LHC bounds on coloured scalars

We analyze the constraints on coloured scalar bosons imposed by the current LHC data at √ s = 13 TeV. Specifically, we consider an additional electroweak doublet of colouroctet scalars, satisfying the principle of Minimal Flavour Violation in order to fulfill the stringent experimental limits on flavour-changing neutral currents. We demonstrate that coloured scalars with masses below 800 GeV are already excluded, provided they are not fermiophobic.


I. INTRODUCTION
The discovery of the Higgs boson by the ATLAS and CMS Collaborations [1,2] can be understood not just like the confirmation of the mechanism for generating the masses of the Standard Model (SM) particles but like the discovery of the first element of a possible more extended scalar sector. Although the current data are compatible with the SM, there remain many unanswered fundamental questions which motivate the existence of new physics (NP). Furthermore, with the available data, there is still plenty of room for extensions of the scalar sector of the SM at the TeV scale.
Electroweak models with extended scalar sectors usually introduce two potentially worrisome problems: dangerous contributions to the mass ratio ρ ≡ M 2 W =ðM 2 Z cos 2 θ W Þ and unsuppressed flavor-changing neutral-current (FCNC) transitions. The phenomenological requirement that ρ − 1 must be very small, ≤ Oð10 −3 Þ, selects SUð2Þ L singlets and doublets as the preferred scalar candidates, while unwanted FCNCs are usually avoided introducing some ad hoc discrete symmetry that forces each type of SM righthanded fermion to couple only to a single scalar doublet. This guarantees that all Yukawa matrices are diagonal in the mass basis and keeps the resulting flavor structure stable under quantum corrections (natural flavor conservation) [3]. However, there are more generic possibilities in order to suppress FCNCs.
The principle of minimal flavor violation (MFV) [4,5] constitutes a much weaker (and general) assumption that also leads to a very effective suppression of FCNCs. It is based on the hypothesis that all Yukawa matrices are proportional to the same flavor structures that break the SUð3Þ Q L ⊗ SUð3Þ u R ⊗ SUð3Þ d R symmetry. In multi-Higgs doublet models, this leads to the flavor alignment of all scalar Yukawa couplings to a given right-handed fermion type [6][7][8]. FCNCs are then absent at tree level, and all flavor-changing phenomena are controlled by the charged-current Cabibbo-Kobayashi-Maskawa (CKM) matrix [9,10].
Manohar and Wise realized that the principle of MFV can only be satisfied by those scalar representations transforming under SUð3Þ C ⊗ SUð2Þ L ⊗ Uð1Þ Y like ð1; 2Þ 1=2 or ð8; 2Þ 1=2 [11]. Therefore, models with additional scalar doublets that are either color singlet or color octet become very interesting candidates for possible extensions of the scalar sector of the SM. The phenomenological implications of having additional color-singlet scalar doublets have been extensively studied [8,. In this work, we focus on the other possibility, scalar extensions in which we add a color-octet electroweak doublet satisfying the principle of MFV, the so-called Manohar and Wise (MW) model [11]. This model is also motivated by the fact that some SUð5Þ and SOð10Þ unification theories predict color-octet electroweak doublets with masses around the electroweak scale [48][49][50][51][52][53].
The MW model has also been widely studied in the literature, and its parameter space has been constrained with theoretical considerations, such as unitarity and vacuum stability [54,55], and phenomenological analyses. The presence of color-octet scalars affects the SM Higgs production and its decay to diphotons [56][57][58][59][60], electroweak precision observables like the oblique parameters [11,61] or R b [62,63], flavor observables such as neutral-meson mixing or the decay B s → l þ l − [64], and the anomalous electric and magnetic dipole moments of the quarks [65].
Being colored particles, the MW scalars could be massively produced at the LHC, provided they are light enough to be kinematically accessible. However, most of the works analyzing the direct production of these scalars have been performed before the huge release of data on heavy-particle searches at the LHC [66][67][68][69][70][71][72]. The most complete analysis on the direct production of these scalars was done in Ref. [70], where the single and pair production of neutral scalars was studied. Slightly stronger limits were obtained in Ref. [72] for the particular case of a neutral colored pseudoscalar, assuming that it decays exclusively into top quarks (top-philic limit). The recent release of high-luminosity experimental searches for massive resonances decaying into heavy quarks [73][74][75][76] makes it now possible to improve those bounds. In this work, we focus on the associated production of neutral and charged scalars with heavy quarks that was not considered in previous analyses. With these searches, we are able to cover some regions of the parameter space that have not been studied before. Furthermore, we also update the experimental data on the single production of neutral scalars decaying to top quarks, at ffiffi ffi s p ¼ 13 TeV. However, the limits obtained from this channel are still not better than those obtained for ffiffi ffi s p ¼ 8 TeV in Ref. [70].
Note that in order to compare our theoretical predictions with the experimental data, we make use of specific analyses performed by the ATLAS and CMS Collaborations. In the case of the single production, we compare our predictions with the experimental limits found for Z 0 particles and Kaluza-Klein gluons and gravitons [73]. For the associated production, we exploit the limits obtained for the type-II two-Higgs-doublet model (2HDM) in Refs. [74][75][76]. We have considered that the kinematics of our production channels and the experimental analyses we compare with are very similar, which is a reasonable assumption.
Although in this work we have focused in the MW model, it is worth to mention that there are other possibilities to suppress the NP contributions to FCNCs, like the next-tominimal-flavor-violation scenario where the NP dominantly couples with the third-generation quarks and is quasialigned with the Yukawa matrices [77]. Phenomenological signatures of colored-octet scalars in this type of models have been studied in Refs. [78,79].
We first provide in Sec. II a very brief description of the MW model and define the relevant parameters. Our phenomenological analyses are detailed in Sec. III, where we analyze the single production of neutral scalars and the associated production of neutral and charged scalars with top quarks. We finally summarize in Sec. IV.

II. THE MW MODEL
As we mentioned in the Introduction, the MW model adds a new scalar field to the SM with the SUð3Þ C ⊗ SUð2Þ L ⊗ Uð1Þ Y quantum numbers ð8; 2Þ 1=2 . Therefore, the scalar sector will be formed by the usual Higgs doublet ϕ¼ðϕ þ ;ϕ 0 Þ T plus an SUð3Þ C -octet field S A ¼ðS A;þ ;S A;0 Þ T . Since they are colorful particles, the new scalars cannot acquire a vacuum expectation value (vev), neither can they mix with the SM Higgs doublet.
The most general potential that can be build with this scalar sector takes the form, where i and j are SUð2Þ L indices, the traces are in color space, and we have used the notation S ¼ S A T A , with T A the generators of the SUð3Þ C group. All potential parameters are real except λ 3 , λ 4 , and λ 5 , but we can choose λ 3 to be real performing a global phase rotation of the S multiplet field. The vev of the SM Higgs doublet, hϕ 0 i ¼ v= ffiffi ffi 2 p , generates a mass splitting among the physical colored scalars, where m S AE is the charged-scalar mass, while m S 0 R and m S 0 I are the masses of the CP-even and CP-odd neutral scalars, respectively.
The interaction of the octet scalars with the gauge bosons is generated by the kinetic term, through the covariant derivative, with G μ ¼ G A μ T A the octet gluon field. The factor of 2 in Eq. (3) generates the correct canonical normalization for the fields.
The last remaining piece is the Yukawa interaction of the color-octet scalar multiplet. The MFV assumption implies that the Yukawa flavor matrices of the S field are proportional to the SM ones, where η D and η U are, in general, complex parameters. The scalar-fermion interactions are then proportional to the corresponding fermion masses since The S 0 R;I G 2 interaction plays an important role in the production and decay of the colored neutral scalars. The dimension-4 Lagrangian does not contain any direct coupling of the neutral octet scalars to two gluon fields. However, this coupling gets generated by quantum effects through scalar and fermion loops. The corresponding vertex can be represented by the dimension-6 gaugeinvariant effective Lagrangian, where G A μν is the gluon strength tensor andG A μν ¼ 1 2 ϵ μναβ G A αβ . The Wilson coefficients F R and F I are easily obtained at the one-loop level. This was first calculated in Ref. [62], neglecting the mass splitting between the scalars. If the splitting is not neglected, these coefficients are given by in the CP-conserving limit, i.e., considering all parameters to be real. In these expressions, we have made use of the functions, I q ðzÞ ¼z½2þð4z−1ÞF ðzÞ; I s ðzÞ ¼−z½1þ2zF ðzÞ; ð8Þ where F ðzÞ ¼ In order to reduce the total number of parameters, we will work in the CP-conserving limit. Furthermore, the couplings m 2 S , λ 1 , λ 2 , and λ 3 are only needed for the determination of the scalar masses; since we only have three different colored scalars, we can then remove 1 degree of freedom (d.o.f.). Moreover, the parameters λ 4 and λ 5 are just relevant for the decay of the CP-even neutral scalar to gluons, which depends on their sum λ 4 þ λ 5 , allowing us to remove another d.o.f. Finally, the four-point interactions of the colored scalars will be irrelevant for this analysis, so we can take λ 6−11 ¼ 0 for simplicity.
With this considerations, we end up with only 6 d.o.f., These parameters must satisfy some theoretical requirements and phenomenological constraints. Perturbative unitarity enforces the modulus of λ 4;5 to be smaller than 13 [54]. Flavor observables [64] and R b [62] put an upper bound on the up-type Yukawa, jη U j < 2, for a chargedscalar mass below 1 TeV. The analogous limits on jη D j are currently very weak, allowing it to go beyond 100. Finally, the mass splittings are constrained by the oblique parameters S, T, and U. For light scalar masses (≤ 350 GeV), Δm 2 S I (Δm 2 S R ) is constrained to be smaller than 60 GeV if λ 2 and λ 3 have the same (opposite) sign [67]. The oblique constraint becomes weaker for larger masses.

III. COLLIDER PHENOMENOLOGY
In our phenomenological analysis, we have considered the processes pp → S 0 R;I → tt, pp → S 0 R;I tt → tttt, and pp → S þt b → tbt b. In order to constrain the model, we will focus on those regions of the parameter space in which the colored scalar particles decay mainly to quarks, i.e., where decay modes into another colored scalar and a weak boson such as S 0 R;I → S þ W − , S þ → S 0 R;I W þ , and S 0 R;I → S 0 I;R Z 0 are suppressed. Of course, these channels are only possible for the heavier scalars decaying to the lighter ones because the opposite possibility is kinematically forbidden. Therefore, when we consider the production of every scalar, we select λ 2 and λ 3 in such a way that this scalar is the lightest one, avoiding then the unwanted decay modes. Hence, the study of the CP-even neutral scalars, the one of the CP-odd neutral scalars, and the one of the charged scalars analyze different regions of λ 2 and λ 3 , but once we combine the results coming from the three searches, we are able to cover all possible values of these parameters.
We have generated collision events with the program MG5_aMC@NLO [80], using first FEYNRULES [81,82] to produce the universal Feynrules output needed to run our model. All calculations in MG5_aMC@NLO were performed at tree level. In the event generation, we have used the particle distribution functions (PDF) set NNPDF23_nlo_ as_0119 of the package LHAPDF_6.1.6 [83]. The center-of-mass energy has been fixed to ffiffi ffi s p ¼ 13 TeV, while the experimental data had an integrated luminosity of 36.1 fb −1 for the ATLAS data [73][74][75] and 137 fb −1 for the CMS data [76]. Furthermore, we have estimated the theoretical uncertainties varying the values of the renormalization and factorization scales by a factor between 1=2 and 2.

A. Single production of neutral scalars
The neutral colored scalars, apart from decaying into top quarks, can also decay to gluons and to bottom quarks. The decay amplitude to gluons depends on the parameters η U and η D for the CP-odd scalars, and also on λ 4;5 for the CP-even scalars; the decay to tt is proportional to η U , and the decay to bb is proportional to η D . Furthermore, as mentioned before, we select the values of λ 2 and λ 3 in such a way that the decay into another colored scalar is kinematically forbidden, i.e., λ 2 þ 2λ 3 < 0 and λ 3 < 0 for the CP-even scalar, and λ 2 − 2λ 3 < 0 and λ 3 > 0 for the CP-odd scalar.
In order to maximize the production cross section of a colored neutral CP-even scalar through gluon fusion, we have chosen λ 4;5 with the opposite sign to η U and η D , so that all contributions to F R in Eq. (7) interfere constructively. Notice that this Wilson coefficient also governs the S 0 R decay into two gluons. Hence, both the single production cross section and the decay amplitude to gluons are higher in this case compared with the alternative sign choice for λ 4;5 . However, as shown in Ref. [70], the scalar contribution is much smaller than the top one, unless jλ 4;5 j is very large. Therefore, one gets finally similar results with both signs. The resulting limits are just slightly better with our choice of sign due to the increase of the cross section; the corresponding increase of the decay width into gluons does not worsen them because the branching ratio of the S 0 R decay into quarks remains always very close to one, for the values of η U considered here. We have varied jλ 4;5 j from 0 to 10, which is almost at the limit of its perturbative unitarity region (jλ 4;5 j ≲ 13) [54]. As it has been previously commented, the parameter η U is strongly constrained by flavor observables and, for masses of the colored scalars smaller than 1 TeV, its absolute value cannot be higher than 2 [62,64], so we go at most up to this value of jη U j.
For the single production, we have taken the experimental data from the ATLAS search for Z 0 bosons and Kaluza-Klein gravitons and gluons, in Ref. [73]. The interference with the SM production amplitudes has not been considered in those searches, although it could have some effects on the signal shape for scalars and pseudoscalars decaying to top quarks [84]. Therefore, some care has to be taken with the limits obtained from these channels, until more direct searches for colored scalars and pseudoscalars are released. Note also that for these searches we have not considered any QCD corrections. In the SM, these corrections enhance the single production of the Higgs boson at the LHC by a factor around 1.5 [85], and we could expect a similar or even larger contribution for our colored scalars. However, there are additional Feynman diagrams contributing to the production of the colored scalars, which have not been calculated yet. Since the missing QCD corrections are expected to increase the production cross section, our limits are then quite conservative. Figure 1 compares the (95% C.L.) experimental limits on the production of heavy particles that decay into top-quark pairs [84] with the calculated production cross section times branching ratio for the CP-even colored scalar, as a function of the scalar mass. The model-dependence of the ATLAS exclusion limits can be appreciated from the broad range of bounds obtained for the different explicit models analyzed: Z 0 bosons, Kaluza-Klein gluons, and Kaluza-Klein gravitons decaying into tt. The production of CP-even scalars is dominated by the gluon-fusion mechanism, which depends on η U , η D , and λ 4;5 . Therefore, this experimental constraint can be easily avoided, taking smallenough values for these parameters. When jη U j is of order one, the branching ratio for the decay S 0 R → tt is almost one, provided m S 0 R > 2m t , for all values of λ 4;5 within its perturbative unitarity region. Owing to the large mass splitting between the top and bottom quarks, the value of η D is also almost irrelevant for this branching ratio, so we have just taken η D ¼ 1 as a representative value. FIG. 1. Cross section times branching ratio for the single production of S 0 R and its decay to tt. In (a) we set η U to 1 and vary λ 4;5 from −10 to 0, and in (b) λ 4;5 ¼ −10 and η U is varied from 1 to 2. The nonsolid lines correspond to the experimental limits from Ref. [73].
It is evident from the figure that this channel does not provide strong constraints on the scalar mass. For η U ¼ 1 and λ 4;5 ¼ −10, we can just infer that m S 0 R should be heavier than 500 GeV, comparing the predicted cross section with the ATLAS bound on the production of Kaluza-Klein gravitons. For smaller values of η U and/or jλ 4;5 j, one does not obtain any useful constraints. Obviously, the missing QCD corrections could not change much the situation, which, furthermore, gets slightly worse with the opposite sign choice for λ 4;5 . In spite of the much higher statistics accumulated at 13 TeV, compared with the data sample analyzed in Ref. [70], the emerging limits are still not better than those extracted from the 8 TeV data and are not competitive with the ones obtained from the associated-production process, which are shown afterwards.
The analogous limits on the single production of the CP-odd scalar S 0 I are shown in Fig. 2. The behavior is very similar to the CP-even case with λ 4;5 ¼ 0, because the production of CP-odd scalars does not depend on λ 4;5 . The production cross section depends only on η U and η D , and it grows with the modulus of these parameters, although the dependence on η D is again extremely weak. The numerical differences between the left (η D ¼ 1) and right (η D ¼ 10) panels can hardly be seen in the figure. For values of jη U j ≤ 2, the predicted signal remains below the experimental limits in the whole range of m S 0 I analyzed. Therefore, this channel does not provide any constraint.

B. Associated production of neutral scalars and top quarks
The production of the colored neutral scalars in association with top-quarks proceeds through the two mechanisms indicated in Fig. 3. Both diagrams contain a single S 0 R;I tt vertex with a coupling η U . However, in the kinematical region of interest, the production through the left mechanism is dominated by an on shell intermediate S 0 R;I particle decaying into tt. Therefore, as long as the branching ratio of this decay is close to one, the dependence with η U will be small (gets canceled by the total decay width contained in the on shell propagator). The second mechanism on the right will then become more relevant for values of jη U j in its higher allowed range. However, for values of jη U j of order 10 −1 or smaller, the left diagram dominates, provided the branching ratio of the decay to top quarks is of order one. Note that when the left diagram dominates, the resulting amplitude does not change under a global rescaling of the parameters η U , η D , and λ 4;5 . 1 The limits on the production of neutral scalars in association with top-quark pairs are shown in Figs. 4 and 5, for the CP-even and CP-odd cases, respectively. As mentioned before, we compare our theoretical predictions for σ · Br of pp → S 0 R;I tt → tttt with the experimental upper bounds obtained for the type-II 2HDM by ATLAS [74] and CMS [76]. 2 Note that the ATLAS analysis covers a higher mass region, reaching masses up to 1 TeV, while the CMS one, although being more restrictive, only applies to m S 0 R;I ≤ 650 GeV. Regarding the QCD corrections, they increase the SM prediction for this process (pp → Htt → tttt) by a factor of 1.2 [86]. We expect them to be more relevant for the colored scalars, because they can directly couple to gluons. However, these corrections have not been calculated yet. Since they are expected to increase the prediction for σ · Br, our limits without QCD corrections are then conservative.
The four panels in Fig. 4 show the predicted production of the CP-even scalar, as a function of m S 0 R , for several representative choices of the relevant parameters. In all cases, we have considered that the CP-even neutral scalar is the lightest colored scalar. Taking η U ¼ 0.1 and η D ≤ 1 (top-left panel), the predictions are well above the experimental bounds in the whole range of m S 0 R , for any value of jλ 4;5 j, resulting in a lower limit of 1 TeV for the mass of the neutral scalar. This limit remains valid for larger values of jη U j because they result in much larger theoretical predictions for σ · Br. The expected signal decreases with increasing values of jη D j, mainly because it increases the branching ratio of the decay to bottom quarks, relaxing the experimental constraint. The decay amplitude to gluons also depends on η D , and on its relative sign with η U and λ 4;5 , but the effect is extremely weak. Indeed, the decay to gluons is also affected by λ 4;5 but, as we discussed before, the effect is small for values of this parameter inside its perturbative unitarity region. Looking at the top panels, we can see that the difference between setting λ 4;5 to 0 or to −10 is very small. For positive values of λ 4;5 , the effect would be even smaller because the decay amplitude to gluons would grow less while varying λ 4;5 . With η U ¼ 0. 1 and FIG. 4. Cross section of the associated production of S 0 R with tt times its branching ratio into tt, as a function of m S 0 R , for representative choices of the parameters. In the top panels, η U ¼ 0.1 and λ 4;5 is varied from 0 to −10, while η D ¼ 1 (left) or η D ¼ 3 (right). The bottom panels correspond to λ 4;5 ¼ −10, with η U varying from 0.1 to 1, and η D ¼ 5 (left) or η D ¼ 10 (right). The experimental bounds are taken from Refs. [74,76]. η D ¼ 3 (top-right panel), one still obtains m S 0 R ≥ 850 GeV, for values of λ 4;5 inside its perturbative unitarity region. As we keep increasing η D (lower panels) the constraints become worse, but even with η D ¼ 10 the mass of the CP-even scalar is constrained to be higher than 1 TeV for η U ¼ 1.
For the analysis of the CP-odd neutral scalar, we have chosen λ 2 and λ 3 in such a way that this particle is the lightest colored scalar, as mentioned before. As illustrated by Fig. 5, the results are similar to the previous ones; although in this case, there is no dependence on λ 4;5 , which makes the analysis simpler. The decay width into gluons can be neglected compared to the decays into top and bottom quarks because it is a loop process with only top and bottom quarks in the loop. The dependence of the resulting limits with η U and η D follows the same trend as before: the higher the value of jη U =η D j, the better the constraint. We find again a lower limit of 1 TeV for values of jη U j as small as 10 −1 when η D ¼ 1 (top-left panel). For η U ¼ 0.1 and η D ¼ 3 (top-right panel), the constraint is m S 0 I ≥ 850 GeV, as in the CP-even case. However, thanks to the CMS data, we can find a better limit on the CP-odd scalar, m S 0 I ≥ 575 GeV, for η D ¼ 10 and η U ¼ 0.1 (bottom-right panel).

C. Production of charged scalars
The analysis of charged-scalar production is much simpler because they cannot decay to gluons; therefore, there is no dependence on λ 4;5 . Now, we restrict ourselves to the region of parameter space with λ 2 þ 2λ 3 > 0 and λ 2 − 2λ 3 > 0, so that the decays S þ → S 0 R;I W þ are kinematically forbidden. Thus, the only possible two-body decays of the charged scalars are into two quarks, and therefore, the decay to heavy quarks will always dominate. Hence, the cross section only depends on η U and η D .
The production mechanisms, indicated in Fig. 6, are analogous to the associated production of a neutral scalar with the obvious changes on the particle charges. Again, the mechanism on the right dominates when jη U j ∼ Oð1Þ, while the left diagram is the relevant one for jη U j of order 10 −1 or smaller. But now the branching ratio of the decay into top and bottom quarks is always close to one, independently of the values of η U and η D , provided the charged scalar is not fermiophobic (η U;D ¼ 0). The only requirement needed for our analysis is that the charged scalar indeed decays before entering the detector, which is the case even when η U and η D are of order 10 −5 . Figure 7 compares our theoretical predictions for the cross section times branching ratio of the process pp → S þt b → tbt b, with the (95% C.L.) experimental bounds, as a function of m S AE . We profit here from the ATLAS analysis performed in the context of the type-II 2HDM model [75]. In the left panel, η D ¼ 1 and η U is varied from 0 to 3, while in the right panel the most pessimistic possibility, η U ¼ 0, is taken and η D is varied between 1 and 20. From this plot, we are able to constrain the mass of the charged scalars to be higher than 800 GeV, even for η U ¼ 0 and for any value of jη D j as small as 10 −5 .
For higher values of jη U j, the constraints become stronger, and for η U ¼ 2, we push this limit up to 900 GeV. Note that this is the maximum value of η U allowed by flavor constraints, for charged-scalar masses up to 1 TeV [62,64], with a 95% C.L.
The production cross section depends very weakly on jη D j. This is easy to understand, since η D couples proportionally to the bottom quark mass while η U brings an m t factor. Even when η U ¼ 0 (right panel), the sensitivity to jη D j is quite mild because the left diagram dominates in this case.

IV. CONCLUSION
In this work, we have found lower limits for the masses of all the colored scalars of the MW model: the CP-even and CP-odd neutrals and the charged colored scalars. Our phenomenological study has been performed in the CP-conserving limit. In order to avoid unwanted decay modes into other colored scalars, for each separate analysis, we have selected the parameters of the scalar potential in such a way that the analyzed scalar was the lightest. Thus, once we find a limit for the mass of this scalar, the others must be necessarily heavier, so the limit applies to all of them. Combining the results obtained for all possible mass splittings, the least restrictive limit is valid for any value of the parameters λ 2 and λ 3 .
We have analyzed the single production of neutral scalars and the associated production of neutral and charged scalars with heavy quarks, using the available LHC data at ffiffi ffi s p ¼ 13 TeV. From the study of the single production, we could not find better constraints than those already obtained in Ref. [70], at lower LHC energies. However, our analysis of the associated production results in relevant limits on the scalar masses, which are significantly better than those extracted in previous works. The associated production of the CP-even neutral scalar depends on λ 4;5 , η U , and η D . Nevertheless, we have found a very small sensitivity to λ 4;5 for values of jη U j higher than 10 −1 . If this is the case, the limits on the scalar mass only depend on the relation between η U and η D . For jη D =η U j ≤ 10, a CP-even scalar is excluded in the whole range explored by ATLAS, which extends up to 1 TeV. The limit gets relaxed for larger values of this ratio. Taking jη D =η U j ∼ 30, one still finds m S 0 R > 850 GeV. However, no limit is found when jη D =η U j > 100. Similar results are obtained for jη U j smaller than 10 −1 , as long as jλ 4;5 =η U j ≤ 100. The associated production of the CP-odd scalar behaves in a similar way, but it does not depend on λ 4;5 , which results in somewhat stronger constraints for some parameter configurations. The limits on m S 0 I are only sensitive to the ratio η D =η U and are valid for values of η U as small as 10 −7 . For jη D =η U j ≤ 10, one finds again that a CP-odd scalar is excluded in the full kinematical range analyzed, i.e., m S 0 I > 1 TeV. The limit gets relaxed to 850 GeV, when jη D =η U j ∼ 30, as in the CP-even case. However, even for jη D =η U j ∼ 100, one still finds a relevant limit of m S 0 I > 575 GeV. Since the charged scalars cannot decay into gluons, the analysis of their associated production is much simpler. We found a quite strong lower bound, m S AE > 800 GeV, provided the charged scalars are not fermiophobic. This constraint remains valid as long as jη D j>10 −5 or jη U j> 10 −7 . The limit becomes stronger with increasing values of jη U j, reaching 900 GeV for jη U j ¼ 2, the maximum value of this parameter allowed by flavor constraints.
Combining all searches, we find an absolute lower bound of 800 GeV on the masses of all colored scalars, which is valid under very mild requirements on the relevant model parameters: jη D =η U j<30, jλ 4;5 =η U j<100, and jη U j > 10 −7 . Note that this limit applies for all possible mass splittings among the scalars. Our analysis puts then a very severe limitation on the hypothetical existence of light colored scalars with masses below the TeV. In practice, they seem to be strongly excluded for all reasonable choices of parameters. The only possibility that remains still viable are fermiophobic scalars with η U;D ¼ 0. The analysis of this extreme case would require a quite different and much more tricky phenomenological approach.
Much larger data samples are going to be accumulated in the forthcoming run 3 and, specially, at the High Luminosity LHC (HL-LHC). Prospective analyses for the pair production of pseudoscalar color-octet particles decaying into the ðttÞðttÞ final state [72] have been presented in Ref. [87]. From these results, we can infer that our 850 GeV limit on the CP-odd scalar mass, when jη D =η U j ∼ 30, could be pushed to 1 TeV. In order to find a rough estimate for the High Energy LHC (HE-LHC), we have generated some events at a center-of-mass energy of 27 TeV in the same channel, finding an exclusion lower limit of 1.4 TeV in the top-philic case. Larger mass scales could of course be reachable at the FCC-hh, since the pair-production cross section would be orders of magnitude larger, for a given scalar mass value [88]. In the singleproduction channel, the achievable sensitivity at the FCC-hh can be roughly estimated from existing simulations of strongly interacting vector particles, which quote limits on σðpp → Z 0 Þ × BrðZ 0 → ttÞ in the 10 −4 to 10 −3 pb range, for Z 0 masses between 10 and 35 TeV [88]. Extrapolating these estimates to lower masses, exclusion limits around 3 to 5 TeV for the octet scalar mass appear to be feasible. More detailed analyses would be needed to reliably assess the discovery potential of future colliders.

ACKNOWLEDGMENTS
We would like to thank A. Peñuelas for her useful comments on the manuscript. We also thank M. Perello and M. Vos for helpful discussions on experimental aspects of the analysis. This work has been supported in part by the Spanish Government and European Regional Development