Cosmological Aspects of Higgs Vacuum Metastability

The current central experimental values of the parameters of the Standard Model give rise to a striking conclusion: metastability of the electroweak vacuum is favoured over absolute stability. A metastable vacuum for the Higgs boson implies that it is possible, and in fact inevitable, that a vacuum decay takes place with catastrophic consequences for the Universe. The metastability of the Higgs vacuum is especially significant for cosmology, because there are many mechanisms that could have triggered the decay of the electroweak vacuum in the early Universe. We present a comprehensive review of the implications from Higgs vacuum metastability for cosmology along with a pedagogical discussion of the related theoretical topics, including renormalization group improvement, quantum field theory in curved spacetime and vacuum decay in field theory.


Introduction
One of the most striking results of the discovery of Higgs boson [1,2] has been that its mass lies in a regime that predicts the current vacuum state to be a false vacuum, that is, there is a lower energy vacuum state available to which the electroweak vacuum can decay into [3,4]. That this was a possibility in the Standard Model (SM) has been known for a long time [5][6][7][8][9][10]. The precise behavior of the Higgs potential is sensitive to the experimental inputs, in particular the physical masses for the Higgs and the top quark and also physics beyond the SM. The current best estimates of the Higgs and top quark masses [11], M h = 125.18 ± 0.16 GeV, M t = 173.1 ± 0.9 GeV, (1.1) place the Standard Model squarely in the metastable region. As in any quantum system, there are three main ways in which the vacuum decay can happen. They are illustrated in Fig. 1. If the system is initially in the false vacuum state, the transition would take place through quantum tunneling. On the other hand, if there is sufficient energy available, for example in a thermal equilibrium state, it may be possible for the system to move classically over the barrier. The third way consists of quantum tunneling from an excited initial state. This is often the dominant process if the temperature is too low for the fully classical process. All three mechanisms can be relevant for the decay of the electroweak vacuum state, and their rates depending on the conditions. In each of them, the transition happens initially locally in a small volume, nucleating a small bubble of the true vacuum. The bubble then starts to expand, reaching the speed of light very quickly, any destroying everything in its way.
If the Universe was infinitely old, even an arbitrarily low vacuum decay rate would be incompatible with our existence. The implications of vacuum metastability can therefore only be considered in the cosmological context, taking into account the finite age and the cosmological history of the Universe. Although the vacuum decay rate is extremely slow in the present day, that was not necessarily the case in the early Universe. High Hubble rates during inflation and high temperatures afterwards could have potentially increased the rate significantly. Therefore the fact that we still observe the Universe in its electroweak vacuum state allows us to place constraints on the cosmological history, for example the reheat temperature and the scale of inflation, and on Standard Model parameters, such as particle masses and the coupling between the Higgs field and spacetime curvature.
In this review we discuss the implications from Higgs vacuum metastability in early Universe cosmology and describe the current state of the literature. We also discuss all the theoretical frameworks, with detailed derivations, that are needed for the final results. In Section 2 we present renormalization group improvement in flat space by using the Yukawa theory as an example before discussing the full SM. Section 3 contains an overview of quantum field theory on curved backgrounds relevant for our purposes, including the modifications to the SM. In Section 4 we go through the various ways vacuum decay can occur. In Section 5 we discuss the connection to cosmology and in Section 6 we present our concluding remarks.
Our sign conventions for the metric and curvature tensors are (−, −, −) in the classification of [12] and throughout we will use units where the reduced Planck constant, the Boltzmann constant and the speed of light are set to unity, ≡ k B ≡ c ≡ 1. The reduced Planck mass is given by Newton's constant as M P ≡ (8πG) −1/2 ≈ 2.435 × 10 18 GeV. (1.2) We will use ϕ for the vacuum expectation value (VEV) of a spectator field (usually the Higgs), φ for the inflaton and Φ for the SM Higgs doublet. The inflaton potential is U (φ) and the Higgs potential V (ϕ). The physical Higgs and top masses read M h and M t .

Example: Yukawa theory
The possibility of quantum corrections destabilizing a classically stable vacuum has been known for quite some time [5,[13][14][15][16][17]. Although our focus will be strictly on the SM, one should keep in mind that the instability that potentially arises in the SM is only a specific example of a more general phenomenon that could manifest in a variety of other theories of elementary particles. For this reason all the essential features of the vacuum instability in the SM can be illustrated with the simple Yukawa theory, which we will now discuss before moving on to the full Standard Model in Section 2.3.
The action containing a massless, quartically self-interacting scalar field ϕ Yukawacoupled to a massless Dirac fermion ψ is Classically, the potential for the scalar field is simply which quite trivially has a well-defined state of lowest energy at the origin. When quantized the potential for the field ϕ becomes modified by quantum corrections which may be investigated within the usual framework of quantum field theory [18]. Importantly, it has been for a long time understood that in some instances predictions in a quantum theory can deviate significantly from those of the classical case. A prime example of such behaviour is radiatively induced symmetry breaking [19].
In the one-loop approximation the result for the quantum corrected or effective potential for the Yukawa model has the form (see, for example, Ref. [20]) with M 2 ϕ (µ) ≡ 3λ(µ)ϕ 2 (µ) ; M 2 ψ (µ) ≡ g 2 (µ)ϕ 2 (µ) . (2.5) In the above we have explicitly denoted the dependence on the renormalization scale µ, which is an arbitrary energy scale, which one needs to choose in order to define the renormalised parameters of the theory. There is also a similar dependence in ϕ(µ) which now refers to the renormalized one-point function of the quantized field, which is related to the bare field via the field renormalization constant [18] ϕ bare = Z(µ)ϕ(µ) . (2.6) In the one-loop effective potential (2.4), the contribution from the fermion ψ comes with a minus sign. For sufficient high values of g, it can overtake the classical contribution and lead to a region with negative potential energy In the limit of large field values ϕ → ∞, one may write the potential as implying that if the potential has a barrier and starts to decrease without bound at high field values [13]. When λ is larger than the critical threshold λ cr the quantum correction approaches +∞ indicating that an arbitrary small deviation from λ cr leads either to +∞ or −∞ at large enough field values.
Hence we have seen that in the Yukawa theory the low-field vacuum will be separated by a barrier from an infinitely deep well on the other side. Even if the barrier is very robust, after a sufficiently long time the system initialized in the classical vacuum must eventually make a transition to the other side of the barrier and evolve towards the state of minimum energy.
A potential unbounded from below is a problematic concept and it is often assumed that, perhaps due to non-perturbative physics invisible to a loop expansion, some mechanism reverses the behavior of the potential at very high energies. This means that the minimum energy is in fact bounded from below, and the effect of the quantum corrections is to generate second local minimum beyond the barrier as depicted in Fig. 1. In theories containing U (1) gauge fields, such as the SM, the reversal of the potential can be shown to happen and the issue of an infinitely deep well does not arise. In the effective theory framework, which arguably is the correct way of viewing the SM, this issue is also not present as one will always encounter a finite scale beyond which the calculation becomes unreliable. Indeed, gravitational corrections are a prime example of a modification that is expected to become significant at large field values.
From a practical point of view, whether or not the potential is infinitely or deep of has a second or more accurately a true minimum beyond the barrier is not important for the generic prediction that the vacuum at the origin should eventually decay if the potential possesses regions with lower energy than at the origin.
However, conclusions based on the behaviour of the perturbative one-loop result (2.4) may be premature. This is because for very large field values the logarithms become nonpertubatively large making the loop expansion invalid: generically one would expect higher powers of the logarithmic contributions in the square brackets of Eq. (2.7) to be generated by higher orders in the expansion, as for example is evident in the results of Ref. [21]. Concretely, for our Yukawa theory (2.1) this requirement means that we can only draw conclusions in the region where 4g 4 64π 2 log g 2 ϕ 2 µ 2 1 and 9λ 2 64π 2 log In principle, the smaller the logarithms the more accurate the result.

Renormalization group improvement
By making use of renormalization group (RG) techniques it is possible to improve the accuracy of an existing perturbative expression such that the issue of large logarithms may be avoided [22][23][24][25].
Demanding that the effective potential (2.4) does not depend on the renormalisation scale µ gives rise to the Callan-Symazik equation [26][27][28] where we have defined the beta functions and the anomalous dimension in the usual manner with γ from the field renormalization constant in (2.6), which has a dependence on the renormalization scale Z ≡ Z(µ). Deriving the beta functions and the anomalous dimension for the Yukawa theory is a well-known calculation (see for example Ref. [23]) and here we simply state the results where for completeness we have included the beta function also for a mass parameter of the scalar field. The beta functions tell us how the values of the renormalised parameters "run", i.e., depend on the scale choice µ. For example, assuming renormalised coupling value g(µ 0 ) at some scale choice µ 0 , one may solve the running of the Yukawa coupling g(µ) from Eq. (2.14), . (2.16) This shows that increasing µ leads to a larger g(µ), and that the coupling g(µ) appears to diverge at scale , (2.17) which is known as the Landau pole [29]. However, well before the Landau pole is reached, the loop expansion ceases to be valid. For more information on the effect of running couplings we refer the reader to more or less any textbook on quantum field theory (for example Refs. [18,30]). Even though the full effective potential V eff (ϕ) has to be independent of the scale choice µ, for any finite-order perturbative result that is only true up to neglected higher-order terms. This means that some scale choices will work better than the others, and by a judicious choice, one can improve the accuracy of the perturbative result. In general, one would choose the scale µ to optimise the perturbative expansion in such a way that the loop corrections are small as indicated in Eq. (2.9). However, for the effective potential (2.7), the loop corrections depend on the field value ϕ. Therefore each given choice of scale would only work well over a relatively narrow range of field values.
To ensure that Eq. (2.9) remains satisfied at any field values, one can take this approach further and make the renormalisation scale a function of the field ϕ, µ = µ * (ϕ), (2.18) so that the expansion is optimised at all field values. This procedure is generically called renormalization group improvement (RGI). 1 This way one can define the renormalization group improved (RGI) effective potential as One should note that in this expression ϕ refers to the field defined at the field-dependent renormalisation scale, ϕ = ϕ(µ * (ϕ)) (for more discussion, see Ref. [20]), and that at any finite order in perturbation theory the resulting function V RGI (ϕ) depends on the choice of the function µ * (ϕ). In principle, one could choose µ * in such a way that the loop correction vanishes exactly. For the one-loop potential (2.7) in the Yukawa theory, this would give where both the couplings g and λ and the field ϕ are renormalised at scale µ exact * (ϕ), and therefore the equation defines the scale µ * (ϕ) implicitly. With this choice, the RGI effective potential V RGI (ϕ) is given simply by the tree-level potential with ϕ-dependent couplings, In more realistic theories it is often impractical to choose µ * (ϕ) that cancels the loop correction exactly [20]. Instead, one chooses some simpler function that keeps the loop correction sufficiently small. The most common choice in the literature is simply Because the loop correction in Eq. (2.7) does not vanish for this scale choice it should still be included in the effective potential. It is nevertheless, fairly common to make the further approximation of dropping it, and writing the tree-level RGI effective potential simply as For weak couplings this is not a good approximation, though. Eq. (2.20) shows that the loop correction vanishes for µ * ≈ gϕ, and therefore a good approximation to RGI effective potential is From this we can see that the use of the tree-level RGI potential (2.23) with the scale choice (2.22) gets the barrier position wrong by a factor of g and the barrier height by a factor of g 4 . Therefore one should either keep the one-loop correction, or use a more accurate scale choice. From the beta function (2.13) for λ, we see that if g 2 λ, λ can become negative at high scales µ. It is conventional to define the instability scale µ Λ as the scale where this happens, If µ Λ < ∞, the effective potential (2.21) becomes negative at high field values, too, implying an instability. Again, the root cause is a negative contribution from the fermions, this time in the beta functions. The solution for the running λ(µ) can be obtained analytically, but is unfortunately quite complicated (see e.g. Ref. [23]). However, it is easy to see that the critical value of the coupling, below which the instability appears, is Close to this critical value one may provide relatively simple analytical results. Suppose we have initial conditions given at some reference scale µ 0 for the running parameters g(µ 0 ) and λ(µ 0 ) the latter of which we parametrize as a fixed value λ cr and a perturbation δλ as By solving Eqs. (2.13)-(2.14) explicitly one may show that λ(µ) has the following expansion (2.28) From Eq. (2.16) it is apparent that, because g(µ) is a monotonically increasing function of µ, the RGI effective potential (2.21) is unbounded from below at large field values, for an arbitrarily small positive perturbation δλ > 0. For comparison, the threshold (2.8) in the unimproved case was λ cr /g 2 = 2/3, somewhat lower than the RGI result (2.26). The above makes apparent a very important generic feature: renormalization group improvement can lead to conclusions that are qualitatively different from the unimproved results. In particular, sizes of couplings deemed as well-behaved and hence giving rise to a stable potential may in fact reveal to result in an instability by the RG improved results. This also implies that close to the critical value the higher loop corrections become quite important as even a small change may tilt the conclusion from stable to unstable, or vice versa. This is also suggested by the fact that the couplings run very gradually and the precise value of the instability scale is very sensitive to small corrections: even a tiny change in the initial values or the running may change µ Λ by several orders of magnitude. These features are illustrated in the example below.
For concreteness, let us consider a numerical example that highlights the importance of renormalization group improvement. Specifically, we choose the Yukawa theory with a Behaviour of the 1-loop RGI effective potential (2.19) (green), the tree-level RGI effective potential (2.23) (blue), and the non-improved result (2.5) (red) with the choices (2.27) at the reference scale µ 0 . The RGI scale choice was µ * (ϕ) = ϕ.
negligible mass parameter and with the initial conditions defined at the renormalisation scale µ 0 as as which from (2.27) can be seen to correspond to a choice that is below the critical value by δλ = 10 −2 . (2.30) Since Eq. (2.29) satisfies λ(µ 0 ) > 2 3 g 2 (µ 0 ) the unimproved effective potential (2.8) implies no instability. This is however not the case after renormalization group improvement as shown in Fig. 2. We must however make sure that the above scale is such that all parameters remain perturbative, in particular for the Yukawa theory we need to check that the g-coupling is sufficiently small. For our parametrization this can be loosely expressed as 2g 2 (µ Λ ) 4π and perturbativity is easily demonstrated with the help of Eq. (2.16). This check is quite important since if g(µ) reaches a large value before µ Λ , it will render the entire derivation inconsistent.
What is also apparent from Fig. 2 that there is a clear difference between the treelevel RGI approximation (2.23) and the full RGI result (2.19), when using the simple non-exact scale choice (2.22). In many applications this would result in a non-negligible inaccuracy, but as shown in Eq. (2.24), it changes the barrier position by a factor O(g) and height by O(g 4 ), which can be important for vacuum stability. This sensitivity to quantum corrections and the choice of µ * comes from the fact that the instability occurs precisely at the point where the tree-level contribution vanishes.

Effective potential in the Standard Model
The SM has a far richer particle content than the simple Yukawa theory of Section 2.1, but the main reason for the possible vacuum instability remains the same: Quantum corrections from the fermions contribute with a minus sign and if significant enough can lead to the formation of regions with lower potential energy than the electroweak vacuum. In the SM the effect is mostly due to the top quark, because it is by far the heaviest and thus has the largest Yukawa coupling. As discussed in Section 4, general field theory principles then dictate that after a sufficiently long time has passed the system should relax into the configuration with the lowest energy resulting in the decay of the electroweak vacuum.
Through increasing experimental accuracy and improved analytic estimates in recent years it has become apparent that the central values for the couplings of the SM allow extrapolation to energy scales close to the Planck scale and that they are in fact incompatible with the situation where the electroweak vacuum would be the state of lowest energy. Some important early works addressing the question of vacuum instability are Refs. [5,[13][14][15][16][17]. The full body of work studying aspects of the vacuum instability is vast (to say the least) and includes Refs. [3, 4, 6-10, 25, 31-87].
The modern high precision era of vacuum instability investigations can be thought to have been initiated by the detailed analyses performed in Refs. [3,4], which presented the first complete next-to-next-to-leading order analysis of the Standard Model Higgs potential and the running couplings.
The current state-of-the-art calculation for the running of Standard Model parameters uses two-loop matching conditions, three-loop RG evolution and pure QCD corrections to four-loop order [79]. The running of the Higgs self-coupling λ is shown in Fig. 3  This depends sensitively on the top and Higgs masses: At 1σ the range is 1.16 × 10 9 GeV < µ Λ < 2.37 × 10 11 GeV, and the case in which λ(µ) is never negative at still included within 3σ uncertainty. Using the three-loop running, and including the one-loop correction in the RGI effective potential with the scale choice µ * (ϕ) = ϕ, the top of the potential barrier lies at ϕ bar = 4.64 × 10 10 GeV, (2.32) and the barrier height is ∆V (ϕ bar ) = V (ϕ bar ) − V (ϕ fv ) = 3.46 × 10 38 GeV 4 = (4.31 × 10 9 GeV) 4 .  of the potential barrier, ϕ bar = 7.70 × 10 9 GeV. Using the unimproved one-loop effective potential with parameters renormalised at the electroweak scale gives as even lower value ϕ bar = 5.78 × 10 4 GeV. This demonstrates that, as discussed in Section 2.2, the use of renormalisation group improvement and the inclusion of at least the one-loop correction in the RGI effective potential are both crucial for accurate results. A slightly more formal issue that must also be kept in mind is that the barrier position ϕ bar is in fact gauge dependent and strictly speaking has limited physical significance [88][89][90][91]. The value of the potential at its extrema are however gauge independent as demanded by the famous Nielsen identity [92]. In the simplest approximation the probability of vacuum decay involves only the values of the potential at the extrema and subtleties involving gauge dependence are evaded. Furthermore, more precise calculations of the rate of vacuum decay, since it is a physical process, can be expected to always be cast into a gauge-invariant form.

Spectator field on a curved background
In the extreme conditions of the early Universe, gravity plays a significant role. In order to investigate the consequences from Higgs metastability we must therefore make use of an approach that incorporates also gravitational effects. This can be achieved in the framework of quantum field theory in a curved spacetime. The study of quantum fields theory on curved backgrounds is hardly a recent endeavour. For a thorough discussion on the subject we refer the reader to the standard textbooks, such as Refs. [93][94][95].
As a representative model we choose an action consisting only of a self-interacting scalar field where the curved background is visible in the metric dependence of integration measure, |g|, the covariant derivative ∇ µ and in the appearance of the non-minimal coupling ξ that connects the field to the scalar curvature of gravity R. The necessity of an operator ∝ Rϕ 2 was discovered already in Refs. [96][97][98], the reasons for which we will elaborate in Section 3.5. It will turn out to be a key ingredient for the implications of the vacuum (in)stability in the early Universe.
Since our discussion assumes a classical curved background with no fluctuations of the metric g µν some effects visible in a complete quantum gravity approach are possibly missed. For energy scales below the Planck threshold and for spectator fields with a negligible effect on the evolution of the background modifications from quantum gravity are expected to be suppressed. For the case of a quasi de Sitter background this was verified in detail in Ref. [99] for the SM Higgs. The reason why quantum gravity is not relevant for a potential SM metastability can be understood from the simple fact that the instability scale (see Section 2.3) is significantly lower than the Planck mass

Homogeneous and isotropic spacetime
From the cosmological point of view it is often sufficient to consider the special case of a homogeneous and isotropic spacetime with the Friedmann-Lemaître-Robertson-Walker (FLRW) line-element given in cosmic time as where a(t) ≡ a is the scale factor describing cosmic acceleration. We will furthermore assume that the energy and pressure densities of the background, ρ and p, are connected via the constant equation of state parameter w as where T µν is the energy-momentum tensor of the background. With the line-element (3.3) the Einstein equation reduces to the the Friedmann equations which allow one to easily find expressions for the Hubble rate H ≡ȧ/a and the scale factor as functions of w For the purposes of this discussion the most important quantity characterizing gravitational effects will be the scalar curvature of gravity R, which may be written as a function of the equation of state parameter and the Hubble rate

Amplified fluctuations
Let us then concentrate on a free quantum theory by setting λ = 0. For this case the action (3.1) leads to the equation of motion + m 2 + ξR φ = 0 , (3.8) whose solutions, as usual, can be expressed as a mode expansion where k is the co-moving momentum and k ≡ |k|. In the above we have also made use of conformal time defined as where the primes denote derivatives with respect to conformal time and we have defined the effective mass Equation (3.11) may be interpreted as that of a harmonic oscillator with a timedependent mass. The crucial point is that for many cosmologically relevant combinations of m, ξ and w the M 2 -contribution is in fact negative. A prime example would be a massless minimally coupled scalar field during cosmological inflation for which m = 0, ξ = 0 and w = −1 giving M 2 = −2H 2 . If M 2 < 0 it is a simple matter to show that the modes with (k/a) 2 + M 2 < 0 contain an exponentially growing branch, which implies that a large field fluctuation can be generated. This effect coming from an imaginary mass-like contribution is sometimes called tachyonic or spinodal instability/amplification [100]. We note that even if no tachyonicity occurs, a large fluctuation can nonetheless be generated if there is a rapid i.e. a non-adiabatic change in M.
A more precise way of understanding the generation of a large fluctuation is by calculating the infrared (IR) portion of the variance i.e. a loop with a low-momentum cut-off Λ IR . This shows that in many situations that can broadly be characterized as having When the theory is not free interactions will via backreaction prevent the generation of arbitrary large fluctuations. In practice one may understand this as the emergence of positive mass-like contributions from the interactions making the field heavy and thus preventing tachyonic or non-adiabatic amplification. The functioning of this mechanism usually allows a significant φ 2 term indicating that quite generally an IR divergence in the free theory implies a large fluctuation when interactions are included.
This rather simple discussion leads to an important implication in regards the vacuum instability problem in the cosmological setting: even if in flat space the decay of a metastable vacuum is enormously unlikely, this may not have been the case during the earlier cosmological epochs when a transition over the barrier can be induced by a large fluctuation generated by the dynamics on a curved background.

Quantum theory in de Sitter space
Even in the simple special case of a de Sitter background it is difficult to perform analytic calculations for an interacting quantum theory. This is mostly due to the non-trivial infrared behavior of quantum fields in de Sitter space [102][103][104][105][106][107][108][109]. A manifestation of this is the lack of a perturbative expansion based on a non-interacting propagator due to the infrared divergence as described in Eq. (3.13). The infrared properties of de Sitter space have attracted significant attention over the years and we refer the interested reader to the review [110] for more information and references.
One popular way forward is to use techniques based on the so-called two-particleirreducible (2PI) diagrams, which are essentially non-perturbative resummations of distinct classes of Feynman diagrams. The 2PI approach is attractive in that it is derivable via first principles from quantum field theory without any approximations. Hence, in principle it can be used up to arbitrary accuracy. Unfortunately, only the leading terms that come by the Hartree approximation are analytically tractable. Applications of 2PI techniques to de Sitter space include Refs. [111][112][113][114][115][116][117].
A non-perturbative framework for calculating quantum correlators in de Sitter space was laid out in Refs. [118,119]. This technique is generally known as the stochastic formalism and is surprisingly straightforward calculationally. It is based on the insight that to a good approximation in de Sitter space one may neglect the quantum nature of the problem and devise a set-up in which the correlators may be calculated from a classical probability distribution P (t, ϕ). If the scalar fieldφ is light, m H, coarse graining over horizon sized patches allows one to approximate its dynamics with a Langevin equatioṅ where V (ϕ) is the classical potential and f (t) is a white noise term satisfying The reason why the 'hat' notation has been dropped from ϕ is that Eq. (3.14) contains only classical stochastic quantities i.e. the quantum features are no longer visible. The Langevin equation (3.14) can be cast in the form of a Fokker-Planck equation for the probability density P (t, ϕ) [119] After a sufficiently long time has passed one would expect that P (t, ϕ) reaches a constant equilibrium distribution. WhenṖ (t, ϕ) = 0, Eq. (3.16) has a simple analytic solution as where N is a normalization factor. As an example, for a theory with only a quartic term V (ϕ) = (λ/4)ϕ 4 , which in many cases is the relevant approximation for the SM Higgs in the early Universe, this results in the equilibrium probability distribution The corresponding field variance becomes This means that the Higgs field develops a non-zero value ϕ ∼ λ −1/4 H, which is sometimes called a condensate [120][121][122][123][124][125][126][127].
The central assumption that leads to the stochastic description is that the effect of the ultraviolet physics on the infrared behaviour can be described as a white noise term in the Langevin equation (3.14). The ultraviolet modes also contribute to the effective potential V (ϕ) in the Fokker-Planck equation (3.16), as was discussed in Section 2.2 in flat space. Especially when investigating the vacuum stability of the SM it is imperative that the quantum corrections are incorporated in the stochastic approach, for example by making use of the RGI effective potential as the input in Eq. (3.16).

Curvature corrections to the effective potential
It is clearly evident from the derivations of Section 3.3 that a scalar field in curved spacetime feels the curvature of the background. It then follows that also the effective potential must receive a contribution from curvature. In order to reliably investigate the implications from the SM metastability in the early Universe these contributions then must be included in a discussion of quantum corrections to the potential.
Deriving the effective potential for a quantized scalar field on a curved background is naturally much more difficult than in flat space: for many backgrounds even the case of a free scalar field admits no closed form solutions for the mode equation [93]. Another complication that arises is that choosing the boundary condition i.e. the specific quantum state in which the effective potential is calculated is far from obvious. This is due to the fact that in curved space the concept of a particle and hence the vacuum state is no longer well-defined globally, but depends on the specific dynamics and perceptions of a given particular observer [143].
However, even on an arbitrary curved background some things remain universal: renormalizability of a quantum field theory imposes the requirement that all quantum states should have coinciding divergences. From this it follows that it is possible to derive an effective potential retaining terms only originating from the very high ultraviolet (UV), which is a contribution that is always present irrespective of the quantum state one is interested in. Such an expression would then allow one to determine all the generated operators and their respective runnings, as RG effects are ultimately the result of UV physics.
Let us once more study the Yukawa theory of Section 2.1 only this time in curved spacetime and without neglecting the mass parameter for the scalar. In curved spacetime the action reads The most convenient way of deriving the effective potential is the Heat Kernel method reviewed in Ref. [144], see also [142]. This approach has been known for a long time, see Refs. [145][146][147][148][149][150] for early work. We will make use of the resummed form of the Heat Kernel expansion derived in Refs. [151,152], which for the action (3.20) gives (for details, see Ref. [20]) where the one-loop quantum corrections from the scalar and the fermion, V respectively, and the curved space effective masses M 2 ϕ and M 2 ψ are now The R µν and R µνα β are the Ricci and Riemann tensors, respectively. We have introduced the absolute values in the logarithms to ensure that the result is never complex. A complex effective potential in flat space can be interpreted as a finite lifetime of the quantum state [153], but this is ultimately an infrared effect and hence not correctly represented in an UV expansion. Therefore, the effective potential in curved space derived with the Heat Kernel expansion correctly represents the local physics and can for example be used to determine the running of parameters in curved space and the possible generation of new operators (see the next section), but in order to answer questions about vacuum decay one needs additional technology, which is discussed in Section 4. What the above clearly shows is that on a curved background a highly non-trivial dependence on the curvature emerges: A curved spacetime leads to the generation of additional operators that couple to the scalar field. Importantly, the non-minimal term ∝ Rϕ 2 directly coupling ϕ to R is not the only one, but terms ∝ R 2 , R µν R µν and R µνδη R µνδη are also unavoidable and they couple to the scalar field via the logarithmic loop contributions. These terms are not necessarily small, for example in de Sitter space with a constant Hubble rate H the various curvature contributions may be written as and in the early Universe the Hubble rate can be several orders of magnitude larger than any mass parameter of the SM. Simply put, since curvature is felt by the scalar field its inclusion in the calculation is vital for making robust predictions because the scale provided by H often is the largest scale of the problem.

Running couplings in curved space
The basic principles laid out in the flat space analysis of Section 2.2 remain unchanged when the background in no longer flat: Demanding a result independent of the renormalization scale µ leads to the Callan-Symanzik equation from which the beta functions may be solved given the anomalous dimension γ. Since γ is a dimensionless number it will receive no contributions from constants associated with the curvature of space such as ξ. Otherwise parameters only visible in the action when the background is curved would nonetheless influence the RG running of, say, λ. Similar arguments imply that all beta functions present in flat space remain unchanged when the background is curved. As one may see from Eqs. (3.21) -(3.23) operators that are not present in the treelevel action (3.20) are generated by the loop correction. This means that even if one renormalizes these terms to zero, the may resurface via RG running. Ultimately, this is the reason behind the non-minimal term ∝ Rϕ 2 already in (3.20). For the same reason in our theory we must include the following purely gravitational action which along with Eqs. (2.12)-(2.15) provide a complete set of RG equations for the Yukawa theory in curved spacetime. A crucial difference to the flat space case arises when implementing renormalization group improvement. In Section 2.2 we exploited the fact that the full quantum result must be independent of the renormalization scale µ in order to optimize the pertubative expansion. Namely, we made the choice (2.18) in order to keep the logarithms small also at large scales. In curved space the logarithms in the loop corrections (3.22) and (3.23) have dependence on the scalar curvature R, and therefore it must be included in the optimization. The one-loop calculation shows that the exact scale choice that would fully cancel the loop correction is not possible across the whole range of field values [20]. Instead, a sensible choice for the optimized scale µ * is a linear combination of ϕ 2 and R i.e.
where the parameters a and b are chosen in such a way that the logarithms remain under control. Equation (3.28) highlights an often neglected effect arising in curved spaces after renormalization group improvement: In a curved background the optimal scale choice depends significantly on curvature. This phenomenon may be characterized as curvature induced running and was recently studied in detail for the full SM in Ref. [20]. In situations where the curvature of the background is significant it can give the dominant contribution to the scale. Considering the metastability of the SM in the early Universe this in fact is often the case as during and after inflation one may have a Hubble rate much larger than the instability scale, H µ Λ .

The Standard Model
The Standard Model particle content can be expressed with the Lagrangian The first three terms in Eq. (3.29) describe the contributions coming from the gauge fields, the fermions and the Higgs doublet Φ whose one point function we write as Φ ≡ ϕ, from now one dropping the hats. The 'GF' and 'GH' are the gauge fixing and ghost Lagrangians, respectively. Here we show explicitly only the Higgs contribution (for the full result see Ref. [20]) with the SM covariant derivative where ∇ µ contains the connection appropriate for Einsteinian gravity, g, g and Y (A a µ and A µ ) are the SU (2) and U (1) gauge couplings (fields) and σ a the Pauli matrices.
As de Sitter space is the most important application of our results here we show the perturbative 1-loop correction for the SM in a spacetime with an equation of state w = −1 i.e. a constant Hubble rate H (see Section 3.2) where the sum is over all degrees of freedom of the SM, which may be found in tables 1 and 2. The masses are defined as and the ζ i are the gauge fixing parameters. The flat space beta functions have of course been known for some time, see for example Refs. [4,25]. The complete set of SM beta functions to 1-loop order for the SM was however first calculated only in Ref. [20]. The 1-loop SM beta functions for couplings associated with gravity coming from the action (3.26), ξ, V Λ , κ, α 1 , α 2 and α 3 can be solved from Table 1. The 1-loop effective potential (3.32) contributions with tree-level couplings to the Higgs. Ψ stands for W ± and Z 0 bosons, the 6 quarks q, the 3 charged leptons l, the Higgs h. The Goldstone bosons are χ W and χ Z and ghosts c W and c Z . The masses may be found in Eq. (3.33).
Eq. (3.32) and read with y i being a Yukawa coupling for a fermion type i. Table 2. Contributions to the effective potential (3.32) with no coupling to the Higgs at tree-level. The Ψ include the photon γ, the 8 gluons g, the 3 neutrinos ν and the respective ghosts c γ and c g .
Much like in the flat space case in Eq. (2.19) we can write the RGI effective potential by choosing an optimized scale µ * (ϕ, R) in such a way that the loop correction is small [20]. In curved space in addition to the Lagrangian from Eq. (3.30) we must include the purely gravitational terms from Eq. (3.26) in addition to the one loop contributions (3.32) giving rise to where µ * generally depends on both ϕ and R, and we have assumed |R| |m 2 h |, which is usually true for the SM Higgs in the early Universe.
When the Hubble rate is above electroweak scales it is quite obvious that the highly non-trivial curvature dependence apparent in Eq. (3.32) and also in Eq. (3.41) with the optimized scale (3.28) cannot be neglected: it is just as, if not more, important as what would have been obtained by using only a flat space derivation. The most obvious difference is the emergence of the direct non-minimal coupling between the Higgs and the scalar curvature R. Due to the curvature dependence of the optimized renormalization scale in curved space (3.28), which can be traced back to the curvature dependence of the one-loop correction (3.32), the generation of the non-minimal coupling in the current cosmological paradigm is unavoidable. It will be sourced by the changing Hubble rate H. Furthermore, as can be read from the beta function (3.34), ξ = 0 is not a fixed point of the RG flow. Depending on the sign of ξR, the non-minimal coupling can have a stabilizing or destabilizing effect, which can be very significant in the early Universe.  In Fig. 4 we illustrate the behavior of effective potential for the full SM in de Sitter space including the one loop quantum correction (3.32). We have chosen to set the renormalised non-minimal coupling ξ to zero at the electroweak scale. We use the field renormalized at the physical top mass as the x-axis. It is clearly evident that in curved space the potential may have drastically different predictions to flat space. As can be read off from the beta function for ξ (3.34), if ξ = 0 at some low scale, it will run to negative values at high scales. Furthermore, since in de Sitter space R = 12H 2 > 0, a negative ξ can prevent the emergence of a potential barrier, even if robustly present on a flat background, as visible in Fig. 4.

Quantum tunnelling and bubble nucleation
The main mechanism behind vacuum decay in the Standard Model is essentially a direct extension of ordinary quantum tunnelling to quantum field theories. In ordinary quantum mechanics, the wave-function for particles trapped by a potential barrier can penetrate the classically forbidden region of the barrier, leading to a non-zero probability to be found on the other side. The transition rate for particles of energy E incident on a barrier described by potential W (x) can be estimated using the WKB method [154], where x 1 , x 2 are the turning points of the potential. As is clear from this expression, the tunnelling rate is suppressed by wide and tall barriers. Although Eq. (4.1) can in principle be evaluated directly, we will follow a different approach that readily generalises to quantum field theories [155,156]. The idea is to use the equation of motion, The region (x 1 , x 2 ) is classically forbidden, since W (x) − E > 0 there. We can apply a trick, however, by analytically continuing time to an imaginary value: τ = it, which gives a Euclidean equation of motion, The most notable feature of these equations is that the potential has effectively been inverted. This means that we can find a classical solution that rolls through the barrier between the turning points x 1 and x 2 . If we can find this solution, it allows us to re-express the integral in Eq. (4.1) as where S E is the Euclidean action corresponding to Eq. (4.3) is a bounce solution of the Euclidean equations of motion satisfying x (τ 1 ) = x (τ 2 ) = 0, and x fv (τ ) is a constant solution, sitting in the false vacuum with energy E. The 'bounce' solution is so named because we see, by energy conservation, that it starts at x 1 , rolls down the inverted potential before 'bouncing' off x 2 and rolling back. By finding this solution and evaluating its action, we can compute the rate for tunnelling through a barrier. This argument generalised straightforwardly to many-body quantum systems, where we use the action (4.6) With more than one degree of freedom, however, there are actually an infinite number of paths that q i (τ ) could take when passing through the barrier, corresponding to an infinite number of solutions. However, since the decay rate is exponentially dependent on the action, Γ ∝ e −S E [q i ] , it is clear that only the solution with smallest Euclidean action will contribute significantly, as this will dominate the decay rate (in other words, the tunnelling takes the 'path of least resistance').
The generalisation from a many body system, q i , to a quantum field theory with scalar field ϕ(x) is then straightforward, The integral here is over flat four-dimensional Euclidean space, and note that the opposite sign of the potential leads to an opposing sign in the equations of motion, Although it is tempting to interpret V (φ) as the potential to be tunneled through, this is only somewhat true. The analogue of W (q i ) in Eq. (4.6) is a functional of the field configuration ϕ(x), given by an integral over three-dimensional space space, where ∇ϕ represents the spatial derivative of the field. In the analogy with quantum mechanics, this term should be considered part of the potential, as its many body equivalent is a nearest-neighbour interaction between adjacent degrees of freedom, q i , q i±1 . This means, in particular, while in quantum mechanics, the particle emerges after tunneling at a point x 2 that has the same potential energy, W (x 1 ) = W (x 2 ), in quantum field theory, the field emerges lower down the potential V . In a field theory, the analogue of x 2 is a field configuration, ϕ(x), given by slicing the bounce solution at its mid-way point. This is a nucleated 'true-vacuum' bubble, whose decay rate is determined by the Euclidean action of the bounce solution, ϕ B . As we will see in Section 4.7, the dominant Euclidean solutions have O(4) symmetry, which means that the bubble nucleates with O(3, 1) symmetry. This causes it to expand at near the speed of light, resulting in the space around a nucleation point being converted to the true vacuum, releasing energy into the bubble wall. Apart from the destruction that this would unleash, and the different masses of fundamental particles in the bubble interior, the result is also gravitational collapse of the bubble [157], making its nucleation in our past light-cone completely incompatible with the trivial observation that the vacuum has not decayed (yet).
In cosmological applications, but also other areas, it is also important to consider the effect of thermally induced fluctuations over the barrier. Brown and Weinberg [156] describe how thermal effects can be included in the above argument. At non-zero temperature, we must integrate over the possible excited states, and the decay exponent which depends on energy, where B(E) is the (energy dependent) difference in Euclidean action between the bounce solution and the excited state of energy E. This integral is dominated by the energy that minimises the exponent βE + B(E), which is easily shown to satisfy where τ 1 , τ 2 are the initial and final values in imaginary time of the (energy dependent) bounce solution. In other words, the bounce solution is periodic in imaginary time, with period controlled by the temperature.
In quantum field theory, the decay rate per unit volume and time of a metastable vacuum decays was first discussed by Coleman [155,158], and is given by is the difference between the Euclidean action of a so called bounce solution ϕ B of the Euclidean (Wick rotated) equations of motion, and the action of the constant solution ϕ fv which sits in the false vacuum. S denotes the second functional derivative of the Euclidean action of a given solution, and det denotes the functional determinant after extracting the four zero-mode fluctuations which correspond to translations of the bounce (these are responsible for the formula giving a decay rate per unit volume). Precise calculations of the pre-factor A in the Standard Model were performed in [6], and involve computing the fluctuations around the bounce solution of all fields that couple to the Higgs. This requires renormalising the loop corrections, and also to avoid double-counting, expanding around the tree-level bounce, rather than the bounce in the loop corrected potential.
In the gravitational case, the prefactor A is harder to compute. The main issue is that it includes both Higgs and gravitational fluctuations, and without a way of renormalizing the resulting graviton loops, the calculation becomes much harder. Various attempts have been made to do this using the fluctuations discussed in Section 4.5 (see Refs. [159][160][161] for example), but a full description, especially for the Standard Model case, is not yet available.
In most cases, it is reasonable to estimate the prefactor A using dimensional analysis. Because A has dimension four, one would expect where µ the characteristic energy scale of the instanton solution. Due to the exponential dependence on the decay exponent, B, this will not lead to large errors, and therefore we will use this result in the absence of more accurate estimates.

Asymptotically flat spacetime at zero temperature
In flat Minkowski space, the bounce solution corresponds to a saddle point of the Euclidean action, with one negative eigenvalue (see Section 4.5). Since Eq. (4.12) depends exponentially on the bounce action, only the lowest action bounce solutions will contribute. In flat space, it is always the case that the lowest action solution has O(4) symmetry [162]. This means that the equations of motion for the bounce can be reduced tö subject to the boundary conditionsφ(0) = 0 and ϕ(r → ∞) → ϕ fv . These ensure that the bounce action is finite and thus gives non-zero contribution to the decay rate. There are always trivial solutions corresponding to the minimal of the potential V (ϕ), but they do not contribute to vacuum decay because they have no negative eigenvalues. For example, in a theory with a constant negative quartic coupling, that is, there exists the Lee-Weinberg or Fubini bounce [163,164]. This is a solution of the form: where the arbitrary parameter r B characterises the size of the bounce (and thus the nucleated bubble). This arbitrary parameter appears in the theory because the potential Eq. (4.17) is conformally invariant, and thus bounces of all scales contribute equally with action In fact, similar bounces contribute approximately in the Standard Model, where the running of the couplings breaks this approximate conformal symmetry, so that bounces of order the scale at which λ is most negative (which is the minimum of the λ(µ) running curve) dominate the decay rate [6]. The complete calculation would also include gravity, and would therefore involve finding the corresponding saddle point of the action where R is the Ricci scalar. The leading gravitational correction to Eq. (4.19) is [165] ∆S gravity = 256π 3 45(r B M P λ) 2 . Another approach is to solve the bounce equations numerically, which makes it possible to use the exact field and Einstein equations and the full effective potential. The difference is a second order correction [165]. Using the tree-level RGI effective potential ( A non-minimal value of the Higgs curvature coupling ξ changes the action and the shape of the bounce solution (and thus the scale that dominates tunneling) [85,[165][166][167][168]. Fig. 5 shows the bounce action B as a function of ξ, computed numerically in Ref. [166]. As the plot shows, the action is smallest near the conformal value ξ = 1/6. For ξ ≈ 1/6, the result agrees well with the perturbative calculation [85]), For comparison, for the same parameters, the numerically computed decay exponent in flat space is [166] B flat = 1805.8, Calculations of the prefactor A show that the decay rate (4.12) is well approximated by [6] Γ ∼ µ 4 min e −B ∼ 10 −716 GeV 4 , (4.26) where the numerical value corresponds to the action (4.22). This agrees with the estimate from dimensional analysis (4.14). Note, however, that the rate is very sensitive to the top quark and Higgs boson masses, and also to higher-dimensional operators [66,78]. The presence of a small black hole can catalyse vacuum decay and make it significantly faster [169][170][171][172][173]. The action of the vacuum decay instanton in the presence of a seed black hole is given by where M seed and M remnant are the masses of the seed black hole and the left over remnant black hole. For black holes of mass M seed 10 5 M P ≈ 1g the vacuum decay rate becomes unsuppressed. This can be interpreted [173,174] as a thermal effect due to the black hole temperature T seed = M 2 P /M seed . The catalysis of vacuum decay does not necessarily rule out cosmological scenarios with primordial black holes, because positive values of nonminimal coupling ξ would suppress the vacuum decay in the presence of a black hole [175].

Non-zero temperature
The presence of a heat bath with non-zero temperature has a significant impact on the vacuum decay rate Γ [36,176]. On one hand, the thermal bath modifies the effective potential of the Higgs field, and on the other hand, as discussed in Section 4.1, it modifies the process itself because it can start from an excited state rather than the vacuum state.
At one-loop level, the finite-temperature effective potential can be written as [36] V where n i and M 2 i are given in Table 1 (taking H = 0). In the high-temperature limit, T M h , this can be approximated by where γ 2 ≈ 1 12 Therefore the thermal fluctuations give rise to a positive contribution to the quadratic term. This raises the height of the potential barrier, and therefore would appear to suppress the decay rate.
At non-zero temperatures the decay process is described by a periodic instanton solution with period β in the Euclidean time direction. In the high-temperature limit, the solution becomes independent of the Euclidean time, and has the interpretation of a classical sphaleron configuration. The instanton action is therefore given by where E sph is the energy of the sphaleron, which is the three-dimensional saddle point configuration analogous to the Coleman bounce (4. 16), and satisfies the equation Using the approximation of constant negative λ, the action is [36] Because γ 1, this is smaller than the zero-temperature action (4.19). Therefore the net effect of the non-zero temperature is to increase the vacuum decay rate compared to the zero-temperature case.
More accurately, the sphaleron solutions have been calculated numerically in Ref. [85,177]. At high temperatures T 10 16 GeV, the action is roughly B(T 10 16 GeV) ∼ 300. (4.34) When the temperature decreases, the action increases, so that B(10 14 GeV) ∼ 400. Salvio et al. [85] obtained fully four-dimensional instanton solutions numerically, without assuming independence on the Euclidean time, and found that the three-dimensional sphaleron solutions have always the lowest action and are therefore the dominant solutions. They also showed that including the two-loop corrections to the quadratic term (4.30) or the one-loop correction to the Higgs kinetic term gives only small correction to the action.
Taking also the prefactor into account, the vacuum decay rate at non-zero temperature is [177,178]

Vacuum decay in de Sitter space
In extending from flat space to curved space, the theorem [162] that guarantees O(4) symmetry of the bounce no longer applies. There is some evidence, however, that in background metrics that do respect this symmetry, O(4) symmetric solutions should still dominate [179]. This would include the special case of particular interest in this review -an inflationary, or de Sitter background. 3 A Wick rotated metric can be placed in a co-ordinate system that makes the O(4) symmetry of the bounce immediately manifest, where χ is a radial variable, dΩ 2 3 is the 3-sphere metric, and a 2 (χ) is a scale factor that physically describes the radius of curvature of a surface at constant χ. The bounce equations of motion then take the form [157] ϕ + 3ȧ aφ − V (ϕ) = 0 (4.37) We will consider the case in which the false vacuum has a positive energy density, V (ϕ fv ) > 0, and therefore non-zero Hubble rate (4.39) The boundary conditions the bounce solution must satisfy require special attention: a(0) = 0 is required because of the definition of a(χ) as a radius of curvature of a surface of constant χ, while we requireφ(0) =φ(χ max ) = 0, where χ max > 0 is defined by a(χ max ) = 0. These boundary conditions avoid the co-ordinate singularities at χ = 0, χ max giving infinite results, but allow for the peculiar property that the bounces are compact, and do not approach the false vacuum anywhere.
One way of understanding this peculiar feature was discussed by Brown and Weinberg [156]. They considered vacuum decay in de Sitter space, specifically the static patch co-ordinates where the metric takes the form where dΩ 2 n−2 is the n − 2-sphere metric (in this case, n = 4). The important feature of these co-ordinates is that they are valid only up to the horizon at r = 1/H. The Euclidean action can then be re-written as where h ij is the remaining spatial metric. Brown and Weinberg interpreted this to mean that tunneling takes place on a compact Euclidean space, with a curved three-dimensional geometry. This compactness condition is reflected in the boundary conditionsφ(0) = ϕ(χ max ), which inevitably produce a compact bounce solution. They observed that the same effect would be seen in considering a spatially curved universe with this same spatial geometry, but with a non-zero temperature, This corresponds to the Gibbons-Hawking temperature of de Sitter space [143], and implies that bounces in de Sitter space may have a thermal interpretation. [180]. This is a constant solution, for which ϕ = ϕ bar sits at the top of the barrier for the entire Euclidean period, and the scale factor is given by

The simplest solution of Eqs. (4.37) and (4.38) is the Hawking-Moss solution
(4.43) Hence χ max = π/H HM . The action difference of Eq. (4.13) is then easily computed analytically to be . (4.44) A particularly important limit is that in which ∆V (ϕ bar ) = V (ϕ bar ) − V (ϕ fv ) V (ϕ fv ). In that case, Eq. (4.44) is approximately where H 2 = V (ϕ fv )/3M 2 P is the background Hubble rate. The prefactor (4.14) in the decay rate can be expected to be at the scale of the Hubble, and therefore the vacuum decay rate due to the Hawking-Moss instanton can be approximated by Eq. (4.45) has a simple thermal interpretation: It is the ratio of the energy required to excite an entire Hubble volume, 4π/3H 3 from the false vacuum to the top of the barrier, divided by the background Gibbons-Hawking temperature (4.42). Therefore it can be understood as Boltzmann suppression in classical statistical physics.
The bounce equations (4.37) and (4.38) also often have Coleman-de Luccia (CdL) instantons, in which the field increases monotonically from ϕ(0) < ϕ bar to ϕ(χ min ) > ϕ bar . For low false vacuum Hubble rates, H µ min , a CdL solution can be found as a perturbative correction to Eq. On the other hand, the Hawking-Moss action (4.44) decreases rapidly as the Hubble rate increases. It crosses below B CdL at Hubble rate [182] H cross = 1.931 × 10 8 GeV.
(4.48) Decay exponent, B  At Hubble rates below this, H > H cross vacuum decay is dominated by the Colemande Luccia instanton, which describes quantum tunneling through the potential barriers, whereas above this, H > H cross , the dominant process is the Hawking-Moss instanton. This is discussed further in Section 4.6.
In addition to the HM and CdL solutions, one may also find oscillating solutions [183][184][185][186], which cross the top of the barrier ϕ bar multiple times between χ = 0 and χ = χ max , and additional CdL-like solutions with higher action [182,184]. The latter were found numerically in the Standard Model in Ref. [182]. Because these solution have a higher action than the HM and CdL solutions, they are highly subdominant as vacuum decay channels. Oscillating solutions also have more than one negative eigenvalues [159,187].

Negative eigenvalues
In order for a stationary point of the action to describe vacuum decay, it has to have precisely one negative eigenvalue. The reason is that the decay rate of a metastable vacuum is determined by the imaginary part of the energy as computed by the effective action [158], and thus only solutions that contribute an imaginary part to the vacuum energy will contribute to metastability.
This requirement comes in via the functional determinant which encodes the quantum corrections to the bounce solution. This functional determinant is given by a product over the eigenvalues for fluctuations around the relevant bounce solution. In flat space, these all satisfy [158] − ∇ µ ∇ µ δϕ + V (ϕ B )δϕ = λδϕ, (4.49) where ϕ B is the solution expanded around. The O(4) symmetric bounce solutions in flat space can be shown to have at least one negative eigenvalue, since they possess zero modes corresponding to translations of the bounce around the space-time. In fact, there must only be one such eigenvalue. Solutions with more negative eigenvalues do not contribute to tunneling rates, because while they are stationary points of the Euclidean action, they are not minima of the barrier penetration integral (4.1) obtained from the WKB approximation [154]. The situation is somewhat different in the gravitational case, however, due to the fact that in addition to the scalar field, we can also consider metric fluctuations about a bounce solution. Several authors have constructed a quadratic action for fluctuations about a bounce in curved space [160,161,187]. This takes the gauge invariant form where and f is a complicated function of a and ϕ which can be found in Refs. [160,161,187]. The analysis of this Lagrangian is complicated, but some conclusions can be drawn. To begin with, it is possible to argue that expanded around a CdL bounce solution, this action always has an infinite number of negative eigenvalues. This is the so called 'negative mode problem' [160,161,187]. The argument, as expressed in Ref. [160], is that we can re-write Q using Eq. (4.38) as (4.52) Note that the bounce always has a point satisfyingȧ = 0, which is the largest value obtained by a(χ). Consequently, there is always a region where Q is negative, so for the l = 0 modes it is possible to construct a negative kinetic term in Eq. (4.50). This means that sufficiently rapidly varying fluctuations will have action unbounded below, so there is an infinite tower of high frequency, rapidly oscillating fluctuations that all have negative eigenvalues. Note that for l = 1 the quadratic Lagrangian is zero (these are the zero-modes associated to translations of the bounce), and for l > 1, both numerator and denominator in Eq. (4.50) are negative, thus the kinetic terms are always positive. Since Q = 1 in flat space (obtained by taking the M P → ∞ limit), it is clear that these 'rapidly oscillating' modes are somehow associated to the gravitational sector. At first, this seems concerning, however, it was pointed out in Ref. [160] that these high frequency oscillations are inherently associated with quantum gravity contributions, and thus may not affect tunneling. If we focus on the 'slowly varying' modes, the structure of these is much more similar to the analogous flat space bounces. The conclusion we should draw then, is that a solution is relevant only if there is a single slowly varying negative eigenvalue.

Hawking-Moss/Coleman-de Luccia transition
As discussed in Section 4.4, there are two types of solutions that contribute to vacuum decay in de Sitter space. The first is the Hawking-Moss solution (4.43), and the second is the Coleman-de Luccia solution, which crosses the barrier once. By considering the negative eigenvalues of the HM solution, one gains insight into which solutions exist and contribute to vacuum decay at a given Hubble rate.
Higher modes will all be positive if and only if This imposes a lower bound on H HM , below which the Hawking-Moss solution has multiple negative eigenvalues. Hence, it cannot contribute to vacuum decay for Hubble rates below the critical threshold [154,156]. An alternative way of expressing this is in terms of a critical Hubble rate. If we define H 2 = V (ϕ fv )/3M 2 P to be the background Hubble rate in the false vacuum, then the condition for Hawking-Moss solutions to contribute to vacuum decay is H > H crit where (4.56) Here, ∆V (ϕ) ≡ V (ϕ) − V (ϕ fv ). However, the second term generally only contributes significantly if the difference in height between the top of the barrier and the false vacuum is comparable to the Planck Mass. For most potentials, only the second derivative at the top of the barrier matters. At low Hubble rates, H < H crit , the Hawking-Moss solution does not contribute to vacuum decay, but on the other hand, a CdL solution is guaranteed to exist [189]. In most potentials, the CdL solution smoothly merges into the Hawking-Moss solution as the Hubble rate approached H crit from below, and and the Hawking-Moss solution becomes relevant [184,190]. Close to the critical Hubble rate, H ∼ H crit , one can define the quantity [190][191][192] which divides potentials into two classes [182,190]. Those with ∆ < 0 are 'typical' potentials, for which the perturbative solution only exists for H < H crit [190], while those with ∆ > 0 only have perturbative solutions for H > H crit . When a perturbative solution exists, its action is given by [190] where ϕ 0 is the true vacuum side value of the bounce (which approaches ϕ HM in the H → H crit limit) and ϕ bar is the top of the barrier.
Hence one can see that if ∆ < 0, a CdL solution with lower action, B CdL < B HM , exists for H < H crit , and approaches the Hawking-Moss solution smoothly as H → H crit , until it vanishes at H crit . At the same point, the second eigenvalue of the HM solution turns positive, and therefore the HM solution starts to contribute to vacuum decay.
On the other hand, if ∆ > 0, which is the case for the Standard Model Higgs potential [182], the perturbative CdL solution exists only for H > H crit . Below H crit , the HM solution has two negative eigenvalues, which means that it does not contribute to vacuum decay. Instead, the relevant solution is the CdL solution, which also has a lower action (see Fig. 6). When the Hubble rate is increased, a second, perturbative CdL solution appears smoothly at H = Hcrit, at the same as the second eigenvalue of the HM solution becomes positive. At H > H crit there are, therefore, at least two distinct CdL solutions, and in fact, numerical calculations indicate that there are at least four [182]. For the parameters used in Fig. 6, the critical Hubble rate is H crit = 1.193 × 10 8 GeV.

Evolution of bubbles after nucleation
The bounce solution ϕ B determines the field configuration to which the vacuum state tunnels [156,158], and therefore sets the initial conditions for its later evolution. It is the equivalent of the second turning point on the true vacuum side, x 2 , appearing in Eq. (4.1). In ordinary quantum mechanics, a particle with energy E emerges on the true vacuum side of the barrier at x 2 (E) after tunnelling. This is related to the bounce solution, which starts at x 1 , rolls until reaching x 2 , and then bounces back to x 1 , thus x 2 represents a slice of the bounce solution half way through.
In complete analogy, the field emerges at a configuration corresponding to a slice half way through the bounce solution (in Euclidean time). In flat space tunnelling, the bounce is ϕ B (χ) where χ 2 = τ 2 + r 2 , and thus touches the false vacuum at τ → ±∞. Hence the mid-way points occurs at τ = 0 and the solution emerges with φ(x, 0) = ϕ B (r). One can then use this as an initial condition at t = 0 for the Lorentzian field equations, However, this is not really necessary, as the O(4) symmetry of the bounce solution carries over into O(3, 1) solution [158], and thus the solution can be read off as ϕ(x, t) = ϕ B ( r 2 − t 2 ) for r > t. (4.60) From this one can see that the bubble wall is moving outwards asymptotically at the speed of light. The inside of the light cone corresponds to an anti-de Sitter spacetime collapsing into a singularity [172,178,193,194]. The situation in de Sitter space is considerably more complicated, but the conclusion is the same [156]. First, de Sitter bounces can be thought of as bounces at finite temperature on a curved spatial background described by constant time slices of the static patch of de Sitter space, The temperature in this case is the Gibbons-Hawking temperature (4.42) of de Sitter space. Bounces at finite temperature β = 1/k B T correspond to periodic bounces in Euclidean space [156], with period τ period = β. In this case, the bounce starts at the false vacuum at τ = −π/H, hits its mid-point at τ = 0, and returns to the false vacuum side at τ = π/H. Thus, the τ = 0 hypersurface describes the final state of the field after tunnelling. Analytic continuation of the metric back to real space can be performed using the approach of [172]. The O(4) symmetric Euclidean metric is of the form (4.62) where in the de Sitter case, Since it is straightforward to analytically continue the flat space metric back to real space via the transformation τ = it, then the same thing can be done with any conformally flat metric, by changing variables toτ ,r such that (4.64) which is achieved by choosing f (χ) such that f (χ) = f /a, f (0) = 0. In the de Sitter case, this means f (χ) = C sin(Hχ) 1 + cos(Hχ) = C tan(Hχ/2), (4.65) where C is an arbitrary constant -we can choose it to be 1. This co-ordinate system is obtained from the O(4) symmetric co-ordinates viã One then transforms back to real space exactly as in flat space, viaτ = it. The co-ordinate χ is then related tot andr via It should be noted thatt,r as defined only cover ther >t portion of de Sitter space. Because the metric is manifestly conformally flat in these co-ordinates, we can see that this corresponds to the portion of de Sitter space outside the light-cone, which lies at r = ±t. Doing this for de Sitter yields the real space metric (4.69) which at first glance, is not obviously de Sitter space. However, the transformation , (4.71) can be readily shown to yield Eq. (4.61), thus this is indeed a valid analytic continuation of the Euclidean 4-sphere back to de Sitter space.
To describe the subsequent evolution of the bubble, it is argued in Ref. [172] that φ(r, t) = φ B (χ(r, t)) matches the symmetry of the O(4) symmetric bounce, just as in flat space, with χ(r, t) defined by Eq. (4.68). As mentioned before, this describes only the evolution of the scalar field outside the light-cone. Forr <t, it is necessary to solve the Euclidean equations directly. That calculation demonstrates explicitly that the formation of a singularity in the negative-potential region is inevitable [172], confirming previous calculations using the thin wall approximation in Ref. [157].
As for the evolution outside the light-cone, it can be seen that, much as in flat space, a point of constant field value ϕ 0 corresponding to χ 0 where ϕ 0 = ϕ(χ 0 ), satisfies r(t) = t 2 + f 2 (χ 0 (φ 0 )), (4.72) which means that it rapidly approaches the speed of light ast → ∞. Thus, just as in flat space, the bubble expands outwards at the speed of light. Even if the bubble wall moves outward at the speed of light, it does not necessarily grow to fill the whole Universe, if it is trapped behind an event horizon. Scenarios in which bubbles of true vacuum form primordial black holes have been discussed [195][196][197][198]. This can happen if inflation ends before the space inside the bubble hits the singularity. When the Universe reheats, thermal corrections (4.28) stabilise the Higgs potential, preventing the collapse. The bubble then collapses into a black hole, and the primordial black holes produced in this way could potentially constitute part or all of the dark matter in the Universe [197]. This scenario requires fine tuning to avoid the singularity or new heavy degrees of freedom that modify the potential at high field values [198]. The same scenario can also produce potentially observable gravitational waves [199].

Cosmological history
For the Universe to be currently in a metastable state rather than in its true ground state, it is not enough that the decay rate is slow today. The Universe also had to somehow end up in the metastable electroweak-scale state, and the decay rate had to be sufficiently slow in the past for the Universe to stay there through the whole history of the Universe. The former requirement depends on the initial conditions of the Universe, which are often assumed to involve Planck-scale field values, and therefore one needs to explain how the Higgs could have relaxed into the electroweak-scale vacuum without getting trapped into the negative-energy true vacuum. The latter condition, the survival of the current metastable state through the history of the Universe, requires that no bubbles of true vacuum were nucleated in our past light cone [178]. This is because, once nucleated, a bubble of true vacuum expands at the speed of light and destroys everything in its way. If even a single bubble had nucleated at any time, anywhere in our past lightcone, it would have already hit us.
To describe the history of the Universe, we approximate it with the FLRW metric (3.3). The scale factor a(t) satisfies the Friedmann equation (3.5) where ρ is the energy density of the Universe. When the dominating energy forms can be described by ideal fluids, one can write an equation of state p = wρ, which relates the pressure p to the energy density ρ through the equation of state parameter w. From the first law of thermodynamics it then follows that the energy density scales with the expansion of the Universe as ρ ∝ a −3(1+w) .

(5.2)
Observations indicate the the Universe currently contains three forms of energy: radiation (w = 1/3), matter (w = 0) and dark energy, which we assume to be a cosmological constant with w = −1. The total energy density can be therefore written as a function of the scale factor as where Ω Λ = 0.69, Ω mat = 0.31 and Ω rad = 5.4 × 10 −5 are the observed energy fractions of cosmological constant, matter and radiation, respectively [11], and ρ 0 tot is the current total energy density. The Universe is therefore currently dominated by dark energy, but in past it was dominated by matter and, at even earlier times, by radiation. Observations also show that in its very early stages, before radiation-dominated epoch, the Universe went through a period of accelerating expansion known as inflation, during which the equation of state was, again, w ≈ −1.
To find the expected number of bubbles in the past lightcone, it is convenient to write the FLRW metric in terms of the conformal time η as in Eq. (3.10). In these coordinates, light satisfies |d r/dη| = 1, so if we denote the current conformal time by η 0 , the comoving radius of our past light cone at conformal time η is r(η) = η 0 − η.
The dependence of the scale factor on the conformal time is determined by the Friedmann equation (3.5), which in terms of the conformal time is The bubble nucleation rate Γ may have been very different in different stages of the early evolution of the Universe. It depends on the curvature of spacetime and temperature, and also potentially on any perturbations or non-equilibrium processes that could catalyse or trigger the decay process and therefore it is function of the scale factor, Γ = Γ(a). This allows us to write an expression for the expected number of bubbles N in our past lightcone (after some initial time η ini ) as If this number is much greater than one, it would be unlikely that our part of the Universe could have survived until today, and therefore our existence requires N 1. (5.7)

Late Universe
Let us first consider the post-inflationary Universe described by the energy density (5.3) and assume that the bubble nucleation rate Γ(a) in the past was at least as high as its current Minkowski space value Γ 0 , i.e., Γ(a) ≥ Γ 0 . In this case the expected number of bubbles is Hence, the constraint on the nucleation rate Γ 0 from the post-inflationary era is on the bounce action. By calculating the nucleation rate Γ 0 , theories can be divided into categories: stable, metastable and unstable. If the rate exceeds the bound (5.9), the Universe would not have survived until the present day, and hence the vacuum is said to be unstable. If the rate is non-zero but satisfies Eq. (5.9), the vacuum would not have decayed by the present time but would decay in the future, and hence it is said to be metastable. Finally, if the decay rate is strictly zero, which is the case when the current vacuum state is the global minimum of the potential, then the vacuum is said to be stable.   [11]. In the green region, the current vacuum is absolutely stable, in the yellow region it satisfies the bound (5.9), and in the red region it is so unstable that it would not have survived until the present day. The instability boundary includes gravitational backreaction [166] and is shown for ξ = 0 and ξ = ±1000. of the non-minimal curvature coupling. The blue dashed line shows the instability bound (5.62) obtained by taking the thermal history of the Universe into account [177] and assuming reheat temperature T RH = 10 16 GeV. Fig. 7 shows the stability diagram of the Standard Model based on Ref. [166] (see Section 4.2 for discussion), in terms of the Higgs mass M h , top mass M t , for three different values of the non-minimal coupling ξ. The ellipses show the 68%, 95%, and 99% contours based on the experimental and theoretical uncertainties in the masses.
It is worth mentioning that one could invoke the anthropic principle to evade the bound (5.9). Even if the expected number of bubbles N is large, there is always a nonzero probability that no bubbles were nucleated. Life can obviously only exist in those parts of the Universe that have no bubble nucleation event in their past light cone, and therefore that is necessarily what we observe, no matter how low the probability is a priori. One can therefore argue that observations do not require N 1. However, the anthropic argument does not rule out bubbles hitting us in the future, and therefore, if the Universe survives for a further period of time, that imposes a bound that is not subject to the anthropic principle. For this, the quantity that matters is the time derivative of the expected number of bubbles, This imposes constraints that are numerically weaker but cannot be avoided by anthropic reasoning. To be concrete, one can carry out an experiment by waiting for a period of time t exp , for example one year. If, at the end of the time period, the experimenter has not been hit by a bubble wall, this gives a constraint For the post-inflationary Universe this is and for t exp = 1yr, one obtains the bound Γ 0 2.9 × 10 10 H 4 0 , or B 520. (5.14) This is weaker than Eq. (5.9), but because of the very strong dependence of Γ 0 on the top and Higgs masses, it does not change the stability constraints on them significantly.

Inflation
Although most of the spacetime volume of our past lightcone comes from the late times, the vacuum decay rate Γ(a) was much higher in the very early Universe. Depending on the cosmological scenario, it can be high enough to violate the bound (5.7), and this can be used to constrain theories. The earliest stage in the evolution of the Universe that we have evidence for is inflation, a period of accelerating expansion, which made the Universe spatially flat, homogeneous and isotropic and also generated the initial seeds for structure formation. In simplest models of inflation, the energy density driving it is in the form of the potential energy V (φ) of a scalar field φ known as the inflaton. The inflaton field is nearly homogeneous, and satisfies the equation of motion φ + 3Hφ + V (φ) = 0. (5.15) During inflation the potential satisfies the slow-roll conditions, These conditions guarantee the existence of a solution in which the first term in Eq. (5.15) is subdominant, and the inflaton field rolls slowly down the potential V (φ). As a consequence, the energy density ρ ≈ V (φ) and the Hubble rate are approximately constant. The Hubble rate during inflation, H inf , is largely unknown. Observationally it is constrained from above by the limits on primordial B-mode polarisation in the cosmic microwave background radiation. This gives an upper bound r < 0.09 on the tensor-toscalar ratio [200], which implies H inf 3.3 × 10 −5 M P ≈ 8.0 × 10 13 GeV at the time when the observable scales left the horizon. In a realistic inflationary model, the Hubble rate decreases with time, and would therefore be lower at the end of inflation. Although there are models in which the Hubble rate is well below the tensor bound, it is generally expected to be close to it, and in the simplest single-field inflation models it even exceeds it. It is therefore considered to be likely that the Hubble rate was significantly higher than the Higgs mass m H ≈ 125 GeV.
The minimal inflationary model is Higgs inflation [201], in which the non-minimal curvature coupling of the Higgs field is large, ξ ∼ −49000 √ λ. This allows it to play the role of the inflaton, without the need for a separate inflaton field. During inflation, the Higgs field has a large value ϕ ∼ M P /|ξ|, which means that the existence of a negativeenergy minimum would appear to pose a problem for the scenario, because if the Higgs field gets trapped there, it would lead to a rapid collapse of the Universe instead of inflation. However, inclusion of higher-dimensional operators and finite temperature effects can avoid this problem [202]. Of course, if the actual top and Higgs masses lie in the stable region (see Fig. 7), no problem arises. Furthermore, if they are just below the stability boundary, the effective Higgs potential would have an inflection point which would allow the scenario known as critical Higgs inflation [203], in which the Higgs field values are significantly lower than in conventional Higgs inflation. In the following our focus will be on the conventional scenario in which the inflaton is a separate field, and therefore we will not discuss Higgs inflation in detail. A thorough and up-to-date review of Higgs inflation, covering also the vacuum stability issues, is given in Ref. [204].
Even in the scenario in which the inflaton is not the Standard Model Higgs field, one could expect on general grounds that the natural initial value for the Higgs field is at the Planck scale ϕ ∼ M P [205]. In that case the existence of a negative-energy true vacuum between the electroweak and Planck scales would appear to be a problem, just like in Higgs inflation. Therefore one either has to assume special initial conditions that guarantee ϕ ϕ bar everywhere, or find a mechanism that allows the Higgs field to roll to small values without gettting trapped in the negative energy true vacuum.
In addition, even if that problem is solved, one still needs to avoid the nucleation bubbles of true vacuum, and hence satisfy the bound (5.7). Approximating inflation with a de Sitter space with constant Hubble rate H inf , the expected number of bubbles (5.6) in our past lightcone originating from inflation is where Γ inf is the vacuum decay rate, and V inf is the volume of the inflationary part of our past light cone. One can write this as 18) where a inf is the scale factor at the end of inflation, H inf is the Hubble rate during inflation, and N tot is the total number of e-foldings of inflation. In principle, if inflation lasted for an infinite amount of time, the volume of the inflationary past light cone would be infinite. In practice, inflation has a finite duration in most models, and the first term usually dominates in Eq. (5.18).
The factor (a inf H inf /a 0 H 0 ) is the ratio of the comoving Hubble lengths today and at the end of inflation. It can be expressed as 19) where N is the number of e-foldings from the moment the largest observable scales left the horizon during inflation, to the end of inflation. It depends somewhat on the cosmological history, but is approximately [206] N ≈ 60 + ln V 1/4 inf 10 16 GeV . (5.20) This means that the spacetime volume of the inflationary past light cone is From Eq. (5.7), one then obtains a bound on the decay rate during inflation In the literature, the vacuum stability during inflation is often discussed in terms of the survival probability P survival , which can be defined either as the fraction of volume that remains in the metastable vacuum at the end of inflation, or as the probability that a given Hubble volume remains in the metastable vacuum until the end of inflation. This is related to N by N ≈ e N (1 − P survival ), (5.23) and therefore the bound (5.9) can be written as One can use the bounds (5.22) or (5.24) to constrain the Hubble rate during inflation H inf and other parameters of the theory. This computation can be done in two ways, either using the instanton calculation of the tunneling rate discussed in Section 4, or using the stochastic Starobinsky-Yokoyama approach discussed in Section 3.4. The instanton calculation includes both quantum tunneling and classical excitation, and it can incorporate interactions and gravitational backreaction at short distances. Because it requires analytic continuation, it only works with constant Hubble rate H inf , but it can still be expected to be a good approximation when the Hubble rate is slowly varying. In contrast, the stochastic approach can describe a time-dependent Hubble rate and gives a more detailed picture of the time evolution, but it includes only the classical excitation process and does not include interactions on sub-Hubble scales.
In the stochastic approach, the dynamics is described by either the Langevin equation (3.14), or by the Fokker-Planck equation (3.16), which gives the time evolution of the one-point probability distribution P (t, ϕ) of the Higgs field ϕ.
If the Higgs field is assumed to vanish initially, ϕ = 0, the probability distribution grows initially as This is obtained by ignoring the Higgs potential V (ϕ), which should be a good approximation at early times. After some time the potential becomes important and starts to limit this growth. If the Hubble rate H is constant, the field approaches asymptotically the equilibrium distribution (3.17), and it is also a good approximation if the Hubble rate is varying sufficiently slowly. Considering the tree-level potential V (ϕ) = λϕ 4 /4 with constant λ > 0, the typical (rms) value of the field is given by Eq. (3.19) as (5.26) where the last expression is for the experimental value of the Higgs self coupling λ ≈ 0.13. If H 10 10 GeV, these field values are beyond the position (2.32) of the maximum of the potential. This means that for such values of the Hubble rate, inflationary fluctuations of the Higgs field would be able to throw the Higgs field over the potential barrier, triggering the vacuum instability [77,122,178,181,195,196,205,[207][208][209][210][211][212][213]. This would place a rough upper bound on the Hubble rate, H 10 10 GeV. (5.27) To make the bound more precise, Espinosa et al. [178] solved the equation for the initial state P (0, ϕ) = δ(ϕ), with the boundary condition P (ϕ bar , t) = 0 to account for the destruction of any Hubble volume where ϕ > ϕ bar . They then defined the survival probability of the vacuum as P survival (t) = ϕ bar −ϕ bar P (h, t). (5.28) Because of the boundary conditions, the survival probability is not conserved but decreases with time, and from the late-time asymptotic decay, The numerical value of this constraint depends on the number of e-foldings N and, in particular, the height of the potential barrier, which is highly dependent on the precise Higgs and top masses. The bound on the ratio H/ϕ bar is much less sensitive to the mass values, and therefore also quote the bounds in units of ϕ bar rather than GeV. To obtain indicative bounds in physical units, one can use the central estimate for ϕ bar in Eq. (2.32). Using N = 60, the bound (5.32) becomes H 0.067ϕ bar . (5.33) The same result be also obtained using the instanton approach [207], which gives the decay rate (4.46),  6 shows that for Hubble rates near ϕ bar , the relevant instanton solution is the Hawking-Moss instanton, whose action (4.45) agrees with the exponent in Eq. (5.31) in the limit where the barrier height is much less than the false vacuum energy. The instanton and Fokker-Planck calculations are therefore in good agreement in this case.
As discussed in Section 4.6, the relevant instanton for lower Hubble rates, H < H cross , is the Coleman-de Luccia solution [182]. However, this is below the bound (5.32) and the Coleman-de Luccia action is very high, B ∼ 1800, so that it gives a negligible decay rate, and therefore this does not change the bound (5.32).
There has been some debate about the correct field value used for the boundary condition (5.28) in the Fokker-Planck calculation. Hook et al. [195] applied the boundary condition P (t, ϕ cl ) = 0 at ϕ = ϕ cl , determined from the condition This condition means that at h > h cl the classical motion of the field due to the potential gradient dominates over the quantum noise. Therefore it allows field trajectories that cross the top of the barrier but return to the metastable side because of the quantum fluctuations. This leads to a slower decay rate in the case of the high Hubble rate, East et al. [194] considered the cutoff point the value ϕ sr , where This is the value above which the Higgs field no longer satisfies the slow roll condition and therefore the stochastic approach fails. The choice of the boundary condition becomes less important when H ϕ bar , and therefore it does not affect the bound Eq. (5.32) very much. By solving the Fokker-Planck equation numerically, the authors obtained the bound H 0.067ϕ bar , (5.39) which coincides numerically with Eq. (5.32).
There are aspects of physics that are not included in the approximations leading to the bound (5.32), and which can therefore provide a way to evade the bound. First, the high spacetime curvature R = 12H 2 during inflation modifies the effective potential both at the tree level through the non-minimal coupling ξ and through the curvature-dependence of the loop corrections. The non-minimal coupling gives rise to an effective curvature-dependent mass term (3.12), If ξ is positive, it increases the potential height between the electroweak and true vacua and helps to stabilise the electroweak vacuum even if the Hubble rate is well above the bound (5.32) [178,208,212]. On the other hand, negative values of ξ make the vacuum less stable. For ξ < 0, Joti et al. [192] obtained the bound The stabilising effects of the non-minimal coupling have also been discussed in Refs. [20,181,193,[213][214][215][216][217] The curvature dependent loop corrections mean that the non-minimal coupling ξ runs with the renormalisation scale, and even it is zero at low energies, it runs to a negative value ξ ≈ −0.03 at the relevant scales for the instability µ Λ ∼ 10 10 GeV [212]. Curvature contributions to the loop corrections to the rest of the effective potential can be approximated using renormalisation group improvement [212], by choosing the renormalisation scale as µ * ≈ H when H ϕ, rather than µ * ≈ ϕ which had been used previously. Using the curvature-dependent renormalisation scale, such as Eq. (3.28), has become the norm in the more recent literature [194,196,218]. Having µ * ∼ H means that for sufficiently high Hubble rates the effective coupling becomes negative, λ(µ * ) < 0, and the potential barrier disappears completely, unless ξ is sufficiently large. Both of these effects, running ξ and the curvature-dependent renormalisation scale, tend to de-stabilise the vacuum. Taking them into account gives the bound [212] ξ 0.06 for H ϕ bar .
The full curvature-dependent effective potential was computed at one-loop order in Ref. [20], and confirms this expectation. The stability bounds as a function of the Hubble rate H and the non-minimal coupling ξ are shown in Fig. 8. For comparison, the bound from particle collider experiments is |ξ| 2.6 × 10 15 [219]. A sufficiently large positive non-minimal coupling ξ can also avoid the Higgs field initial condition problem. It was found in Ref. [217] that if ξ H/10 −4 M P , (5.43) the positive curvature contribution to the effective potential allows the Higgs field to roll from Planck-scale values to its electroweak minimum during inflation without getting trapped into the negative-energy true vacuum. The bound (5.32) also does not take into account any direct coupling between the Higgs and the inflaton field φ. Although a direct coupling is not radiatively generated, in general it is possible and the precise form it would have and its effects on vacuum stability depend on the detailed of the inflaton sector. The simplest example is a coupling of the form λ φh φ 2 h 2 in chaotic inflation with a quadratic potential. During inflation, the inflaton field has a high value φ M P , and therefore the coupling produces an effective mass term for the Higgs field, Coupling values λ φh 10 −6 would not spoil the flatness of the inflaton potential [205,220], and if λ φh 10 −10 , it would stabilise the vacuum during inflation and allow the Higgs field to roll to its current small field values even if starts from a Planck-scale value at the beginning of inflation [205,209,220]. This coupling has also been discussed in Refs. [213].
Considering the non-minimal curvature coupling ξ and the direct Higgs-inflaton coupling λ φh together, Ref. [221] finds the constraint 10 −10 λ φh + 10 −10 ξ 10 −6 , (5. 45) in the quadratic chaotic inflation model. Other forms of the Higgs-inflaton coupling have been considered in Refs. [195,211,222,223]. There are also other effects that could potentially stabilise the vacuum state during inflation. Non-zero temperature T 6 × 10 13 GeV during inflation [209], moduli fields [224], coupling to a spectator scalar field [225], or top quark production [218] could all generate an effective stabilising term in the effective potential.

Reheating
The end of inflation can be defined as the point at which the Universe no longer undergoes accelerated expansion, which occurs when w = −1/3. This marks the beginning of the so-called reheating phase during which the energy density stored as potential energy gets converted into the hot thermal plasma of the Big Bang. If the acceleration is sourced by a slowly rolling inflaton φ, during reheating the slow-roll conditions seize to hold and the inflaton will begin a phase where its (average) kinetic energy is comparable to its potential energy. This usually manifests as coherent oscillations around the minimum of the potential. Reheating is said to be completed when the energy density of the hot Big Bang overtakes that of the inflaton sector, which often proceeds via direct couplings allowing the inflaton to decay into SM constituents. It is however worth pointing out that it is perfectly possible to have successful reheating without any couplings between the inflaton and the SM sector, for examples of such models see Refs. [226][227][228][229].
An inflaton field coherently oscillating around the minimum of its potential may source a very potent non-perturbative amplification of quantum modes, which takes place during the early stages of reheating and is hence often referred to as preheating [230,231]. If a phase of preheating occurs, it does not lead to the completion of reheating as the created particles tend to shut off any non-perturbative behaviour through backreaction and a perturbative decay channel is often required to ensure the complete decay of the inflaton.
From the point of view of a possible vacuum destabilization, preheating is a crucial epoch because vacuum decay is potentially induced by a large amplification of the Higgs field [232]. It is important to note that at the time of preheating, the Universe has not yet reheated to a high temperature, and therefore the thermal effects discussed in Section 4.3 cannot stabilise the vacuum state.
Let us proceed to consider the familiar Lagrangian appropriate for the Higgs doublet in curved space (3.30). We consider Hubble rates well above the electroweak scale, H M h , and therefore we can we neglect the tree-level mass parameter, and use the action We also assume a single-field model of inflation with a canonical kinetic term and the potential U (φ). The inflaton φ is taken to dominate the energy density of the Universe completely and because of this the Higgs field may be considered as a subdominant spectator that can be neglected in the Einstein equation. Using then in the Friedmann equations (3.5), we can solve for the Ricci scalar R After inflation ends, the inflaton field φ rolls down its potential, and initially oscillates coherently about its minimum φ min , until it eventually decays. We assume that the inflaton potential vanishes at the minimum, U (φ min ), as is usually the case. We can see from Eq. (5.48) that during every oscillation, when φ ≈ φ min , the Ricci scalar becomes negative, R < 0. This, in turn, means that the non-minimal term ∼ ξRϕ 2 gives rise to a tachyonic mass term (3.12) for the Higgs field. As already discussed in Section 3.3, this gives rise to significant excitation of the field. The fact that the non-minimal term can lead to extremely efficient particle creation during preheating was first discussed in Refs. [233,234].
Particle creation from a periodically tachyonic effective mass was analyzed in detail in Ref. [235] where it was named tachyonic resonance. It is much more extreme than the resonant effects usually taking place during preheating. Hence a dangerous fluctuation of the Higgs field can be generated during a single oscillation of the inflaton.
For concreteness, we now focus on the case of a quadratic inflaton potential Although as a complete model of inflation, this not compatible with observations [236], it approximates the shape of the potential around the minimum in general single-field models. The behaviour of the inflaton field during its coherent oscillations can be approximately written as φ = φ 0 (t) cos(mt) where φ 0 is a slowly changing amplitude φ 0 (t) = √ 6H(t)M P /m [230,231]. We will focus only on a very brief time period immediately after inflation, when no thermalization has yet taken place. In cosmic time the properly normalized mode is obtained from Eq. (3.9) as f (η) → f (t)/ √ a giving the mode equation By using the Friedmann equations (3.5) in this approximation the mode equation can be cast in the Mathieu form [232,233]

51)
Making use of the analysis in Ref. [235] we can derive an analytical result for the occupation number of the Higgs field n k after the first oscillation where Ω 2 ≡ −ω 2 . The ω 2 is the term in the square brackets in Eq. (5.51) and ∆z covers the time period when ω 2 < 0. Including only the IR modes k < aH, neglecting the expansion of space, the self-interaction and furthermore assuming ξ 1 we can estimate the generated Higgs fluctuations at horizon scale, φ 2 aH , after the first oscillation of the inflaton as [232] . Note that a large and positive ξ gives rise to a destabilizing effect after inflation. This is opposite to what happens during inflation when it suppresses fluctuations by effectively making the field heavy (see Section 5.3).
In general, once a significant particle density is produced it tends to work against any further particle production [231]. For the Higgs the main backreaction comes from the selfinteraction term, which contributes to the effective mass (3.12), along with the curvature terms visible in (5.50), as very similarly as we derived in the 1-loop approximation for a scalar singlet in Eq. (3.24) of Section 3.5. In order for tachyonic particle creation to take place one must have ξ|R| 6λ φ 2 for ξ 1. However, in Section (3.6) it was shown that the Hubble rate contributes to the RG scale through curvature induced running (see Eq. (3.28)). If H ϕ bar , the fourpoint coupling is negative, implying that the backreaction in fact enhances the instability, and will not suppress tachyonic particle creation even if a large variance is generated. Backreaction also arises from the gravitational disturbance of the generated particle density. In order to reach this threshold one must create enough particles such that their energy density approaches 3H 2 M 2 P . The relevance of gravitation backreaction we can estimate from the approximate energy density for the Higgs [232] ρ Higgs ≈ 24ξH 2 φ 2 + 6λ φ 2 2 .
(5.56) Figure 9. Regions where the Higgs fluctuations ∆ϕ generated after inflation by a single inflaton oscillation are greater than the barrier height ϕ bar for a model with no direct couplings between the Higgs and the inflaton. The Higgs fluctuations are given by Eq. (5.53) while taking into account backreaction effects from self-interactions and gravity. The amplitude of the inflaton at the end of inflation is assumed to satisfy φ 0 = 0.3M P and we have used the value ϕ bar = 4.64 × 10 10 GeV from Section 2.3. Regions below ξ = 5 have been cut according to the bound obtained in Ref. [238].
When ρ Higgs ∼ 3H 2 M 2 P the Higgs starts to influence the dynamics of spacetime requiring a non-linear analysis. Below we will assume that when the gravitational backreaction threshold is reached the particle production will seize.
More detailed calculations of the process have been carried out using linearised approximation [215] and lattice field theory simulations [237,238]. The most detailed analysis, carried out in Ref. [238], used the tree-level RGI effective potential with three-loop running, and considered different top quark masses. The main conclusion was that the instability is triggered with high probability for ξ 4 − 5, for a top quark mass M t ≈ 173.3GeV. This implies an upper bound on ξ in the context of quadratic chaotic inflation.
The regions where a dangerously large fluctuation of the Higgs is generated after a single oscillation of the inflaton are shown in Fig. 9. The rather complicated shapes are the result of the interplay of the variance (5.53) and the constraints coming from self interactions and gravitational backreaction. In Fig. 9 we have assumed that the amplitude of the inflaton at the end of inflation satisfies φ 0 = 0.3M P . While this is true for the quadratic (chaotic) model of inflation, it is not true generically. Since Eq. (5.53) is exponentially dependent on φ 0 the predictions are very sensitive to the specifics of inflation. Similarly, the duration of reheating plays a crucial role and for prolonged reheating a possible instability may be further enhanced. The derivation of Eq. (5.53) is based on the adiabatic approximation [235], which can be shown to break down for small ξ [239]. Furthermore, very little particle creation is expected when close to the (approximately) conformally invariant point ξ = 1/6. For these reasons and the lattice results of Ref. [238] we have conservatively cut out regions with ξ ≤ 5.
As discussed in the previous section, the set-up in Eq. (5.46) assuming that the inflaton is decoupled from the SM is in many ways the minimal one. Couplings between the inflaton and the SM sector may of course be introduced or even required by a specific reheating model. Vacuum stability during preheating in models with no non-minimal coupling but with direct couplings between the inflaton and the Higgs was investigated in Refs. [220,240,241]. In particular, in Refs. [240] it was shown that in some cases vacuum decay during preheating may take place also for low-scale inflation.
In Refs. [215,221,237,242] both the non-minimal coupling and direct Higgs-inflaton couplings were considered. In a sense in this case the Higgs fluctuations are sourced in a complicated manner by the interplay of tachyonic resonance [235] and the (usual) parametric resonance [231]. For the precise coupling ranges where significant particle production takes place and possible implications for instability, see Ref. [221]. We also point out that particle creation resulting from the non-adiabatic change in the background curvature when inflation ends, already shown in Ref. [243], can be enough to probe the unstable region of the effective Higgs potential [232].

Hot Big Bang
After reheating, the Universe entered a thermal radiation-dominated state, in which vacuum decay rate can be approximated by the thermal rate (4.35) at the relevant temperature, and the Hubble rate was given by the equation where g * (T ) is the effective number of degrees of freedom and has the value g * (T ) = 106.75 in the Standard Model at high temperatures. Using Eq. (5.6) one can write the expected number of true vacuum bubbles in our past light cone from this era as [85,178] d N d ln T = 4π 3 In practice, reheating is not instantaneous, and there may have been a period when the Standard Model degrees of freedom were in thermal equilibrium but were not the dominant energy component. In the scenario in which the inflaton field decays slowly and dominates the energy density of the Universe for an extented period, the maximum temperature is [10,177,178]  However, this has only a small effect on the numerical bounds [177]. In Fig. 7, the blue dashed line shows the bound (5.62) calculated with a reheat temperature T RH ∼ 10 16 GeV. As can be seen, the inclusion of the thermal history of the Universe reduces the allowed mass range compared with the zero-temperature bounds. The central experimental values are still allowed, but the instability boundary lies within two standard deviations from them.

Concluding Remarks
The current experimental data show that with a high likelihood, the electroweak vacuum state of the Standard Model is metastable. Even though the vacuum state could be stabilised by new physics beyond the Standard Model, and even in the Standard Model parameters corresponding to a stable vacuum are still allowed by experimental errors, it is important to study the implications of the possible metastability. That allows one to understand whether the metastability is compatible with observations, and if so, what constraints it places on the parameters of the theory.
If the electroweak vacuum really is metastable, then bubbles of the true, negativeenergy vacuum can be nucleated by quantum tunneling or classical excitation, as discussed in Section 4. Once a bubble has formed, it expands at the speed of light, destroying everything in its way. This clearly has not happened yet in our part of the Universe, which means there has not been a single bubble nucleation event in our whole past light cone. In Section 5, we showed how the likelihood of this can be estimated by computing the nucleation rate and integrating it over the past light cone. Because the past light cone includes all of the different cosmological eras, and the nucleation rate and its dependence on theory parameters is different in each era, this provides a rich set of constraints on both the cosmological history and on the Standard Model parameters.
In this review, we have focussed on four different cosmological eras: inflation, preheating, hot radiation-dominated phase, and the late Universe. Vacuum stability in the late Universe one obtains constraints on the Higgs and top masses, and they are made tighter by considering the hot radiation-dominated phase, as summarised in Fig. 7. Survival of the vacuum through inflation and the subsequent preheating phase constrains the Hubble rate during inflation and the Higgs-curvature coupling ξ (Figs. 8 and 9), as well as other aspects of inflationary models. A demonstration of the power of these considerations is that for quadratic chaotic inflation, the non-minimal coupling is constrained to be within the range 0.06 ξ 5, which is 15 orders of magnitude stronger than the experimental bounds from the Large Hadron Collider [219]. Cosmological vacuum decay has a unique connection to gravity via the early Universe, which opens up an observational window to particle physics well beyond what colliders can achieve.
In this work we have reviewed the, already rather significant, body of work investigating the cosmological consequences of the SM Higgs possessing a metastable potential. We have also discussed the relevant theoretical frameworks required for such studies. The multidisciplinary nature of the problem is perhaps one of the reasons behind the ongoing significant interest as particle physics, quantum field theory and gravity all play a prominent role. Although the specifics of the theory behind early Universe dynamics are not currently known what has become quite apparent is that a metastable Higgs potential generically leads to non-trivial constraints, which are completely invisible to colliders.
On the other hand, despite the large number of existing studies, much remains to be explored. For example, at the moment very few works exist that go beyond the simple quadratic model of inflation. This is equally true for the inflationary and reheating epochs. There is also a great deal of scope for improving calculation techniques in order to obtain more precise and robust constraints, for example by going beyond the semiclassical approximation or fully including gravitational effects. The work on the cosmological aspects of Higgs vacuum metastability is only starting.