Constraining the Natural MSSM through tunneling to color-breaking vacua at zero and non-zero temperature

We re-evaluate the constraints on the parameter space of the minimal supersymmetric standard model from tunneling to charge- and/or color-breaking minima, taking into account thermal corrections. We pay particular attention to the region known as the Natural MSSM, where the masses of the scalar partners of the top quarks are within an order of magnitude or so of the electroweak scale. These constraints arise from the interaction between these scalar tops and the Higgs fields, which allows the possibility of parameter points having deep charge- and color-breaking true vacua. In addition to requiring that our electro-weak-symmetry-breaking, yet QCD- and electromagnetism-preserving vacuum has a sufficiently long lifetime at zero temperature, also demanding stability against thermal tunneling further restricts the allowed parameter space.


Introduction
The mechanism of spontaneous symmetry breaking by the vacuum expectation value for a scalar field is an essential component of the standard model of particle physics (SM) [1][2][3], which has proven itself to be an accurate description of Nature all the way to the tera-electron-Volt scale. The discovery of the bosonic resonance at 125 GeV at the Large Hadron Collider (LHC) [4,5] is consistent with the Higgs boson of the spontaneous symmetry breaking of the SM, leading one to take the issue of minimizing the scalar potential seriously.
The minimal supersymmetric extension of the SM (the MSSM) has a much more complex scalar potential by merit of there being many more scalar fields (partners for each SM fermion as well as a second Higgs SU (2) L doublet) which interact with the Higgs fields. The large effect of the extra loops on the mass of the Higgs boson along with the non-observation of supersymmetric partners thus far has led to the pragmatic region of the MSSM parameter space known as the Natural MSSM [6][7][8]. This is the region where the masses of all the partners are very large but for those with the largest contributions to the Higgs mass [9][10][11][12][13][14], which should have masses not very far above the electroweak scale so that there is little finely tuned cancellation between loop contributions to the minimization conditions, and thus is in some sense natural [15][16][17][18]. Thus the stopst (scalar partners of the top quarks) should have TeV-scale soft supersymmetry-breaking parameters while all others are assumed to have very large masses. The partners of the bottom quarks and tau leptons could also be in the TeV-scale, but in this letter we consider only stops, noting that our algorithm is trivially generalizable and is already implemented in the public code Vevacious [19].
While the interaction between stops and the Higgs fields allow the mass of the Higgs boson to reach 125 GeV in the MSSM, it also leads to the possibility of the scalar potential having undesired minima apart from the desired symmetry-breaking (DSB) vacuum, where only the neutral components of the Higgs doublets get non-zero VEVs. Even though a parameter point may be chosen where the scalar potential has a minimum where the stops do not have non-zero VEVs, there is no guarantee that this is the global minimum: there may be deeper charge-and color-breaking (CCB) minima to which the Universe may tunnel [20][21][22][23][24][25][26][27][28][29][30]. However, even if the DSB vacuum is only metastable, the parameter point is still acceptable if the expected tunneling time is of the order of the age of the known Universe [31][32][33]. Also, given the convincing success of the Big Bang theory, acceptable parameter points with metastable DSB vacua should also have a high probability of surviving tunneling to the true CCB vacua through thermal fluctuations.
In section 2 we lay out the algorithm by which we compute whether a parameter point is excluded by the DSB vacuum having a very low probability of surviving to the present day either by a high probability of critical bubbles of true vacuum forming through quantum fluctuations in our past light-cone at zero temperature, or by such bubbles forming through thermal fluctuations during the period when the Universe was at sufficiently high temperature. In section 3, we show how much of the parameter space is excluded by such conditions, and compare this to previous work. Finally we conclude in section 4.

Parameter point selection and stability evaluation
We categorize the stability or metastability of a parameter point by a multi-stage process. First, a consistent set of Lagrangian parameters at a fixed renormalization scale is generated by SPheno [34,35], such that the MSSM physics at the DSB vacuum is consistent with the SM inputs (m Z , G F , etc.), and these parameters are stored in a file in the SLHA format which is passed to Vevacious, using a model file automatically generated by SARAH [36][37][38][39][40]; for consistency of input, the version of SPheno was also generated by SARAH. Vevacious is a publicly-available code [19] that then prepares the minimization conditions for the tree-level potential as input for the publicly-available binary HOM4PS2 [41] that finds all possible solutions to the particular minimization conditions of the parameter point. These are then used by Vevacious as starting points for gradient-based minimization by MINUIT [42] through PyMinuit [43] to minimize the full one-loop potential with thermal corrections at a given temperature. If a minimum deeper than the DSB vacuum is found, the probability of tunneling out of the false DSB vacuum is then calculated through CosmoTransitions [44]. For a full discussion of the calculation of the bounce action and its conversion to a tunneling time from a false vacuum to a true vacuum, we refer the reader to the Vevacious manual [19], the CosmoTransitions manual [44], and the seminal papers on tunneling out of false vacua [45,46].
If a parameter point is found to have a deeper CCB minimum, we label it as metastable, otherwise we label it stable 2 . We then divide the metastable points into shortlived points which would tunnel out of the false DSB vacuum in three giga-years or less (corresponding to a survival probability of lasting 13.8 Gy of one per-cent or less), and the rest as long-lived. Finally, we divide the long-lived points into thermally excluded, by having a probability of the DSB vacuum surviving thermal fluctuations of one per-cent or less, or allowed, by having a survival probability of greater than one per-cent, as described in more detail in the following subsection.

Thermal corrections
Since the temperature of the Universe has been negligible for most of its existence, it is quite reasonable to calculate the tunneling time assuming that the four-dimensional bounce action S 4 is the dominant contribution to the decay width of the false vacuum. However, for sufficiently high temperatures, the dominant contribution may come from solitons that are O(3) cylindrical in Euclidean space rather than O(4) spherical [47].
If the thermal contribution dominates, the expression for the decay width per unit volume Γ/V at a temperature T changes accordingly: where A is a quantity of energy dimension four, which is related to the ratio of eigenfunctions of the determinants of the action's second functional derivative, and S 3 (T ) is the bounce action integrated over three dimensions rather than four, with the integration over time simply replaced by division by temperature because of the constant value along the Euclidean time direction. The leading thermal corrections to the potential are at one loop, and given by where the sum is over degrees of freedom: bosons as sets of real scalars, fermions as sets of Weyl fermions, and with J + for a real bosonic degree of freedom and J − for a Weyl fermion (note that we incorporate the negative sign into the definition of J − in contrast to Ref. [48]). The probability P (T i , T f ) of not tunneling between the time when the Universe is at temperature T i and when it is at temperature

Evaluating the survival probability
Even the numerical evaluation of the action is computationally intense and while one could attempt to numerically integrate eq. (4), this is impractical for more than a handful of parameter points. Hence we exclude parameter points based on an upper bound on the survival probability under some approximations, which requires S 3 (T ) to be evaluated only once.
Firstly, the factor A(T ) is taken to be T 4 , as the evaluation of the eigenfunctions of the determinant is so hard that they are usually estimated on dimensional grounds anyway, which is justified as the exponent of the action is much more important [49]. Any deviation would effectively contribute ln(AT −4 ) to S 3 (T )/T , and S 3 (T )/T is ∼ 240 for survival probabilities that are not extremely close to zero or one.
Secondly, we assume that the Universe is radiation dominated during its evolution from T i to T f and that entropy is approximately conserved between T i and today, as it is appropriate for the MSSM. Entropy conservation implies that V (T 0 )/V (T ) = s(T )/s(T 0 ), where s is the entropy density and T 0 = 2.73 K is the temperature of the Universe today. Using the relation for dt/dT during radiation domination, we can replace in eq. (4) The tunneling is assumed to be dominated at a temperature above that at which the DSB vacuum evaporates, so all degrees of freedom of the SM are taken to be relativistic, while non-SM particles are assumed to be still non-relativistic at this temperature. This is because if the dimensionful terms such as soft SUSY-breaking terms are of the order of some scale Q, the CCB minimum should be deeper than the DSB vacuum by about Q 4 and effective thermal contributions to the masses of about T are likely to make tunneling impossible by T Q. Hence T should be less than the typical masses of the non-SM particles. Thus g * s (T ) ≡ g * (T ) = 106.75, entirely due to the SM particles.
Putting it all together, we take Thirdly, as the evaluation of S 3 (T ) is very costly in CPU time, we assume that S 3 (T ) is a monotonically increasing function of T . As the magnitudes of the field values increase along the path from the false vacuum to the CCB vacuum, the masses of the degrees of freedom increase (barring occasional cancellations). Hence the thermal contributions lower the effective potential less near the CCB vacuum than near the false vacuum, hence increasing T leads to the absolute height of the energy barrier decreasing but the barrier height relative to the false vacuum, which is the important quantity, increases, and thus Given this, and all that remains is to find the optimal T = T opt to maximize this quantity to find an upper bound on the survival probability P (T i = T opt , T f = 0) for the DSB vacuum. Hence if we can choose T opt before attempting to calculate S 3 (T opt ), we only need make one evaluation of S 3 (T opt ).
The evaluation of the three-dimensional bounce action along a straight path in "field space" from the false vacuum to the true vacuum 3 , denoted S straight 3 (T ), is much quicker to calculate than searching for the optimal path, so for each parameter point S straight 3 (T ) was calculated for a set of temperatures between the temperature at which the DSB vacuum evaporates and the critical temperature T crit at which tunneling to the CCB minimum becomes impossible, then was fitted as (T crit −T ) −2 times a polynomial in T , since the action should diverge as (T crit − T ) −2 as T approaches T crit [49,50]. This fitted function was then numerically minimized to estimate the value of T = T opt which minimizes P (T i = T opt , T f = 0), which was then used to evaluate the right-hand side of eq. (9), taken as the upper bound on the survival probability of the false vacuum.
The estimated optimal T opt was then used to evaluate S 3 (T opt ) properly, along the correct tunneling path (not the straight path) between the CCB vacuum at temperature T opt (found by gradient-based minimization of the full one-loop thermal potential starting from the minimum at T = 0) and the "DSB vacuum at T opt ", which is where gradient-based minimization starting from the position of the DSB vacuum at T = 0 ends up: above the evaporation temperature, this should be the field origin, and indeed was for each parameter point, also demonstrating that the field origin is a true minimum of the potential at T = T opt . 3 The full set of equations of motion of the critical bubble are not solved by this path [44], but would be solved by adding a term to the effective potential raising the energy barrier away from this path in the appropriate way. A critical bubble of the unmodified potential must then have an action less than the action for a critical bubble for the modified potential.
The above procedure has been incorporated into version 1.1 of Vevacious and has been made public for download from HepForge.

Range of validity
As discussed in Ref. [19], one should not trust a fixedorder loop expansion for VEVs very much larger than the renormalization scale Q. Likewise, thermal tunneling dominated at temperatures T Q might not be very accurate. One would hope that the incorporation of running parameters and leading logarithmic corrections to the thermal contributions [48] would stabilize the results acceptably. While we are working on extending Vevacious to include these enhancements, the results presented here are based purely on the one-loop effective potential with running parameters evaluated at a fixed Q. However, for every single one of our parameter points, the VEVs of the CCB minima were within a factor of a few of Q and the thermal tunneling was also dominated by T Q. Hence the logarithms associated with higher orders are not large, and the one-loop expansion of the thermal potential is valid throughout the entire field space considered.
We note that exclusion based on thermal tunneling is dependent on the thermal history of the Universe: if combined with a model where the Universe is never hot enough to allow tunneling at the optimal T , the parameter point is still valid. Indeed, given appropriate initial conditions, consistency with big bang nucleosynthesis requires reheat temperatures only above a few MeV (see e.g. [51]). The exclusions presented here are nonetheless important since in the most commonly hypothesized cosmologies, T ∼ 10 5 GeV is already considered very low [52][53][54][55][56][57].
Finally, we do not address the question of whether there are additional CCB minima at extremely large VEVs 10 16 GeV which can only be reliably calculated with current methods using running parameters and even then only under restricted circumstances [58], nor do we consider the effects of inflation and re-heating [52].

Parameter scan
While spontaneous symmetry breaking in the SM is triggered by a negative mass-squared term in the Lagrangian for the Higgs field 4 , it is neither a necessary nor sufficient condition for any scalar field in a multi-scalar theory to develop a non-zero VEV [60]. In particular, a positive mass-squared for the stop fields does not preclude a parameter point from having a CCB minimum, especially if the trilinear couplings T U 33 = Y t A t and Y t µ for H utLtR and H dtLtR respectively are large compared to the square roots of the soft SUSY-breaking mass-squareds m 2 Q33 and m 2 U 33 .
Parameter Range tan β 5 -60 m 2 Q33 500 2 GeV 2 -1500 2 GeV 2 m 2 U 33 500 2 GeV 2 -1500 2 GeV 2 µ 100 GeV -500 GeV T U 33 -3000 GeV -3000 GeV Table 1: Parameter ranges used in the scan. The soft SUSY-breaking mass-squared parameter for the SU (2) L doublet squarks is given by m 2 Q and that of the SU (2) L singlet up-type squarks by m 2 U . All mass-squared matrices for the scalar partners of SM fermions were diagonal, and all diagonal entries but those shown above were set to 1500 2 GeV 2 . The soft SUSY-breaking mass terms for the U (1) Y , SU (2) L , and SU (3)c gauginos were 100 GeV, 300 GeV, and 1000 GeV, respectively. The soft SUSY-breaking coefficient for the trilinear Hut LtR interaction T U 33 is often written as At × Yt; all other soft SUSY-breaking trilinear terms were set to zero. Finally, the mass of the pseudoscalar Higgs boson was set to 1000 GeV. The renormalization scale for each parameter point was the mean of the physicalt masses at the DSB vacuum.
Given then that we are investigating the Natural MSSM and restricting ourselves to the possibility of tunneling to minima witht VEVs, we choose the region in parameter space described by table 1. The large value of the pseudoscalar Higgs mass places the scan firmly in the decoupling regime of the MSSM Higgs sector [61]. To ensure that scalar partners other than the stops are not relevant to the analysis, we set them to have large masses-squared and zero soft SUSY-breaking trilinear interactions. Since the gluino also can have a non-negligible contribution to the mass of the lightest scalar Higgs, we chose to keep it at 1000 GeV and took masses for the other gauginos roughly according to a typical hierarchy that is expected from unification of the gauge forces [62]. Our parameter scan thus largely overlaps with that of Ref. [63].

Comparison in methodology to previous works
Much early work in the area of tunneling to CCB minima in the MSSM focused on analytic expressions derived from the tree-level potential to determine whether there would be a CCB global minimum [21,22,24,25,27,64], though it has been known for some time that such expressions are neither necessary nor sufficient [26,28], and only general outlines of algorithms could be given [29]. It has also been known for some time that they gave no hint as to whether the tunneling time out of the DSB false vacuum could be phenomenologically acceptable [23,33].
The algorithm used by Vevacious improves upon these by finding all the minima of the tree-level potential, not just those that may lie on special lines in field space, as well as incorporating loop corrections, which, despite various claims in the literature [29,65], are important [66][67][68].
One may note the overlap in objective with the works of Refs. [69] and [63]: CCB minima with stop VEVs are searched for in a similar parameter space, and metastable points are categorized as acceptably long-lived or not based on tunneling times calculated by CosmoTransitions. The major improvement over these works is that we also ex- clude points based on low probabilities to survive thermal fluctuations when the Universe was at a temperature of the order of 1 TeV. However, we also note that we improve upon the zero-temperature results of these works in two significant ways: the first is that our use of the homotopy continuation method guarantees that we find all the minima of the tree-level potential, as opposed to a random seeding of the field space followed by gradient minimization in Ref. [69], which obviously cannot guarantee that the random seeding did not miss a CCB minimum, or a brute-force four-dimensional grid scan in Ref. [63], which may miss minima just beyond the range of the grid. The second way is that we use the full one-loop effective potential rather than the tree-level potential. Though one would hope that the loop corrections do not significantly alter the tree-level conclusions, it is not always the case, and the tree-level results can be rather sensitive to the renormalization scale chosen for the running parameters, while the loop corrections stabilize the dependence on the scale [68].

Constraining the parameter space of the Natural MSSM
Our primary result is that a large proportion of the parameter space where the Higgs boson mass is even slightly compatible with the measurement of 125 GeV [4,5] is ruled out by thermal tunneling even though the tunneling time at zero temperature is much longer than the observed age of the Universe. This is presented in figure 1, where the parameter points of our five-dimensional scan are projected onto a two-dimensional plane with the axes being the mass of the lightest Higgs scalar m h and the ratio X t /M S , where X t = A t − µ cot β, and M S is the square root of the product of the tree-levelt masses evaluated at the DSB minimum as this should keep higher order corrections small [70]. One would expect this ratio to be correlated with the probability of tunneling out of the DSB vacuum, as a combination of the comparisons mentioned in section 2.2. (Though tunneling was evaluated at one loop, the value of m h for each parameter point is based on a full diagrammatic one-loop calculation including the effects of the external momenta [71] and, in addition, the known two-loop corrections are included [72][73][74][75].) However, while increasing |X t |/M S is correlated with decreasing stability in some sense, in the phenomenologically interesting region where m h > 123 GeV, it fails to discriminate effectively between acceptable points with stable or high survival probability DSB vacua and those with low survival probabilities for their DSB vacua. In fact, no projection of our scan onto a plane in terms of simple combinations of the input parameters showed any clear discriminatory power, and thus we conclude that a full calculation is inevitably necessary.

Comparison to previous results
Even though it was derived under the assumption that the Yukawa coupling is much smaller than the gauge couplings, which is obviously wrong for the top sector, and even though it has been known to be neither necessary nor sufficient [26,28], the condition has been used in place of a proper analysis as a check that parameter points have stable DSB vacua. It has been demonstrated numerically that it is neither necessary nor sufficient, nor meaningfully correlated with long-/shortlived metastable vacua [68], but for completeness we show how our results are if we exclude points which fail the condition in figure 2. Coincidentally, the condition happens Figure 2: Results without displaying those points which would be excluded by condition (10). The color coding is the same as in figure 1, with green points contrasted against the others in the top picture, purple in the middle, and blue in the bottom picture. We note that this condition excludes all points with short-lived DSB vacua at zero temperature.
to exclude all the points with DSB vacua that are shortlived at zero temperature, but it both unnecessarily excludes stable and acceptably long-lived metastable points at larger |X t | and fails to exclude most of the points which are excluded by thermal tunneling with m h > 123 GeV. An attempt to account for acceptably long-lived DSB vacua by empirically fitting coefficients [33] led to the following condition: but not even one of the points in our scan was excluded by this condition, hence we consider it irrelevant. Hence we stress again that one should not rely on analytic conditions which are derived using simplifying assumptions or which are based on ostensible patterns found in a particular numerical analysis. For a serious check of the stability of the scalar potential, a full-fledged numerical evaluation for each point is usually inevitable. However, while the typical running time per metastable point (those plotted in red, purple, or blue) per CPU core with Vevacious 1.1 was 10-30 minutes, one can easily change a setting so that it will evaluate whether a parameter point is stable or metastable (green or not) within seconds, which for example should be sufficient for purposes of finding conservatively acceptable parameter regions.
If we ignore thermal tunneling, our results qualitatively agree with Refs. [69] and [63]. Though the parameter space overlap with Ref. [69] is not as great, we largely agree with the ratios of X t to M S where the CCB minima become deeper than the DSB minima and where the tunneling time becomes unacceptably short. Likewise, we agree with the ratios one can read off the figures in Ref. [63], but note that the values of m h therein are inconsistent with our calculation (using SPheno), the calculation in Ref. [69] (using SuSpect), and the results in Ref. [76] (using SuSpect and SoftSUSY) 5 . The Higgs masses calculated by SoftSUSY, SPheno, and SuSpect are within the theoretical uncertainity of 2-4 GeV using two-loop corrections. In contrast, the differences between these codes and FeynHiggs [78,79] are usually larger because of the different renormalization scheme: the difference for figure 1 is a steady increase in m h with increasing |X t /M S |, with good agreement for |X t /M S | = 0, to an increase of about 3 GeV for X t /M S = −2.5 and an increase of about 5 GeV for X t /M S = +2.5.

Conclusion
We have presented an exploration of what regions of the Natural MSSM parameter space can be excluded by demanding at least a one per-cent survival probability for the vacuum with the desired symmetry breaking against tunneling to charge-and color-breaking vacua at non-zero temperatures. In order to do so, we extended the feature set of Vevacious to include the functionality to exclude parameter points based on thermal tunneling, which we have made publicly available: Vevacious 1.1 is available for download from HepForge: http://www.hepforge.org/downloads/vevacious. Stability against thermal tunneling is a relevant constraint, especially in the parameter space of the MSSM where the mass of the lightest Higgs boson is consistent with observations. While exclusion based on zero-temperature tunneling can also exclude regions of the parameter space, points that have sufficiently long lifetimes at zero temperature may have very low probability to avoid ending up in a CCB vacuum by the time the temperature drops to a negligible value. Unfortunately, the dependence on the Lagrangian parameters is not simple, and a full analysis of any given parameter point seems necessary, though straightforward given the availability of Vevacious.
We have also showed that results on metastability based on previous tree-level analyses are not significantly affected by the zero-temperature one-loop corrections, as opposed to the effects at finite-temperature.