Effects of discontinuities of the derivatives of the inflaton potential

We study the effects of a class of features of the inflaton potential, corresponding to discontinuities in its derivatives. We perform fully numerical calculations and derive analytical approximations for the curvature perturbations spectrum and the bispectrum which are in good agreement with the numerical results. The spectrum of primordial perturbations has oscillations around the scale k0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$k_0$$\end{document} which leaves the horizon at the time τ0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau _0$$\end{document} when the feature occurs, with the amplitude and phase of the oscillations determined by the size and the order of the discontinuity. The large scale bispectrum in the squeezed and equilateral limits have a very similar form and are linearly suppressed. Both in the squeezed and the equilateral small scale limit the bispectrum has an oscillatory behavior whose phase depends on the parameters determining the discontinuity, and whose amplitude is inversely proportional to the scale. Given the generality of this class of features they could be used to model or classify phenomenologically different types of non-Gaussian features encountered in observational data such as the cosmic microwave background radiation or large scale structure.


Introduction
In the last few decades the outstanding advances in observational cosmology have allowed for the first time to test theoretical cosmological models [1][2][3][4]. Among the most important sources of cosmological observational data we can mention the Sloan Digital Sky Survey (SDSS), the Wilkinson Microwave Anisotropy Probe (WMAP), and the Planck mission, and other ground-based and sub-orbital experiments [5,6]. According to the standard cosmological model the cosmic microwave background (CMB) radiation consists of photons which decoupled from the primordial plasma at the time when protons and electrons combined to form neutral light atoms. Although this radiation is extremely isotropic a e-mail: aer@physics.uoc.gr there are small fluctuations in the temperature of the order of T /T ∼ 10 −5 . Since the CMB radiation was emitted at a redshift of about 1100 it provides a unique window on the early universe [7][8][9].
Inflation theory [10] explains the anisotropies of the CMB temperature as the consequence of primordial curvature perturbations whose statistical properties can be described by the n-points correlation functions. If the perturbations followed a perfectly Gaussian distribution the two points correlation function would be enough, but even the most recent observations are compatible with some non-Gaussianity corresponding to f local NL = 2.5 ± 5.7 and f equil NL = −16 ± 70 [8,11], motivating the theoretical study of the conditions which could have generated it. Some recent developments in the study of models which could generate non-Gaussianity and in the detection can be found for example in [12][13][14][15].
The theoretical study of the effects of features of the inflaton potential was started in the seminal works of Starobinsky [16], and once CMB observational data became available it was shown that features can be used to model the glitches of the power spectrum [17,18]. Some other interesting studies and reviews in this area can be found for example in [19][20][21][22]. In this paper we focus on the effects of features of the inflaton potential on the primordial curvature perturbations, considering a class corresponding to a discontinuity in the derivatives of the potential. Our model is a generalization of other features which have been studied earlier such as the Starobinsky model or the mass step [23]. These kinds of features could have arisen through different mechanisms such as for example particle production [24], or phase transitions [25], but in this paper we study their effects from a purely phenomenological point of view, without investigating their fundamental origin.
There is also an important observational motivation for studying this kind of potentials: recent analyses of CMB observations based on cubic Hermite interpolating polynomials (PCHIP) for the primordial curvature perturbations spec-trum [26] have in fact shown some evidence for a feature around the wave number k = 0.002 Mpc −1 , which is in good qualitative agreement with the results of our calculations for some of the potentials we consider.
The paper is organized as follows: first we define the features, then we give both a numerical and analytical solution for the background, and finally we provide both numerical and analytical calculations of the spectrum and the bispectrum, giving details of the squeeze and equilateral limit and showing the effects of varying the different parameters defining the feature, i.e., its amplitude and the order n of the discontinuous derivatives.

Inflation
We consider inflationary models with a single scalar field and a standard kinetic term according to the action [27,28] where M Pl = (8π G) −1/2 is the reduced Planck mass. Varying the action with respect to the metric tensor and the scalar field we get Friedmann equation and the equation of motion of the inflaton, where dots and ∂ φ indicate derivatives with respect to time and scalar field, respectively, and H is the Hubble parameter. We adopt the following definitions of the slow-roll parameters:

The model
We consider a single scalar field φ with potential where where V a and V b are different in order to ensure the continuity of the potential at φ 0 .
The value of φ 0 determines the scale at which the effects of the feature appear in the power spectrum and the bispectrum of curvature perturbations, and as such it is a free parameter, which can be fixed phenomenologically based on experimental data. It in fact determines the value of conformal time when φ(τ 0 ) = φ 0 , and, as will be shown in the following sections, the features in the spectrum and bispectrum appear around the scale k 0 = −1/τ 0 which is leaving the horizon at that time.
The potential has a discontinuity in the derivatives at φ 0 , which is the cause of the temporary slow-roll regime violation which produce the non-Gaussian features. The continuity condition for the potential at τ 0 gives In this paper we study potentials dominated by the vacuum energies V b and V a before and after the feature. The potential in Eq. (5) is similar to the one studied in [29,30], but it only coincides with it in the special case (φ 0 = 0, p = 2). Another important difference is that inflation in our model is driven by the dominating vacuum energy term V a .

Analytic solution of the background equations
The Friedmann equation and the equation of motion for the inflaton in terms of conformal time τ take the form where primes indicate derivatives with respect to conformal time.
Since V (φ) is dominated by the vacuum energy we can use the de Sitter approximation, in which H is a constant and the scale factor is given by Before the feature the equation of motion of the inflaton is which has the solution where The slow-roll regime corresponds to φ − b = 0. After the feature the equation of motion of the inflaton becomes In order to find an analytical solution we can expand the last term in Eq. (14) to second order in conformal time around τ 0 From now on quantities evaluated at τ 0 are denoted by the subscript 0. From Eq. (12) we have an analytical expression for the first and second derivative of the field at τ 0 , where we have only kept terms linear in λ + because during the slow regime higher order terms can be safely neglected according to Assuming the same slow-roll regime condition we can expand Eq. (15) to linear order in λ + and get which admits an analytical solution of the form where The constants of integration φ ± a are determined by imposing the continuity conditions for φ and φ at τ 0 , which give We can also find an analytical approximation for the slowroll parameters after the feature by substituting Eq. (20) in Eq. (4),

Numerical solution of the background equations
The background evolution can be obtained by solving the system of coupled differential equations for a(τ ) and φ(τ ), or alternatively H (τ ) and φ(τ ). In the numerical integration we chose the following value for the different parameters defining the model: This choice of the parameters is made in order to satisfy the Planck normalization on small scales. If the term λ φ is of the order of the vacuum energy V b , then from Eq. (2) we have implying that in this case the de Sitter approximation used to obtain an analytical solution in the previous section is not valid, as shown in Fig. 1, and the numerical integration is necessary to obtain reliable results. As stated previously, we will focus on vacuum energy dominated models for which the de Sitter approximation is valid, so Fig. 1 is given only to show the limits of its validity, but in all the cases we consider it turns out to be quite accurate as shown in the figures comparing analytical results, based on the de Sitter approximation, and to numerical results, which take into account the small variation of the Hubble parameter. We adopt a system of units in which M Pl = 1. As can be seen in the left plot the de Sitter approximation is not valid at late times when λφ n > V 0 , so that a full numerical integration of the background equations is required, while from the right plot we can see that when λφ n < V 0 the de Sitter approximation is quite accurate The potential as a function of the field is shown in Fig. 2 for different types of features. The effects of the features on the scalar field φ and on the slow-roll parameters are shown in Figs. 3, 4, 5 and 6. When |λ| is kept constant, larger values of n tend to produce larger variations of both and η, while when n is kept constant, larger values of |λ| tend to produce larger variations of both and η. As can be seen in Figs. 3 and 4 the numerical and analytical solution for the scalar field are in good agreement. The analytical approximation is also good for the slow-roll parameters as shown in Figs. 5 and 6. We can conclude that Eq. (20) is a good approximation for the background solution within the limits of validity of the assumptions used to derive it, and we will use it in the following sections to calculate curvature perturbations. The analytical solution we derived around the time τ 0 , when the feature occurs, should be accurate as long as the de Sitter approximation is valid and |λ + | 1.

Spectrum of curvature perturbations
In order to study curvature perturbations we need to expand perturbatively the action respect to the background FRLW solution [31,32]. We adopt the comoving gauge, in which there is no fluctuation in the scalar field, δφ = 0. The second and third order actions are, respectively, and δL/δζ | 1 is the variation of the quadratic action with respect to ζ [31]. The Lagrange equations for the second order action give the equation for the curvature perturbations ζ ∂ ∂t The Fourier transform of the above equation, using conformal time, gives where z ≡ a √ 2 and k is the comoving wave number. It is convenient to define the variable [23] in terms of which Eq. (33) takes the form As can be seen in Fig. 7, small scale modes, which are sub-horizon at time τ 0 , are affected by the feature. Modes that had left the horizon at that time are unaffected, since they were already frozen. In Fig. 8 the power spectrum of  primordial curvature perturbations P ζ is plotted for different types of features.

Analytical approximation for curvature perturbations
At the time of the feature there is a discontinuity in φ which implies that z contains a Dirac delta function [24]. We can evaluate the discontinuity in z /z by integrating the Dirac delta function around the feature time, Before the feature we assume the Bunch-Davies vacuum, [33] v(τ, The curvature perturbations modes after the feature are approximated as a linear combination of the positive and negative frequency modes before the feature [16,19] according to where are the Bogoliubov coefficients, and τ k = −1/k is the horizon crossing time for the mode k. The coefficients α, β are determined by imposing the continuity of the modes and their derivative at the feature time [24]. In Figs. 9 and 10 we show the comparison between the numerical results for the mode function and the analytical approximation for small scales. The parameters used for the feature are n = 4 and λ = 3.9 × 10 −19 . In Fig. 9 we show the real and imaginary part of the dependence of the mode functions on the scale and for a particular time, namely, after ten e-folds after the feature. In Fig. 10 we show the evolution of the real and imaginary parts of the mode function at a particular scale 100k 0 .

Analytical approximation for the spectrum
The two-point function is where the power spectrum of curvature perturbations is defined as where τ e is the time at which inflation ends. After substituting in the above definition the analytical approximations obtained in the previous sections we get This generalizes the result obtained in [16] for a potential with a discontinuous first derivative to the more general case considered in this paper, corresponding to Eq. (36). In Fig. 8 we compare the analytical expression for the power spectrum given by Eq. (42) with the numerical results obtained by integrating numerically both the background and the perturbation equations. The analytical result is quite accurate and it improves the results obtained in [23] because we use the analytical approximation for the perturbations modes also for modes which were superhorizon slightly before τ 0 , improving substantially the agreement with numerical results. This is due to the well-known fact that modes are not completely frozen at τ k = −1/k but keep evolving for a few e-folds thereafter, so that also scales slightly greater than k 0 are mildly affected by the features.
It should be noted that Eq. (42) we obtained is not depending on the slow-roll approximation, since it is derived directly from the definition of the power spectrum in Eq. (41), and this explains why it is in so well in agreement with fully numerical calculations despite the temporary violation of slow-roll regime produced by the features. The in the denominator of Eq. (42) comes in fact from the analytical solution for the perturbation modes, which is not based on any slow-roll expansion because it comes from z ≡ a √ 2 in Eqs. (33) and (35) which are valid at any order in slow roll.
As can be seen in Fig. 8 negative values of λ correspond to a suppression of the spectrum on large and intermediate scales, while positive values produce a suppression on small scales. While the analysis of observational data goes beyond the scope of this paper, we can see in Fig. 8 that an appropriate choice of parameters gives spectra in good qualitative agreement with the features recently found when parameterizing the free primordial power spectrum with a piecewise cubic Hermite interpolating polynomial [26].

Calculation of the bispectrum
The Fourier transform of the three-point correlation function [8,11], also known as the bispectrum B ζ , is given by and it should vanish if the curvature perturbations are Gaussian [9,34]. Therefore, deviations from non-Gaussianity can be determined by measuring the implications of a nonvanishing B ζ on the CMB radiation. Following the procedure in Ref. [23], after a field redefinition, the third order action can be written From this action the interaction Hamiltonian can be written in terms of conformal time as The three-point correlation function is given by [9,31] |ζ(τ e , and after substitution the expression for the bispectrum B ζ [9,35] is where is the imaginary part and we evaluate the integral from τ 0 to τ e , where τ e is some time sufficiently long after the Hubble crossing horizon, when the modes are frozen [14,15,34,36,37]. It is common to study non-Gaussianity using the parameter f NL defined by where Replacing P ζ in Eq. (48) we obtain f NL in terms of our dimensionless definition of the spectrum P ζ (k), In this paper we will study non-Gaussianity using a different quantity defined as where k * is a pivot scale at which the power spectrum is evaluated which corresponds approximately to the scale of normalization of the spectrum, i.e. P ζ (k * ) ≈ 2.2×10 −9 . Our definition of F NL reduces to f NL in the equilateral limit if the spectrum is approximately scale invariant, but in general f NL and F NL are different, and for example in the squeezed limit they are not the same. For this reason they cannot be compared directly but F NL still provides useful information as regards the non-Gaussian behavior of B ζ . In Figs

Analytical approximation for the bispectrum
In order to obtain an analytical approximation for the bispectrum we use Eq. (38) for curvature perturbations and Eq. (25) for slow-roll parameters. This implies that also the different approximations are used in different cases as explained in more detail in the following sections. All the results presented in this section should be considered taking into account the existence of a cut-off scale beyond which the Heaviside approximation of a smooth transition is not valid as discussed in more detail in [23,24]. We provide analytical expressions for the bispectrum for different types of features, i.e., different values of n and λ, in the squeezed and equilateral limit for both large and small scales. The analytical results are shown in dashed black lines in Figs. 11, 12, 13 and 14, where it can be seen that the approximations for the bispectrum are in good agreement with the numerical results.

Large scales
In the large scale limit when k i < k 0 , i = 1, 2, 3, the curvature perturbations modes are frozen in the time interval of interest, since there is no time evolution for τ > τ 0 > τ k i . Thus in Eq. (47) all the modes functions can be evaluated at τ 0 and pulled out of the integrals while the terms ζ * (τ, k i ) can be set to zero. Following this approximation we get    k, k, k) is plotted for small scales. On the left λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green). On the right n is constant, n = 3, while λ = 6.0 × 10 −19 (blue), λ = 1.2 × 10 −18 (red), and λ = 2.4 × 10 −18 (green). The dashed black lines correspond to the analytical approximation Fig. 15 On the top the large scales equilateral limit F NL (k, k, k) and squeezed limit F NL (k 0 /500, k, k) are plotted respectively on the right and left. On the bottom the small scales equilateral limit F NL (k, k, k) and squeezed limit F NL (k, 1000k 0 , 1000k 0 ) are plotted respectively on the right and left. The dashed lines are the analytical approximations. All the plots are for n = 3, λ = −8 × 10 −20 (orange) and λ = −6 × 10 −20 (cyan). This choice of parameters is the same used in Fig. 16, for models able to account for the observed large scale suppression where we have used the approximations for the slow-roll parameters in Eq. (25) in the integration. Now we use the analytical approximations for the perturbation to obtain at large scales In the squeezed limit with k 1 k 2 = k 3 ≡ k and k < k 0 this expression reduces to the following analytical formula: As shown in Fig. 11 the analytical approximation is in good agreement with the numerical results. Here and in any other approximation for the F NL as defined in Eq. (51) we use Eq. (42) for the spectrum P ζ . In the large scale equilateral limit, when k 1 = k 2 = k 3 ≡ k k 0 , Eq. (53) becomes The numerical result and Eq. (55) are in good agreement, as shown in Fig. 12. In Figs. 11 and 12 we have evaluated both the numerical and the analytical expressions at time τ e corresponding approximately to 10 e-folds after the feature [14,15,34,36]. As can be seen in Figs. 11 and 12 the large scale bispectrum in the squeezed and equilateral limits have a very similar form and are linearly suppressed.

Small scales
In the small scale limit, when k i > k 0 , it is convenient to re-write the expression for F NL as where In the above expressions we have used Eq. (25) to derive an approximation for η a 2 as All the integrals we need to compute have a similar form, so it is useful to re-write Eqs. (57) and (58) as where we have defined The above integrals can be computed analytically in terms of functions and are given in detail in "Appendix".
It is now possible to obtain a fully analytical template in the squeeze limit, when k 2 = k 3 , and k 2 k 1 > k 0 , In the equilateral limit, when k ≡ k 1 = k 2 = k 3 and k > k 0 , instead we have Numerical results and the analytical templates are in good agreement both in the squeezed and equilateral limits as shown in Figs. 13 and 14. In the squeezed and equilateral small scale limits the bispectrum has an oscillatory behavior whose phase and amplitude depend on the value of the parameters n and λ as can be seen in Figs. 13 and 14. The amplitude is inversely proportional to the scale as shown in Figs. 13 and 14. As previously observed, all the results derived can be trusted only up to cut-off scales beyond which the Heaviside approximation is not valid, as discussed in more detail in [23,24]. The same applies to other similar models previously studied, such as the well-known Starobinsky model [16].

Behavior of the small scale bispectrum
As seen in Figs. 13 and 14 both the equilateral and the squeezed limit small scale bispectrum do not behave in the same way as the spectrum and the slow parameters respect to the variation of n. To clarify this we can write the bispectrum B ζ in Eq. (47) as where First of all we can see from Figs. 17 and 20 that the dominant contribution to the bispectrum comes from the term (B 1 ) (B 2 ), which in fact behaves in the same way as the bispectrum with respect to the variation of n (Figs. 18, 19, 21,   Fig. 17 The real part of B 1 is plotted on the left, and the imaginary part on the right, for the small scale squeezed limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) Fig. 18 The real part of B 2 is plotted on the left, and the imaginary part on the right, for the small scale squeezed limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) Fig. 19 The products (B 1 ) (B 2 ) and (B 1 ) (B 2 ) are plotted for the small scale squeezed limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) Fig. 20 The real part of B 1 is plotted on the left, and the imaginary part on the right, for the small scale equilateral limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) 22). The two terms (B 2 ) and (B 1 ) behave like the spectrum, i.e., are larger for larger values of n, but, since (B 1 ) is negative, their product (B 1 ) (B 2 ), and consequently the bispectrum which is dominated by it, behaves in the opposite way, i.e. it is decreasing when n in increasing. The effect is not noticeable in the case of n = 2/3, because in this case (B 2 ) is very close to zero, while it is clear for n = 3 and n = 4.

Conclusions
We have studied the effects of a general type of features produced by discontinuities of the derivatives of the potential. We found that each different type of feature has distinctive effects on the spectrum and bispectrum of the curvature perturbations which depend both on the order n and on the amplitude λ of the discontinuity. The spectrum of primordial Fig. 21 The real part of B 2 is plotted on the left, and the imaginary part on the right, for the small scale equilateral limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) Fig. 22 The products (B 1 ) (B 2 ) and (B 1 ) (B 2 ) are plotted for the small scale squeezed limit. The parameter λ is constant, λ = 3.9 × 10 −19 , while n = 2/3 (blue), n = 3 (red), and n = 4 (green) curvature perturbations shows oscillations around the scale k 0 , which leaves the horizon at the time τ 0 when the feature occurs, with amplitude and phase determined by the parameters n and λ.
Both in the squeezed and equilateral small scale limit the bispectrum has an oscillatory behavior whose phase depends on the parameters determining the discontinuity, and whose amplitude is inversely proportional to the scale. The large scale bispectrum in the squeezed and equilateral limits have a very similar form and are linearly suppressed.
The analytical approximation for the spectrum is in good agreement with the numerical results, and improves substantially the accuracy for large scales respect to previous studies. The analytical approximations for the bispectrum are in good agreement with numerical calculations at large scales in both the squeeze and the equilateral limit. At small scales we found an analytical template which is in very good agreement with the numerical calculations both in the squeezed and the equilateral limit, and it is able to account for both the oscillations and the amplitude of the bispectrum.
The type of feature we have studied generalizes previous models such as the Starobinsky model or the mass step [23], providing a general framework to classify and model phenomenologically non-Gaussian features in CMB obser-vations or in large scale structure survey data. In the future it would be interesting to find the parameters which better fit different non-Gaussian features in observational data and to investigate what more fundamental physical mechanism, such as phase transitions for example, could actually produce these features.

(A.2)
In order to simplify the calculation we fix (τ ) = (τ 0 ) only in the analytical approximation for perturbations modes in Eq. (39), while we keep (τ ), η(τ ) as functions of conformal time when they appear explicitly in the integrand.
After some rather cumbersome calculations the final result can be written in this general form