Critical Higgs Inflation and Second Order Gravitational Wave Signatures

The self coupling $\lambda$ of the Higgs boson in the Standard Model may show critical behavior, i.e. the Higgs potential may have a point at an energy scale $\sim 10^{17-18}$ GeV where both the first and second derivatives (almost) vanish. Since $\lambda$ is very small in this region, the Higgs boson can serve as inflaton even if its nonminimal coupling to the curvature scalar is only ${\cal O}(10)$, thereby alleviating concerns about the perturbative unitarity of the theory. We find that just before the Higgs as inflaton enters the flat region of the potential the usual slow--roll conditions are violated. This leads to ``overshooting'' behavior, which in turn strongly enhances scalar curvature perturbations because of the excitation of entropic (non--adiabatic) perturbations. For appropriate choice of the free parameters these large density perturbations occur at length scales relevant for the formation of primordial black holes. Even if these perturbations are not quite large enough to trigger copious black hole formation, they source second order tensor perturbations, i.e. primordial gravitational waves; the corresponding energy density can be detected by the proposed space-based gravitational wave detectors DECIGO and BBO.


Introduction
Inflation is a beautiful paradigm for the evolution of the very early universe: it not only solves the problems of standard cosmology [1,2], but also generates the initial seeds for the formation of large structures via quantum fluctuations. The simplest inflationary models feature a single scalar field that slowly "rolls down" a rather flat potential ("slow-roll" inflation). The energy density during inflation is then dominated by the potential, which leads to an approximately exponential expansion of the universe. Often a separate "inflaton" field is introduced for this purpose, but it would obviously be more economical to instead use the single scalar Higgs field Φ of the Standard Model (SM) of particle physics as inflaton. At TeV energies the Higgs self coupling λ is O(0.1); a coupling of this size leads to a rather steep potential, which needs to be "flattened" by a large non-minimal coupling to the Ricci scalar, ξΦ † ΦR [3]. This yields [3][4][5][6] n s 0.97, r 0.0034 in agreement with observation, but for a coupling ξ ∼ 10 4 . 1 Such a large coupling may violate perturbative unitarity [7][8][9]. This issue has been discussed at length in the literature [10][11][12][13][14][15][16]. Several authors have also attempted to resolve this problem by considering Higgs inflation with an additional field beyond the SM [17][18][19][20][21][22][23]; however, these models lack the simplicity of the original suggestion.
On the other hand, at the large field values where Higgs inflation may have occurred, the value of λ is expected to be quite different than at the weak scale. This difference is described by renormalization group equations (RGE). At the one-loop order λ is driven to larger values by Higgs self-interactions (i.e. the λ 2 term contributes with positive sign in the RGE) and by electroweak gauge interactions, but is reduced by Yukawa interactions, the by far most important one being that of the top. The evolution of the top Yukawa coupling in turn is also affected by QCD interactions. For the measured values of the mass of the physical Higgs boson (which determines λ at the electroweak energy scale), the top mass (which determines the top Yukawa coupling at the electroweak scale), and of the electroweak and strong gauge couplings, the two-loop RGE indicate that λ shows critical behavior at energy scale ∼ 10 17−18 GeV, i.e. λ and its first derivative, described by its beta function, can both be very small [25]. This implies that the potential becomes very flat in this region [26]. This can give rise to "critical Higgs inflation" (CHI). 2 Since λ is small, one only needs a nonminimal coupling ξ ∼ O(10) [27][28][29]; see also [30][31][32][33] for recent investigations concerning CHI and [34] for a comprehensive review of Higgs inflation.
In addition to reproducing the measured CMB power spectrum accurately, recently some other cosmological implications of the CHI scenario have been investigated. In particular ref. [35] showed that curvature perturbations are greatly enhanced at the length scales that leave the horizon when the inflaton field enters the very flat region of the potential; this might even lead to the formation of a cosmologically significant abundance of primordial black holes (PBH). 3 In fact, in the simplest approximation the strength of the density perturbations is inversely proportional to the first derivative V of the inflaton potential V . It is thus tempting to associate a spike in the spectrum of density perturbations with an "ultraslow roll" (USR) phase while the inflaton traverses a very flat piece of the potential. We will see in Section 2 that this is not the whole story: the largest enhancement actually does not happen in the USR phase, but during a transitionary "overshooting" stage just before USR where the inflaton potential has a sizable curvature, so that the slow-roll conditions are violated.
The enhanced scalar curvature perturbations are expected to source tensor perturbations at second order, as investigated in recent papers [37][38][39][40]. In this paper, we extend the analysis in [35] and investigate the second order gravitational wave (GW) signatures arising from large scalar curvature perturbations in the CHI scenario. We show that such signatures can be detected by several proposed space based GW experiments. The calculation of the PBH density is fraught with considerable uncertainty [41,42]. Our result indicates that, at least for CHI inflation, an inflationary GW signal should be detectable in all cases that could conceivably lead to sizable PBH production, i.e. a failure to detect the latter in future experiments would exclude the possibility that PBHs contribute significantly to cosmological dark matter.
The reminder of this paper is organized as follows. We first revisit the curvature perturbations under adiabatic and non-adiabatic conditions in Section 2; in particular, we show that for non-adiabatic conditions the amplitude density perturbations does not necessarily remain constant after horizon crossing, in contrast to the usual slow-roll (SR) treatment. In Section 3 we set up the CHI scenario. Using the general Mukhanov-Sasaki formalism to compute the power spectrum, we show that the standard SR approximation to calculate the power spectrum fails when the inflaton enters an overshooting phase, even if we use "Hubble" rather than "potential" SR parameters. In Section 4 we discuss second order GW signatures induced by the scalar density perturbations. Finally we summarize our results in Section 5. In appendix A we give a quick review of the Mukhanov-Sasaki equation and its analytical solution for a (quasi) de Sitter spacetime.

Adiabatic and Entropic Perturbation
In SR single field inflation the quantum fluctuations are adiabatic [56]. As a result the perturbations of all inflaton field dependent quantities X i share the same phase trajectory [57]: where i, j donate different observables,X i represents the average (background) of X i , δX i ≡ X i −X i its perturbation, andẊ the time derivative. In particular, using X 1 = p (pressure) and X 2 = ρ (energy density), adiabatic perturbation satisfy Energy density and pressure are defined via the energy-momentum tensor T ν µ , with ρ = T 0 0 and p = − 3 i=1 T i i /3. However in some cases the perturbation may not be adiabatic, for example when there are multiple fields interacting with the inflaton [56] or while the universe undergoes a non-SR inflationary phase, see subsection 2.2. Thus more generally the pressure perturbation can be decomposed into an adiabatic part and an entropic (i.e. non-adiabatic) one: i.e. δp en ≡ δp−δp ad . We will show in the next subsection that this distinction has significant impact on the evolution of curvature perturbations.

Evolution of Curvature Perturbations in SR, USR and Overshooting Phases
In order to relate inflationary predictions and cosmic microwave background (CMB) measurements the gauge invariant scalar quantity called curvature perturbation is usually introduced; it is defined by [58,59] Here H denotes the Hubble parameter, and Ψ(t, x) is a scalar function of coordinates. Physically ζ represents the spatial curvature of hypersurfaces with uniform energy density [60].
Since the power of the two-point correlation function of ζ is related to the CMB temperature anisotropies, one has to compute the power spectrum of ζ for a given inflationary model. This is usually done in Fourier space, where ζ(t, k) is associated to perturbations at a comoving length scale 1/k with k = |k|. As we will review below, under SR conditions the power spectrum can be computed when some mode k crosses the horizon, since ζ is frozen at super-horizon scale, i.e. it remains constant once k aH where a is the (dimensionless) scale factor in the Friedman-Robertson-Walker (FRW) metric. However whenever the universe deviates from SR expansion, we must in general take the super-horizon evolution of ζ into account; unfortunately this considerably complicates the accurate computation of the power spectrum. Using energy conservation it can be shown that the evolution of ζ is given by [60] Here δp en is the non-adiabatic component of the pressure perturbation, and Π is defined as where Ψ B is a Bardeen potential [58] which does not depend on k. We thus see that at super horizon scales, i.e. for k aH, the second term in eq.(2.5) can be neglected. If in addition the perturbations are adiabatic, i.e. if δp en can be neglected, then ζ is conserved on super-horizon scales. We will see that this is true if the universe follows a SR expansion, i.e. for a quasi-de Sitter spacetime.
Since we define δX = X(φ) − X(φ), δX can include terms that are of higher order in the field perturbation δφ. We assume that only a single (inflaton) field φ has sizable perturbations. Moreover, we are interested only in super-horizon modes where all gradient terms can be neglected. The energy density and pressure are thus given by The same equations also describeρ andp, which φ →φ on the right-hand side. Superhorizon size non-adiabatic pressure perturbations are then given by: (2.8) Using Taylor expansion for the potential up to second order 4 of δφ , we obtain (2.9) For a strict de Sitter spacetime V has to be constant, i.e. all derivatives of V vanish. It is easy to see that δp en = 0 in this case. However, during realistic SR inflation, |φ| V φ , so that eq.(2.9) reduces to (2.10) In order to see that δp en is indeed very small during SR, consider the equation of motion for δφ: [57,61] δφ In momentum space this becomes In analyses of inflationary dynamics it is often useful to trade the time for the number of e-folds N via dN = Hdt; eq.(2.12) then becomes (2.13) At super-horizon scales (k aH) the third term in eq.(2.13) can be neglected. Moreover, during SR the total energy density is dominated by the potential energy, so that 5 H 2 1 3 V . Finally, we introduce the second potential SR parameter η V ≡ V V , which has to be small during SR inflation. Eq.(2.13) can then be written as: For constant η V with |η V | 1 the solution of eq.(2.14) is given by where the constants C 1,2 are of order H/(2π), which determines the size of |δφ| due to quantum fluctuations during SR inflation. This 5 We set the reduced Planck scale MP 2.4 · 10 18 GeV to 1 in the following.
denotes the first potential SR parameter. Thus we seeφ is also rather small and nearly constant Inserting these estimates in eq.(2.10) and using eqs.(2.5) and (2.7) we find for the time evolution of the curvature perturbations at super-horizon scales: (2.16) Note that H V ∼ P ζ ∼ |ζ|, where P ζ is the power in the perturbation. At the length scales probed by the CMB, P ζ ∼ 10 −9 is very small. More importantly, the last line in eq.(2.16) shows that the time derivative of ζ is suppressed by the SR parameter η V once the exponentially decaying part of δφ k can be ignored. This completes our argument that in the SR regime ζ is (nearly) constant once the mode crosses out of the horizon.
However whenever the universe deviates from SR expansion,ζ may no longer be negligible even in the super-horizon regime. Of special interest to us is the situation where the acceleration term is much larger than the derivative of the potential, i.e. |φ| |V φ |. In this case eq.(2.9) becomes (2.17) In the last step we have neglected V φ relative to V φ δφ.
Our assumption |φ| |V φ | is equivalent to having the second (Hubble) SR parameter η H ≡ −φ/(Hφ) ≈ 3, which is manifestly not smaller than unity, i.e. the SR conditions are violated. Some recent papers state that this scenario corresponds to an USR phase. We disagree with this interpretation. The expression "ultra-slow" roll implies that the inflaton field evolves even more slowly than during SR, which happens when the potential becomes very flat, in which case the SR parameters should also be small. In other words, the spacetime during USR should be even more de Sitter like than that during SR, so the perturbations should be even more adiabatic than during SR, and the evolution of ζ at super-horizon scales should be even more suppressed. Hence the SR approximation for the power spectrum, where it is computed at horizon crossing, should work even better in a true USR phase, rather than breaking down.
So the phase with |φ| |V φ | cannot correspond to USR, but to an intermediate transition "overshooting" stage (also mentioned in [51]) between SR to USR, where the curvature of the potential is sizable but the first derivative is already rather small. We will see below that critical Higgs inflation can indeed lead to a situation where η H 3 for several e-folds of inflation.
In order to get a first qualitative understanding of such an "overshooting" stage, we insert the final result of eq.(2.17) into eq.(2.5): (2.18) Note thatφ + 3Hφ = 0 impliesφ ∝ e −3N so that 1/φ 2 ∝ e 6N . As a result the derivative dζ/dN grows exponentially during this overshooting region. Since neither δφ nor the curvature V are (approximately) constant during this overshooting epoch, eq.(2.18) is not so well suited for a quantitative treatment of the evolution of the density perturbations; this can be done using the Mukhanov-Sasaki equation, as will be described in the next section. However, we can already conclude that |φ| |V φ | implies that the curvature perturbation is not frozen at the super-horizon scales, and even increases significantly if the potential has a large positive curvature V . Of course, the usual SR treatment of approximating the final power spectrum by its value at horizon crossing will then no longer work. We consider eq.(2.18) and its consequences to be one of the central results of this paper, which is applicable whenever an overshooting epoch occurs during the evolution of the inflaton field. In the next chapter we will explore the quantitative consequences for the case of critical Higgs inflation.
Before concluding this section we briefly discuss the evolution of the perturbations after inflation ends. During matter domination the pressure is by definition negligible. During radiation domination, p ρ/3 holds locally, which again implies δp en = 0. Hence density perturbations remain frozen on super-horizon scales after inflation.

Critical Higgs Inflation
In this section we discuss critical Higgs inflation, with emphasis on the enhancement of density perturbations associated with an overshooting region. We first describe the basic setup in the Jordan and Einstein frames. In the second subsection we analyze the inflationary dynamics in the Einstein frame and connect it to CMB observables. In sec. 3.3 we show that SR conditions are violated in the overshooting region, just before the inflaton enters the very flat part of the potential. We then review the Mukhanov-Sasaki formalism, which we use in sec. 3.4 for a detailed numerical investigation.

Formalism
Starting point of the analysis is the action in the Jordan frame (in Planckian units, where In the second line we have introduced the function f (h) = 1 2 (1+ξh 2 ). The crucial observation [35] is that for realistic values of the relevant SM parameters, the running Higgs self coupling λ attains a minimum at scale µ. Near this minimum it can then be expanded as: The running non-minimal coupling ξ to the Ricci scalar is also expanded around scale µ: since ξ does not have an extremum at scale µ, the leading energy dependence is described by a term linear, rather than quadratic, in ln(h/µ).
While the matter part of the Jordan frame action has its canonical form, this is not true for the gravitational part, unless |ξ(h)h 2 | 1. In order to use standard results for the inflationary dynamics we transform to the Einstein frame, where gravity is described by the well-known Einstein-Hilbert action and the inflaton is described by a canonically normalized field χ. To that end we first utilize a conformal transformation to the Einstein frame:g Then we use a field redefinition to obtain a canonical kinetic term [62]; it is defined by: (3.5) After these transformations the action becomes While the gravitational part as well as the kinetic energy term in the action now have the standard form, the inflationary potential has become more complicated: It is convenient to introduce the quantities The inflaton potential can then be written as Note that for nonminimal coupling ξ = 0 the potential approaches a constant as x → ∞; it is this "flattening" which allows inflation. Consistency with the CMB observables (see below) and with current measurements of SM parameters can be obtained for parameter values in the ranges [35] . In order to compare our calculations, especially the power spectrum, with those in [35] based on the SR approximation, we mainly work with their representative set of parameters: The inflaton potential for these values of the parameters is shown in Fig. 1. It features an inflection point 6 at x = x c = 0.784. Again following ref. [35], we introduce one more free parameter β in order to study slight deviations from a perfect inflection point: where a(x c , c) and b(x c , c) are the values of the parameters which lead to V (x c ) = V (x c ) = 0. This is of interest since the inflaton field can linger near a true inflection point for a very large number of e-folds. This modification can give a slight slope to the ultra-flat region. Of course, the shape of the overshooting region, in particular V , will also be slightly modified: the larger the slope in the ultra-flat region is, the smaller V will be in the overshooting regime. We will use β in the range 10 −5 to 10 −4 .

Parameters of the CMB Power Spectrum
The inflaton dynamics in the Einstein frame is given by the Klein-Gordon equation in curved space:χ Using the relation between the number of e-folds and time, dN = Hdt, we can rewrite eq.(3.11) as [48,63] The two Hubble SR parameters are defined as and (3.14) SR inflation requires H 1 and |η H | 1. We have seen in sec. 2.2 that under the SR approximation, curvature perturbations are basically frozen at super-horizon scales. The power spectrum is therefore usually computed at horizon crossing, defined by k = aH, and can be analytically given by [64] (see appendix A for details): where N cross denotes the number of e-folds at horizon crossing. 7 The scale dependence of P ζ is usually parameterized as a power law, with spectral index n s given by The deviation from an exact power law is described by the "running" of the spectral index, parameterized through the quantity α: The final CMB observable of phenomenological interest is the tensor to scalar ratio r, i.e. the perturbations in tensor modes (which can be probed via the polarization of the CMB) normalized to the scalar perturbations. To leading order in SR parameters, Now our task is to solve eq.(3.12), from which the parameters of the CMB power spectrum can be computed. We find it more convenient to calculate the evolution of x, rather than the canonically normalized field χ; this is because we have an explicit expression for V (x), see eq.(3.8), and thus also for and V (x). By using eq.(3.5), eq.(3.12) and x = h/µ, we find the differential equation for x is: is too complicated to solve analytically. For a numerical solution we have to choose initial values for x and dx/dN at some N = N 0 . The initial value of x should evidently be above the field values where the CMB scales cross the horizon, so that our solution covers the entire range of scales probed by the CMB and other cosmological observations. On the other hand, it would be wasteful to choose x(N 0 ) to be much larger than the field values probed by the CMB, since this earlier evolution leaves no observable traces anyway. In practice we have used x(0) = 8.5. At these large field values the potential is very flat; if the initial kinetic energy of the field is not very large, the field evolution will therefore quickly approach the SR solution. 8 The initial choice of dx/dN is therefore largely irrelevant; we chose dx/dN = −0.21, which corresponds to assuming the SR solution already at N = N 0 = 0.
With these initial conditions, eq.(3.19) can be then solved numerically. Once x(N ) is known, the evolution of the canonical field χ can be obtained by integrating eq.(3.5): The constant of integration χ(0) can be fixed by using the fact that g(x) → 1 for x → 0; the natural choice is thus χ h for ξh 1, which corresponds to χ(0) = 6.94. Fig. 2 shows the evolution of x as well as the canonically normalized field χ with N for our standard set of parameters (3.9) with β = 10 −5 . We see that the field at first gradually decreases with increasing N . The evolution of χ becomes quite fast at N 30, signaling a break-down of SR. However, from N 35 both x and χ remain nearly constant for more than 30 e-folds; this is when the inflaton traverses the very flat part of the potential around the pseudo-critical point. Inflation ends when the inflaton leaves this region.
Once the dynamics of the inflaton field is known, the parameters of the CMB power spectrum can be computed. Using our standard parameter set (3.9) and β = 10 −5 we find that the CMB "pivot scale" k = 0.05 Mpc −1 crosses out of the horizon at N end −N CMB ≈ 68.
The numerical values of the CMB parameters at this scale are which is consistent with the Planck 2018 results [65] at the with 2σ level. The large value of N end − N CMB is to a large extent due to the USR phase. This number of e-folds can be reduced by increasing β, which increases the slope of the potential near the pseudo-critical point. For example, using β = 10 −4 , we find the same predictions as given by eq.

Slow-roll Violation
For our standard set of parameters, CMB scales first crossed out of the horizon during a SR phase, i.e. the SR approximation works very well for the predictions collected in eq.(3.21). However, Fig. 2 also shows that the canonically normalized inflaton field χ moves rather fast around N = 35. In this subsection we show that the SR conditions are indeed violated in this "overshooting" region. The dependence of the potential and its first and second derivatives, both with respect to x and with respect to χ, are plotted as function of N in Fig. 3a. The first derivatives remain positive and fairly small throughout. The second derivatives are initially small and negative, but increase in size as the inflaton field approaches the overshooting region, where the second derivative changes very rapidly from large negative to large positive values; in the region around the pseudo-critical point the second derivatives are again very small.
When discussing non-adiabatic pressure perturbations in Sec. 2, we had assumed that the second time derivative of the inflaton field is much larger in magnitude than the slope of the potential. Fig. 3b shows that this is indeed the case for some range of N around 35. In this case eq.(3.11) becomesχ In SR, both H and |η H | should be (much) smaller than 1; the result (3.23) clearly violates this. This is further illustrated in Fig. 4, which shows the evolution of the SR parameters with N for the same set of parameters. We show both the "Hubble" SR parameters defined in eqs. (3.13) and (3.14) and their "potential" analogues, defined via SR implies that H V and η H η V ; we see that in our case this holds for N < 30 as well as 9 for N > 38. H and V always remain significantly smaller than 1, but vary rapidly, and differ markedly, in the overshooting region.
We have shown in subsection 2.2 that entropic perturbation can be excited if the SR conditions are violated, in which case the curvature perturbation are no longer conserved at super horizon scale. Hence the usual estimate (3.15) of the power spectrum can no longer be justified for modes that first crossed out of the horizon near the overshooting region. In the next subsection we instead use the Mukhanov-Sasaki (MS) formalism to compute the power spectrum numerically.

Mukhanov-Sasaki Formalism
Our numerical treatment of the evolution of the curvature perturbations is based on the MS equation; a quick derivation of this equation and its analytical solution in the quasi-de Sitter limit is reviewed in appendix A. It is usually written in terms of the Mukhanov variable v k ≡ −zζ k , (3.25) where z is defined by In these variables, the MS equation reads: where τ denotes the conformal time, i.e. dτ = dt a . Rewriting eq.(3.27) using the number of e-folds N instead of the conformal time gives [48] Under SR conditions the curvature perturbation ζ k is frozen at super-horizon scales, hence the power spectrum of ζ k is usually computed at horizon crossing. However, we have seen in the previous subsection that the SR approximation fails in the overshooting regime. In order to account for the evolution of the curvature perturbation also at super-horizon scales the power spectrum should be computed at the end of inflation: This can usually only be done numerically. Recall also that super-horizon perturbations are frozen after inflation, as we showed at the very end of Sec. 2.
In order to solve eq.(3.28), initial conditions are needed. We follow the usual procedure, which assumes the Bunch-Davies vacuum at very early times [66]: Since v k is a complex quantity, 10 in practice it is more convenient to solve for its real and imaginary parts separately. To this end one can rewrite the initial condition eq.(3.30) as [48]: Here N i is the "initial" point where we start the numerical integration of the MS equation.
In principle the Bunch-Davis initial conditions (3.30) should be imposed at τ → −∞, which also corresponds to N → −∞ if the CMB pivot scale crossed the horizon at N 0, as we assumed in the last three figures. Physically this does not make much sense, since we don't know how many e-folds of inflation happened before that time. Moreover, the ansatz (3.30) remains a very good approximation of the exact solution of the MS equation as long as the perturbation is well within the horizon, i.e. for k 2 a 2 H 2 . Let the mode k cross the horizon at N = N k,cross . In practice it is then usually sufficient to use ∆N ≡ N k,cross − N i = 2 ∼ 3. We have checked that reducing N i , which costs a lot of CPU time since the term k 2 /(a 2 H 2 ) in eq.(3.28) grows ∝ e 2∆N requiring correspondingly reduced step sizes to attain numerical convergence, does not change the final result appreciably. However, we will see below that choosing too small a value for ∆N can lead to inaccuracies.

Power Spectrum for Critical Higgs Inflation
We now apply the MS formalism to CHI. In order to compute the power spectrum we have to integrate eq.(3.28) with the initial conditions eq.(3.30) till the end of inflation, and then plug the solution into eq.(3.29).
The result for our standard parameter set with β = 10 −5 is shown in Fig. 5. We see that the SR approximation fails badly for modes crossing the horizon in the vicinity of the overshooting region. In particular, the SR approximation gets both the location and the depth of the dip in the power spectrum wrong by more than one order of magnitude. The approximation (3.15) using H underestimates the maximum of P ζ by only a factor of about 1.4, but gets the location k max of the true maximum off by an order of magnitude, and underestimates the power at k max by about five orders of magnitude. Using the approximation (3.15) but replacing H by V , as is done in much of the older literature on inflation, actually gets k max approximately right, but overestimates the power at this scale by more than two orders of magnitude. In contrast, the SR approximation works well both for the large scales probed by the CMB and for the much smaller scales that cross out of the horizon in the USR regime after the end of the overshooting epoch. The power at these small scales exceeds that at CMB scales by roughly five orders of magnitude.
In order to better understand the red curve in Fig. 5, in Fig. 6 we show the evolution of P ζ with N for four representative values of k. These results have been obtained by While the results of Fig. 6 have been obtained from eq.(3.28), the qualitative behavior is more easily understood by combining eqs. (3.27) and (3.25), which yields the equivalent differential equation This equation again has to be satisfied by both the real and imaginary parts of ζ k . For sub-horizon modes, where k 2 a 2 H 2 , the last term in eq.(3.33) dominates; this by itself leads to an oscillatory behavior of ζ k , with amplitude increasing ∝ e N/2 and with exponentially decreasing oscillation frequency. For SR conditions, H , |η H | 1, the second term in eq.(3.33) is a damping term, which reduces the amplitude of the oscillations ∝ e −3N/2 . Altogether this yields P ζ ∝ e −2N , which explains the initial steep decline in all four cases depicted in Fig. 6.
Of course, the term ∝ k 2 in eq.(3.33) decreases ∝ e −2N , due to the exponential growth of a(N ); by definition the coefficient multiplying ζ k in this term equals 1 at N = N cross . Moreover, the SR conditions are badly violated in the overshooting region. The evolution of the power depends on the where N cross lies relative to the overshooting region. To see this, let us discuss the four cases depicted in Fig. 6 one by one. Fig. 6a: Here we chose k = 2.5 × 10 10 Mpc −1 , so that horizon crossing takes place at N = 29.7, where the SR conditions still hold. As shown in Fig. 6a, the power spectrum for this mode first approaches a constant after horizon crossing. Here the last term in eq.(3.33) is negligible. As long as the coefficient of the second term is close to +3, the absolute value of the first derivative of ζ k keeps decreasing exponentially with increasing N ; this corresponds to an overdamped oscillator. The solution for this range of N can thus be written as ζ k (N ) = C 1 + C 2 e −3(N −Ncross) , where C 1 , C 2 are two constants determined by the initial conditions 11 . Let N SR denote the number of e-folds which ζ k undergoes in the SR regime after horizon crossing, but before overshooting; then this epoch suppresses the first derivative of ζ k by a factor e −3N SR . Since the derivative of ζ k is small, ζ k itself is basically constant. This solution is no longer valid in the overshooting region, where η H ≈ 3 while H remains rather small, so that (3 + H − 2η H ) ≈ −3, i.e. the second term in eq.(3.33) changes sign relative to the SR epoch. This means that now the first derivative to ζ k begins to grow exponentially in magnitude, however without changing sign. At the end of this epoch one thus has ζ k (N ) = C 3 + C 4 e +3N OS , where N OS is the total "length" of the overshooting epoch, i.e. the number of e-folds during which η H ≈ 3. 12 By matching the two solutions for ζ k at the point where the overshooting epoch begins, one finds C 1 = C 3 and C 4 = −C 2 e −3N SR . So after overshooting ends, the value of ζ k is approximately given by ζ k = C 1 − C 2 e 3(N OS −N SR ) . Since afterwards the SR conditions hold again, ζ k remains approximately constant, i.e. this result still holds at the end of inflation.  The overshooting region will therefore only have significant impact on the final power for modes that crossed out of the horizon not much more than N OS e-folds before its onset. From the left frame of Fig. 7 we read off that for our numerical example overshooting starts at N ≈ 34, with N OS ≈ 3.5. For the case considered in Fig. 6a N OS and N SR are comparable. For much smaller co-moving wave number k, N SR N OS , so that the effect of the overshooting region on the final power is not significant. This explains why the standard SR approximation works for k < 10 10 Mpc −1 in Fig. 5.
The detailed evolution of the real and imaginary parts of ζ k is shown in Fig. 8. Note that the overshooting region has significant impact on ζ k itself (as opposed to its derivative) only beginning at N ≈ 37, where the exponential growth of the modulus of the derivative compensated its exponential suppression between horizon crossing and the onset of the overshooting epoch. Since overshooting already ends at N ≈ 38, its total effect is still moderate for this value of k. Notice, however, that the second flat region lies well below the first one, which corresponds to the prediction of the usual analytical SR estimates. This is because in the overdamped oscillator phase just after horizon crossing, the first derivative of ζ k always has the opposite sign as ζ k itself, for both the real and imaginary part. The exponential decrease of the modulus of the derivatives will therefore decrease |ζ k |, and thus P ζ . We will come back to this point shortly. Fig. 6b: for co-moving k = 10 11 Mpc −1 , P ζ nearly vanishes for a value of N during the overshooting epoch. 14 Now the nominal horizon crossing at N = 31.2 occurs just before the onset of the overshooting epoch, which means we cannot always assume k aH when we discuss the evolution of the curvature perturbation around the overshooting regime. In the following discussion we denote the coefficient of the second and third terms in eq.(3.33) by A and B, respectively; their dependence on N is depicted in Fig. 7. After horizon crossing, the evolution of ζ k undergoes four stages, which are shown in Fig. 9: Fig. 9a: in this region, A is always positive and therefore acts as a friction 13 Actually after the overshooting dynamics ends, the inflaton enters the USR phase, where the matter perturbation is even more adiabatic compared to that in a SR.
14 To the best of our knowledge, a similar behavior as shown in Fig. 6b was first explored in [67] and recently was mentioned in [49,54]. term, and B decreases exponentially. Since we are already beyond horizon crossing, B 2 −A 2 /4 < 0, i.e. eq.(3.33) approximately describes an over-damped oscillator. This means that ζ k does not oscillate any more, but decreases rather slowly. This explains the first, short flat region in Fig. 6b. In this region the second derivative of ζ k can be neglected, thus the curvature perturbation satisfies dζ k dN ≈ B A ζ k . This approximation ceases to hold somewhat before the value N 0 where A turns to zero, i.e. where the overshooting region starts; recall that in our case N 0 ≈ 34.
• 34 < N < 37, Fig. 9b: for N ≥ N 0 , B ≤ 0.01 (see Fig. 7) has become essentially negligible, while A changes from positive to negative hence acts as a driving term. As discussed above this leads to an exponential increase of the first derivative of ζ k . Since the epoch of exponentially decreasing first derivatives is considerably shorter than for k = 2.5 · 10 10 Mpc −1 , ζ k itself now begins to vary appreciably already at N ≈ 36.
Remarkably, shortly thereafter both the real and imaginary parts cross the zero point nearly at the same time, leading to ζ k → 0. This can be understood from the approximate solution for ζ k in this range of N : , we see that the real and imaginary parts of ζ k (N ) go through zero at the same point, which is also the origin of the very sharp minimum depicted in Fig. 6b. Note that for the previous case, k = 2.5 · 10 10 Mpc −1 , the overshooting epoch ended before ζ k reached zero; the sharp minimum of the final power spectrum depicted in Fig. 5 corresponds to that value of k where the overshooting epoch lasts just long enough to drive ζ k to zero, and then ends. In the case at hand instead |ζ k | again increases exponentially beyond the zero crossing.
• 37 < N < 38, Fig. 9c: this is the stage just after the overshooting epoch, i.e. A is again positive so that the modulus of the first derivative of ζ k is decreasing exponentially again. For a while |ζ k | keeps increasing, albeit more slowly than before. Of course, B ≤ 10 −5 is now completely negligible.
• 38 < N < 70, Fig. 9d: the overshooting phase has ended and the universe comes back to (U)SR inflation, therefore curvature perturbations are frozen again at the super horizon scale and matter perturbation evolve adiabatically again. This corresponds the second flat region shown in Fig. 6b. For this value of k the second flat region is already higher than the first one, thus the analytical SR approximation underestimates the final power.  Nevertheless for N < 34 eq.(3.33) again describes a damped oscillator, the amplitude of the oscillation decreasing ∝ e −N . This is also shown in Fig. 10a. However, at N = N 0 34 the second term in eq.(3.33) changes sign, eventually reaching −3 as shown in Fig. 7. For N 0 < N < N cross eq.(3.33) therefore leads to oscillations whose amplitude grows ∝ e 2N . For N > N cross the approximately exponential growth continues for a while, this time for |ζ k | itself which no longer oscillates. As before, for N > 37 the derivative of ζ k begins to decrease exponentially in magnitude, which leads to ζ k itself becoming essentially constant for N > 38, after the end of the overshooting epoch. Not surprisingly, for this mode the analytical SR estimate for the final power also fails. Fig. 6d: for the mode with co-moving wave number k = 2 × 10 14 Mpc −1 , horizon crossing takes place at N = 39.5, i.e. during USR well after the end of the overshooting epoch. We again see an (initially very rapid) oscillation whose amplitude first drops ∝ e −N and then increases ∝ e 2N once N > 34, i.e. in the overshooting regime. Since for this value of k the overshooting epoch ends before horizon crossing, for N > 38 the function ζ k again undergoes a few oscillations with exponentially decreasing amplitude, before settling into an overdamped oscillator mode, i.e. approaching a constant.
In this case the analytical SR approximation for the final power actually works quite well. On the one hand this may not be surprising, since the arguments of Sec. 2 imply that perturbations are now adiabatic at super-horizon scales. On the other hand, it may be surprising that one still gets the correct result by imposing the initial conditions (3.30) just a couple of e-folds before horizon crossing. This implies that these initial conditions capture the dynamics of the MS equation even in the overshooting regime, as long as the mode is still (deep) inside the horizon.
In fact, eq.(3.30) shows that |v k | is simply a constant (independent of N ) for subhorizon modes. The dynamics is therefore entirely captured by the factor 1/z which relates ζ k to v k , see eq.(3.26). This contains a factor 1/a ∝ e −N , which dominates the N −dependence in the (U)SR regime where the SR parameter H is approximately constant (and small). However, eq.(3.12) shows that in the overshooting regime, where the term containing the derivative of the potential can be neglected, , so that altogether |ζ k | ∝ 1/(a √ H ) ∝ e +2N , as we had inferred from the MS equation. It should be noted that this agreement between the analytical treatment and the solution of the MS equation crucially hinges on using the "Hubble" version of the SR parameters; we saw in Fig. 4 that H and V are very different in the overshooting epoch.
For modes with even larger k, the situation is very similar to the case with k = 2 × 10 14 Mpc −1 , i.e. the analytical SR approximation for the final power agrees with the numerical result, as shown in Fig. 5, as long as the modes cross out of the horizon (well) before the end of inflation. The approximation fails again for modes with very large k which cross the horizon near the end of inflation where SR again fails; however, we know of no way to probe those modes observationally.
From the above discussion it is easy to understand that the peak in the power spectrum shown in Fig. 5 occurs for the mode which crosses the horizon just at the beginning of the overshooting regime. In this case, N SR = 0, and the exponential increase of |dζ k /dN | is maximized. This also greatly enhances the final value of |ζ k |, leading to the maximum in the spectrum. We find the maximal scalar power spectrum is P ζ ≈ 1.1 × 10 −4 for k = 1.1 · 10 13 Mpc −1 . According to ref. [35], density fluctuations of this size are large enough to lead to significant formation of primordial black holes, which might even constitute a sizable fraction of all dark matter. In the next section we will investigate another cosmological consequence of such a large curvature perturbation, namely the amplification of primordial gravitational wave signatures due to second order effects.
Before closing this section we comment on possibilities to increase the power spectrum even further. According to ref. [54], quantum diffusion effects can in principle further enhance the power; however, we checked that in our case one always has |χ| H 2 /(2π), which indicates that quantum diffusion does not change the evolution of the inflaton field significantly. Moreover, refs. [48,68] argue that a shallow local minimum of the inflaton potential can also enhance the power spectrum; this agrees with our finding in eq.(2.18), because the curvature of the potential V is maximal near a local minimum. However, the inflaton might get stuck in a local minimum, in which case inflation would never end. In contrast, the scenario we presented leads to a well-behaved inflationary epoch, in agreement with current observations.

Second Order Gravitational Wave Signatures
As well known, primordial perturbations of the inflaton field source primordial gravitation waves. Usually the strength of the GW signal is estimated in linear order in perturbations; for SR inflation, this leads to the famous prediction r = 16 V , where r is the tensor-to-scalar ratio. However, in some cases effects that are second order in the curvature perturbations can also contribute significantly to the primordial GW signal [69][70][71]. As has recently been emphasized in [37], which analyses a polynomial potential with an inflection point, this occurs in particular when an overshooting regime enhances the power spectrum. In the following analysis, we mainly follow the formalism given in [37,71].
In order to obtain the current gravitational wave density, we have to solve eq.(4.1) with source given by eq.(4.2). To that end we'll apply the Green's function method of ref. [71]. Rewriting eq.(4.1) with v := ah, we get The solution of eq.(4.1) can then be written as where g is the Green's function for eq.(4.3), which satisfies: Once the tensor perturbations are known, we can further compute the contribution of these primordial gravitational waves to the total energy budget of the universe, which can be written as [71]: .
The result is shown in Fig. 11, which also shows the sensitivity of several planned gravitational wave detectors. We saw at the end of the last Section that the maximum of the power spectrum is at k = 1.1 × 10 13 Mpc −1 , which corresponds to frequency f = 0.017 Hz. This is near the frequency of maximal sensitivity of the upcoming space mission LISA, which may just barely be able to detect this signal if the parameter β = 10 −5 (red curve), while for β = 10 −4 (green) the signal is below the foreseen LISA sensitivity. Recall that the CMB predictions for both values of β are consistent with latest Planck measurements, see Sec. 3.2. Since a larger β makes the potential less flat in the USR region and reduces V in the overshooting region, it is expected that the corresponding curvature perturbation is less enhanced compared to that with a smaller β according to eq.(2.18). This explains why the peak of GW signatures with β = 10 −4 is lower. However, the second generation space missions DECIGO and BBO should easily detect this signal even for β = 10 −4 . Fig. 11 also shows that the peak of the CHI signal lies at frequencies that are too large for the pulsar timing arrays even after SKA comes on-line. The size of the signal is well below the sensitivity of advanced LIGO, and even below that of the planned Einstein Telescope (ET).

Summary and Conclusions
In this paper, we have revisited critical Higgs inflation, carefully computing the power spectrum as well as the gravitational wave signatures induced by second order effects.
In Sec. 2 we analyzed the evolution of curvature perturbations under (ultra-)slow roll as well as overshooting conditions in general terms. In the former, the second derivative of the inflaton field with respect to time can be neglected in the equation of motion; we showed that the matter perturbations are adiabatic in this case, which further implies that they are frozen at super-horizon scales. This allows one to calculate the final power spectrum (at the end of inflation), which seeds all observed structures in the universe, by simply computing the power spectrum at horizon crossing. We emphasize that this also holds for ultra-slow roll (USR), which in our model describes the epoch when the inflaton field is near the (almost) saddle point of the potential. Here the deviations from adiabacity are even smaller than in the SR case, so that the usual approximate treatment is even more accurate.
In contrast, when the inflaton enters an overshooting phase where the acceleration |χ| is much larger than the derivative of the potential |V |, we showed that entropic perturbations are excited, which are not conserved at super-horizon scales. In this case, the standard SR approximation to calculate the power spectrum is expected to break down. To our knowledge this is the first time that the significance of entropic perturbations has been discussed in this context. Our eq.(2.18) shows that the enhancement of the perturbations after horizon crossing but during the overshooting epoch will increase for larger curvature V of the potential. This can be very useful for inflationary model building if one wants to strongly enhance the power spectrum, e.g. in order to produce primordial black holes.
In Sec. 3 we illustrated these general results by analyzing the CHI scenario in detail. For judiciously chosen parameters, an overshooting epoch appears between the SR and USR eras. During the overshooting stage the Hubble SR parameters vary rapidly, which implies the universe deviates significantly from SR evolution. As a result the usual analytical approximation to compute the power spectrum fails. We instead solved the Mukhanov-Sasaki equation numerically to compute the power spectrum at the end of inflation. We find that the modes which cross the horizon just before or during the overshooting epoch are greatly enhanced. The power spectrum can reach values of order ∼ 10 −4 for k ∼ 10 13 Mpc −1 ; this is to be compared to values of order 10 −9 at the (much smaller) k−values probed by the CMB anisotropies.
In the course of this discussion we found a version of the MS equation very useful which holds for the k−space perturbation ζ k directly, rather than for the related quantity v k which is usually employed, see eq. (3.33). This allowed us to understand the numerical results in detail: why SR approximation agrees with MS formalism for modes that cross out of the horizon well before (small k) or well after (large k) the overshooting epoch; why there is a sharp minimum in the power spectrum, for scales that cross the horizon a few e-folds before the overshooting region; and where the maximum of the spectrum lies. These findings are generic and can also be applied to explain the numerical power spectrum results for other inflation models featuring a near-inflection point, for example [37,50,54], where detailed explanations concerning the numerical results are not given.
Finally we analyzed the second order GW signatures induced by the enhanced scalar perturbations. The strength of this signal is proportional to the square of the scalar power spectrum. The peak of the latter at co-moving wave number of order 10 13 Mpc −1 corresponds to a peak of the GW signal at a frequency of 0.017 Hz. We find that for our choices of parameters, the GW signal should remain detectable up to frequency of order 1 Hz by two planned second-generation space based GW experiments, DECIGO and BBO. Detection of this signal is a firm prediction, if the power spectrum is enhanced to the level that might allow significant production of PBHs. This statement holds also in other models of inflation proposed recently [36,37,[43][44][45][46][47][48][49][50]. Hence if future GW experiments fail to detect this signal, one could conclude that no significant PBH formation occurred immediately after inflation.
In this paper we did not consider effects due to non-Gaussianity. Since PBHs only form in regions with large overdensity, their formation rate can be greatly enhanced if there are significant non-Gaussian tails in the distribution function of the density perturbations [83,84]. It should be noted that the calculation of the PBH formation rate is in any case somewhat uncertain; however, the second order GW signal computed in the Gaussian approximation should be detectable by second generation space missions for the entire range of density perturbations that could plausibly lead to sizable PBH formation, unless non-Gaussianities are quite large for the relevant modes. In this context it is important to note that according to a recent analysis [85] primordial non-Gaussianities will also enhance the second order GW signal itself. We leave a detailed investigation of the impact on non-Gaussianities on CHI inflation for future work.
Unfortunately this definition is not unique, since eq.(A.11) only fixes the product v kâk , i.e. the vacuum state defined by eq.(A.21) depends on the form of mode function v k . In other words, eq.(A.21) defines the vacuum state uniquely only once v k has been fixed.
To that end we consider the limit τ → −∞, such that |kτ | 1 or k aH, where k = |k|; this corresponds to perturbations with wavelength much smaller than the Hubble horizon. In this limit the mode function v k behaves like a massless field in Minkowski spacetime, since the z term in eq.(A.10) can be neglected compared to k 2 : v k + k 2 v k = 0 . Eq.(A.24) fixes the mode function v k and thus also the vacuum state (up to some angledependent phase factor, which is not physically relevant); this is usually referred to as the Bunch-Davies vacuum.

A.3 Analytical Solution in Quasi-de Sitter Spacetime
We now describe the analytical solution of the MS equation in the limit where the Hubble parameter H is nearly constant. This also means that the Hubble SR parameter H is small and nearly constant, thus the time derivative of H can be neglected; these conditions are met during (U)SR inflation. Using eq.(A.4), we then obtain: The general analytical solution of this equation is given by [60] v where α and β are integration constants. The initial conditions in eq.(A.24) imply α = 1 and β = 0, which leads to the Bunch-Davies mode functions v k = e −ikτ √ 2k 1 − i kτ . (A.28)

A.4 Power Spectrum in Quasi-de Sitter Spacetime
Having solved the Mukhanov-Sasaki equation, we can compute the power spectrum of the field, δφ k = a −1v k : 0|δφ k (τ ) δφ k (τ )|0 = (2π) 3 δ(k + k ) |v k (τ )| 2 a 2 = (2π) 3 δ(k + k ) H 2 2k 3 (1 + k 2 τ 2 ) , where we have used eq(A.28) as well as the expression for the scale factor a(τ ) = − 1 Hτ which holds for constant H, i.e. during (U)SR inflation. On super-horizon scales, |kτ | 1 or equivalently k aH, eq.(A.29) becomes 0|δφ k (τ ) δφ k (τ )|0 → (2π) 3 δ(k + k ) H 2 2k 3 , (A. 30) or in a dimensionless version (recall that we are using Planckian units where M p = 1): (A.31) Eq.(A.31) also implies δφ 2 k = H/(2π), which is the frequently used formula for the quantum fluctuations of light fields (with mass smaller than H) during SR inflation. As shown in Section 2.2, during SR inflation curvature perturbations are frozen at superhorizon scale, thus the power spectrum can be computed at the horizon crossing, i.e. for a(t k )H(t k ) = k [60]: The corresponding dimensionless power spectrum is Eq.(A.33) is widely used in the literature when discussing SR inflation, where H is small and its variation with time can be neglected. Moreover, during SR inflation, the energy is mainly dominated by the potential, thus we have H 2 = V 3 . Using in addition the SR solution for the equation of motion of the inflaton field,φ = −V 3H , allows us to rewrite H as: V is usually called the potential SR parameter. This leads to another frequently used formula for the power spectrum: (A. 35) However, for non SR inflation -in particular, during the overshooting epoch which we have explored in this paper -H changes rapidly and z /z = 2/τ 2 . There eqs.(A.33), (A.34) and eq.(A.35) are no longer valid; note in particular that eq.(A.35) predicts a diverging power spectrum at a true saddle point where V , and hence V , vanishes. In this case we must solve the Mukhanov-Sasaki eq.(A.10) numerically in order to reliably estimate the power spectrum at the end of inflation.