# Non-Gaussianity from Particle Production during Inflation

Authors: Neil Barnaby
Journal: Advances in Astronomy
 (2010)

---

## Abstract

In a variety of models the motion of the inflaton may trigger the production of some non-inflaton particles during inflation, for example via parametric resonance or a phase transition. Such models have attracted interest recently for a variety of reasons, including the possibility of slowing the motion of the inflaton on a steep potential. In this review we show that interactions between the produced particles and the inflaton condensate can lead to a qualitatively new mechanism for generating cosmological fluctuations from inflation. We illustrate this effect using a simple prototype modelg2(ϕ−ϕ0)2χ2 for the interaction between the inflaton, ϕ, and iso-inflaton, χ. Such interactions are quite natural in a variety of inflation models from supersymmetry and string theory. Using both lattice field theory and analytical calculations, we study the production of χ particles and their subsequent rescatterings off the condensate ϕ(t), which generates bremsstrahlung radiation of light inflaton fluctuations δϕ. This mechanism leads to observable features in the primordial power spectrum. We derive observational constraints
on such features and discuss their implications for popular models of inflation. Inflationary particle production also leads to a very novel kind of nongaussian signature which may be observable in future missions.


---

## Body

## 1. Introduction

In recent years the inflationary paradigm has become a cornerstone of modern cosmology. In the simplest scenario the observed cosmological perturbations are seeded by the quantum vacuum fluctuations of the inflaton field [1–5]. This mechanism predicts a nearly scale invariant spectrum of adiabatic primordial fluctuations, consistent with recent observational data [6]. In addition to this standard mechanism, there are also several alternatives for generating cosmological perturbations from inflation; examples include modulated fluctuations [7–10] and the curvaton mechanism [11]. These various scenarios all lead to similar predictions for the power spectrum. On the other hand, nongaussian statistics (such as the bispectrum) provide a powerful tool to observationally discriminate between different mechanisms for generating the curvature perturbation. In this paper, which is based on [12–15], we will present a qualitatively new mechanism for generating cosmological perturbations during inflation. We discuss in detail the predictions of this new scenario for both the spectrum and nongaussianity of the primordial curvature fluctuations, showing how this new mechanism may be observationally distinguished from previous approaches.


### 1.1. Non-Gaussianity from Inflation

The possibility to discriminate between various inflationary scenarios has led to a recent surge of interest in computing and measuring nongaussian statistics. Although single field, slow roll models are known to produce negligible nongaussianity [16–18], there are now a variety of scenarios available in the literature which may predict an observable signature. Departures from gaussianity are often parametrized in the following form:(1)ζ(x)=ζg(x)+35fNL[ζg2(x)-〈ζg2(x)〉],
where ζ(x) is the primordial curvature perturbation, ζg(x) is a Gaussian random field, and fNL characterizes the degree of nongaussianity. The ansatz (1) is known as the “local” form of nongaussianity.Although the local ansatz (1) has received significant attention, it is certainly not the only well-motivated model for a nongaussian curvature perturbation. For example, the nongaussian part of ζ(x) need not be correlated with the gaussian part. Consider a primordial curvature perturbation of the form(2)ζ(x)=ζg(x)+FNL[χg(x)],
where FNL is some nonlinear (not necessarily quadratic) function, and χg(x) is a gaussian field which is uncorrelated with ζg(x). Both (1) and (2) are local in position space; however, these two types of nongaussianity will have very different observational implications. The uncorrelated ansatz (2) for the primordial curvature perturbation can arise, for example, in models with preheating into light fields [19–21]. (See also [22, 23] for more discussion of nongaussianity from preheating and [24, 25] for another model where nongaussianity is generated at the end of inflation.)A useful quantity to consider is thebispectrum, B(k1,k2,k3), which is the 3-point correlation function of the Fourier transform of the primordial curvature perturbation(3)〈ζk1ζk2ζk3〉=(2π)3δ(k1+k3+k3)B(ki),
where ki≡|ki|. The delta function appearing in (3) reflects translational invariance and ensures that B(ki) depends on three momenta ki which form a triangle: k1+k2+k3=0. Rotational invariance implies that B(ki) is symmetric in its arguments.If we assume the ansatz (1) for the primordial curvature perturbation, then B(ki) has a very particular dependence on momenta; it peaks in the squeezed limit where one of the wavenumbers is much smaller than the remaining two (e.g., k1≪k2,k3). Such a bispectrum is referred to as having a squeezed shape. However, other shapes of bispectrum are worth considering. A bispectrum is referred to as “equilateral” if it peaks when k1=k2=k3 and “flattened” if it peaks when one of the wavenumbers is half the size of the remaining two (e.g., 2k1=k2=k3).Without assuming any specific form for the primordial curvature perturbation, such as (1) or (2), one may characterize an arbitrary bispectrum (3) by specifying its shape, running, and size [26]. As discussed above, the shape refers to the configuration of triangle on which B(ki) is maximal (squeeze, equilateral, or flattened). The running of the bispectrum refers to how the magnitude of B(ki) depends on the overall size of the triangle. For example, in the case of scale invariant fluctuations, the bispectrum must scale as B(λk1,λk2,λk3)=λ-6B(k1,k2,k3). Finally, the overall size of the bispectrum is often quantified by evaluating the magnitude of B(ki) on some fixed equilateral triangle. However, the skewness of the probability density function (defined later) might provide a better measure of the size of nongaussianity.Different types of nongaussian signatures are correlated with properties of the underlying inflation model. Let us first consider some examples with small running (by “small” here we refer to any model where the running of the bispectrum is proportional to slow-variation parameters or arises due to loop effects. This doesnot necessarily mean that such running cannot lead to interesting observational signatures, see; [27–30]).(1)
A large bispectrum of local shape, along with iso-curvature effects, is associated with models where multiple fields are light (or otherwise dynamically important) during inflation. Examples include the curvaton mechanism [31–34] or models with turning points along the inflationary trajectory [35–39]. The observational bound on local type nongaussianity, coming from the WMAP7 [6] data, is -10<fNLlocal<74 [40] at 95% confidence level. When combined with large scale structure (LSS) data the bound becomes somewhat stronger: -1<fNLlocal<65 [41].(2)
A large local bispectrumwithout any isocurvature fluctuations can only be produced by nonlocal inflation models [42–44]. For any single-field inflation model described by a local low-energy effective field theory, the results of [45] imply that the ratio of the 3-point correlation function to the square of the 2-point function must be of order of the spectral tilt, in the squeezed limit. Hence, it has been argued that a large squeezed bispectrum must be associated with the presence of multiple light degrees of freedom and hence iso-curvature effects. However, in [42–44] it was shown that single field nonlocal inflation models can produce a large squeezed bispectrum in the regime where the underlying scale of nonlocality is much larger than the Hubble scale during inflation. Such constructions evade the no-go theorem of [45] precisely because they violate the usual assumption of cluster decomposition. Moreover, models of this type are not subsumed by the general analysis of [46] since nonlocal field theories with infinitely many derivatives cannot be obtained in the regime of low-energy effective field theory. It is nevertheless sensible to study such constructions since they may be derived from ultra-violet (UV) complete frameworks, such as string field theory or p-adic string theory. See [47–49] for details concerning the underlying consistency of nonlocal field theories, and see [50] for a succinct review of nonlocal cosmology.(3)
A large equilateral bispectrum is typically associated with a small sound speed for the inflaton perturbations [26], such as in Dirac-Born-Infeld (DBI) inflation models [51, 52]. However, such a signature may also be obtained in multifield gelaton [53] or trapped inflation [54] models. The observational bound on equilateral type nongaussianity is -125<fNLequil<435 at 95% confidence level [55].(4)
A large flattened bispectrum is associated with nonvacuum initial conditions [26, 56–58]. (To our knowledge there is no explicit computation of the observational bound on flattened nongaussianity. In [56], a template (the enfolded model) was proposed. The analysis of [55] is sufficiently general to study this shape; however, they do not explicitly place bounds on fNLflat but instead constrain an alternative shape (the orthogonal model) which is a superposition of flattened and equilateral shapes.)If we relax the assumption that the bispectrum is close to scale invariant, then a much richer variety of nongaussian signatures is possible. For example, in models with sharp steps in the inflaton potential [59, 60] the bispectrum is large only for triangles with a particular characteristic size. We will refer to such a signature as a localized nongaussian feature. Localized nongaussianities are not well constrained by current observation but may be observable in future missions.Given the significant role that nongaussianity may play in discriminating between different models of the early universe, it is of crucial importance to explore and classify all possible consistent signatures for the bispectrum and other nongaussian statistics. Indeed, in this paper we will describe a new kind of signature—uncorrelated nongaussian features—which is predicted in a variety of simple and well-motivated models of inflation, but which has nevertheless been overlooked in previous literature.


### 1.2. Inflationary Particle Production

Recently, a new mechanism for generating cosmological perturbations during inflation was proposed [12]. This new mechanism, dubbed infra-red (IR) cascading, is qualitatively different from previous proposals (such as the curvaton or modulated fluctuations) in that it does not rely on the quantum vacuum fluctuations of some light scalar fields during inflation. Rather, the scenario involves the production of massive iso-curvature particles during inflation. These subsequently rescatter off the slow-roll condensate to generate bremsstrahlung radiation of light inflaton fluctuations (which induce curvature perturbations and temperature anisotropies in the usual manner). IR cascading can also be distinguished from previous mechanisms from the observational perspective: this new mechanism leads to novel features in both the spectrum and bispectrum.In principle, IR cascading may occur in any model where non-inflaton (iso-curvature) particles are produced during inflation. Models of this type have attracted considerable interest recently; examples have been studied where particle production occurs via parametric resonance [12, 13, 54, 61–66], as a result of a phase transition [19, 20, 67–74] or otherwise [75]. Recent interest in inflationary particle production has been stimulated by various considerations.(1)
Particle production arises naturally in a number of microscopically realistic models of inflation, including examples from string theory [54] and supersymmetric (SUSY) field theory [76]. In particular, inflationary particle production is a generic feature of open string inflation models [13], such as brane/axion monodromy [77–79].(2)
The energetic cost of producing particles during inflation has a dissipative effect on the dynamics of the inflaton. Particle production may therefore slow the motion of the inflaton, even on a steep potential. This gives rise to a new inflationary mechanism, calledtrapped inflation [54, 80, 81], which may circumvent some of the fine-tuning problems associated with standard slow-roll inflation. See [54] for an explicit string theory realization of trapped inflation and [81] for a generalization to higher-dimensional moduli spaces and enhanced symmetry loci. The idea of using dissipative dynamics to slow the motion of the inflaton is qualitatively similar to warm inflation [82] and also to the variant of natural inflation [83, 84] proposed recently by Anber and Sorbo [75].(3)
Observable features in the primordial power spectrum, generated by particle production and IR cascading, offer a novel example of the non-decoupling of high scale physics in the Cosmic Microwave Background (CMB) [61, 85–87]. In the most interesting examples, the produced particles are extremely massive for (almost) the entire history of the universe; however, their effect cannot be integrated out due to the nonadiabatic time dependence of the iso-inflaton mode functions during particle production. In [61] particle production during large field inflation was proposed as a possible probe of Planck-scale physics.In this paper we study in detail the impact of particle production and IR cascading on the observable primordial curvature perturbations. In order to illustrate the basic physics we focus on a very simple and general prototype model where the inflaton,ϕ, and iso-inflaton, χ, fields interact via the coupling(4)ℒint=-g22(ϕ-ϕ0)2χ2.
We expect, however, that our results will generalize in a straightforward way to more complicated models, such as fermion iso-inflaton fields or gauged interactions, wherein the physics of particle production and rescattering is essentially the same. Our result may also have implications for inflationary phase transitions, because spinodal decomposition can be interpreted as a kind of particle production, and similar bilinear interactions will induce rescattering effects.Scalar field interactions of type (4) have also been studied recently in connection with nonequilibrium Quantum Field Theory (QFT) [88–90], in particular with applications to the theory of preheating after inflation [91–95] and also moduli trapping [80, 81] at enhanced symmetry points. Although our focus is on particle production during inflation (as opposed to during preheating, after inflation) some of our results nevertheless have implications for preheating, moduli trapping, and also non-equilibrium QFT more generally. For example, in [12] analytical and numerical studies of rescattering and IR cascading during inflation made it possible to observe, for the first time, the dynamical approach to the turbulent scaling regime that was discovered in [96, 97].Particle production during inflation in model (4) leads to observable features in the primordial power spectrum, P(k). A number of recent studies have found evidence for localized features in P(k) that are incompatible with the simplest power-law model P(k)~kns-1 [62, 73, 98–110]. Although these observed features may simply be statistical anomalies (see, e.g., [111]), there remains the tantalizing possibility that they represent some new physics beyond the simplest slow-roll model. Upcoming polarization data may play an important role in distinguishing these possibilities [73]. In the meantime, it is interesting to determine the extent to which such features may be explained by a simple and well-motivated model such as (4). Moreover, because (4) is a complete microscopic model (as opposed to a phenomenological modification of the power spectrum) it is possible to predict a host of correlated observables, such as features in the scalar bispectrum and tensor power spectrum. Hence, it should be possible to robustly rule out (or confirm) the possibility that some massive iso-curvature particles were produced during inflation.If detected, features from particle production and IR cascading will provide a rare and powerful new window into the microphysics driving inflation. This scenario opens up the possibility of learning some details about how the inflaton couples to other particles in nature, as opposed to simply reconstructing the inflaton potential along the slow-roll trajectory. Moreover, due to the non-decoupling discussed above, features from particle production and IR cascading may probe new (beyond the standard model) physics at extraordinarily high energy scales.The outline of this paper is as follows. In Section2 we provide a brief, qualitative overview of the dynamics of particle production and IR cascading in model (4). In Section 3 we study in detail this same dynamics using fully nonlinear lattice field theory simulations. In Section 4 we provide an analytical theory of particle production and IR cascading in an expanding universe. A complimentary analytical analysis, using second-order cosmological perturbation theory, is provided in Section 5. In Section 6 we consider the observational constraints on inflationary particle production using a variety of data sets. In Section 7 we provide several explicit microscopic realizations of our scenario and study the implications of our observational constraints on models of string theory inflation, in particular brane monodromy. In Section 8 we quantify and characterize the nongaussianity generated by particle production and IR cascading. Finally, in Section 9, we conclude and discuss possible future directions.


## 1.1. Non-Gaussianity from Inflation

The possibility to discriminate between various inflationary scenarios has led to a recent surge of interest in computing and measuring nongaussian statistics. Although single field, slow roll models are known to produce negligible nongaussianity [16–18], there are now a variety of scenarios available in the literature which may predict an observable signature. Departures from gaussianity are often parametrized in the following form:(1)ζ(x)=ζg(x)+35fNL[ζg2(x)-〈ζg2(x)〉],
where ζ(x) is the primordial curvature perturbation, ζg(x) is a Gaussian random field, and fNL characterizes the degree of nongaussianity. The ansatz (1) is known as the “local” form of nongaussianity.Although the local ansatz (1) has received significant attention, it is certainly not the only well-motivated model for a nongaussian curvature perturbation. For example, the nongaussian part of ζ(x) need not be correlated with the gaussian part. Consider a primordial curvature perturbation of the form(2)ζ(x)=ζg(x)+FNL[χg(x)],
where FNL is some nonlinear (not necessarily quadratic) function, and χg(x) is a gaussian field which is uncorrelated with ζg(x). Both (1) and (2) are local in position space; however, these two types of nongaussianity will have very different observational implications. The uncorrelated ansatz (2) for the primordial curvature perturbation can arise, for example, in models with preheating into light fields [19–21]. (See also [22, 23] for more discussion of nongaussianity from preheating and [24, 25] for another model where nongaussianity is generated at the end of inflation.)A useful quantity to consider is thebispectrum, B(k1,k2,k3), which is the 3-point correlation function of the Fourier transform of the primordial curvature perturbation(3)〈ζk1ζk2ζk3〉=(2π)3δ(k1+k3+k3)B(ki),
where ki≡|ki|. The delta function appearing in (3) reflects translational invariance and ensures that B(ki) depends on three momenta ki which form a triangle: k1+k2+k3=0. Rotational invariance implies that B(ki) is symmetric in its arguments.If we assume the ansatz (1) for the primordial curvature perturbation, then B(ki) has a very particular dependence on momenta; it peaks in the squeezed limit where one of the wavenumbers is much smaller than the remaining two (e.g., k1≪k2,k3). Such a bispectrum is referred to as having a squeezed shape. However, other shapes of bispectrum are worth considering. A bispectrum is referred to as “equilateral” if it peaks when k1=k2=k3 and “flattened” if it peaks when one of the wavenumbers is half the size of the remaining two (e.g., 2k1=k2=k3).Without assuming any specific form for the primordial curvature perturbation, such as (1) or (2), one may characterize an arbitrary bispectrum (3) by specifying its shape, running, and size [26]. As discussed above, the shape refers to the configuration of triangle on which B(ki) is maximal (squeeze, equilateral, or flattened). The running of the bispectrum refers to how the magnitude of B(ki) depends on the overall size of the triangle. For example, in the case of scale invariant fluctuations, the bispectrum must scale as B(λk1,λk2,λk3)=λ-6B(k1,k2,k3). Finally, the overall size of the bispectrum is often quantified by evaluating the magnitude of B(ki) on some fixed equilateral triangle. However, the skewness of the probability density function (defined later) might provide a better measure of the size of nongaussianity.Different types of nongaussian signatures are correlated with properties of the underlying inflation model. Let us first consider some examples with small running (by “small” here we refer to any model where the running of the bispectrum is proportional to slow-variation parameters or arises due to loop effects. This doesnot necessarily mean that such running cannot lead to interesting observational signatures, see; [27–30]).(1)
A large bispectrum of local shape, along with iso-curvature effects, is associated with models where multiple fields are light (or otherwise dynamically important) during inflation. Examples include the curvaton mechanism [31–34] or models with turning points along the inflationary trajectory [35–39]. The observational bound on local type nongaussianity, coming from the WMAP7 [6] data, is -10<fNLlocal<74 [40] at 95% confidence level. When combined with large scale structure (LSS) data the bound becomes somewhat stronger: -1<fNLlocal<65 [41].(2)
A large local bispectrumwithout any isocurvature fluctuations can only be produced by nonlocal inflation models [42–44]. For any single-field inflation model described by a local low-energy effective field theory, the results of [45] imply that the ratio of the 3-point correlation function to the square of the 2-point function must be of order of the spectral tilt, in the squeezed limit. Hence, it has been argued that a large squeezed bispectrum must be associated with the presence of multiple light degrees of freedom and hence iso-curvature effects. However, in [42–44] it was shown that single field nonlocal inflation models can produce a large squeezed bispectrum in the regime where the underlying scale of nonlocality is much larger than the Hubble scale during inflation. Such constructions evade the no-go theorem of [45] precisely because they violate the usual assumption of cluster decomposition. Moreover, models of this type are not subsumed by the general analysis of [46] since nonlocal field theories with infinitely many derivatives cannot be obtained in the regime of low-energy effective field theory. It is nevertheless sensible to study such constructions since they may be derived from ultra-violet (UV) complete frameworks, such as string field theory or p-adic string theory. See [47–49] for details concerning the underlying consistency of nonlocal field theories, and see [50] for a succinct review of nonlocal cosmology.(3)
A large equilateral bispectrum is typically associated with a small sound speed for the inflaton perturbations [26], such as in Dirac-Born-Infeld (DBI) inflation models [51, 52]. However, such a signature may also be obtained in multifield gelaton [53] or trapped inflation [54] models. The observational bound on equilateral type nongaussianity is -125<fNLequil<435 at 95% confidence level [55].(4)
A large flattened bispectrum is associated with nonvacuum initial conditions [26, 56–58]. (To our knowledge there is no explicit computation of the observational bound on flattened nongaussianity. In [56], a template (the enfolded model) was proposed. The analysis of [55] is sufficiently general to study this shape; however, they do not explicitly place bounds on fNLflat but instead constrain an alternative shape (the orthogonal model) which is a superposition of flattened and equilateral shapes.)If we relax the assumption that the bispectrum is close to scale invariant, then a much richer variety of nongaussian signatures is possible. For example, in models with sharp steps in the inflaton potential [59, 60] the bispectrum is large only for triangles with a particular characteristic size. We will refer to such a signature as a localized nongaussian feature. Localized nongaussianities are not well constrained by current observation but may be observable in future missions.Given the significant role that nongaussianity may play in discriminating between different models of the early universe, it is of crucial importance to explore and classify all possible consistent signatures for the bispectrum and other nongaussian statistics. Indeed, in this paper we will describe a new kind of signature—uncorrelated nongaussian features—which is predicted in a variety of simple and well-motivated models of inflation, but which has nevertheless been overlooked in previous literature.


## 1.2. Inflationary Particle Production

Recently, a new mechanism for generating cosmological perturbations during inflation was proposed [12]. This new mechanism, dubbed infra-red (IR) cascading, is qualitatively different from previous proposals (such as the curvaton or modulated fluctuations) in that it does not rely on the quantum vacuum fluctuations of some light scalar fields during inflation. Rather, the scenario involves the production of massive iso-curvature particles during inflation. These subsequently rescatter off the slow-roll condensate to generate bremsstrahlung radiation of light inflaton fluctuations (which induce curvature perturbations and temperature anisotropies in the usual manner). IR cascading can also be distinguished from previous mechanisms from the observational perspective: this new mechanism leads to novel features in both the spectrum and bispectrum.In principle, IR cascading may occur in any model where non-inflaton (iso-curvature) particles are produced during inflation. Models of this type have attracted considerable interest recently; examples have been studied where particle production occurs via parametric resonance [12, 13, 54, 61–66], as a result of a phase transition [19, 20, 67–74] or otherwise [75]. Recent interest in inflationary particle production has been stimulated by various considerations.(1)
Particle production arises naturally in a number of microscopically realistic models of inflation, including examples from string theory [54] and supersymmetric (SUSY) field theory [76]. In particular, inflationary particle production is a generic feature of open string inflation models [13], such as brane/axion monodromy [77–79].(2)
The energetic cost of producing particles during inflation has a dissipative effect on the dynamics of the inflaton. Particle production may therefore slow the motion of the inflaton, even on a steep potential. This gives rise to a new inflationary mechanism, calledtrapped inflation [54, 80, 81], which may circumvent some of the fine-tuning problems associated with standard slow-roll inflation. See [54] for an explicit string theory realization of trapped inflation and [81] for a generalization to higher-dimensional moduli spaces and enhanced symmetry loci. The idea of using dissipative dynamics to slow the motion of the inflaton is qualitatively similar to warm inflation [82] and also to the variant of natural inflation [83, 84] proposed recently by Anber and Sorbo [75].(3)
Observable features in the primordial power spectrum, generated by particle production and IR cascading, offer a novel example of the non-decoupling of high scale physics in the Cosmic Microwave Background (CMB) [61, 85–87]. In the most interesting examples, the produced particles are extremely massive for (almost) the entire history of the universe; however, their effect cannot be integrated out due to the nonadiabatic time dependence of the iso-inflaton mode functions during particle production. In [61] particle production during large field inflation was proposed as a possible probe of Planck-scale physics.In this paper we study in detail the impact of particle production and IR cascading on the observable primordial curvature perturbations. In order to illustrate the basic physics we focus on a very simple and general prototype model where the inflaton,ϕ, and iso-inflaton, χ, fields interact via the coupling(4)ℒint=-g22(ϕ-ϕ0)2χ2.
We expect, however, that our results will generalize in a straightforward way to more complicated models, such as fermion iso-inflaton fields or gauged interactions, wherein the physics of particle production and rescattering is essentially the same. Our result may also have implications for inflationary phase transitions, because spinodal decomposition can be interpreted as a kind of particle production, and similar bilinear interactions will induce rescattering effects.Scalar field interactions of type (4) have also been studied recently in connection with nonequilibrium Quantum Field Theory (QFT) [88–90], in particular with applications to the theory of preheating after inflation [91–95] and also moduli trapping [80, 81] at enhanced symmetry points. Although our focus is on particle production during inflation (as opposed to during preheating, after inflation) some of our results nevertheless have implications for preheating, moduli trapping, and also non-equilibrium QFT more generally. For example, in [12] analytical and numerical studies of rescattering and IR cascading during inflation made it possible to observe, for the first time, the dynamical approach to the turbulent scaling regime that was discovered in [96, 97].Particle production during inflation in model (4) leads to observable features in the primordial power spectrum, P(k). A number of recent studies have found evidence for localized features in P(k) that are incompatible with the simplest power-law model P(k)~kns-1 [62, 73, 98–110]. Although these observed features may simply be statistical anomalies (see, e.g., [111]), there remains the tantalizing possibility that they represent some new physics beyond the simplest slow-roll model. Upcoming polarization data may play an important role in distinguishing these possibilities [73]. In the meantime, it is interesting to determine the extent to which such features may be explained by a simple and well-motivated model such as (4). Moreover, because (4) is a complete microscopic model (as opposed to a phenomenological modification of the power spectrum) it is possible to predict a host of correlated observables, such as features in the scalar bispectrum and tensor power spectrum. Hence, it should be possible to robustly rule out (or confirm) the possibility that some massive iso-curvature particles were produced during inflation.If detected, features from particle production and IR cascading will provide a rare and powerful new window into the microphysics driving inflation. This scenario opens up the possibility of learning some details about how the inflaton couples to other particles in nature, as opposed to simply reconstructing the inflaton potential along the slow-roll trajectory. Moreover, due to the non-decoupling discussed above, features from particle production and IR cascading may probe new (beyond the standard model) physics at extraordinarily high energy scales.The outline of this paper is as follows. In Section2 we provide a brief, qualitative overview of the dynamics of particle production and IR cascading in model (4). In Section 3 we study in detail this same dynamics using fully nonlinear lattice field theory simulations. In Section 4 we provide an analytical theory of particle production and IR cascading in an expanding universe. A complimentary analytical analysis, using second-order cosmological perturbation theory, is provided in Section 5. In Section 6 we consider the observational constraints on inflationary particle production using a variety of data sets. In Section 7 we provide several explicit microscopic realizations of our scenario and study the implications of our observational constraints on models of string theory inflation, in particular brane monodromy. In Section 8 we quantify and characterize the nongaussianity generated by particle production and IR cascading. Finally, in Section 9, we conclude and discuss possible future directions.


## 2. Overview and Summary of the Mechanism

In this section we provide a brief overview of the dynamics of particle production and IR cascading in model (4) and also summarize the resulting observational signatures. In the remainder of this paper we will flesh out the details of this mechanism with analytical and numerical calculations.We consider the following model:(5)S=∫d4x-g[Mp22R-12(∂ϕ)2-V(ϕ)-12(∂χ)2-g22(ϕ-ϕ0)2χ2],
where R is the Ricci curvature constructed from the metric gμν, ϕ is the inflaton field, and χ is the iso-inflaton. As usual, we assume a flat FRW space-time with scale factor a(t)(6)ds2≡gμνdxμdxν=-dt2+a2(t)dx2
and employ the reduced Planck mass Mp≅2.43×1018GeV. We leave the potential V(ϕ) driving inflation unspecified for assuming that it is sufficiently flat in the usual sense, that is, ϵ≪1, |η|≪1 where(7)ϵ≡Mp22(V′V)2,η≡Mp2V′′V
are the usual slow-roll parameters.Note that one might wish to supplement (5) by its supersymmetric completion in order to protect the flatness of the inflaton potential from large radiative corrections coming from loops of the χ field. We expect that our results will carry over in a straightforward way to SUSY models and also to more complicated scenarios such as higher spin iso-inflaton fields and (possibly) inflationary phase transitions.The coupling(g2/2)(ϕ-ϕ0)2χ2 in (5) is introduced to ensure that the iso-inflaton field can become instantaneously massless at some point ϕ=ϕ0 along the inflaton trajectory (which we assume occurs during the observable range of e-foldings of inflation). At this moment χ particles will be produced by quantum effects.Let us first consider the homogeneous dynamics of the inflaton field,ϕ(t). Near the point ϕ=ϕ0 we can generically expand(8)ϕ(t)≅ϕ0+vt,
where v≡ϕ̇(0), and we have arbitrarily set the origin of time so that t=0 corresponds to the moment when ϕ=ϕ0. (We are, of course, assuming that ϕ̇(0)≠0.) The interaction (4) induces an effective (time varying) mass for the χ particles of the form(9)mχ2=g2(ϕ-ϕ0)2≅k⋆4t2,
where we have defined the characteristic scale(10)k⋆=g|v|.
It is straightforward to verify that the simple expression (9) will be a good approximation for (H|t|)-1≲𝒪(ϵ,η) which, in most models, will be true for the entire observable 60 e-foldings of inflation.Note that, without needing to specify the background inflationary potentialV(ϕ), we can write the ratio k⋆/H as(11)k⋆H=g2π𝒫ζ1/2,
where 𝒫ζ1/2=5×10-5 is the usual amplitude of the vacuum fluctuations from inflation. In this work we assume that k⋆>H which is easily satisfied for reasonable values of the coupling g2>10-7. In particular, for g2~0.1 we have k⋆/H~30.The scenario we have in mind is the following. Inflation starts at some field valueϕ>ϕ0 and the inflaton rolls toward the point ϕ=ϕ0. Initially, the iso-inflaton field is extremely massive mχ≫H and hence it stays pinned in the vacuum, χ=0, and does not contribute to superhorizon curvature fluctuations. Eventually, at t=0, the inflaton rolls through the point ϕ=ϕ0 where mχ=0 and χ particles are produced. To describe this burst of particle production one must solve for the following equation for the χparticle mode functions in an expanding universe:(12)χ̈k+3Hχ̇k+[k2a2+k⋆4t2]χk=0.
Equations of this type are well-studied in the context of preheating after inflation [92] and moduli trapping [80]. The initial conditions for (12) should be chosen to ensure that the q-number field χ is in the adiabatic vacuum in the asymptotic past (see Sections 3 and 4 for more details). In the regime k⋆>H particle production is fast compared to the expansion time and one can solve (12) very accurately for the occupation number of the created χ particles(13)nk=e-πk2/k⋆2.
Very quickly after the moment t=0, within a time Δt~k⋆-1≪H-1, these produced χ particles become nonrelativistic (mχ>H), and their number density starts to dilute as a-3.Following the initial burst of particle production there are two distinct physical effects which take place. First, the energetic cost of producing the gas of massive out-of-equilibriumχ particles drains energy from the inflaton condensate, forcing ϕ̇ to drop abruptly. This velocity dip is the result of the backreaction of the produced χ fluctuations on homogeneous condensate ϕ(t). The second physical effect is that the produced massive χ particles rescatter off the condensate via the diagram in Figure 1 and emit bremsstrahlung radiation of light inflaton fluctuations (particles).Figure 1
Rescattering diagram.Backreaction and rescattering leave distinct imprints in the observable cosmological perturbations. Let us first discuss the impact of backreaction. In Figure2 we plot the velocity dip resulting from the backreaction of the produced χ particle on the homogeneous inflaton condensate ϕ(t). From this figure we see that the quantity ϕ̈/(Hϕ̇) becomes large in the dip. This violation of slow roll is a transient effect; at late times the produced χ particles become extremely massive and their number density dilutes as a-3.Figure 2
|ϕ̇|/(Mpm) plotted against mt for g2=0.1 (where m=V,ϕϕ is the effective inflaton mass). Time t=0 corresponds to the moment when ϕ=ϕ0 and χ particles are produced copiously. The solid red line is the lattice field theory result taking into account the full dynamics of rescattering and IR cascading while the dashed blue line is the result of a mean-field theory treatment which ignores rescattering [64]. The dot-dashed black line is the inflationary trajectory in the absence of particle creation.One can understand the temporary slowing down of the inflaton from an analytical perspective. Backreaction is taken into account using the mean-field equation(14)ϕ̈+3Hϕ̇+V,ϕ+g2(ϕ-ϕ0)〈χ2〉=0,
where the vacuum average is computed following [80, 92] (15)〈χ2〉≅nχa-3g|ϕ-ϕ0|.
In (14) we have implicitly assumed that the usual Coleman-Weinberg corrections to the inflaton potential have already been absorbed into V(ϕ), hence the vacuum average 〈χ2〉 should include only the effects of nonadiabatic particle production. (Here nχ=∫(d3k/(2π)3)nk~k⋆3 is the total number density of produced χ particles, and the factor a-3 reflects the usual volume dilution of non-relativistic matter.) In Figure 2 we have plotted the solution of (14) along with the exact result obtained from lattice field theory simulations, illustrating the accuracy of this simple treatment.Using the mean-field approach, one finds that the transient violation of slow roll leads to a “ringing pattern” (damped oscillations) in the power spectrumPϕ(k)=(k3/2π2)|δϕk|2 of inflaton fluctuations [64]. This ringing pattern is localized around wavenumbers which left the horizon at the moment when particle production occurred. The effect is very much analogous to Fresnel diffraction at a sharp edge.The second physical effect, rescattering, was considered for the first time in the context of inflationary particle production in [12]. Figure 1 illustrates the dominant process: bremsstrahlung emission of long-wavelength δϕ fluctuations from rescattering of the produced χ particles off the condensate. The time scale for such processes is set by the microscopic scale, k⋆-1, and is thus very short compared to the expansion time, H-1. Moreover, the production of inflaton fluctuations δϕ deep in the infrared (IR) is extremely energetically inexpensive, since the inflaton is very nearly massless. The combination of the short time scale for rescattering and the energetic cheapness of radiating IR δϕ leads to a rapid build-up of power in long wavelength inflaton modes: IR cascading. This effect leads to a bump-like feature in the power spectrum of inflaton fluctuations, very different from the ringing pattern associated with backreaction. The bump-like feature from rescattering dominates over the ringing pattern from backreaction for all values of parameters.In [12] model (5) was studied using lattice field theory simulations, without neglecting any physical processes (that is to say that full nonlinear structure of the theory, including backreaction and rescattering effects, was accounted for consistently). However, this same dynamics can be understood analytically by solving the equation for the inflaton fluctuations δϕ in the approximation that all interactions are neglected, except for the diagram in Figure 1. The appropriate equation is(16)δϕ̈+3Hδϕ̇-∇⃗2a2δϕ+V,ϕϕδϕ≅-g2[ϕ(t)-ϕ0]χ2.
See [14] for a detailed analytical theory. The solution of (16) may be split into two parts: the solution of the homogeneous equation and the particular solution which is due to the source term. The former simply corresponds to the usual scale invariant quantum vacuum fluctuations from inflation. The particular solution, on the other hand, corresponds to inflaton fluctuations generated by rescattering. The abrupt growth of χ inhomogeneities at t=0 sources the particular solution and generates inflaton fluctuations which subsequently cross the horizon and freeze in.As mentioned earlier, rescattering generates a bump-like contribution to the primordial power spectrum of the curvature perturbations. To good approximation this may be described by a simple semi analytic fitting function(17)P(k)=As(kk0)ns-1+AIR(πe3)3/2(kkIR)3e-(π/2)(k/kIR)2,
where the first term corresponds to the usual vacuum fluctuations from inflation (with amplitude As and spectral index ns) while the second term corresponds to the bump-like feature from particle production and IR cascading. The amplitude of this feature (AIR) depends on g2 while the location (kIR) depends on ϕ0.In [13] the simple fitting function (17) was used to place observational constraints on inflationary particle production using a variety of cosmological data sets. Current data are consistent with rather large spectral distortions of the type (17). Features as large as 𝒪(10%) of the usual scale-invariant fluctuations from inflation are allowed, in the case that kIR falls within the range of scales relevant for CMB experiments. Such a feature corresponds to a realistic coupling g2~0.01. Even larger values of g2 are allowed if the feature is localized on smaller scales. In Figure 3 we have illustrated the primordial power spectrum in model (5) for a representative choice of parameters. We also plot the CMB angular Temperature-Temperature (TT) power spectrum for the same parameters.(a) shows a sample bump in the power spectrum with amplitudeAIR=2.5×10-10 which corresponds to a coupling g2~0.01. The feature is located at kIR=0.01Mpc-1. This example represents a distortion of 𝒪(10%) as compared to the usual vacuum fluctuations and is consistent with the data at 2σ. (b) shows the CMB angular TT power spectrum for this example, illustrating that the distortion shows up mostly in the first peak.
(a)(b)The prototype model (5) may be realized microscopically in a variety of different particle physics frameworks. In particular, particle production is a rather generic feature of open string inflation models [13] where the inflaton, ϕ, has a geometrical interpretation as the position of some mobile D-brane. In this context the iso-inflaton, χ, corresponds to a low-lying open string excitation which is stretched between the mobile inflationary brane and any other (spectator) branes which inhabit the compactification volume. If the inflationary and spectator branes become coincident during inflation, then the symmetry of the system is enhanced [80] and some low-lying stretched string states will become instantaneously massless, mimicking interaction (4) (see also [54]). An explicit realization of this scenario is provided by brane/axion monodromy models [77–79]. Our observational constraints on inflationary particle production may be used to place bounds on parameters of the underlying string model [13].The bump-like feature inP(k), illustrated in Figure 3, must be associated with a nongaussian feature in the bispectrum [12, 14]. Indeed, it is evident already from inspection of (16) that the inflaton fluctuations generated by rescattering are significantly nongaussian; the particular solution of (16) is bi-linear in the gaussian field χ. The nongaussian signature from IR cascading is rather novel. The nongaussian part of ζ is uncorrelated with the gaussian part. Moreover, the bispectrum B(ki) is very far from scale invariant; it peaks strongly for triangles with a characteristic size ~kIR, corresponding to the location of the bump in the power spectrum (17). The shape of the bispectrum therefore depends sensitively on the size of the triangle and is not well described by any of the templates that have been proposed in the literature to date.The magnitude of this new kind of nongaussianity may be quite large. To quantify the effect it is useful to introduce the probability density function (PDF),P(ζ), which is the probability that the curvature perturbation has a fluctuation of size ζ. If we define the central moments of the PDF as(18)〈ζn〉=∫ζnP(ζ)dζ,
then a useful measure of nongaussianity is the dimensionless skewness of the PDF, defined by(19)Ŝ3≡〈ζ3〉c〈ζ2〉3/2,
where the subscript c indicates that only the connected part of the correlator should be included. The skewness Ŝ3 encodes information about the bispectrum B(ki) integrated over all size and shape configurations and thus provides a meaningful single number to compare the nongaussianity of inflation models which may have very different shapes or running [27].If we chooseg2~0.01 (which is compatible with observation for all values of ϕ0), then model (5) produces the same value of Ŝ3 as a local model (1) with |fNL|~102. This large value suggests that nongaussianity from particle production during inflation may be observable in future missions.Depending on model parameters, the nongaussian features predicted by model (5) may lead to a rich variety of observable consequences for the CMB or Large Scale Structure (LSS). The phenomenology of this model is quite different from other constructions that have been proposed to obtain large nongaussianity from inflation. However, the underlying microscopic description (5) is extremely simple and, indeed, rather generic from the low-energy perspective. Explicit realizations of interaction (4) have been obtained from string theory and SUSY. Moreover, in order to obtain an observable signature it was not necessary to fine-tune the inflationary trajectory or appeal to re-summation of an infinite series of high-dimension operators.


## 3. Numerical Study of Rescattering and IR Cascading

### 3.1. HLattice Simulations

In this section we study numerically the creation ofδϕ fluctuations by rescattering of the produced χ particles off the condensate ϕ(t) in model (5). To this end, we have written a new lattice field theory code, HLattice [112], for simulating the interactions of scalar fields in a cosmological setting. HLattice can be used to simulate the dynamics of any number of interacting scalar fields with arbitrary scalar potential and metric on field space [113]. We solve the Klein-Gordon equations for the scalar field dynamics in an expanding FRW space-time and also solve the Friedmann equation self-consistently for the scale factor, a(t). Since the production of long wavelength δϕ modes is so energetically inexpensive, a major requirement for successfully capturing this effect is respecting energy conservation to very high accuracy. HLattice conserves energy with an accuracy of order ~10-8, as compared to 10-3–10-5, which has been obtained using previous codes such as DEFROST [114] or LATTICEASY [115]. A minimum accuracy of order 10-4 is required for the problem at hand.The box size of our5123 simulations corresponds to a comoving scale which is initially (20/2π)~3 times the horizon size H-1, while k⋆≅60gH. We run our simulations for roughly 3 e-foldings from the initial moment t=0 when the χ particles are produced although a single e-folding would have been sufficient to capture the effect. For the sake of illustration, we have chosen the standard chaotic inflation potential V=m2ϕ2/2 with m=10-68πMp for our numerical analysis. However, our results do not depend sensitively on the choice of background inflation model. (The model independence of our result arises simply because all the dynamics of rescattering and IR cascading occurs within a single e-folding from the moment when ϕ=ϕ0. Over such a short time it will always be a good approximation to expand ϕ(t)≅ϕ0+vt. Hence the dependence on the background dynamics arises only through v=ϕ̇(0) which is determined by the Hubble scale and the observed amplitude of curvature perturbations. This claim of model independence is born out by explicit analytical calculations in the next section.) We have considered both ϕ0=28πMp and ϕ0=3.28πMp and also three different values of the coupling constant: g2=0.01,0.1,1. As expected, the coupling g2 determines the magnitude of the effect while ϕ0 simply shifts the location of the power spectrum feature. For this choice of inflationary potential, the choice ϕ0=3.28πMp corresponds to putting the feature on scale slightly smaller than today’s horizon. On the other hand, ϕ0=28πMp corresponds to placing the feature on scales much smaller than those probed by the CMB (we considered this case in order to be able to directly contrast our results with [64]).In order to capture the quantum production ofχ particles using classical lattice simulations, we start our numerical evolution very shortly after particle production has occurred, when the χk modes are nearly adiabatic, but before any significant inflaton fluctuations have been produced. In practice, this corresponds to initializing the simulation at a time tinitial=𝒪(k⋆-1). The initial conditions for the modes χk(t) are given by the usual Bogoliubov computation [80, 92]. These are chosen to reproduce the occupation number nk=e-πk2/k⋆2, while ensuring that the source term for the δϕ fluctuations is turned on smoothly at the initial time. As long as the initial conditions are chosen appropriately, our results are not sensitive to the choice of tinitial.At the initial time, the occupation numbers in the inflaton and iso-inflaton fluctuations are small. However, very quickly the massiveχ particles are diluted away by the expansion of the universe, and the occupation number of the produced IR δϕ fluctuations grows large compared to unity. Thus, classical lattice field theory simulations are sufficient to capture the late-time dynamics. (In the next section we will provide a quantum mechanical treatment of the dynamics of particle production and IR cascading during inflation, which will serve as an a posteriori justification for our classical lattice calculation.)Our approach is very similar to the methodology that has been employed successfully in studies of preheating after inflation for many years [114, 115]. In that case the initial fluctuations of the fields are chosen to reproduce the exact behaviour of the quantum correlation functions. The occupation numbers of the fields are small at the initial time. However, these grow rapidly as a result of the preheating instability, and classical simulations are sufficient to capture the late-time dynamics.


### 3.2. Numerical Results

We have studied the fully nonlinear dynamics ofχ particle production and the subsequent interactions of the produced χ with the inflaton field in model (5), as described above. We are interested in the power spectrum of the inflaton fluctuations(20)Pϕ(k)=k32π2|δϕk|2.
This contains a contribution coming from the usual quantum vacuum fluctuations from inflation that is close to the usual power-law form kns-1 on large scales. Such a contribution would be present even in the absence of particle production and is not particularly interesting for us. In order to isolate the effects of rescattering we have subtracted off this component in Figures 4, 5, and 6. In all cases we have normalized Pϕ to the amplitude of the usual vacuum fluctuations from inflation, H2/(2π)2.Figure 4
The power spectrum of inflaton modes induced by rescattering (normalized to the usual vacuum fluctuations) as a function ofln(k/k⋆), plotted for three representative time steps in the evolution, showing the cascading of power into the IR. For each time step we plot the analytical result (the solid line) and the data points obtained using lattice field theory simulations (diamonds). The time steps correspond to the following values of the scale factor: a=1.03,1.04,2.20 (where a=1 at the moment when ϕ=ϕ0). By this time the amplitude of fluctuations is saturated due to the expansion of the universe. The vertical lines show the range of scales from our lattice simulation.Figure 5
The dependence of the power spectrumPϕ on the coupling g2. The three curves correspond to Pϕ for g2=0.01,0.1,1, evaluated at a fixed value of the scale factor, a=2.20. We see that even for small values of g2 the inflaton modes induced by rescattering constitute a significant fraction of the usual vacuum fluctuations after only a single e-folding.Figure 6
The power spectrum of inflaton modes induced by rescattering (normalized to the usual vacuum fluctuations) as a function ofln(k/k⋆), plotted for three representative time steps in the late-time evolution. This figure illustrates the final stages of IR cascading; we see the peak of the bump-like feature slide to k~e-3k⋆, at which point the associated mode functions δϕk have crossed the horizon and become frozen. At later times in the evolution the peak of the feature and also the IR tail (~k3) remain fixed. Modes associated with the UV end of the spectrum are still inside the horizon and continue to evolve as δϕk~a-1, which explains the damping of the k>e-2k⋆ part of the spectrum. For each time step we plot the analytical result (the solid line) and the data points obtained using lattice field theory simulations (diamonds).Figure4 shows time evolution of the power in the inflaton fluctuations generated by rescattering, for three different time steps early in the evolution. This figure illustrates how multiple rescatterings lead to a dynamical cascading of power into the IR. To illustrate the magnitude of this effect, the horizontal yellow line corresponds to the amplitude of the usual vacuum fluctuations from inflation. For g2≳0.06, the fluctuations from rescattering come to dominate over the vacuum fluctuations within a single e-folding. In Figure 5 we illustrate how the magnitude of the spectral distortion depends on the coupling, g2. (The apparent change in the location of the feature for different values of g2 arises because we are plotting the power spectrum as a function of ln(k/k⋆) and k⋆ depends on g.)At late times, the IR portion of the power spectrum illustrated in Figure4 will remain fixed since the modes associated with these scales have gone outside the horizon and become frozen. On the other hand, the UV portion of this curve corresponds to modes that are still inside the horizon, hence we expect δϕk~a-1 and the UV tail of the power spectrum should damp as a-2, due to the Hubble expansion. We observe precisely this behaviour in our lattice field theory simulations, and this is illustrated in Figure 6, which displays the dynamics of IR cascading over a much longer time scale.Within a fewe-foldings from the time of particle production, the entire bump-like feature from IR cascading becomes frozen outside the horizon. At this point the fluctuations have become classical, large-scale adiabatic density perturbations and are observable in the present epoch (presuming that ϕ=ϕ0 occurs during the observable range of e-foldings). In Figure 3 we have illustrated this bump-like feature in both the primordial power spectrum and angular TT spectrum, for a representative choice of parameters.


### 3.3. Backreaction Effects

As discussed previously, the production ofχ fluctuations at t=0 backreacts on the homogeneous ϕ(t) causing a transient violation of slow roll. We can study this backreaction numerically, by averaging the inhomogeneous field ϕ̇(t,x) over the simulation box. The result is plotted in Figure 2. We have also plotted the analytical solution of the mean field (14), showing that this agrees with the exact numerical result.The dynamics illustrated in Figure2 is easy to understand physically. The production of χ particles at t=0 drains kinetic energy from the condensate and hence ϕ̇ must decrease abruptly. However, within a few e-foldings of the moment t=0, the produced iso-inflaton particles become extremely massive and are diluted by the expansion as a-3. At late times the inflaton velocity ϕ̇ must tend to the slow-roll value. Notice that the velocity ϕ̇ including backreaction effects is not changed significantly, as compared to the usual slow-roll result. This illustrates the energetic cheapness of particle production and IR cascading in model (5).The transient violation of slow roll illustrated in Figure2 is expected to induce a ringing pattern in the vacuum fluctuations from inflation [64]. This effect is accounted for automatically in our HLattice simulations. However, we would like to disentangle the effect of backreaction on the cosmological fluctuations from the effect of rescattering. This will be useful in order to compare the relative importance of different physical processes and also to guide our analytical efforts in the next section. To this end, we consider the evolution of the curvature perturbation on comoving hypersurfaces, ℛ. In linear theory the equation for the Fourier modes ℛk is well known(21)ℛk′′+2z′zℛk′+k2ℛk=0.
Here the prime denotes derivatives with respect to conformal time τ=∫(dt/a) and z≡aϕ̇/H. Equation (21) is only strictly valid in the absence of entropy perturbations. However, in our case the χ field is extremely massive mχ2≫H2 for nearly the entire duration of inflation, hence one may expect that direct iso-curvature contributions to ℛ are small. We have solved (21) numerically. In order to take backreaction effects into account we compute the dynamics of z(t)=a(t)ϕ̇(t)/H(t) by averaging over our HLattice simulation box. Next, we solve (21) given this background evolution and compute the power spectrum(22)Pℛ(k)=k32π2|ℛk|2.
The result is very close to the usual power-law form kns-1, with small superposed oscillations resulting from the transient violation of slow roll; see Figure 7. In order to make the ringing pattern more visible, we have subtracted off the usual (nearly) scale-invariant result which would be obtained in the absence of particle production. For comparison, we also plot the bump-like feature from rescattering and IR cascading. This latter contribution was obtained using the results for Pϕ(k) from the previous subsection and the naive formula ℛ~(H/ϕ̇)δϕ (so that Pℛ~(2ϵMp2)-1Pϕ).Figure 7
A comparison of curvature fluctuations from different physical effects. The dashed blue line is the usual (nearly) scale invariant vacuum fluctuations from inflation. The red solid line is the bump-like feature induced by rescattering and IR cascading. The dotted blue line is the ringing pattern resulting from the momentary slowing-down of the inflaton (computed using the mean field approach of [64]). The vertical lines show aH at the beginning of particle production and after ~3e-foldings. This figure clearly illustrates the dominance of IR cascading over backreaction effects. For illustration we have taken g2=0.1, but the dominance is generic for all values of the coupling.From Figure7 we see that IR cascading has a much more significant impact on the observable curvature fluctuations than does backreaction. Indeed, for g2=0.1 the transient violation of slow roll yields an order 10-2 correction to the vacuum fluctuations while the correction from IR cascading is of order 101. This dominance is generic for all values of the coupling. Thus, in developing an analytical theory of particle production during inflation, it is a very good approximation to completely ignore backreaction effects.


## 3.1. HLattice Simulations

In this section we study numerically the creation ofδϕ fluctuations by rescattering of the produced χ particles off the condensate ϕ(t) in model (5). To this end, we have written a new lattice field theory code, HLattice [112], for simulating the interactions of scalar fields in a cosmological setting. HLattice can be used to simulate the dynamics of any number of interacting scalar fields with arbitrary scalar potential and metric on field space [113]. We solve the Klein-Gordon equations for the scalar field dynamics in an expanding FRW space-time and also solve the Friedmann equation self-consistently for the scale factor, a(t). Since the production of long wavelength δϕ modes is so energetically inexpensive, a major requirement for successfully capturing this effect is respecting energy conservation to very high accuracy. HLattice conserves energy with an accuracy of order ~10-8, as compared to 10-3–10-5, which has been obtained using previous codes such as DEFROST [114] or LATTICEASY [115]. A minimum accuracy of order 10-4 is required for the problem at hand.The box size of our5123 simulations corresponds to a comoving scale which is initially (20/2π)~3 times the horizon size H-1, while k⋆≅60gH. We run our simulations for roughly 3 e-foldings from the initial moment t=0 when the χ particles are produced although a single e-folding would have been sufficient to capture the effect. For the sake of illustration, we have chosen the standard chaotic inflation potential V=m2ϕ2/2 with m=10-68πMp for our numerical analysis. However, our results do not depend sensitively on the choice of background inflation model. (The model independence of our result arises simply because all the dynamics of rescattering and IR cascading occurs within a single e-folding from the moment when ϕ=ϕ0. Over such a short time it will always be a good approximation to expand ϕ(t)≅ϕ0+vt. Hence the dependence on the background dynamics arises only through v=ϕ̇(0) which is determined by the Hubble scale and the observed amplitude of curvature perturbations. This claim of model independence is born out by explicit analytical calculations in the next section.) We have considered both ϕ0=28πMp and ϕ0=3.28πMp and also three different values of the coupling constant: g2=0.01,0.1,1. As expected, the coupling g2 determines the magnitude of the effect while ϕ0 simply shifts the location of the power spectrum feature. For this choice of inflationary potential, the choice ϕ0=3.28πMp corresponds to putting the feature on scale slightly smaller than today’s horizon. On the other hand, ϕ0=28πMp corresponds to placing the feature on scales much smaller than those probed by the CMB (we considered this case in order to be able to directly contrast our results with [64]).In order to capture the quantum production ofχ particles using classical lattice simulations, we start our numerical evolution very shortly after particle production has occurred, when the χk modes are nearly adiabatic, but before any significant inflaton fluctuations have been produced. In practice, this corresponds to initializing the simulation at a time tinitial=𝒪(k⋆-1). The initial conditions for the modes χk(t) are given by the usual Bogoliubov computation [80, 92]. These are chosen to reproduce the occupation number nk=e-πk2/k⋆2, while ensuring that the source term for the δϕ fluctuations is turned on smoothly at the initial time. As long as the initial conditions are chosen appropriately, our results are not sensitive to the choice of tinitial.At the initial time, the occupation numbers in the inflaton and iso-inflaton fluctuations are small. However, very quickly the massiveχ particles are diluted away by the expansion of the universe, and the occupation number of the produced IR δϕ fluctuations grows large compared to unity. Thus, classical lattice field theory simulations are sufficient to capture the late-time dynamics. (In the next section we will provide a quantum mechanical treatment of the dynamics of particle production and IR cascading during inflation, which will serve as an a posteriori justification for our classical lattice calculation.)Our approach is very similar to the methodology that has been employed successfully in studies of preheating after inflation for many years [114, 115]. In that case the initial fluctuations of the fields are chosen to reproduce the exact behaviour of the quantum correlation functions. The occupation numbers of the fields are small at the initial time. However, these grow rapidly as a result of the preheating instability, and classical simulations are sufficient to capture the late-time dynamics.


## 3.2. Numerical Results

We have studied the fully nonlinear dynamics ofχ particle production and the subsequent interactions of the produced χ with the inflaton field in model (5), as described above. We are interested in the power spectrum of the inflaton fluctuations(20)Pϕ(k)=k32π2|δϕk|2.
This contains a contribution coming from the usual quantum vacuum fluctuations from inflation that is close to the usual power-law form kns-1 on large scales. Such a contribution would be present even in the absence of particle production and is not particularly interesting for us. In order to isolate the effects of rescattering we have subtracted off this component in Figures 4, 5, and 6. In all cases we have normalized Pϕ to the amplitude of the usual vacuum fluctuations from inflation, H2/(2π)2.Figure 4
The power spectrum of inflaton modes induced by rescattering (normalized to the usual vacuum fluctuations) as a function ofln(k/k⋆), plotted for three representative time steps in the evolution, showing the cascading of power into the IR. For each time step we plot the analytical result (the solid line) and the data points obtained using lattice field theory simulations (diamonds). The time steps correspond to the following values of the scale factor: a=1.03,1.04,2.20 (where a=1 at the moment when ϕ=ϕ0). By this time the amplitude of fluctuations is saturated due to the expansion of the universe. The vertical lines show the range of scales from our lattice simulation.Figure 5
The dependence of the power spectrumPϕ on the coupling g2. The three curves correspond to Pϕ for g2=0.01,0.1,1, evaluated at a fixed value of the scale factor, a=2.20. We see that even for small values of g2 the inflaton modes induced by rescattering constitute a significant fraction of the usual vacuum fluctuations after only a single e-folding.Figure 6
The power spectrum of inflaton modes induced by rescattering (normalized to the usual vacuum fluctuations) as a function ofln(k/k⋆), plotted for three representative time steps in the late-time evolution. This figure illustrates the final stages of IR cascading; we see the peak of the bump-like feature slide to k~e-3k⋆, at which point the associated mode functions δϕk have crossed the horizon and become frozen. At later times in the evolution the peak of the feature and also the IR tail (~k3) remain fixed. Modes associated with the UV end of the spectrum are still inside the horizon and continue to evolve as δϕk~a-1, which explains the damping of the k>e-2k⋆ part of the spectrum. For each time step we plot the analytical result (the solid line) and the data points obtained using lattice field theory simulations (diamonds).Figure4 shows time evolution of the power in the inflaton fluctuations generated by rescattering, for three different time steps early in the evolution. This figure illustrates how multiple rescatterings lead to a dynamical cascading of power into the IR. To illustrate the magnitude of this effect, the horizontal yellow line corresponds to the amplitude of the usual vacuum fluctuations from inflation. For g2≳0.06, the fluctuations from rescattering come to dominate over the vacuum fluctuations within a single e-folding. In Figure 5 we illustrate how the magnitude of the spectral distortion depends on the coupling, g2. (The apparent change in the location of the feature for different values of g2 arises because we are plotting the power spectrum as a function of ln(k/k⋆) and k⋆ depends on g.)At late times, the IR portion of the power spectrum illustrated in Figure4 will remain fixed since the modes associated with these scales have gone outside the horizon and become frozen. On the other hand, the UV portion of this curve corresponds to modes that are still inside the horizon, hence we expect δϕk~a-1 and the UV tail of the power spectrum should damp as a-2, due to the Hubble expansion. We observe precisely this behaviour in our lattice field theory simulations, and this is illustrated in Figure 6, which displays the dynamics of IR cascading over a much longer time scale.Within a fewe-foldings from the time of particle production, the entire bump-like feature from IR cascading becomes frozen outside the horizon. At this point the fluctuations have become classical, large-scale adiabatic density perturbations and are observable in the present epoch (presuming that ϕ=ϕ0 occurs during the observable range of e-foldings). In Figure 3 we have illustrated this bump-like feature in both the primordial power spectrum and angular TT spectrum, for a representative choice of parameters.


## 3.3. Backreaction Effects

As discussed previously, the production ofχ fluctuations at t=0 backreacts on the homogeneous ϕ(t) causing a transient violation of slow roll. We can study this backreaction numerically, by averaging the inhomogeneous field ϕ̇(t,x) over the simulation box. The result is plotted in Figure 2. We have also plotted the analytical solution of the mean field (14), showing that this agrees with the exact numerical result.The dynamics illustrated in Figure2 is easy to understand physically. The production of χ particles at t=0 drains kinetic energy from the condensate and hence ϕ̇ must decrease abruptly. However, within a few e-foldings of the moment t=0, the produced iso-inflaton particles become extremely massive and are diluted by the expansion as a-3. At late times the inflaton velocity ϕ̇ must tend to the slow-roll value. Notice that the velocity ϕ̇ including backreaction effects is not changed significantly, as compared to the usual slow-roll result. This illustrates the energetic cheapness of particle production and IR cascading in model (5).The transient violation of slow roll illustrated in Figure2 is expected to induce a ringing pattern in the vacuum fluctuations from inflation [64]. This effect is accounted for automatically in our HLattice simulations. However, we would like to disentangle the effect of backreaction on the cosmological fluctuations from the effect of rescattering. This will be useful in order to compare the relative importance of different physical processes and also to guide our analytical efforts in the next section. To this end, we consider the evolution of the curvature perturbation on comoving hypersurfaces, ℛ. In linear theory the equation for the Fourier modes ℛk is well known(21)ℛk′′+2z′zℛk′+k2ℛk=0.
Here the prime denotes derivatives with respect to conformal time τ=∫(dt/a) and z≡aϕ̇/H. Equation (21) is only strictly valid in the absence of entropy perturbations. However, in our case the χ field is extremely massive mχ2≫H2 for nearly the entire duration of inflation, hence one may expect that direct iso-curvature contributions to ℛ are small. We have solved (21) numerically. In order to take backreaction effects into account we compute the dynamics of z(t)=a(t)ϕ̇(t)/H(t) by averaging over our HLattice simulation box. Next, we solve (21) given this background evolution and compute the power spectrum(22)Pℛ(k)=k32π2|ℛk|2.
The result is very close to the usual power-law form kns-1, with small superposed oscillations resulting from the transient violation of slow roll; see Figure 7. In order to make the ringing pattern more visible, we have subtracted off the usual (nearly) scale-invariant result which would be obtained in the absence of particle production. For comparison, we also plot the bump-like feature from rescattering and IR cascading. This latter contribution was obtained using the results for Pϕ(k) from the previous subsection and the naive formula ℛ~(H/ϕ̇)δϕ (so that Pℛ~(2ϵMp2)-1Pϕ).Figure 7
A comparison of curvature fluctuations from different physical effects. The dashed blue line is the usual (nearly) scale invariant vacuum fluctuations from inflation. The red solid line is the bump-like feature induced by rescattering and IR cascading. The dotted blue line is the ringing pattern resulting from the momentary slowing-down of the inflaton (computed using the mean field approach of [64]). The vertical lines show aH at the beginning of particle production and after ~3e-foldings. This figure clearly illustrates the dominance of IR cascading over backreaction effects. For illustration we have taken g2=0.1, but the dominance is generic for all values of the coupling.From Figure7 we see that IR cascading has a much more significant impact on the observable curvature fluctuations than does backreaction. Indeed, for g2=0.1 the transient violation of slow roll yields an order 10-2 correction to the vacuum fluctuations while the correction from IR cascading is of order 101. This dominance is generic for all values of the coupling. Thus, in developing an analytical theory of particle production during inflation, it is a very good approximation to completely ignore backreaction effects.


## 4. Analytical Formalism

In the last section, we have studied particle production, rescattering, and IR cascading using nonlinear lattice field theory simulations. In this section we will develop a detailed analytical theory, in order to understand those results from a physical perspective. These results were first presented in [14]. We consider, again, model (5). The equations of motion that we wish to solve are(23)-ϕ+V′(ϕ)+g2(ϕ-ϕ0)χ2=0,(24)-χ+g2(ϕ-ϕ0)2χ=0,
where =gμν∇μ∇ν is the covariant d'Alembertian. It will be useful to work with conformal time τ, related to cosmic time t via adτ=dt. In terms of conformal time the metric takes the form(25)ds2=-dt2+a2(t)dx·dx=a2(τ)[-dτ2+dx·dx].
We denote derivatives with respect to cosmic time as ḟ≡∂tf and with respect to conformal time as f′≡∂τf. The Hubble parameter H=ȧ/a has conformal time analogue ℋ=a′/a. For an inflationary (quasi-de Sitter) phase (H≅const) one has(26)a=-1Hτ11-ϵ,ℋ=-1τ11-ϵ
to leading order in the slow roll parameter ϵ≪1.As discussed in Section2, the motion of the homogeneous inflaton ϕ(t) leads to the production of a gas of χ particles at the moment t=0 when ϕ=ϕ0. The first step in our analytical computation is to describe this burst of particle production in an expanding universe. Following the initial burst, both backreaction and rescattering effects take place. Our formalism will focus on the latter effect, which is much more important, and we provide only a cursory treatment of backreaction.


### 4.1. Particle Production in an Expanding Universe

The first step in our scenario is the quantum mechanical production ofχ particles due to the motion of ϕ. To understand this effect we must solve the equation for the χ fluctuations in the rolling inflaton background. Approximating ϕ≅ϕ0+vt (24) gives(27)χ̈+3Hχ̇-∇⃗2a2χ+k⋆4t2χ=0,
where k⋆≡g|v|. We remind the reader that k⋆≫H for reasonable values of the coupling; see (11).The flat space analogue of (27) is very well understood from studies of broad-band parametric resonance during preheating [92] and also moduli trapping at enhanced symmetry points [80]. One does not expect this treatment to differ significantly in our case since both the time scale for particle production Δt and the characteristic wavelength of the produced fluctuations λ are small compared to the Hubble scale: Δt~λ-1~k⋆-1≪H-1. Hence, we expect that the occupation number of produced χ particles will not differ significantly from the flat-space result (13), at least on scales k≳H. Furthermore, notice that the χ field is extremely massive for most of the inflation(28)mχ2H2≅k⋆4t2H2.
Since k⋆≫H, it follows that mχ2≫H2, except in a tiny interval H|Δt|~(H/k⋆)2 which amounts to roughly 10-3  e-foldings for g2~0.1. Therefore, we do not expect any significant fluctuations of χ to be produced on superhorizon scales k≲H.Let us now consider the solutions of (27). We work with conformal time τ and write the Fourier transform of the quantum field χ as(29)χ(τ,x)=∫d3k(2π)3/2ξkχ(τ)a(τ)eik·x.
Note the explicit factor of a-1 in (29) which is introduced to give ξkχ a canonical kinetic term. The q-number valued Fourier transform ξkχ(τ) can be written as(30)ξkχ(τ)=akχk(τ)+a-k†χk⋆(τ),
where the annihilation/creation operators satisfy the usual commutation relation(31)[ak,ak'†]=δ(3)(k-k′),
and the c-number valued mode functions χk obey the following oscillator-like equation:(32)χk′′(τ)+ωk2(τ)χk(τ)=0.
The time-dependent frequency is(33)ωk2(τ)=k2+a2mχ2(τ)-a′′a≅k2+1τ2[k⋆4H2t2(τ)-2],
where mχ2(τ)=g2(ϕ-ϕ0)2≅k⋆4t2(τ) is the time-dependent effective mass of the χ particles, and(34)t(τ)=1Hln(-1Hτ)
is the usual cosmic time variable. We have arbitrarily set the origin of conformal time so that τ=-1/H corresponds to the moment when ϕ=ϕ0.In Figure8(a) we have plotted a representative solution of (32) in order to illustrate the qualitative behaviour of the modes χk. In Figure 8(b) we plot the occupation number nk of particles with momentum k, defined as the energy of the mode (1/2)|χk′|2+(1/2)ωk2|χk|2 divided by the energy ωk of each particle. Explicitly, we define(35)nk=ωk2[|χk′|2ωk2+|χk|2]-12
where the term -1/2 comes from extracting the zero-point energy of the linear harmonic oscillator (see [92] for a review). From Figure 8(a) we see that, near the massless point t=0, the fluctuations χk get a “kick”, and from Figure 8(b) we see that the occupation number nk jumps abruptly at this same moment.(a) illustrates the time dependence of the solutionsχk of (32) for a representative choice of parameters. The oscillatory behaviour at early times represents the adiabatic initial condition. At t=0 the effective frequency (33) varies nonadiabatically and the fluctuations get a “kick”. (b) plots the occupation number (35) for the same mode. No particles are present in the adiabatic in-going regime. This figure illustrates how the violations of adiabaticity at t=0 lead to production of χ particles.
(a)(b)Let us now try to understand analytically the behaviour of the solutions of (32). At early times t≪-k⋆-1, the frequency ωk varies adiabatically(36)|ωk′ωk2|≪1.
In this in-going adiabatic regime the modes χk are not excited and the solution of (32) is well described by the adiabatic solution χk(τ)=fk(τ) where(37)fk(τ)≡12ωk(τ)exp[-i∫τdτ′ωk(τ)].
We have normalized (37) to be pure positive frequency so that the state of the iso-inflaton field at early times corresponds to the adiabatic vacuum with no χ particles. (Inserting (37) into (35) one finds nk=0 for the adiabatic solution, as expected.)The adiabatic solution (37) ceases to be a good approximation very close to the moment when ϕ=ϕ0, that is at times |t|≲k⋆-1. In this regime the adiabaticity condition (36) is violated for modes with wavenumber H≲k≲k⋆ and χ particles within this momentum band are produced. During the non-adiabatic regime we can still represent the solutions of (32) in terms of the functions fk(τ) as(38)χk(τ)=αk(τ)fk(τ)+βk(τ)fk⋆(τ).
This expression affords a solution of (32) provided the time-dependent Bogoliubov coefficients obey the following set of coupled equations:(39)αk′(τ)=ωk′(τ)2ωk(τ)exp[+2i∫τdτ′ωk(τ′)]βk(τ),βk′(τ)=ωk′(τ)2ωk(τ)exp[-2i∫τdτ′ωk(τ′)]αk(τ).
The Bogoliubov coefficients are normalized as |αk|2-|βk|2=1, and the assumption that no χ particles are present in the asymptotic past (this assumption is justified since any initial excitation of χ would have been damped out exponentially fast by the expansion of the universe) fixes the initial conditions αk=1, βk=0 for t→-∞. This is known as the adiabatic initial condition.From the structure of (39) it is clear that violation of condition (36) near t=0 leads to a rapid growth in the |βk| coefficient. The time variation of βk can be interpreted as a corresponding growth in the occupation number of the χ particles(40)nk=|βk|2.At late times (t≳k⋆-1) adiabaticity is restored and the growth of nk=|βk|2 must saturate. By inspection of (39) we can see that the Bogoliubov coefficients must tend to constant values in the out-going adiabatic regime. Therefore, within less than an e-folding from the moment of particle production the solution χk of (32) can be represented as a simple superposition of positive frequency fk modes and negative frequency fk⋆ modes. Our goal now is to derive an analytical expression for the modes χk which is valid in this out-going adiabatic region.Let us first study the adiabatic solutionfk(τ). If we focus on the interesting region of phase space, H≲k≲k⋆, then the adiabatic solution (37) is very well approximated by(41)fk(τ)≅1a1/2k⋆2t(τ)e-(i/2)k⋆2t2(τ),
where t(τ) is defined by (34). It is interesting to note that (41) is identical to the analogous flat-space result [12], except for the factor of a-1/2. Taking into account also the explicit factor of a-1 in our definition of the Fourier transform (29), we recover the expected large-scale behaviour for a massive field in de Sitter space, that is, χ~a-3/2. This dependence on the scale factor is easy to understand physically, it simply reflects the volume dilution of nonrelativistic particles: ρχ~mχ2χ2~a-3.Next, we seek an expression for the Bogoliubov coefficientsαk, βk in the out-going adiabatic regime t≳k⋆-1. From (39) it is clear that the value of the Bogoliubov coefficients at late times can depend only on dynamics during the interval |t|≲k⋆-1 where adiabaticity condition (36) is violated. This interval is tiny compared to the expansion time, and we are justified in treating a(τ) as a constant during this phase. Hence, it follows that the flat space computation of the Bogoliubov coefficients [80, 92] must apply, at least for scales k≳H. To a very good approximation we therefore have the well-known result(42)αk≅1+e-πk2/k⋆2,(43)βk≅-ie-πk2/(2k⋆2)
in the out-going adiabatic regime. Equation (43) gives the usual expression (13) for the co-moving occupation number of particles produced by a singe burst of broad-band parametric resonance:(44)nk=|βk|2=e-πk2/k⋆2.Finally, we arrive at an expression for the out-going adiabaticχ modes which is accurate for interesting scales k⋆≲k≲H. Putting together results (41) and (38) along with well-known expressions (42) and (43) we arrive at(45)χk(τ)≅1+e-πk2/k⋆21a1/2k⋆2t(τ)e-(i/2)k⋆2t2(τ)-ie-πk2/(2k⋆2)1a1/2k⋆2t(τ)e+(i/2)k⋆2t2(τ),
valid for t≳k⋆-1. Equation (45) is the main result of this subsection. We will now justify that this expression is quite sufficient for our purposes.For modes deep in the UV,k≳k⋆, our expression (45), is not accurate. (Expression (41) for the adiabatic modes fk is not valid at high momenta where ωk≅k.) However, such high momentum particles are not produced, condition (36) is always satisfied for k≫k⋆. Note that the absence of particle production deep in the UV is built into our expression (45): as k→∞ this function tends to the vacuum solution χk→fk.Our expression (45) is also not valid deep in the IR, for modes k<H. To justify this neglect requires somewhat more care. Notice that, even very far from the massless point, t=0, long wavelength modes k≪H should not be thought of as particle like. The large-scale mode functions are not oscillatory but rather damp exponentially fast as χ~a-3/2. Hence, even if we started with some super-horizon fluctuations of χ at the beginning of inflation, these would be suppressed by an exponentially small factor before the time when particle production occurs. Any super-horizon fluctuation generated near t=0 would need to be exponentially huge to overcome this damping. However, resonant particle production during inflation does not lead to exponential growth of mode functions. (In this regard our scenario is very different from preheating at the end of inflation. In the latter case the inflaton passes many times through the massless point mχ=0, and there are, correspondingly, many bursts of particle production. After many oscillations of the inflaton field, the χ particle occupation numbers build up to become exponentially large, and, averaged over many oscillations of the background, the χ mode functions grow exponentially. However, in our case there is only a single burst of particle production at t=0. The resulting occupation number (13) is always less than unity, and the solutions of (32) never display exponential growth.)To verify explicitly that there is no significant effect for super-horizon fluctuations let us consider solving (27), neglecting gradient terms. The equation we wish to solve, then, is(46)∂t2(a3/2χ)+[k⋆4t2-94H2](a3/2χ)=0.
(For simplicity we take ϵ=0 for this paragraph; however, this has no effect on our results.) The solution of this equation may be written in terms of parabolic cylinder functions Dν(z) as(47)χ(t,x)~1a3/2(C1D-1/2+(9H2/8k⋆2)i[(1+i)k⋆t]+C2D-1/2-(9H2/8k⋆2)i[(-1+i)k⋆t]).
For our purposes the precise values of the coefficients C1, C2 are not important. Rather, it suffices to note that for k⋆|t|≳1 function (47) behaves as(48)χ(t,x)~|t|-1/2e-3Ht/2×[oscillatory].
This explicit large-scale asymptotics confirms our previous claims that the super-horizon fluctuations of χ damp to zero exponentially fast, as a-3/2~e-3Ht/2. As discussed previously, this damping is easy to understand in terms of the volume dilution of non-relativistic particles. We can also understand the power-law damping that appears in (48) from a physical perspective. The properly normalized modes behave as a3/2χ~ωk-1/2 while on large scales we have ωk~|mχ|~|t|. Hence, the late-time damping factor t-1/2 which appears in (48) reflects the fact that the χ particles become ever more massive as ϕ rolls away from the point ϕ0.Finally, it is straightforward to see that function (47) does not display any exponential growth near t=0. Hence, we conclude that there is no significant generation of super-horizon χ fluctuations due to particle production. (This is strictly true only in the linearized theory. It is possible that χ particles are generated by nonlinear effects such as rescattering. However, even such second-order χ fluctuations will be extremely massive compared to the Hubble scale and must therefore suffer exponential damping a-3/2 on large scales.)In this subsection we have seen that the quantum production ofχ particles in an expanding universe proceeds very much as it does in flat space. This is reasonable since particle production occurs on a time scale short compared to the expansion time and involves modes which are inside the horizon at the time of production.


### 4.2. Inflaton Fluctuations

In Section4.1 we studied the quantum production of χ particles which occurs when ϕ rolls past the massless point ϕ=ϕ0. Subsequently, there are two distinct physical processes which take place: backreaction and rescattering. As we have argued in Section 3, the former effect has a negligible impact of the observable spectrum of cosmological perturbations. Hence, we will not study this effect in any detail (see [61, 63, 64] for analytical calculations). Instead we provide a cursory treatment of backreaction in Appendix A, in order to clear up some common misconceptions.In this subsection we study the rescattering of producedχ particle off the inflaton condensate. The dominant process to consider is the diagram illustrated in Figure 1, corresponding to bremsstrahlung emission of δϕ fluctuations (particles) in the background of the external field. (There is also a subdominant process of the type χχ→δϕδϕ which is phase-space suppressed.) Taking into account only the rescattering diagram illustrated in Figure 1 is equivalent to solving the following equation for the q-number inflaton fluctuation:(49)[∂t2+3H∂t-∇⃗2a2+m2]δϕ=-g2[ϕ(t)-ϕ0]χ2,
where we have introduced the notation m2≡V,ϕϕ for the inflaton effective mass. (Note that we are not assuming a background potential of the form m2ϕ2/2, only that V,ϕϕ≠0 in the vicinity of the point ϕ=ϕ0.)Equation (49) may be derived by noting that (5) gives an interaction of the form g2(ϕ-ϕ0)δϕχ2 between the inflaton and iso-inflaton, in the background of the external field ϕ(t). Equivalently, one may construct this equation by a straightforward iterative solution of (23).We work in conformal time and define the q-number Fourier transformξkϕ(τ) of the inflaton fluctuation analogously to (29):(50)δϕ(τ,x)=∫d3k(2π)3/2ξkϕ(τ)a(τ)eik·x.
(To avoid potential confusion we again draw the attention of the reader to the explicit factor a-1 in our convention for the Fourier transform.) The equation of motion (49) now takes the form(51)[∂τ2+k2+a2m2-a′′a]ξkϕ(τ)=-gk⋆2a(τ)t(τ)∫d3k′(2π)3/2ξk′χξk-k′χ(τ).
The solution of (51) consists of two parts: the solution of the homogeneous equation and the particular solution which is due to the source. The former corresponds, physically, to the usual vacuum fluctuations from inflation. On the other hand, the particular solution corresponds physically to the secondary inflaton modes which are generated by rescattering.


### 4.3. Homogeneous Solution and Green Function

We consider first the homogeneous solution of (51). Since the homogeneous solution is a Gaussian field, we may expand the q-number Fourier transform in terms of annihilation/creation operators bk, bk† and c-number mode functions ϕk(τ) as(52)ξkϕ(τ)=bkϕk(τ)+b-k†ϕk⋆(τ).
Here the inflaton annihilation/creation operators bk, bk† obey(53)[bk,bk′†]=δ(3)(k-k′)
and commute with the annihilation/creation operators of the χ-field:(54)[ak,bk′]=[ak,bk′†]=0.Using (26) and (7) it is straightforward to see that the homogeneous inflaton mode functions obey the following equation:(55)∂τ2ϕk+[k2-1τ2(ν2-14)]ϕk=0,
where we have defined(56)ν≅32-η+ϵ.
The properly normalized mode function solutions are well known and may be written in terms of the Hankel function of the first kind as(57)ϕk(τ)=π2-τHν(1)(-kτ).
This solution corresponds to the usual quantum vacuum fluctuations of the inflaton field during inflation.In passing, let us compute the power spectrum of the quantum vacuum fluctuations from inflation. Using the solutions (57) we have(58)Pϕvac(k)=k32π2|ϕka|2≅H2(2π)2(kaH)ns-1
on large scales k≪aH. The explicit factor of a-2 in (58) appears to cancel the a-1 in our definition of the Fourier transform (50). The spectral index is(59)ns-1=3-2ν≅2η-2ϵ
using (56).Given the solution (57) of the homogeneous equation, it is now trivial to construct the retarded Green function for (51). This may be written in terms of the free theory mode functions (57) as(60)Gk(τ-τ′)=iΘ(τ-τ′)[ϕk(τ)ϕk⋆(τ′)-ϕk⋆(τ)ϕk(τ′)]=iπ4Θ(τ-τ′)ττ′[Hν(1)(-kτ)Hν(1)(-kτ′)⋆-Hν(1)(-kτ)⋆Hν(1)(-kτ′)].


### 4.4. Particular Solution: Rescattering Effects

We now consider the particular solution of (51). This is readily constructed using the Green function (60) as(61)ξkϕ(τ)=-gk⋆2(2π)3/2∫dτ′d3k′Gk(τ-τ′)a(τ′)t(τ′)ξk′χξk-k′χ(τ′).
Notice that this particular solution is statistically independent of the homogeneous solution (52). In other words, the particular solution can be expanded in terms of the annihilation/creation operators ak,ak† associated with the χ field whereas the homogeneous solution is written in terms of the annihilation/creation operators bk,bk† associated with the inflaton vacuum fluctuations. These two sets of operators commute with one another.We will ultimately be interested in computing then-point correlation functions of the particular solution (62). For example, carefully carrying out the Wick contractions, the connected contribution to the 2-point function is(62)〈ξk1ϕξk2ϕ(τ)〉=2g2k⋆4(2π)3δ(3)(k1+k2)×∫dτ′dτ′′a(τ′)a(τ′′)t(τ′)t(τ′′)Gk1(τ-τ′)Gk2(τ-τ′′)×∫d3k′χk1-k′(τ′)χk1-k′⋆(τ′′)χk′(τ′)χk′⋆(τ′′).
The power spectrum of δϕ fluctuations generated by rescattering is then defined in terms of the 2-point function in the usual manner(63)〈ξkϕ(t)ξk′ϕ(τ)〉≡δ(3)(k+k′)2π2k3a(τ)2Pϕresc.
(The explicit factor of a2 in definition (63) appears to cancel the factor of a-1 in our convention for Fourier transforms (50).)The total power spectrum is simply the sum of the contribution from the vacuum fluctuations (58) and the contribution from rescattering (63):(64)Pϕ(k)=Pϕvac(k)+Pϕresc(k).
There are no cross-terms, owing to the fac that ak and bk commute.


### 4.5. Renormalization

We now wish to evaluate the 2-point correlator (62). In principle, this is straightforward: first substitute result (45) for the χk modes and result (60) for the Green function into (62), next evaluate the integrals. However, there is a subtlety. The resulting power spectrum is formally infinite. Moreover, the 2-point correlation function (62) receives contributions from two distinct effects. There is a contribution from particle production, which we are interested in. However, there is also a contribution coming from quantum vacuum fluctuations of the χ field interacting nonlinearly with the inflaton. The latter contribution would be present even in the absence of particle production, when αk=1, βk=0.In order to isolate the effects of particle production on the inflaton fluctuations, we would like to subtract off the contribution to the 2-point correlation function (62) which is coming from the quantum vacuum fluctuations of χ. This subtraction also has the effect of rendering the power spectrum (63) finite, since it extracts the usual UV divergent contribution associated with the Minkowski-space vacuum fluctuations.As a step towards renormalizing the 2-point correlation function of inflaton fluctuations from rescattering (62), let us first consider the simpler problem of renormalizing the 2-point function of the gaussian field χ. We defined the renormalized 2-point function in momentum space as follows:(65)〈ξk1χ(t1)ξk2χ(t2)〉ren=〈ξk1χ(t1)ξk2χ(t2)〉-〈ξk1χ(t1)ξk2χ(t2)〉in.
In (65) the quantity 〈ξk1χ(t1)ξk2χ(t2)〉in is the contribution which would be present even in the absence of particle production, computed by simply taking solution (38) with αk=1, βk=0. Explicitly, we have(66)〈ξk1χ(t1)ξk2χ(t2)〉in=δ(3)(k1+k2)fk1(t1)fk2⋆(t2),
where fk are the adiabatic solutions (37).To see the impact of this subtraction, let us consider the renormalized variance for the iso-inflaton field,〈χ2〉. Employing prescription (65) we have(67)〈χ2(τ,x)〉ren=∫d3k(2π)3a2(τ)[|χk(τ)|2-12ωk(τ)]=〈χ2(τ,x)〉-δM,
where δM is the contribution from the Coleman-Weinberg potential. This proves that our prescription reproduces the scheme advocated in [80]. The renormalized variance (67) is finite and may be computed explicitly using our solutions (45). We find(68)〈χ2(t,x)〉ren≅nχa-3g|ϕ-ϕ0|,
where(69)nχ≡∫d3k(2π)3nk~k⋆3
is the total co-moving number density of produced χ particles. Result (68) was employed in [64] to quantify the effect of backreaction on the inflaton condensate in the mean field treatment (14). Hence, the renormalization scheme (65) was implicit in that calculation also.At the level of the 2-point function, our renormalization scheme is tantamount to assuming that Coleman-Weinberg corrections are already absorbed into the definition of the inflaton potential,V(ϕ). In general, such corrections might steepen V(ϕ) and spoil slow-roll inflation. Here, we assume that this problem has already been dealt with, either by fine-tuning the bare inflaton potential or else by including extended SUSY (which can minimize dangerous corrections). See also [80] for a related discussion.Having established a scheme for renormalizing the 2-point function of the gaussian fieldχ, it is now straightforward to consider higher-order correlation functions. We simply rewrite the 4-point function as a product of 2-point functions using Wick's theorem. Next, each Wick contraction is renormalized as (65). Applying this prescription to (62) amounts to(70)〈ξk1ϕ(τ)ξk2ϕ(τ)〉ren=2g2k⋆4(2π)3δ(3)(k1+k2)×∫dτ′dτ′′t(τ′)t(τ′′)a(τ′)a(τ′′)Gk1(τ-τ′)Gk2(τ-τ′′)×∫d3k′[χk1-k′(τ′)χk1-k′⋆(τ′′)-fk1-k′(τ′)fk1-k′⋆(τ′′)]×[χk′(τ′)χk′⋆(τ′′)-fk′(t′)fk′⋆(τ′′)],
where fk(t) are the adiabatic solutions defined in (37).


### 4.6. Power Spectrum

We are now in a position to compute the renormalized power spectrum of inflation fluctuations generated by rescattering,Pϕresc(k). We renormalize the 2-point correlator of the inflaton fluctuations generated by rescatter according to (70) and extract the power spectrum by comparison to (63). We have relegated the technical details to Appendix B and here we simply state the final result:(71)Pϕresc(k)=g2k3k⋆16π5[e-πk2/(2k⋆2)22(I2(k,τ)2+|I1(k,τ)|2)+[e-πk2/(4k⋆2)+122e-3πk2/(8k⋆2)]×(I2(k,τ)2-Re[I1(k,τ)])+[8233e-πk2/(3k⋆2)+4255e-3πk2/(5k⋆2)]×g2k3k⋆16π5Im[I1(k,τ)I2(k,τ)]],
where the functions I1, I2 are the curved space generalization of the characteristic integrals defined in [12]. Explicitly we have (72)I1(k,τ)=1a(τ)∫dτ′Gk(τ-τ′)eik⋆2t2(τ′),I2(k,τ)=1a(τ)∫dτ′Gk(τ-τ′).
The characteristic integral I2 can be evaluated analytically; however, the resulting expression is not particularly enlightening. Evaluation of the integral I1 requires numerical methods. More details are in Appendix B. Equation (71) is the main result of this section.


### 4.7. Comparison to Lattice Field Theory Simulations

In Section3 the results of our analytical formalism were plotted alongside the output of fully nonlinear HLattice simulations. It is evident from Figures 4, 5, and 6 that the agreement between these approaches is extremely good, even very late into the evolution and in the regime g2~1. The consistency of perturbative quantum field theory analytics and nonlinear classical lattice simulations provides a highly nontrivial check on our calculation.


### 4.8. The Bispectrum

So far, we have shown how to compute analytically the power spectrum generated by particle production, rescattering, and IR cascading during inflation in model (5). We found that IR cascading leads to a bump-like contribution to the primordial power spectrum of the inflaton fluctuations. However, this same dynamics must also have a nontrivial impact on nongaussian statistics, such as the bispectrum. Indeed, it is already evident from our previous analysis that the inflaton fluctuations generated by rescattering may be significantly nongaussian. From expression (62) we see that the particular solution (due to rescattering) which is bi-linear is the gaussian field χ.We define the bispectrum of the inflaton field fluctuations in terms of the three-point correlation function as(73)〈ξk1ϕξk2ϕξk3ϕ(τ)〉=(2π)3a3(τ)δ(k1+k2+k3)Bϕ(ki).
The factor a3 appears in (73) to cancel the explicit factors of a-1 in our convention (50) for the Fourier transform. It is well known that the nongaussianity associated with the usual quantum vacuum fluctuations of the inflaton is negligible [16–18]; therefore, when evaluating the bispectrum (73) we consider only the particular solution (62) which is due to rescattering. Carefully carrying out the Wick contractions, we find the following result for the renormalized 3-point function:(74)〈ξk1ϕξk2ϕξk3ϕ(τ)〉ren=4g3k⋆6(2π)9/2δ(k1+k2+k3)∏i=13∫dτit(τi)a(τi)Gki(τ-τi)×∫d3p[χk1-p(τ1)χk1-p⋆(τ2)-fk1-p(τ1)fk1-p⋆(τ2)]×[χk3+p(τ2)χk3+p⋆(τ3)-fk3+p(τ2)fk3+p⋆(τ3)]×[χp(τ1)χp⋆(τ3)-fp(τ1)fp⋆(τ3)]+(k2⟷k3),
where the modes χk are defined by (38) and fk are the adiabatic solutions (37). On the last line of (74) we have labeled schematically terms which are identical to the preceding three lines, only with k2 and k3 interchanged. One may verify that this expression is symmetric under interchange of the momenta ki by changing dummy variables of integration.It is now straightforward (but tedious) to plug expressions (41) and (45) into (74) and evaluate the integrals. This computation is tractable analytically because the time and phase-space integrals decouple. The bispectrum is then extracted by comparison to (73). This computation is carried out in detail in [14] where we have shown that Bϕ(ki) peaks only over triangles with a characteristic size, corresponding to the location of the bump in the power spectrum. This result is easy to understand on physical grounds, all the dynamics of rescattering and IR cascading take place over a very short time near the moment t=0. Hence, the effect of this dynamics on the primordial fluctuations must be limited to scales leaving the horizon near the time when particle production occurs.We will provide a cursory discussion of the bispectrumBϕ(ki) in Section 8 when we discuss nongaussianity from particle production during inflation.


### 4.9. Inclusion of a Bare Iso-Inflaton Mass

In passing, it may be interesting to consider how the analysis of this section is modified in the case that our prototype action (5) is supplemented by a bare mass term for the iso-inflaton field of the form δℒ=-(1/2)μ2χ2. Thus, in place of (5) suppose that we consider the model(75)ℒ=-12(∂ϕ)2-V(ϕ)-12(∂χ)2-12μ2χ2-g22(ϕ-ϕ0)2χ2.
Now the χ particles do not become massless at the point ϕ=ϕ0, but rather the effective mass-squared(76)mχ2=μ2+g2(ϕ-ϕ0)2
reaches a minimum value μ2 (which we assume to be positive). Such a correction may arise due to a variety of effects and will reduce the impact of particle production and IR cascading on the observable cosmological fluctuations.Let us briefly consider how the additional bare mass term in (75) alters the dynamics of particle production. The iso-inflaton fluctuations now obey the equation(77)χ̈+3Hχ̇-∇⃗2a2χ+[μ2+k⋆4t2]χ=0
rather than (27). This equation was solved in [80] in the regime where particle production is fast compared to the expansion time. (In the opposite regime, which corresponds to a fine-tuned coupling g2≪10-7, the iso-inflaton will be light for a significant fraction of inflation. In that case the theory (5) must be considered as a multifield inflation model and one can no longer consistently assume 〈χ〉=0. In other words, relaxing the assumption of fast particle production significantly changes the scenario under consideration, and we do not pursue this possibility any further.) The occupation number of produced χ particles is(78)nk=re-πk2/k⋆2
which differs from our previous result (13) by the suppression factor(79)r≡e-πμ2/k⋆2≤1.
Therefore, the effect of the inclusion of a bare mass for the iso-inflaton is to suppress the number density of produced χ particles by an amount r. This suppression reflects the reduced phase space of produced particles: the adiabaticity condition |ωk′/ωk2|≪1 is violated only for modes with k<k⋆2-μ2.The reduction ofnk translates into a suppression for the n-point correlation functions of the iso-inflaton. For example, the renormalized variance 〈χ2〉ren∝∫d3knk is suppressed by a factor of r. The power spectrum of inflaton fluctuations generated by rescattering, Pϕresc, is proportional to the 4-point correlator of χ and hence picks up a suppression factor of order r2. Similarly, the bispectrum Bϕ is proportional to the 6-point correlator of χ and must be reduced by a factor of order r3. The condition(80)μ2≪k⋆2
is equivalent to r≅1 and ensures that the addition of a bare iso-inflaton mass will have a negligible impact on any observable.For models obtained from string theory or supergravity (SUGRA), it is natural to expectμ of the order of the Hubble scale during inflation [116–118]. In the context of SUGRA, the finite energy density driving inflation breaks SUSY and induces soft scalar potentials with curvature of order Vsoft′′~μ2~H2 [117]. In the case of string theory, many scalars are conformally coupled to gravity [118] through an interaction of the form δℒ=-(1/12)Rχ2 where the Ricci scalar is R~H2 during inflation. More generally, any nonminimal coupling δℒ=-(ξ/2)Rχ2 between gravity and the iso-inflaton will induce a contribution of order H to the effective mass of χ, as long as ξ=𝒪(1). In all models where μ2~H2 condition (80) is satisfied for reasonable values of the coupling g2>10-7, see; (11). Thus, we expect that corrections of the form δℒ=-(1/2)μ2χ2 will not alter our results in a wide variety of well-motivated models.


## 4.1. Particle Production in an Expanding Universe

The first step in our scenario is the quantum mechanical production ofχ particles due to the motion of ϕ. To understand this effect we must solve the equation for the χ fluctuations in the rolling inflaton background. Approximating ϕ≅ϕ0+vt (24) gives(27)χ̈+3Hχ̇-∇⃗2a2χ+k⋆4t2χ=0,
where k⋆≡g|v|. We remind the reader that k⋆≫H for reasonable values of the coupling; see (11).The flat space analogue of (27) is very well understood from studies of broad-band parametric resonance during preheating [92] and also moduli trapping at enhanced symmetry points [80]. One does not expect this treatment to differ significantly in our case since both the time scale for particle production Δt and the characteristic wavelength of the produced fluctuations λ are small compared to the Hubble scale: Δt~λ-1~k⋆-1≪H-1. Hence, we expect that the occupation number of produced χ particles will not differ significantly from the flat-space result (13), at least on scales k≳H. Furthermore, notice that the χ field is extremely massive for most of the inflation(28)mχ2H2≅k⋆4t2H2.
Since k⋆≫H, it follows that mχ2≫H2, except in a tiny interval H|Δt|~(H/k⋆)2 which amounts to roughly 10-3  e-foldings for g2~0.1. Therefore, we do not expect any significant fluctuations of χ to be produced on superhorizon scales k≲H.Let us now consider the solutions of (27). We work with conformal time τ and write the Fourier transform of the quantum field χ as(29)χ(τ,x)=∫d3k(2π)3/2ξkχ(τ)a(τ)eik·x.
Note the explicit factor of a-1 in (29) which is introduced to give ξkχ a canonical kinetic term. The q-number valued Fourier transform ξkχ(τ) can be written as(30)ξkχ(τ)=akχk(τ)+a-k†χk⋆(τ),
where the annihilation/creation operators satisfy the usual commutation relation(31)[ak,ak'†]=δ(3)(k-k′),
and the c-number valued mode functions χk obey the following oscillator-like equation:(32)χk′′(τ)+ωk2(τ)χk(τ)=0.
The time-dependent frequency is(33)ωk2(τ)=k2+a2mχ2(τ)-a′′a≅k2+1τ2[k⋆4H2t2(τ)-2],
where mχ2(τ)=g2(ϕ-ϕ0)2≅k⋆4t2(τ) is the time-dependent effective mass of the χ particles, and(34)t(τ)=1Hln(-1Hτ)
is the usual cosmic time variable. We have arbitrarily set the origin of conformal time so that τ=-1/H corresponds to the moment when ϕ=ϕ0.In Figure8(a) we have plotted a representative solution of (32) in order to illustrate the qualitative behaviour of the modes χk. In Figure 8(b) we plot the occupation number nk of particles with momentum k, defined as the energy of the mode (1/2)|χk′|2+(1/2)ωk2|χk|2 divided by the energy ωk of each particle. Explicitly, we define(35)nk=ωk2[|χk′|2ωk2+|χk|2]-12
where the term -1/2 comes from extracting the zero-point energy of the linear harmonic oscillator (see [92] for a review). From Figure 8(a) we see that, near the massless point t=0, the fluctuations χk get a “kick”, and from Figure 8(b) we see that the occupation number nk jumps abruptly at this same moment.(a) illustrates the time dependence of the solutionsχk of (32) for a representative choice of parameters. The oscillatory behaviour at early times represents the adiabatic initial condition. At t=0 the effective frequency (33) varies nonadiabatically and the fluctuations get a “kick”. (b) plots the occupation number (35) for the same mode. No particles are present in the adiabatic in-going regime. This figure illustrates how the violations of adiabaticity at t=0 lead to production of χ particles.
(a)(b)Let us now try to understand analytically the behaviour of the solutions of (32). At early times t≪-k⋆-1, the frequency ωk varies adiabatically(36)|ωk′ωk2|≪1.
In this in-going adiabatic regime the modes χk are not excited and the solution of (32) is well described by the adiabatic solution χk(τ)=fk(τ) where(37)fk(τ)≡12ωk(τ)exp[-i∫τdτ′ωk(τ)].
We have normalized (37) to be pure positive frequency so that the state of the iso-inflaton field at early times corresponds to the adiabatic vacuum with no χ particles. (Inserting (37) into (35) one finds nk=0 for the adiabatic solution, as expected.)The adiabatic solution (37) ceases to be a good approximation very close to the moment when ϕ=ϕ0, that is at times |t|≲k⋆-1. In this regime the adiabaticity condition (36) is violated for modes with wavenumber H≲k≲k⋆ and χ particles within this momentum band are produced. During the non-adiabatic regime we can still represent the solutions of (32) in terms of the functions fk(τ) as(38)χk(τ)=αk(τ)fk(τ)+βk(τ)fk⋆(τ).
This expression affords a solution of (32) provided the time-dependent Bogoliubov coefficients obey the following set of coupled equations:(39)αk′(τ)=ωk′(τ)2ωk(τ)exp[+2i∫τdτ′ωk(τ′)]βk(τ),βk′(τ)=ωk′(τ)2ωk(τ)exp[-2i∫τdτ′ωk(τ′)]αk(τ).
The Bogoliubov coefficients are normalized as |αk|2-|βk|2=1, and the assumption that no χ particles are present in the asymptotic past (this assumption is justified since any initial excitation of χ would have been damped out exponentially fast by the expansion of the universe) fixes the initial conditions αk=1, βk=0 for t→-∞. This is known as the adiabatic initial condition.From the structure of (39) it is clear that violation of condition (36) near t=0 leads to a rapid growth in the |βk| coefficient. The time variation of βk can be interpreted as a corresponding growth in the occupation number of the χ particles(40)nk=|βk|2.At late times (t≳k⋆-1) adiabaticity is restored and the growth of nk=|βk|2 must saturate. By inspection of (39) we can see that the Bogoliubov coefficients must tend to constant values in the out-going adiabatic regime. Therefore, within less than an e-folding from the moment of particle production the solution χk of (32) can be represented as a simple superposition of positive frequency fk modes and negative frequency fk⋆ modes. Our goal now is to derive an analytical expression for the modes χk which is valid in this out-going adiabatic region.Let us first study the adiabatic solutionfk(τ). If we focus on the interesting region of phase space, H≲k≲k⋆, then the adiabatic solution (37) is very well approximated by(41)fk(τ)≅1a1/2k⋆2t(τ)e-(i/2)k⋆2t2(τ),
where t(τ) is defined by (34). It is interesting to note that (41) is identical to the analogous flat-space result [12], except for the factor of a-1/2. Taking into account also the explicit factor of a-1 in our definition of the Fourier transform (29), we recover the expected large-scale behaviour for a massive field in de Sitter space, that is, χ~a-3/2. This dependence on the scale factor is easy to understand physically, it simply reflects the volume dilution of nonrelativistic particles: ρχ~mχ2χ2~a-3.Next, we seek an expression for the Bogoliubov coefficientsαk, βk in the out-going adiabatic regime t≳k⋆-1. From (39) it is clear that the value of the Bogoliubov coefficients at late times can depend only on dynamics during the interval |t|≲k⋆-1 where adiabaticity condition (36) is violated. This interval is tiny compared to the expansion time, and we are justified in treating a(τ) as a constant during this phase. Hence, it follows that the flat space computation of the Bogoliubov coefficients [80, 92] must apply, at least for scales k≳H. To a very good approximation we therefore have the well-known result(42)αk≅1+e-πk2/k⋆2,(43)βk≅-ie-πk2/(2k⋆2)
in the out-going adiabatic regime. Equation (43) gives the usual expression (13) for the co-moving occupation number of particles produced by a singe burst of broad-band parametric resonance:(44)nk=|βk|2=e-πk2/k⋆2.Finally, we arrive at an expression for the out-going adiabaticχ modes which is accurate for interesting scales k⋆≲k≲H. Putting together results (41) and (38) along with well-known expressions (42) and (43) we arrive at(45)χk(τ)≅1+e-πk2/k⋆21a1/2k⋆2t(τ)e-(i/2)k⋆2t2(τ)-ie-πk2/(2k⋆2)1a1/2k⋆2t(τ)e+(i/2)k⋆2t2(τ),
valid for t≳k⋆-1. Equation (45) is the main result of this subsection. We will now justify that this expression is quite sufficient for our purposes.For modes deep in the UV,k≳k⋆, our expression (45), is not accurate. (Expression (41) for the adiabatic modes fk is not valid at high momenta where ωk≅k.) However, such high momentum particles are not produced, condition (36) is always satisfied for k≫k⋆. Note that the absence of particle production deep in the UV is built into our expression (45): as k→∞ this function tends to the vacuum solution χk→fk.Our expression (45) is also not valid deep in the IR, for modes k<H. To justify this neglect requires somewhat more care. Notice that, even very far from the massless point, t=0, long wavelength modes k≪H should not be thought of as particle like. The large-scale mode functions are not oscillatory but rather damp exponentially fast as χ~a-3/2. Hence, even if we started with some super-horizon fluctuations of χ at the beginning of inflation, these would be suppressed by an exponentially small factor before the time when particle production occurs. Any super-horizon fluctuation generated near t=0 would need to be exponentially huge to overcome this damping. However, resonant particle production during inflation does not lead to exponential growth of mode functions. (In this regard our scenario is very different from preheating at the end of inflation. In the latter case the inflaton passes many times through the massless point mχ=0, and there are, correspondingly, many bursts of particle production. After many oscillations of the inflaton field, the χ particle occupation numbers build up to become exponentially large, and, averaged over many oscillations of the background, the χ mode functions grow exponentially. However, in our case there is only a single burst of particle production at t=0. The resulting occupation number (13) is always less than unity, and the solutions of (32) never display exponential growth.)To verify explicitly that there is no significant effect for super-horizon fluctuations let us consider solving (27), neglecting gradient terms. The equation we wish to solve, then, is(46)∂t2(a3/2χ)+[k⋆4t2-94H2](a3/2χ)=0.
(For simplicity we take ϵ=0 for this paragraph; however, this has no effect on our results.) The solution of this equation may be written in terms of parabolic cylinder functions Dν(z) as(47)χ(t,x)~1a3/2(C1D-1/2+(9H2/8k⋆2)i[(1+i)k⋆t]+C2D-1/2-(9H2/8k⋆2)i[(-1+i)k⋆t]).
For our purposes the precise values of the coefficients C1, C2 are not important. Rather, it suffices to note that for k⋆|t|≳1 function (47) behaves as(48)χ(t,x)~|t|-1/2e-3Ht/2×[oscillatory].
This explicit large-scale asymptotics confirms our previous claims that the super-horizon fluctuations of χ damp to zero exponentially fast, as a-3/2~e-3Ht/2. As discussed previously, this damping is easy to understand in terms of the volume dilution of non-relativistic particles. We can also understand the power-law damping that appears in (48) from a physical perspective. The properly normalized modes behave as a3/2χ~ωk-1/2 while on large scales we have ωk~|mχ|~|t|. Hence, the late-time damping factor t-1/2 which appears in (48) reflects the fact that the χ particles become ever more massive as ϕ rolls away from the point ϕ0.Finally, it is straightforward to see that function (47) does not display any exponential growth near t=0. Hence, we conclude that there is no significant generation of super-horizon χ fluctuations due to particle production. (This is strictly true only in the linearized theory. It is possible that χ particles are generated by nonlinear effects such as rescattering. However, even such second-order χ fluctuations will be extremely massive compared to the Hubble scale and must therefore suffer exponential damping a-3/2 on large scales.)In this subsection we have seen that the quantum production ofχ particles in an expanding universe proceeds very much as it does in flat space. This is reasonable since particle production occurs on a time scale short compared to the expansion time and involves modes which are inside the horizon at the time of production.


## 4.2. Inflaton Fluctuations

In Section4.1 we studied the quantum production of χ particles which occurs when ϕ rolls past the massless point ϕ=ϕ0. Subsequently, there are two distinct physical processes which take place: backreaction and rescattering. As we have argued in Section 3, the former effect has a negligible impact of the observable spectrum of cosmological perturbations. Hence, we will not study this effect in any detail (see [61, 63, 64] for analytical calculations). Instead we provide a cursory treatment of backreaction in Appendix A, in order to clear up some common misconceptions.In this subsection we study the rescattering of producedχ particle off the inflaton condensate. The dominant process to consider is the diagram illustrated in Figure 1, corresponding to bremsstrahlung emission of δϕ fluctuations (particles) in the background of the external field. (There is also a subdominant process of the type χχ→δϕδϕ which is phase-space suppressed.) Taking into account only the rescattering diagram illustrated in Figure 1 is equivalent to solving the following equation for the q-number inflaton fluctuation:(49)[∂t2+3H∂t-∇⃗2a2+m2]δϕ=-g2[ϕ(t)-ϕ0]χ2,
where we have introduced the notation m2≡V,ϕϕ for the inflaton effective mass. (Note that we are not assuming a background potential of the form m2ϕ2/2, only that V,ϕϕ≠0 in the vicinity of the point ϕ=ϕ0.)Equation (49) may be derived by noting that (5) gives an interaction of the form g2(ϕ-ϕ0)δϕχ2 between the inflaton and iso-inflaton, in the background of the external field ϕ(t). Equivalently, one may construct this equation by a straightforward iterative solution of (23).We work in conformal time and define the q-number Fourier transformξkϕ(τ) of the inflaton fluctuation analogously to (29):(50)δϕ(τ,x)=∫d3k(2π)3/2ξkϕ(τ)a(τ)eik·x.
(To avoid potential confusion we again draw the attention of the reader to the explicit factor a-1 in our convention for the Fourier transform.) The equation of motion (49) now takes the form(51)[∂τ2+k2+a2m2-a′′a]ξkϕ(τ)=-gk⋆2a(τ)t(τ)∫d3k′(2π)3/2ξk′χξk-k′χ(τ).
The solution of (51) consists of two parts: the solution of the homogeneous equation and the particular solution which is due to the source. The former corresponds, physically, to the usual vacuum fluctuations from inflation. On the other hand, the particular solution corresponds physically to the secondary inflaton modes which are generated by rescattering.


## 4.3. Homogeneous Solution and Green Function

We consider first the homogeneous solution of (51). Since the homogeneous solution is a Gaussian field, we may expand the q-number Fourier transform in terms of annihilation/creation operators bk, bk† and c-number mode functions ϕk(τ) as(52)ξkϕ(τ)=bkϕk(τ)+b-k†ϕk⋆(τ).
Here the inflaton annihilation/creation operators bk, bk† obey(53)[bk,bk′†]=δ(3)(k-k′)
and commute with the annihilation/creation operators of the χ-field:(54)[ak,bk′]=[ak,bk′†]=0.Using (26) and (7) it is straightforward to see that the homogeneous inflaton mode functions obey the following equation:(55)∂τ2ϕk+[k2-1τ2(ν2-14)]ϕk=0,
where we have defined(56)ν≅32-η+ϵ.
The properly normalized mode function solutions are well known and may be written in terms of the Hankel function of the first kind as(57)ϕk(τ)=π2-τHν(1)(-kτ).
This solution corresponds to the usual quantum vacuum fluctuations of the inflaton field during inflation.In passing, let us compute the power spectrum of the quantum vacuum fluctuations from inflation. Using the solutions (57) we have(58)Pϕvac(k)=k32π2|ϕka|2≅H2(2π)2(kaH)ns-1
on large scales k≪aH. The explicit factor of a-2 in (58) appears to cancel the a-1 in our definition of the Fourier transform (50). The spectral index is(59)ns-1=3-2ν≅2η-2ϵ
using (56).Given the solution (57) of the homogeneous equation, it is now trivial to construct the retarded Green function for (51). This may be written in terms of the free theory mode functions (57) as(60)Gk(τ-τ′)=iΘ(τ-τ′)[ϕk(τ)ϕk⋆(τ′)-ϕk⋆(τ)ϕk(τ′)]=iπ4Θ(τ-τ′)ττ′[Hν(1)(-kτ)Hν(1)(-kτ′)⋆-Hν(1)(-kτ)⋆Hν(1)(-kτ′)].


## 4.4. Particular Solution: Rescattering Effects

We now consider the particular solution of (51). This is readily constructed using the Green function (60) as(61)ξkϕ(τ)=-gk⋆2(2π)3/2∫dτ′d3k′Gk(τ-τ′)a(τ′)t(τ′)ξk′χξk-k′χ(τ′).
Notice that this particular solution is statistically independent of the homogeneous solution (52). In other words, the particular solution can be expanded in terms of the annihilation/creation operators ak,ak† associated with the χ field whereas the homogeneous solution is written in terms of the annihilation/creation operators bk,bk† associated with the inflaton vacuum fluctuations. These two sets of operators commute with one another.We will ultimately be interested in computing then-point correlation functions of the particular solution (62). For example, carefully carrying out the Wick contractions, the connected contribution to the 2-point function is(62)〈ξk1ϕξk2ϕ(τ)〉=2g2k⋆4(2π)3δ(3)(k1+k2)×∫dτ′dτ′′a(τ′)a(τ′′)t(τ′)t(τ′′)Gk1(τ-τ′)Gk2(τ-τ′′)×∫d3k′χk1-k′(τ′)χk1-k′⋆(τ′′)χk′(τ′)χk′⋆(τ′′).
The power spectrum of δϕ fluctuations generated by rescattering is then defined in terms of the 2-point function in the usual manner(63)〈ξkϕ(t)ξk′ϕ(τ)〉≡δ(3)(k+k′)2π2k3a(τ)2Pϕresc.
(The explicit factor of a2 in definition (63) appears to cancel the factor of a-1 in our convention for Fourier transforms (50).)The total power spectrum is simply the sum of the contribution from the vacuum fluctuations (58) and the contribution from rescattering (63):(64)Pϕ(k)=Pϕvac(k)+Pϕresc(k).
There are no cross-terms, owing to the fac that ak and bk commute.


## 4.5. Renormalization

We now wish to evaluate the 2-point correlator (62). In principle, this is straightforward: first substitute result (45) for the χk modes and result (60) for the Green function into (62), next evaluate the integrals. However, there is a subtlety. The resulting power spectrum is formally infinite. Moreover, the 2-point correlation function (62) receives contributions from two distinct effects. There is a contribution from particle production, which we are interested in. However, there is also a contribution coming from quantum vacuum fluctuations of the χ field interacting nonlinearly with the inflaton. The latter contribution would be present even in the absence of particle production, when αk=1, βk=0.In order to isolate the effects of particle production on the inflaton fluctuations, we would like to subtract off the contribution to the 2-point correlation function (62) which is coming from the quantum vacuum fluctuations of χ. This subtraction also has the effect of rendering the power spectrum (63) finite, since it extracts the usual UV divergent contribution associated with the Minkowski-space vacuum fluctuations.As a step towards renormalizing the 2-point correlation function of inflaton fluctuations from rescattering (62), let us first consider the simpler problem of renormalizing the 2-point function of the gaussian field χ. We defined the renormalized 2-point function in momentum space as follows:(65)〈ξk1χ(t1)ξk2χ(t2)〉ren=〈ξk1χ(t1)ξk2χ(t2)〉-〈ξk1χ(t1)ξk2χ(t2)〉in.
In (65) the quantity 〈ξk1χ(t1)ξk2χ(t2)〉in is the contribution which would be present even in the absence of particle production, computed by simply taking solution (38) with αk=1, βk=0. Explicitly, we have(66)〈ξk1χ(t1)ξk2χ(t2)〉in=δ(3)(k1+k2)fk1(t1)fk2⋆(t2),
where fk are the adiabatic solutions (37).To see the impact of this subtraction, let us consider the renormalized variance for the iso-inflaton field,〈χ2〉. Employing prescription (65) we have(67)〈χ2(τ,x)〉ren=∫d3k(2π)3a2(τ)[|χk(τ)|2-12ωk(τ)]=〈χ2(τ,x)〉-δM,
where δM is the contribution from the Coleman-Weinberg potential. This proves that our prescription reproduces the scheme advocated in [80]. The renormalized variance (67) is finite and may be computed explicitly using our solutions (45). We find(68)〈χ2(t,x)〉ren≅nχa-3g|ϕ-ϕ0|,
where(69)nχ≡∫d3k(2π)3nk~k⋆3
is the total co-moving number density of produced χ particles. Result (68) was employed in [64] to quantify the effect of backreaction on the inflaton condensate in the mean field treatment (14). Hence, the renormalization scheme (65) was implicit in that calculation also.At the level of the 2-point function, our renormalization scheme is tantamount to assuming that Coleman-Weinberg corrections are already absorbed into the definition of the inflaton potential,V(ϕ). In general, such corrections might steepen V(ϕ) and spoil slow-roll inflation. Here, we assume that this problem has already been dealt with, either by fine-tuning the bare inflaton potential or else by including extended SUSY (which can minimize dangerous corrections). See also [80] for a related discussion.Having established a scheme for renormalizing the 2-point function of the gaussian fieldχ, it is now straightforward to consider higher-order correlation functions. We simply rewrite the 4-point function as a product of 2-point functions using Wick's theorem. Next, each Wick contraction is renormalized as (65). Applying this prescription to (62) amounts to(70)〈ξk1ϕ(τ)ξk2ϕ(τ)〉ren=2g2k⋆4(2π)3δ(3)(k1+k2)×∫dτ′dτ′′t(τ′)t(τ′′)a(τ′)a(τ′′)Gk1(τ-τ′)Gk2(τ-τ′′)×∫d3k′[χk1-k′(τ′)χk1-k′⋆(τ′′)-fk1-k′(τ′)fk1-k′⋆(τ′′)]×[χk′(τ′)χk′⋆(τ′′)-fk′(t′)fk′⋆(τ′′)],
where fk(t) are the adiabatic solutions defined in (37).


## 4.6. Power Spectrum

We are now in a position to compute the renormalized power spectrum of inflation fluctuations generated by rescattering,Pϕresc(k). We renormalize the 2-point correlator of the inflaton fluctuations generated by rescatter according to (70) and extract the power spectrum by comparison to (63). We have relegated the technical details to Appendix B and here we simply state the final result:(71)Pϕresc(k)=g2k3k⋆16π5[e-πk2/(2k⋆2)22(I2(k,τ)2+|I1(k,τ)|2)+[e-πk2/(4k⋆2)+122e-3πk2/(8k⋆2)]×(I2(k,τ)2-Re[I1(k,τ)])+[8233e-πk2/(3k⋆2)+4255e-3πk2/(5k⋆2)]×g2k3k⋆16π5Im[I1(k,τ)I2(k,τ)]],
where the functions I1, I2 are the curved space generalization of the characteristic integrals defined in [12]. Explicitly we have (72)I1(k,τ)=1a(τ)∫dτ′Gk(τ-τ′)eik⋆2t2(τ′),I2(k,τ)=1a(τ)∫dτ′Gk(τ-τ′).
The characteristic integral I2 can be evaluated analytically; however, the resulting expression is not particularly enlightening. Evaluation of the integral I1 requires numerical methods. More details are in Appendix B. Equation (71) is the main result of this section.


## 4.7. Comparison to Lattice Field Theory Simulations

In Section3 the results of our analytical formalism were plotted alongside the output of fully nonlinear HLattice simulations. It is evident from Figures 4, 5, and 6 that the agreement between these approaches is extremely good, even very late into the evolution and in the regime g2~1. The consistency of perturbative quantum field theory analytics and nonlinear classical lattice simulations provides a highly nontrivial check on our calculation.


## 4.8. The Bispectrum

So far, we have shown how to compute analytically the power spectrum generated by particle production, rescattering, and IR cascading during inflation in model (5). We found that IR cascading leads to a bump-like contribution to the primordial power spectrum of the inflaton fluctuations. However, this same dynamics must also have a nontrivial impact on nongaussian statistics, such as the bispectrum. Indeed, it is already evident from our previous analysis that the inflaton fluctuations generated by rescattering may be significantly nongaussian. From expression (62) we see that the particular solution (due to rescattering) which is bi-linear is the gaussian field χ.We define the bispectrum of the inflaton field fluctuations in terms of the three-point correlation function as(73)〈ξk1ϕξk2ϕξk3ϕ(τ)〉=(2π)3a3(τ)δ(k1+k2+k3)Bϕ(ki).
The factor a3 appears in (73) to cancel the explicit factors of a-1 in our convention (50) for the Fourier transform. It is well known that the nongaussianity associated with the usual quantum vacuum fluctuations of the inflaton is negligible [16–18]; therefore, when evaluating the bispectrum (73) we consider only the particular solution (62) which is due to rescattering. Carefully carrying out the Wick contractions, we find the following result for the renormalized 3-point function:(74)〈ξk1ϕξk2ϕξk3ϕ(τ)〉ren=4g3k⋆6(2π)9/2δ(k1+k2+k3)∏i=13∫dτit(τi)a(τi)Gki(τ-τi)×∫d3p[χk1-p(τ1)χk1-p⋆(τ2)-fk1-p(τ1)fk1-p⋆(τ2)]×[χk3+p(τ2)χk3+p⋆(τ3)-fk3+p(τ2)fk3+p⋆(τ3)]×[χp(τ1)χp⋆(τ3)-fp(τ1)fp⋆(τ3)]+(k2⟷k3),
where the modes χk are defined by (38) and fk are the adiabatic solutions (37). On the last line of (74) we have labeled schematically terms which are identical to the preceding three lines, only with k2 and k3 interchanged. One may verify that this expression is symmetric under interchange of the momenta ki by changing dummy variables of integration.It is now straightforward (but tedious) to plug expressions (41) and (45) into (74) and evaluate the integrals. This computation is tractable analytically because the time and phase-space integrals decouple. The bispectrum is then extracted by comparison to (73). This computation is carried out in detail in [14] where we have shown that Bϕ(ki) peaks only over triangles with a characteristic size, corresponding to the location of the bump in the power spectrum. This result is easy to understand on physical grounds, all the dynamics of rescattering and IR cascading take place over a very short time near the moment t=0. Hence, the effect of this dynamics on the primordial fluctuations must be limited to scales leaving the horizon near the time when particle production occurs.We will provide a cursory discussion of the bispectrumBϕ(ki) in Section 8 when we discuss nongaussianity from particle production during inflation.


## 4.9. Inclusion of a Bare Iso-Inflaton Mass

In passing, it may be interesting to consider how the analysis of this section is modified in the case that our prototype action (5) is supplemented by a bare mass term for the iso-inflaton field of the form δℒ=-(1/2)μ2χ2. Thus, in place of (5) suppose that we consider the model(75)ℒ=-12(∂ϕ)2-V(ϕ)-12(∂χ)2-12μ2χ2-g22(ϕ-ϕ0)2χ2.
Now the χ particles do not become massless at the point ϕ=ϕ0, but rather the effective mass-squared(76)mχ2=μ2+g2(ϕ-ϕ0)2
reaches a minimum value μ2 (which we assume to be positive). Such a correction may arise due to a variety of effects and will reduce the impact of particle production and IR cascading on the observable cosmological fluctuations.Let us briefly consider how the additional bare mass term in (75) alters the dynamics of particle production. The iso-inflaton fluctuations now obey the equation(77)χ̈+3Hχ̇-∇⃗2a2χ+[μ2+k⋆4t2]χ=0
rather than (27). This equation was solved in [80] in the regime where particle production is fast compared to the expansion time. (In the opposite regime, which corresponds to a fine-tuned coupling g2≪10-7, the iso-inflaton will be light for a significant fraction of inflation. In that case the theory (5) must be considered as a multifield inflation model and one can no longer consistently assume 〈χ〉=0. In other words, relaxing the assumption of fast particle production significantly changes the scenario under consideration, and we do not pursue this possibility any further.) The occupation number of produced χ particles is(78)nk=re-πk2/k⋆2
which differs from our previous result (13) by the suppression factor(79)r≡e-πμ2/k⋆2≤1.
Therefore, the effect of the inclusion of a bare mass for the iso-inflaton is to suppress the number density of produced χ particles by an amount r. This suppression reflects the reduced phase space of produced particles: the adiabaticity condition |ωk′/ωk2|≪1 is violated only for modes with k<k⋆2-μ2.The reduction ofnk translates into a suppression for the n-point correlation functions of the iso-inflaton. For example, the renormalized variance 〈χ2〉ren∝∫d3knk is suppressed by a factor of r. The power spectrum of inflaton fluctuations generated by rescattering, Pϕresc, is proportional to the 4-point correlator of χ and hence picks up a suppression factor of order r2. Similarly, the bispectrum Bϕ is proportional to the 6-point correlator of χ and must be reduced by a factor of order r3. The condition(80)μ2≪k⋆2
is equivalent to r≅1 and ensures that the addition of a bare iso-inflaton mass will have a negligible impact on any observable.For models obtained from string theory or supergravity (SUGRA), it is natural to expectμ of the order of the Hubble scale during inflation [116–118]. In the context of SUGRA, the finite energy density driving inflation breaks SUSY and induces soft scalar potentials with curvature of order Vsoft′′~μ2~H2 [117]. In the case of string theory, many scalars are conformally coupled to gravity [118] through an interaction of the form δℒ=-(1/12)Rχ2 where the Ricci scalar is R~H2 during inflation. More generally, any nonminimal coupling δℒ=-(ξ/2)Rχ2 between gravity and the iso-inflaton will induce a contribution of order H to the effective mass of χ, as long as ξ=𝒪(1). In all models where μ2~H2 condition (80) is satisfied for reasonable values of the coupling g2>10-7, see; (11). Thus, we expect that corrections of the form δℒ=-(1/2)μ2χ2 will not alter our results in a wide variety of well-motivated models.


## 5. Cosmological Perturbation Theory

In Section4 we developed an analytical theory of particle production and IR cascading during inflation which is in very good agreement with nonlinear lattice field theory simulations. However, this formalism suffers from a neglect of metric perturbations, and, consequently, we were unable to rigorously discuss the gauge invariant curvature perturbation ζ. (This variable is related to the quantity ℛ, defined in Section 3.3 as ζ≅-ℛ on large scales.) Hence, the reader may be concerned about gauge ambiguities in our results. In this section we address such concerns, showing that metric perturbations may be incorporated in a straightforward manner and that their consistent inclusion does not change our results in any significant way. We will do so by showing explicitly that, with appropriate choice of gauge, (49) and (27) for the fluctuations of the inflaton and iso-inflaton still hold, to first approximation. We will also go beyond our previous analysis by explicitly showing that in this same gauge the spectrum of the curvature fluctuations, Pζ, is trivially related to the spectrum of inflaton fluctuations, Pϕ.To render the analysis tractable we would like to take full advantage of the results derived in the last section. To do so, we employ the Seery et al. formalism for working directly with the field equations [119] and make considerable use of results derived by Malik in [120, 121]. (Note that our notations differ somewhat from those employed by Malik. The reader is therefore urged to take care in comparing our formulae.)We expand the inflaton and iso-inflaton fields up to second order in perturbation theory as(81)ϕ(τ,x)=ϕ(τ)+δ1ϕ(τ,x)+12δ2ϕ(τ,x),χ(τ,x)=δ1χ(τ,x)+12δ2χ(τ,x).
The perturbations are defined to average to zero 〈δnϕ〉=〈δnχ〉=0 so that 〈ϕ(t,x)〉=ϕ(t) and 〈χ(t,x)〉=0. (The condition 〈χ〉=0 is ensured by the fact that mχ≫H for nearly the entire duration of inflation.)We employ the flat slicing and threading throughout this section. With this gauge choice the perturbed metric takes the form(82)g00=-a2(1+2ψ1+ψ2),g0i=a2∂i[B1+12B2],gij=a2δij,
so that spatial hypersurfaces are flat. Note also that in this gauge the field perturbations δnϕ, δnχ coincide with the Sasaki-Mukhanov variables [122] at both first and second order.This perturbative approach, of course, neglects the momentary slowdown of the inflaton background due to backreaction. However, we have already shown in Section3.3 that backreaction has a tiny impact on the observable cosmological perturbations (see also [12]).


### 5.1. Gaussian Perturbations

In [120] Malik has derived closed-form evolution equations for the field perturbations δnϕ, δnχ at both first order (n=1) and, second  (n=2) order in perturbation theory. Let us first study the Gaussian perturbations. The closed-form Klein-Gordon equation for δ1ϕ derived in [120] can be written as(83)δ1ϕ′′+2ℋδ1ϕ′-∇⃗2δ1ϕ+[a2m2-3(ϕ′Mp)2]δ1ϕ=0.
Following our previous analysis we expand the first-order perturbation in terms of annihilation/creation operators as(84)δ1ϕ(t,x)=∫d3k(2π)3/2[bkδ1ϕk(τ)a(τ)eik·x+h.c.],
where h.c. denotes the Hermitian conjugate of the preceding term, and we draw the attention of the reader to the the explicit factor of a-1 in our definition of the Fourier transform. Working to leading order in slow-roll parameters we have(85)δ1ϕk′′+[k2+1τ2(-2+3η-9ϵ)]δ1ϕk=0.
This equation coincides exactly with (55), and the properly normalized solutions again take the form of (57). The only difference is that the order of the Hankel function, ν, is now given by(86)ν≅32-η+3ϵ
rather than by (56). The power spectrum of the Gaussian fluctuations is, again, given by (58). The correction to the order of the Hankel function ν translates into a correction to the spectral index: instead of (59) we now have(87)ns-1=2η-6ϵ
which is precisely the standard result [123].Thus, as far as the quantum vacuum fluctuations of the inflaton are concerned, the only impact of consistently including metric perturbations is an𝒪(ϵ) correction to the spectral index ns.Let us now turn our attention to the first-order fluctuations of the iso-inflaton. The closed-form Klein-Gordon equation forδ1χ derived in [120] can be written as(88)δ1χ′′+2ℋδ1χ′-∇⃗2δ1χ+a2k⋆4t2(τ)δ1ϕ=0.
This coincides exactly with (27), which we have already solved. The fact that linear perturbations of χ do not couple to the metric fluctuations follows from the condition 〈χ〉=0.


### 5.2. Non-Gaussian Perturbations

Now let us consider the second-order perturbation equations. The closed-form Klein-Gordon equation forδ2ϕ derived in [120] can be written as(89)δ2ϕ′′+2ℋδ2ϕ′-∇⃗2δ2ϕ+[a2m2-3(ϕ′Mp)2]δ2ϕ=J(τ,x).
As usual, the left-hand side is identical to the first-order equation (83) while the source term J is constructed from a bi-linear combination of the first order quantities δ1ϕ and δ1χ. In order to solve (89) we require explicit expressions for the Green function Gk and the source term J. The Green function is trivial for the case at hand; it is still given by our previous result (56), provided that one takes into account the fact that the order of the Hankel functions ν is now given by (86) rather than (56). In other words, the Green function for the non-Gaussian perturbations (89) differs from the result obtained neglecting metric perturbations only by 𝒪(ϵ) corrections.Next, we would like to consider the source term,J, appearing in (89). Schematically, we can split the source into contributions bi-linear in the Gaussian inflaton fluctuation δ1ϕ and contributions bi-linear in the iso-inflaton δ1χ:(90)J=Jϕ+Jχ.
The contribution Jϕ would be present even in the absence of the iso-inflaton. These correspond, physically, to the usual nongaussian corrections to the inflaton vacuum fluctuations coming from self-interactions. This contribution to the source is well studied in the literature and is known to contribute negligibly to the bispectrum [119]. Thus, in what follows, we will ignore Jϕ.On the other hand, the contributionJχ appearing in (90) depends only on the iso-inflaton fluctuations δ1χ. This contribution can be understood, physically, as generating nongaussian inflaton fluctuations δ2ϕ by rescattering of the produced χ particles off the condensate. Hence, the contribution Jχ may source large nongaussianity and is most interesting for us. It is straightforward to compute Jχ explicitly for our model using the general results of [120]. We find(91)Jχ=-2a2g2(ϕ-ϕ0)(δ1χ)2±2ϵMp[-a2g2(ϕ-ϕ0)2(δ1χ)2-12(∇⃗δ1χ)2+∇-2(∂i(δ1χ)∇⃗2∂i(δ1χ)+∇⃗2(δ1χ)∇⃗2(δ1χ)+δ1χ′∇⃗2δ1χ+(∇⃗δ1χ′)2)],
where the upper sign is for ϕ′>0 and the lower sign is for ϕ′<0. Notice that the contributions to Jχ on the third and fourth line of (91) contain the inverse spatial Laplacian ∇-2 and are thus nonlocal. These terms all contain at least as many gradients as inverse gradients, and hence the large-scale limit is well defined. In [124] it is was argued that these terms nearly always contribute negligibly to the curvature perturbation on large scales.Let us now examine the structure of the iso-inflaton sourceJχ, (91). The first line of (91) goes like a2g2(ϕ-ϕ0)(δ1χ)2. This coincides exactly with the source term in (49) which was already studied in Section 4. On the other hand, the terms on the second, third and fourth lines of (91) are new. These represent corrections to IR cascading which arise due to the consistent inclusion of metric perturbations. We will now argue that these “extra” terms are negligible as compared to the first line. If we denote the energy density in gaussian iso-inflaton fluctuations as ρχ~mχ2(δ1χ)2 then, by inspection, we see that the first line of (91) is parametrically of order ρχ/|ϕ-ϕ0| while the remaining terms are of order ϵρχ/Mp. Hence, we expect the first term to dominate for the field values ϕ≅ϕ0 which are relevant for IR cascading. This suggests that the dominant contribution to Jχ is the term which we have already taken into account in Section 4.Let us now make this argument more quantitative. Inspection reveals that the only “new” contribution to (91) which has any chance of competing with the “old” term a2g2(ϕ-ϕ0)(δ1χ)2 is the one proportional to ϵa2g2(ϕ-ϕ0)2(δ1χ)2/Mp (the first term on the second line). This new correction has the possibility of becoming significant because it grows after particle production, as ϕ rolls away from ϕ0. This growth, which reflects the fact that the energy density in the χ particles increases as they become more massive, cannot persist indefinitely. Within a few e-foldings of particle production the iso-inflaton source term must behave as Jχ~a-3, corresponding to the volume dilution of non-relativistic particles. Hence, in order to justify the analysis of Section 4 we must check that the term (92)Jnew~ϵMpa2g2(ϕ-ϕ0)2(δ1χ)2
does not dominate over the term which we have already considered(93)Jold~a2g2(ϕ-ϕ0)(δ1χ)2
during the relevant time HΔt=𝒪(1) after particle production. It is straightforward to show that(94)JoldJnew~Mpϵ1ϕ-ϕ0~MpHϕ̇ϵ1N~1ϵ1N,
where N=Ht is the number of e-foldings elapsed from particle production to the time when IR cascading has completed. Hence, N=𝒪(1), and we conclude that the second, third, and fourth lines of (91) are (at least) slow-roll suppressed as compared to the first line.In summary, we have shown that consistent inclusion of metric perturbations yields corrections to the inflaton fluctuationsδϕ which fall into two classes.(1)
Slow-roll suppressed corrections to the inflaton vacuum fluctuationsδ1ϕ (these amount to changing the definition of ν in solution (57)). These corrections have two physical effects. First, they yield an 𝒪(ϵ) correction to the spectral index. Second, they modify the propagator Gk by an 𝒪(ϵ) correction.(2)
Corrections to the sourceJ for the nongaussian inflaton perturbation δ2ϕ. These corrections are the second, third, and fourth lines of (91)) which, as we have seen, are slow-roll suppressed.It should be clear that neither of these corrections alters our previous analysis in any significant way.


### 5.3. Correlators

So far, we have shown that a consistent inclusion of metric perturbations does not significantly alter our previous results for the field perturbations. Specifically,δ1χ is identical to our previous solution of (27) for the iso-inflaton while δ1ϕ coincides with the homogeneous solution of (49), up to slow-roll corrections. At second order in perturbation theory, we have seen that(95)δ2ϕ=∫d4x′G(x-x′)Jχ(x′)+𝒪[(δ1ϕ)2].
To a leading order in slow roll, Jχ≅-2a2g2(ϕ-ϕ0)(δ1χ)2, and the first term coincides with our previous result for the particular solution of (49). The terms of order (δ1ϕ)2 represent nongaussian corrections to the vacuum fluctuations from inflation (coming from self-interactions of δϕ and the nonlinearity of gravity). These would be present even in the absence of particle production and are known to have a negligible impact on the spectrum and bispectrum [119].We are ultimately interested in the connectedn-point correlation functions of δϕ. For example, the 2-point function 〈(δϕ)2〉 gets a contribution of the form 〈(δ1ϕ)2〉 which gives the usual nearly scale-invariant large-scale power spectrum from inflation. The cross-term 〈δ1ϕδ2ϕ〉 is of order 〈(δ1ϕ)4〉 and represents a negligible “loop” correction to the scale-invariant spectrum from inflation. (The cross-term does not involve the iso-inflaton since δ1ϕ and δ1χ are statistically independent.) Finally, there is a contribution 〈(δ2ϕ)2〉 which involves terms of order 〈χ4〉 coming from rescattering and terms of order 〈(δ1ϕ)4〉 which represent (more) loop corrections to the scale-invariant spectrum from inflation. Thus, we can schematically write(96)Pϕ(k)=Pϕvac(k)[1+(loops)]+Pϕresc(k).
Here Pϕvac~kns-1 is the usual nearly scale-invariant spectrum from inflation and Pϕresc is the bump-like contribution from rescattering and IR cascading which we have studied in the previous section. The “loop” corrections to Pkvac(k) have been studied in detail in the literature (see, e.g., [125, 126]) and are known to be negligible in most models.We can also make a similar schematic decomposition of the bispectrum by considering the structure of the 3-point correlator〈(δϕ)3〉. Following our previous line of reasoning, it is clear that the dominant contribution comes from rescattering and is of order 〈χ6〉. The terms involving 〈(δ1ϕ)3〉, on the other hand, represent the usual nongaussianity generated during single field slow-roll inflation and are known to be small [119].


### 5.4. The Curvature Perturbation

Ultimately one wishes to compute not the field perturbations,δnϕ, δnχ, but rather the gauge invariant curvature fluctuation, ζ. We expand this in perturbation theory in the usual manner:(97)ζ=ζ1+12ζ2.
In [121] Malik has derived expressions for the large scale curvature perturbation in terms of the Sasaki-Mukhanov variables at both first order and second order in perturbation theory. We remind the reader that in the flat slicing (which we employ) the Sasaki-Mukhanov variable for each field simply coincides with the field perturbation (i.e.,Qϕ=δϕ and Qχ=δχ).At first order in perturbation theory the iso-inflaton does not contribute to the curvature perturbation (since〈χ〉=0), and we have(98)ζ1=-ℋϕ′δ1ϕ.At second order in perturbation theory the expression for the curvature perturbation is more involved. Using the results of [121] and working to leading order in slow-roll parameters we find(99)ζ2≅-ℋϕ′[δ2ϕ-δ2ϕ′3ℋ]+13(ϕ′)2[(δ1χ′)2+a2g2v2t2(τ)(δ1χ)2]+13(ϕ′)2[(δ1ϕ′)2+a2m2(δ1ϕ)2]+2(ℋδ1ϕϕ′)2.
Let us discuss the various contributions to this equation. The third line contributes to the nongaussianity of the vacuum fluctuations during inflation. (The term 2ζ12 on the last line (which would seem to contribute 𝒪(1) to the nonlinearity parameter fNL) stems from using the Malik and Wands [127] definition of ζ2. It can be related to the definition of Lyth and Rodríguez [128] (which also agrees with Maldacena) as ζ2=ζ2LR+2ζ12. Thus, in the absence of particle production (99) would predict fNL~𝒪(ϵ,η), in agreement with [16–18]. See also [19, 20] for a related discussion.) These terms are known to be negligible [16–18, 119].Next, we consider the second line of (99). This represents the direct contribution of the gaussian fluctuations δ1χ to the curvature perturbation. This contribution is tiny since the χ particles are extremely massive for nearly the entire duration of inflation, and hence δ1χ~a-3/2 (see also [54] for a related discussion). The smallness of this contribution to ζ can be understood physically by noting that the super-horizon iso-curvature fluctuations in our model are negligible.Finally, let us consider the contribution on the first line of (99). This contribution is the most interesting. To make contact with observations we must compute the curvature perturbation at late times and on large scales. In Section 4 we have already shown that δnϕ is constant on large scales and at late times for both n=1 and n=2. This is the expected result: the curvature fluctuations are frozen far outside the horizon and in the absence of entropy perturbations [129]. (Note that, in some cases, the curvature fluctuations may evolve significantly after horizon exit [130, 131]. This is a concern in models where there are significant violations of slow roll. In Section 3.3 we have already shown that the transient violation of slow roll has a negligible effect on the curvature fluctuations in our model. Hence, the result ζn~δnϕ~const far outside the horizon is consistent with previous studies.) Hence δ2ϕ′ is completely negligible and the first term on the first line of (99) must dominate over the second term. We conclude that, at late times and on large scales, the second order curvature perturbation is very well approximated by (100)ζ2≅-ℋϕ′δ2ϕ+⋯.In summary, we have shown that the power spectrum of curvature fluctuations from inflation in model (5) is trivially related to the power spectrum of inflaton fluctuations(101)Pζ(k)≅H2ϕ̇2Pϕ(k)=12ϵMp2Pϕ(k)
at both first order and second order in cosmological perturbation theory. This relation is valid at late times and for scales far outside the horizon. Curvature spectrum (101) may be written as(102)Pζ(k)=Pζvac(k)[1+(loops)]+Pζresc(k).
The power spectrum of the inflaton vacuum fluctuations agrees with the usual result obtained in linear theory [123](103)Pζvac(k)≅H28π2ϵMp2(kaH)2η-6ϵ.
In (102) we have schematically labeled the corrections arising from the second line of (99) and the source Jϕ as “loop”. These are nongaussian corrections to the inflaton vacuum fluctuations arising from self-interactions of the inflaton and also the nonlinearity of gravity. Such corrections are negligible. The most interesting contribution to the power spectrum (102) is due to rescattering, Pζresc(k). This quantity is proportional to our previous result (71).In passing, notice that the bispectrumBϕ (defined by (73)) of inflaton fluctuations will differ from the bispectrum B of the curvature fluctuations (defined by (3)) only by a simple rescaling:(104)B(ki)≅-(Hϕ̇)3Bϕ(ki)=-1(2ϵ)3/2Mp3Bϕ(ki).
The dominant contribution to Bϕ comes from rescattering effects 〈δ2ϕ3〉~〈δ1χ6〉.The analysis of this section justifies our neglect of metric fluctuations in Section4.


## 5.1. Gaussian Perturbations

In [120] Malik has derived closed-form evolution equations for the field perturbations δnϕ, δnχ at both first order (n=1) and, second  (n=2) order in perturbation theory. Let us first study the Gaussian perturbations. The closed-form Klein-Gordon equation for δ1ϕ derived in [120] can be written as(83)δ1ϕ′′+2ℋδ1ϕ′-∇⃗2δ1ϕ+[a2m2-3(ϕ′Mp)2]δ1ϕ=0.
Following our previous analysis we expand the first-order perturbation in terms of annihilation/creation operators as(84)δ1ϕ(t,x)=∫d3k(2π)3/2[bkδ1ϕk(τ)a(τ)eik·x+h.c.],
where h.c. denotes the Hermitian conjugate of the preceding term, and we draw the attention of the reader to the the explicit factor of a-1 in our definition of the Fourier transform. Working to leading order in slow-roll parameters we have(85)δ1ϕk′′+[k2+1τ2(-2+3η-9ϵ)]δ1ϕk=0.
This equation coincides exactly with (55), and the properly normalized solutions again take the form of (57). The only difference is that the order of the Hankel function, ν, is now given by(86)ν≅32-η+3ϵ
rather than by (56). The power spectrum of the Gaussian fluctuations is, again, given by (58). The correction to the order of the Hankel function ν translates into a correction to the spectral index: instead of (59) we now have(87)ns-1=2η-6ϵ
which is precisely the standard result [123].Thus, as far as the quantum vacuum fluctuations of the inflaton are concerned, the only impact of consistently including metric perturbations is an𝒪(ϵ) correction to the spectral index ns.Let us now turn our attention to the first-order fluctuations of the iso-inflaton. The closed-form Klein-Gordon equation forδ1χ derived in [120] can be written as(88)δ1χ′′+2ℋδ1χ′-∇⃗2δ1χ+a2k⋆4t2(τ)δ1ϕ=0.
This coincides exactly with (27), which we have already solved. The fact that linear perturbations of χ do not couple to the metric fluctuations follows from the condition 〈χ〉=0.


## 5.2. Non-Gaussian Perturbations

Now let us consider the second-order perturbation equations. The closed-form Klein-Gordon equation forδ2ϕ derived in [120] can be written as(89)δ2ϕ′′+2ℋδ2ϕ′-∇⃗2δ2ϕ+[a2m2-3(ϕ′Mp)2]δ2ϕ=J(τ,x).
As usual, the left-hand side is identical to the first-order equation (83) while the source term J is constructed from a bi-linear combination of the first order quantities δ1ϕ and δ1χ. In order to solve (89) we require explicit expressions for the Green function Gk and the source term J. The Green function is trivial for the case at hand; it is still given by our previous result (56), provided that one takes into account the fact that the order of the Hankel functions ν is now given by (86) rather than (56). In other words, the Green function for the non-Gaussian perturbations (89) differs from the result obtained neglecting metric perturbations only by 𝒪(ϵ) corrections.Next, we would like to consider the source term,J, appearing in (89). Schematically, we can split the source into contributions bi-linear in the Gaussian inflaton fluctuation δ1ϕ and contributions bi-linear in the iso-inflaton δ1χ:(90)J=Jϕ+Jχ.
The contribution Jϕ would be present even in the absence of the iso-inflaton. These correspond, physically, to the usual nongaussian corrections to the inflaton vacuum fluctuations coming from self-interactions. This contribution to the source is well studied in the literature and is known to contribute negligibly to the bispectrum [119]. Thus, in what follows, we will ignore Jϕ.On the other hand, the contributionJχ appearing in (90) depends only on the iso-inflaton fluctuations δ1χ. This contribution can be understood, physically, as generating nongaussian inflaton fluctuations δ2ϕ by rescattering of the produced χ particles off the condensate. Hence, the contribution Jχ may source large nongaussianity and is most interesting for us. It is straightforward to compute Jχ explicitly for our model using the general results of [120]. We find(91)Jχ=-2a2g2(ϕ-ϕ0)(δ1χ)2±2ϵMp[-a2g2(ϕ-ϕ0)2(δ1χ)2-12(∇⃗δ1χ)2+∇-2(∂i(δ1χ)∇⃗2∂i(δ1χ)+∇⃗2(δ1χ)∇⃗2(δ1χ)+δ1χ′∇⃗2δ1χ+(∇⃗δ1χ′)2)],
where the upper sign is for ϕ′>0 and the lower sign is for ϕ′<0. Notice that the contributions to Jχ on the third and fourth line of (91) contain the inverse spatial Laplacian ∇-2 and are thus nonlocal. These terms all contain at least as many gradients as inverse gradients, and hence the large-scale limit is well defined. In [124] it is was argued that these terms nearly always contribute negligibly to the curvature perturbation on large scales.Let us now examine the structure of the iso-inflaton sourceJχ, (91). The first line of (91) goes like a2g2(ϕ-ϕ0)(δ1χ)2. This coincides exactly with the source term in (49) which was already studied in Section 4. On the other hand, the terms on the second, third and fourth lines of (91) are new. These represent corrections to IR cascading which arise due to the consistent inclusion of metric perturbations. We will now argue that these “extra” terms are negligible as compared to the first line. If we denote the energy density in gaussian iso-inflaton fluctuations as ρχ~mχ2(δ1χ)2 then, by inspection, we see that the first line of (91) is parametrically of order ρχ/|ϕ-ϕ0| while the remaining terms are of order ϵρχ/Mp. Hence, we expect the first term to dominate for the field values ϕ≅ϕ0 which are relevant for IR cascading. This suggests that the dominant contribution to Jχ is the term which we have already taken into account in Section 4.Let us now make this argument more quantitative. Inspection reveals that the only “new” contribution to (91) which has any chance of competing with the “old” term a2g2(ϕ-ϕ0)(δ1χ)2 is the one proportional to ϵa2g2(ϕ-ϕ0)2(δ1χ)2/Mp (the first term on the second line). This new correction has the possibility of becoming significant because it grows after particle production, as ϕ rolls away from ϕ0. This growth, which reflects the fact that the energy density in the χ particles increases as they become more massive, cannot persist indefinitely. Within a few e-foldings of particle production the iso-inflaton source term must behave as Jχ~a-3, corresponding to the volume dilution of non-relativistic particles. Hence, in order to justify the analysis of Section 4 we must check that the term (92)Jnew~ϵMpa2g2(ϕ-ϕ0)2(δ1χ)2
does not dominate over the term which we have already considered(93)Jold~a2g2(ϕ-ϕ0)(δ1χ)2
during the relevant time HΔt=𝒪(1) after particle production. It is straightforward to show that(94)JoldJnew~Mpϵ1ϕ-ϕ0~MpHϕ̇ϵ1N~1ϵ1N,
where N=Ht is the number of e-foldings elapsed from particle production to the time when IR cascading has completed. Hence, N=𝒪(1), and we conclude that the second, third, and fourth lines of (91) are (at least) slow-roll suppressed as compared to the first line.In summary, we have shown that consistent inclusion of metric perturbations yields corrections to the inflaton fluctuationsδϕ which fall into two classes.(1)
Slow-roll suppressed corrections to the inflaton vacuum fluctuationsδ1ϕ (these amount to changing the definition of ν in solution (57)). These corrections have two physical effects. First, they yield an 𝒪(ϵ) correction to the spectral index. Second, they modify the propagator Gk by an 𝒪(ϵ) correction.(2)
Corrections to the sourceJ for the nongaussian inflaton perturbation δ2ϕ. These corrections are the second, third, and fourth lines of (91)) which, as we have seen, are slow-roll suppressed.It should be clear that neither of these corrections alters our previous analysis in any significant way.


## 5.3. Correlators

So far, we have shown that a consistent inclusion of metric perturbations does not significantly alter our previous results for the field perturbations. Specifically,δ1χ is identical to our previous solution of (27) for the iso-inflaton while δ1ϕ coincides with the homogeneous solution of (49), up to slow-roll corrections. At second order in perturbation theory, we have seen that(95)δ2ϕ=∫d4x′G(x-x′)Jχ(x′)+𝒪[(δ1ϕ)2].
To a leading order in slow roll, Jχ≅-2a2g2(ϕ-ϕ0)(δ1χ)2, and the first term coincides with our previous result for the particular solution of (49). The terms of order (δ1ϕ)2 represent nongaussian corrections to the vacuum fluctuations from inflation (coming from self-interactions of δϕ and the nonlinearity of gravity). These would be present even in the absence of particle production and are known to have a negligible impact on the spectrum and bispectrum [119].We are ultimately interested in the connectedn-point correlation functions of δϕ. For example, the 2-point function 〈(δϕ)2〉 gets a contribution of the form 〈(δ1ϕ)2〉 which gives the usual nearly scale-invariant large-scale power spectrum from inflation. The cross-term 〈δ1ϕδ2ϕ〉 is of order 〈(δ1ϕ)4〉 and represents a negligible “loop” correction to the scale-invariant spectrum from inflation. (The cross-term does not involve the iso-inflaton since δ1ϕ and δ1χ are statistically independent.) Finally, there is a contribution 〈(δ2ϕ)2〉 which involves terms of order 〈χ4〉 coming from rescattering and terms of order 〈(δ1ϕ)4〉 which represent (more) loop corrections to the scale-invariant spectrum from inflation. Thus, we can schematically write(96)Pϕ(k)=Pϕvac(k)[1+(loops)]+Pϕresc(k).
Here Pϕvac~kns-1 is the usual nearly scale-invariant spectrum from inflation and Pϕresc is the bump-like contribution from rescattering and IR cascading which we have studied in the previous section. The “loop” corrections to Pkvac(k) have been studied in detail in the literature (see, e.g., [125, 126]) and are known to be negligible in most models.We can also make a similar schematic decomposition of the bispectrum by considering the structure of the 3-point correlator〈(δϕ)3〉. Following our previous line of reasoning, it is clear that the dominant contribution comes from rescattering and is of order 〈χ6〉. The terms involving 〈(δ1ϕ)3〉, on the other hand, represent the usual nongaussianity generated during single field slow-roll inflation and are known to be small [119].


## 5.4. The Curvature Perturbation

Ultimately one wishes to compute not the field perturbations,δnϕ, δnχ, but rather the gauge invariant curvature fluctuation, ζ. We expand this in perturbation theory in the usual manner:(97)ζ=ζ1+12ζ2.
In [121] Malik has derived expressions for the large scale curvature perturbation in terms of the Sasaki-Mukhanov variables at both first order and second order in perturbation theory. We remind the reader that in the flat slicing (which we employ) the Sasaki-Mukhanov variable for each field simply coincides with the field perturbation (i.e.,Qϕ=δϕ and Qχ=δχ).At first order in perturbation theory the iso-inflaton does not contribute to the curvature perturbation (since〈χ〉=0), and we have(98)ζ1=-ℋϕ′δ1ϕ.At second order in perturbation theory the expression for the curvature perturbation is more involved. Using the results of [121] and working to leading order in slow-roll parameters we find(99)ζ2≅-ℋϕ′[δ2ϕ-δ2ϕ′3ℋ]+13(ϕ′)2[(δ1χ′)2+a2g2v2t2(τ)(δ1χ)2]+13(ϕ′)2[(δ1ϕ′)2+a2m2(δ1ϕ)2]+2(ℋδ1ϕϕ′)2.
Let us discuss the various contributions to this equation. The third line contributes to the nongaussianity of the vacuum fluctuations during inflation. (The term 2ζ12 on the last line (which would seem to contribute 𝒪(1) to the nonlinearity parameter fNL) stems from using the Malik and Wands [127] definition of ζ2. It can be related to the definition of Lyth and Rodríguez [128] (which also agrees with Maldacena) as ζ2=ζ2LR+2ζ12. Thus, in the absence of particle production (99) would predict fNL~𝒪(ϵ,η), in agreement with [16–18]. See also [19, 20] for a related discussion.) These terms are known to be negligible [16–18, 119].Next, we consider the second line of (99). This represents the direct contribution of the gaussian fluctuations δ1χ to the curvature perturbation. This contribution is tiny since the χ particles are extremely massive for nearly the entire duration of inflation, and hence δ1χ~a-3/2 (see also [54] for a related discussion). The smallness of this contribution to ζ can be understood physically by noting that the super-horizon iso-curvature fluctuations in our model are negligible.Finally, let us consider the contribution on the first line of (99). This contribution is the most interesting. To make contact with observations we must compute the curvature perturbation at late times and on large scales. In Section 4 we have already shown that δnϕ is constant on large scales and at late times for both n=1 and n=2. This is the expected result: the curvature fluctuations are frozen far outside the horizon and in the absence of entropy perturbations [129]. (Note that, in some cases, the curvature fluctuations may evolve significantly after horizon exit [130, 131]. This is a concern in models where there are significant violations of slow roll. In Section 3.3 we have already shown that the transient violation of slow roll has a negligible effect on the curvature fluctuations in our model. Hence, the result ζn~δnϕ~const far outside the horizon is consistent with previous studies.) Hence δ2ϕ′ is completely negligible and the first term on the first line of (99) must dominate over the second term. We conclude that, at late times and on large scales, the second order curvature perturbation is very well approximated by (100)ζ2≅-ℋϕ′δ2ϕ+⋯.In summary, we have shown that the power spectrum of curvature fluctuations from inflation in model (5) is trivially related to the power spectrum of inflaton fluctuations(101)Pζ(k)≅H2ϕ̇2Pϕ(k)=12ϵMp2Pϕ(k)
at both first order and second order in cosmological perturbation theory. This relation is valid at late times and for scales far outside the horizon. Curvature spectrum (101) may be written as(102)Pζ(k)=Pζvac(k)[1+(loops)]+Pζresc(k).
The power spectrum of the inflaton vacuum fluctuations agrees with the usual result obtained in linear theory [123](103)Pζvac(k)≅H28π2ϵMp2(kaH)2η-6ϵ.
In (102) we have schematically labeled the corrections arising from the second line of (99) and the source Jϕ as “loop”. These are nongaussian corrections to the inflaton vacuum fluctuations arising from self-interactions of the inflaton and also the nonlinearity of gravity. Such corrections are negligible. The most interesting contribution to the power spectrum (102) is due to rescattering, Pζresc(k). This quantity is proportional to our previous result (71).In passing, notice that the bispectrumBϕ (defined by (73)) of inflaton fluctuations will differ from the bispectrum B of the curvature fluctuations (defined by (3)) only by a simple rescaling:(104)B(ki)≅-(Hϕ̇)3Bϕ(ki)=-1(2ϵ)3/2Mp3Bϕ(ki).
The dominant contribution to Bϕ comes from rescattering effects 〈δ2ϕ3〉~〈δ1χ6〉.The analysis of this section justifies our neglect of metric fluctuations in Section4.


## 6. Observational Constraints on Particle Production during Inflation

In the previous sections of this paper, we have demonstrated that particle production and IR cascading in model (4) leads to a bump-like feature in the primordial power spectrum. We now aim to determine the extent to which such a spectral distortion is compatible with current cosmological data. The results of this section have first appeared in [13].One motivation for this study is to determine whether features generated by particle production during inflation can explain some of the anomalies in the observed primordial power spectrum,P(k). A number of different studies have hinted at the possible presence of some localized features in the power spectrum [62, 73, 98–110], which are not compatible with the simplest power-law P(k)~kns-1 model. Although such glitches may simply be statistical anomalies [111], there is also the tantalizing possibility that they represent a signature of primordial physics beyond the simplest slow-roll inflation scenario. Forthcoming polarization data may play a crucial role in distinguishing between these possibilities [73]. However, in the meantime, it is interesting to determine the extent to which such features may be explained by microscopically realistic inflation models, such as (5).To answer this question we provide a simple semianalytic fitting function that accurately captures the shape of the feature generated by particle production and IR cascading during inflation. Next, we confront this modified power spectrum with a variety of observational data sets. We find no evidence for a detection; however, we note that observations are consistent with relatively large spectral distortions of the type predicted by model (4). If the feature is located on scales relevant for Cosmic Microwave Background (CMB) experiments, then its amplitude may be as large as 𝒪(10%) of the usual scale-invariant fluctuations, corresponding to g2~0.01. Our results translate into a ϕ0-dependent bound on the coupling g2 which is crucial in order to determine whether the nongaussian signal associated with particle production and IR cascading is detectable in future missions.We also consider the more complicated features which result from multiple bursts of particle production and IR cascading. Such features are a prediction of a number of string theory inflation models, including brane/axion monodromy [77–79]. For appropriate choice of the spacing between the features, we find that the constraint on g2 in this scenario is even weaker than the single-bump case.


### 6.1. A Simple Parametrization of the Power Spectrum

In [12] it was shown that particle production and IR cascading during inflation in model (4) generates a bump-like contribution to the primordial power spectrum. As shown in Figure 9, this feature can be fit with a very simple function Pbump~k3e-πk2/(2k⋆2). The bump-like contribution from IR cascading is complimentary to the usual (nearly) scale-invariant contribution to the primordial power spectrum Pvac~kns-1 coming from the quantum vacuum fluctuations of the inflaton. The total, observable, power spectrum in model (4) is simply the superposition of these two contributions: P(k)~kns-1+k3e-πk2/(2k⋆2) (see (64)). This simple formula can be motivated from analytical considerations [12] and provides a good fit to lattice field theory simulations near the peak of the feature and also in the IR tail. (This fitting formula does not capture the small oscillatory structure in the UV tail of the feature (see Figure 9) which does not concern us since that region is not phenomenologically interesting.)Figure 9
The bump-like features generated by IR cascading. We plot the feature power spectrumPbump(k) obtained from fully nonlinear lattice field theory simulations (the red points) and also the result of an analytical calculation (the dashed blue curve) using the formalism described in [14]. We also superpose the fitting function ~k3e-πk2/(2k⋆2) (the solid black curve) to illustrate the accuracy of this simple formula.It is straightforward to generalize this discussion to allow for multiple bursts of particle production during inflation. Suppose there are multiple pointsϕ=ϕi (i=1,…,n) along the inflationary trajectory where new degrees of freedom χi become massless:(105)ℒint=-∑i=0ngi22(ϕ-ϕi)χi2.
For each instant ti when ϕ=ϕi there will be an associated burst of particle production and subsequent rescattering of the produced massive χi off the condensate ϕ(t). Each of these events proceeds as described above and leads to a new bump-like contribution to the power spectrum. These features simply superpose owing to that fact that each field χi is statistically independent (so that the cross-terms involving χiχj with i≠j in the computation of the two-point function must vanish). Thus, we arrive at the following parametrization of the primordial power spectrum in models with particle production during inflation:(106)P(k)=As(kk0)ns-1+∑i=1nAi(πe3)3/2(kki)3e-(π/2)(k/ki)2,
where As is the amplitude of the usual nearly scale-invariant vacuum fluctuations from inflation, and k0 is the pivot, which we choose to be k0=0.002Mpc-1 following [132]. The constants Ai depend on the couplings gi2 and measure of the size of the features from particle production. We have normalized these amplitudes so that the power in the ith bump, measured at the peak of the feature, is given by Ai. The location of each feature, ki, is related to the number of e-foldings N from the end of inflation to the time when the ith burst of particle production occurs: roughly ln(ki/H)~Ni where N=Ni at the moment when ϕ=ϕi. From a purely phenomenological perspective the locations ki are completely arbitrary.We compare (106) to lattice field theory simulations in order to determine the amplitude Ai in terms of model parameters. We find(107)Ai≅10-6gi15/4.
The power-law form of this relation was determined by inspection of numerical results; however, it can also be motivated by analytical considerations. Assuming standard chaotic inflation V=m2ϕ2/2, we have tested this formula for g2=1,0.1,0.01, taking both ϕ0=28πMp and ϕ0=3.28πMp. We found agreement up to factors order unity in all cases.Theoretical consistency of our calculation of the shape of the feature bounds the coupling as10-7≲gi2≲1 [12]. Hence, power spectrum (106) can be obtained from sensible microphysics only when 10-20≲Ai≲10-6. This constraint still allows for a huge range of observational possibilities: near the upper bound the feature is considerably larger than the vacuum fluctuations while near the lower bound the feature is completely undetectable.Note that for each bump in (106) the IR tail Pbump→k3 as k→0 is similar to the feature considered by Hoi et al. in [108], consistent with causality arguments about the generation of curvature perturbations by local physics.


### 6.2. Data Sets and Analysis

The primordial power spectrum for our model is parametrized as (106). Our aim is to derive observational constraints on the various model parameters As, ns, ki, and Ai using CMB, galaxy power spectrum, and weak lensing data. To this end we use the cosmoMC package [133] to run Markov Chain Monte Carlo (MCMC) calculations to determine the likelihood of the cosmological parameters, including our new parameters Ai and ki. We employ the following data sets.Cosmic Microwave Background (CMB)
Our complete CMB data sets include WMAP-5yr [132, 134], BOOMERANG [135–137], ACBAR [138–141], CBI [142–145], VSA [146], DASI [147, 148], and MAXIMA [149]. We have included the Sunyaev-Zeldovic (SZ) secondary anisotropy [150, 151] for WMAP-5yr, ACBAR, and CBI data sets. The SZ template is obtained from hydrodynamical simulation [152]. Also included for theoretical calculation of CMB power spectra is the CMB lensing contribution.Type Ia Supernova (SN)
We employ the Union Supernova Ia data (307 SN Ia samples) from The Supernova Cosmology Project [153].Large Scale Structure (LSS)
The 2dF Galaxy Redshift Survey (2dFGRS) data [154] and Sloan Digital Sky Survey (SDSS) Luminous Red Galaxy (LRG) data release 4 [155] are utilized.Note that we have used the likelihood code based on the non-linear modeling by Tegmark et al.  [155] (marginalizing the bias b and the Q parameter). However with a large bump in the linear power spectrum, this naive treatment may be not sufficient to characterize the non-linear response to the feature on small scales. Ideally, this should be obtained from N-body simulations; however, such a study is beyond the scope of the current work.There are several other caveats on our results in the high-k regime. First, we assume linear bias for the galaxies, which may not be entirely safe at sufficiently small scales. Moreover, sharp features in the matter power spectrum can cause sharp features in the bias as a function of k.Keeping in mind these caveats our constraints on small scalesk≳0.1Mpc-1 should be taken with a grain of salt and considered as accurate only up to factors order unity.Weak Lensing (WL)
Five WL data sets are used in this paper. The effective survey areaAeff and galaxy number density neff of each survey are listed in Table 1.Table 1
Weak lensing data sets.
Data setsAeffneff(deg2)(arcmin-2)COSMOS [156, 157]1.640CFHTLS-wide [158, 159]2212GaBODS [160, 161]1312.5RCS [160, 161]538VIRMOS-DESCART [159, 162]8.515For COSMOS data we use the CosmoMC plug-in written by Lesgourgues [157], modified to do numerical marginalization on three nuisance parameters in the original code.For the other four weak lensing data sets we use the likelihood given in [163]. To calculate the likelihood we have written a CosmoMC plug-in code, with simplified marginalization on the parameters of galaxy number density function n(z). More details about this plug-in can be found in [164].As for the LSS data, for small scalesk≳0.1Mpc-1 there is the caveat that the nonlinear evolution of the power spectrum in the presence of bump-like distortions may not be treated accurately.


### 6.3. Observational Constraints: A Single Burst of Particle Production

We now present our results for the observational constraints on particle production during inflation, assuming two different scenarios.The minimal scenario to consider is a single burst of particle production during inflation, which corresponds to takingn=1 in (105). The power spectrum is given by (106) with n=1, and, with some abuse of notation, we denote k1≡kIR and A1≡AIR. The prior we have used for AIR is 0≤AIR≤25×10-10 and for kIR is -9.5≤ln[k/Mpc-1]≤1. The former condition ensures that the bump-like feature from IR cascading does not dominate over the observed scale-invariant fluctuations while the latter is necessary in order to have the feature in the observable range of scales. In Figure 10 we plot the marginalized posterior likelihood for the new parameters AIR and kIR describing the magnitude and location of the bump while in Table 2 we give the best fit values for the remaining (vanilla) cosmological parameters.Table 2
Constraints on the standard “vanilla” cosmological parameters for the single-bump model. All errors are 95.4% confidence level.
Ωbh20.0227-0.0010+0.0010Ωch20.1122-0.0044+0.0050θ1.0424-0.0043+0.0042τ0.08-0.03+0.03ns0.956-0.024+0.024ln[1010As]3.206-0.068+0.074ASZ1.62-0.74+0.76Ωm0.264-0.022+0.026σ80.807-0.030+0.034zre10.5-2.7+2.5H071.6-2.4+2.3Figure 10
Marginalized posterior likelihood contours for the parametersAIR and kIR (the magnitude and position of the feature, resp.) in the single-bump model. Black and grey regions correspond to parameter values allowed at 95.4% and 99.7% confidence levels, respectively. At small scales, to the right of the dashed vertical line, our results should be taken with a grain of salt since the nonlinear evolution of the power spectrum may not be modeled correctly in the presence of bump-like distortions.For very large scales≲Gpc-1, the data do not contain much information (due to cosmic variance), and hence the constraint on any modification of the power spectrum is weak. In this region the spectral distortion may be larger than 100% of the usual scale-invariant fluctuations and couplings g2 of order unity are allowed. For smaller scales k≳Gpc-1 the constraints are stronger, and we have, very roughly, AIR/As≲0.1 corresponding to g2≲0.01. For very small scales, k≳0.1Mpc-1, our constraints should be taken with a grain of salt since the nonlinear evolution of the power spectrum may not be modeled correctly in the presence of bump-like distortions. At small scales nonlinear effects tend to wipe out features of this type (see, e.g., [165]), and hence observational constraints for k≳0.1Mpc-1 may be weaker than what is presented in Figure 10. Note that in most of this nonlinear regime we find essentially no constraint on AIR, which is consistent with what would be expected in a more comprehensive treatment.The IR cascading bump in the primordial power spectrum will be accompanied by a corresponding nongaussian feature in the bispectrum, that will be discussed in more detail in the next section. From the perspective of potentially observing this signal it is most interesting if this feature is located on scales probed by CMB experiments. (There is also the fascinating possibility that the nongaussianity from IR cascading could show up in the large-scale structure as in [27, 166–168]. We leave a detailed discussion to future studies.) To get some intuition into what kinds of features in the CMB scales are still allowed by the data we focus on an example with AIR=2.5×10-10 which, using (107), corresponds to a reasonable coupling value g2~0.01. We take the bump to be located at kIR=0.01Mpc-1 and fix the remaining model parameters to As=2.44×10-9, ns=0.97 (which are compatible with the usual values). This sample bump in the power spectrum is illustrated in Figure 3(a) and is consistent with the data at 2σ. In Figure 3(b) we plot the associated angular CMB TT spectrum. This example represents a surprisingly large spectral distortion: the total power in the feature as compared to the scale invariant vacuum fluctuations is Pbump/Pvac~0.1, evaluated at the peak of the bump. We discuss the nongaussianity associated with this feature.


### 6.4. Observational Constraints: Multiple Bursts of Particle Production

Next, we consider a slightly more complicated scenario: multiple bursts of particle production leading many localized features in the power spectrum. For simplicity we assume that all bumps have the same magnitudeAi≡AIR, and we further suppose a fixed number of e-foldings δN between each burst of particle production. This implies that the location of the ith bump will be given by ki=e(i-1)Δk1 where Δ is a model parameter controlling the density of features. We take the number of bursts, n, to be sufficiently large that the series of features extends over the whole observable range. In the next section we will see that these assumptions are not restrictive and that many well-motivated models are consistent with this simple setup.Our multibump model, then, has three parameters:AIR, k1, and Δ. We take the prior on the amplitude to be AIR≤25×10-10 as in Section 6.3. If the features are very widely spaced, Δ≳1, then the constraint on each bump will obviously be identical to the results for the single-bump case presented in Section 6.3. Hence the most interesting case to consider is Δ≲1 so that the bumps are partially overlapping. Our prior for the density of features is therefore 0≤Δ≤1. Finally, the location of the first bump will be a historical accident in realistic models, hence we marginalize over all possible values of k1 and present our constraints and 2d likelihood plots in the space of AIR and Δ. This marginalized likelihood plot is presented in Figure 11. In Table 3 we present the best-fit values for the vanilla cosmological parameters.Table 3
Constraints on the standard “vanilla” cosmological parameters for the multiple-bump model. All error bars are 95.4% confidence level.
Ωbh20.0227-0.0009+0.0009Ωch20.1126-0.0044+0.0049θ1.0424-0.0043+0.0039τ0.078-0.026+0.031ns0.93-0.17+0.04ln[1010As]2.8-0.9+0.4ASZ1.60-0.76+0.77Ωm0.265-0.021+0.026σ80.807-0.030+0.034zre10.3-2.5+2.6H071.4-2.4+2.2Figure 11
Marginalized posterior likelihood contours for the parametersAIR and Δ (the feature amplitude and spacing, resp.) of the multiple-bump model. Black and grey regions correspond to values allowed at 95.4% and 99.7% confidence levels, respectively.From the likelihood plot, Figure11, there is evidently a preferred value of the feature spacing, roughly Δ~0.75, for which the constraints are the weakest. This can be understood as follows. For very high density Δ→0 the localized features from IR cascading smear together, and the total power spectrum (106) is P(k)~As(k/k0)ns-1+C where the size of the constant deformation scales linearly with the density of features: C∝Δ-1. Therefore, the upper bound on the amplitude AIR should scale linearly with Δ. Indeed, this linear trend is very evident from Figure 11 in the small-Δ regime. This linear behaviour must break down at some point since as the features become infinitely widely spaced the constraint on AIR must go to zero. This explains the bump in the likelihood plot, Figure 11, near Δ~0.75.In passing, notice that the behaviourP(k)~As(k/k0)ns-1+C for Δ≪1 also explains why the best-fit As in Table 3 is somewhat lower than the standard value and why the spectral tilt ns-1 is somewhat more red.To get some intuition for the kinds of multibump distortions that are allowed by the data, we consider an example withAIR=1×10-9, Δ=0.75 and fix the vanilla parameters to As=1.04×10-9, ns=0.93. This choice of parameters is consistent with the data at 2σ and corresponds to a reasonable coupling g2~0.02. In Figure 12 we plot the primordial power spectrum P(k) and also the CMB TT angular power spectrum for this example.(a) shows a sample multiple-bump distortion with amplitudeAIR=1×10-9 which corresponds to g2~0.02. The feature spacing is Δ=0.75. We take the vanilla parameters to be As=1.04×10-9, ns=0.93 so that the scale of inflation is slightly lower than in the standard scenario and the spectral tilt is slightly redder. (b) shows the CMB angular TT power spectrum for this example.
(a)(b)


## 6.1. A Simple Parametrization of the Power Spectrum

In [12] it was shown that particle production and IR cascading during inflation in model (4) generates a bump-like contribution to the primordial power spectrum. As shown in Figure 9, this feature can be fit with a very simple function Pbump~k3e-πk2/(2k⋆2). The bump-like contribution from IR cascading is complimentary to the usual (nearly) scale-invariant contribution to the primordial power spectrum Pvac~kns-1 coming from the quantum vacuum fluctuations of the inflaton. The total, observable, power spectrum in model (4) is simply the superposition of these two contributions: P(k)~kns-1+k3e-πk2/(2k⋆2) (see (64)). This simple formula can be motivated from analytical considerations [12] and provides a good fit to lattice field theory simulations near the peak of the feature and also in the IR tail. (This fitting formula does not capture the small oscillatory structure in the UV tail of the feature (see Figure 9) which does not concern us since that region is not phenomenologically interesting.)Figure 9
The bump-like features generated by IR cascading. We plot the feature power spectrumPbump(k) obtained from fully nonlinear lattice field theory simulations (the red points) and also the result of an analytical calculation (the dashed blue curve) using the formalism described in [14]. We also superpose the fitting function ~k3e-πk2/(2k⋆2) (the solid black curve) to illustrate the accuracy of this simple formula.It is straightforward to generalize this discussion to allow for multiple bursts of particle production during inflation. Suppose there are multiple pointsϕ=ϕi (i=1,…,n) along the inflationary trajectory where new degrees of freedom χi become massless:(105)ℒint=-∑i=0ngi22(ϕ-ϕi)χi2.
For each instant ti when ϕ=ϕi there will be an associated burst of particle production and subsequent rescattering of the produced massive χi off the condensate ϕ(t). Each of these events proceeds as described above and leads to a new bump-like contribution to the power spectrum. These features simply superpose owing to that fact that each field χi is statistically independent (so that the cross-terms involving χiχj with i≠j in the computation of the two-point function must vanish). Thus, we arrive at the following parametrization of the primordial power spectrum in models with particle production during inflation:(106)P(k)=As(kk0)ns-1+∑i=1nAi(πe3)3/2(kki)3e-(π/2)(k/ki)2,
where As is the amplitude of the usual nearly scale-invariant vacuum fluctuations from inflation, and k0 is the pivot, which we choose to be k0=0.002Mpc-1 following [132]. The constants Ai depend on the couplings gi2 and measure of the size of the features from particle production. We have normalized these amplitudes so that the power in the ith bump, measured at the peak of the feature, is given by Ai. The location of each feature, ki, is related to the number of e-foldings N from the end of inflation to the time when the ith burst of particle production occurs: roughly ln(ki/H)~Ni where N=Ni at the moment when ϕ=ϕi. From a purely phenomenological perspective the locations ki are completely arbitrary.We compare (106) to lattice field theory simulations in order to determine the amplitude Ai in terms of model parameters. We find(107)Ai≅10-6gi15/4.
The power-law form of this relation was determined by inspection of numerical results; however, it can also be motivated by analytical considerations. Assuming standard chaotic inflation V=m2ϕ2/2, we have tested this formula for g2=1,0.1,0.01, taking both ϕ0=28πMp and ϕ0=3.28πMp. We found agreement up to factors order unity in all cases.Theoretical consistency of our calculation of the shape of the feature bounds the coupling as10-7≲gi2≲1 [12]. Hence, power spectrum (106) can be obtained from sensible microphysics only when 10-20≲Ai≲10-6. This constraint still allows for a huge range of observational possibilities: near the upper bound the feature is considerably larger than the vacuum fluctuations while near the lower bound the feature is completely undetectable.Note that for each bump in (106) the IR tail Pbump→k3 as k→0 is similar to the feature considered by Hoi et al. in [108], consistent with causality arguments about the generation of curvature perturbations by local physics.


## 6.2. Data Sets and Analysis

The primordial power spectrum for our model is parametrized as (106). Our aim is to derive observational constraints on the various model parameters As, ns, ki, and Ai using CMB, galaxy power spectrum, and weak lensing data. To this end we use the cosmoMC package [133] to run Markov Chain Monte Carlo (MCMC) calculations to determine the likelihood of the cosmological parameters, including our new parameters Ai and ki. We employ the following data sets.Cosmic Microwave Background (CMB)
Our complete CMB data sets include WMAP-5yr [132, 134], BOOMERANG [135–137], ACBAR [138–141], CBI [142–145], VSA [146], DASI [147, 148], and MAXIMA [149]. We have included the Sunyaev-Zeldovic (SZ) secondary anisotropy [150, 151] for WMAP-5yr, ACBAR, and CBI data sets. The SZ template is obtained from hydrodynamical simulation [152]. Also included for theoretical calculation of CMB power spectra is the CMB lensing contribution.Type Ia Supernova (SN)
We employ the Union Supernova Ia data (307 SN Ia samples) from The Supernova Cosmology Project [153].Large Scale Structure (LSS)
The 2dF Galaxy Redshift Survey (2dFGRS) data [154] and Sloan Digital Sky Survey (SDSS) Luminous Red Galaxy (LRG) data release 4 [155] are utilized.Note that we have used the likelihood code based on the non-linear modeling by Tegmark et al.  [155] (marginalizing the bias b and the Q parameter). However with a large bump in the linear power spectrum, this naive treatment may be not sufficient to characterize the non-linear response to the feature on small scales. Ideally, this should be obtained from N-body simulations; however, such a study is beyond the scope of the current work.There are several other caveats on our results in the high-k regime. First, we assume linear bias for the galaxies, which may not be entirely safe at sufficiently small scales. Moreover, sharp features in the matter power spectrum can cause sharp features in the bias as a function of k.Keeping in mind these caveats our constraints on small scalesk≳0.1Mpc-1 should be taken with a grain of salt and considered as accurate only up to factors order unity.Weak Lensing (WL)
Five WL data sets are used in this paper. The effective survey areaAeff and galaxy number density neff of each survey are listed in Table 1.Table 1
Weak lensing data sets.
Data setsAeffneff(deg2)(arcmin-2)COSMOS [156, 157]1.640CFHTLS-wide [158, 159]2212GaBODS [160, 161]1312.5RCS [160, 161]538VIRMOS-DESCART [159, 162]8.515For COSMOS data we use the CosmoMC plug-in written by Lesgourgues [157], modified to do numerical marginalization on three nuisance parameters in the original code.For the other four weak lensing data sets we use the likelihood given in [163]. To calculate the likelihood we have written a CosmoMC plug-in code, with simplified marginalization on the parameters of galaxy number density function n(z). More details about this plug-in can be found in [164].As for the LSS data, for small scalesk≳0.1Mpc-1 there is the caveat that the nonlinear evolution of the power spectrum in the presence of bump-like distortions may not be treated accurately.


## 6.3. Observational Constraints: A Single Burst of Particle Production

We now present our results for the observational constraints on particle production during inflation, assuming two different scenarios.The minimal scenario to consider is a single burst of particle production during inflation, which corresponds to takingn=1 in (105). The power spectrum is given by (106) with n=1, and, with some abuse of notation, we denote k1≡kIR and A1≡AIR. The prior we have used for AIR is 0≤AIR≤25×10-10 and for kIR is -9.5≤ln[k/Mpc-1]≤1. The former condition ensures that the bump-like feature from IR cascading does not dominate over the observed scale-invariant fluctuations while the latter is necessary in order to have the feature in the observable range of scales. In Figure 10 we plot the marginalized posterior likelihood for the new parameters AIR and kIR describing the magnitude and location of the bump while in Table 2 we give the best fit values for the remaining (vanilla) cosmological parameters.Table 2
Constraints on the standard “vanilla” cosmological parameters for the single-bump model. All errors are 95.4% confidence level.
Ωbh20.0227-0.0010+0.0010Ωch20.1122-0.0044+0.0050θ1.0424-0.0043+0.0042τ0.08-0.03+0.03ns0.956-0.024+0.024ln[1010As]3.206-0.068+0.074ASZ1.62-0.74+0.76Ωm0.264-0.022+0.026σ80.807-0.030+0.034zre10.5-2.7+2.5H071.6-2.4+2.3Figure 10
Marginalized posterior likelihood contours for the parametersAIR and kIR (the magnitude and position of the feature, resp.) in the single-bump model. Black and grey regions correspond to parameter values allowed at 95.4% and 99.7% confidence levels, respectively. At small scales, to the right of the dashed vertical line, our results should be taken with a grain of salt since the nonlinear evolution of the power spectrum may not be modeled correctly in the presence of bump-like distortions.For very large scales≲Gpc-1, the data do not contain much information (due to cosmic variance), and hence the constraint on any modification of the power spectrum is weak. In this region the spectral distortion may be larger than 100% of the usual scale-invariant fluctuations and couplings g2 of order unity are allowed. For smaller scales k≳Gpc-1 the constraints are stronger, and we have, very roughly, AIR/As≲0.1 corresponding to g2≲0.01. For very small scales, k≳0.1Mpc-1, our constraints should be taken with a grain of salt since the nonlinear evolution of the power spectrum may not be modeled correctly in the presence of bump-like distortions. At small scales nonlinear effects tend to wipe out features of this type (see, e.g., [165]), and hence observational constraints for k≳0.1Mpc-1 may be weaker than what is presented in Figure 10. Note that in most of this nonlinear regime we find essentially no constraint on AIR, which is consistent with what would be expected in a more comprehensive treatment.The IR cascading bump in the primordial power spectrum will be accompanied by a corresponding nongaussian feature in the bispectrum, that will be discussed in more detail in the next section. From the perspective of potentially observing this signal it is most interesting if this feature is located on scales probed by CMB experiments. (There is also the fascinating possibility that the nongaussianity from IR cascading could show up in the large-scale structure as in [27, 166–168]. We leave a detailed discussion to future studies.) To get some intuition into what kinds of features in the CMB scales are still allowed by the data we focus on an example with AIR=2.5×10-10 which, using (107), corresponds to a reasonable coupling value g2~0.01. We take the bump to be located at kIR=0.01Mpc-1 and fix the remaining model parameters to As=2.44×10-9, ns=0.97 (which are compatible with the usual values). This sample bump in the power spectrum is illustrated in Figure 3(a) and is consistent with the data at 2σ. In Figure 3(b) we plot the associated angular CMB TT spectrum. This example represents a surprisingly large spectral distortion: the total power in the feature as compared to the scale invariant vacuum fluctuations is Pbump/Pvac~0.1, evaluated at the peak of the bump. We discuss the nongaussianity associated with this feature.


## 6.4. Observational Constraints: Multiple Bursts of Particle Production

Next, we consider a slightly more complicated scenario: multiple bursts of particle production leading many localized features in the power spectrum. For simplicity we assume that all bumps have the same magnitudeAi≡AIR, and we further suppose a fixed number of e-foldings δN between each burst of particle production. This implies that the location of the ith bump will be given by ki=e(i-1)Δk1 where Δ is a model parameter controlling the density of features. We take the number of bursts, n, to be sufficiently large that the series of features extends over the whole observable range. In the next section we will see that these assumptions are not restrictive and that many well-motivated models are consistent with this simple setup.Our multibump model, then, has three parameters:AIR, k1, and Δ. We take the prior on the amplitude to be AIR≤25×10-10 as in Section 6.3. If the features are very widely spaced, Δ≳1, then the constraint on each bump will obviously be identical to the results for the single-bump case presented in Section 6.3. Hence the most interesting case to consider is Δ≲1 so that the bumps are partially overlapping. Our prior for the density of features is therefore 0≤Δ≤1. Finally, the location of the first bump will be a historical accident in realistic models, hence we marginalize over all possible values of k1 and present our constraints and 2d likelihood plots in the space of AIR and Δ. This marginalized likelihood plot is presented in Figure 11. In Table 3 we present the best-fit values for the vanilla cosmological parameters.Table 3
Constraints on the standard “vanilla” cosmological parameters for the multiple-bump model. All error bars are 95.4% confidence level.
Ωbh20.0227-0.0009+0.0009Ωch20.1126-0.0044+0.0049θ1.0424-0.0043+0.0039τ0.078-0.026+0.031ns0.93-0.17+0.04ln[1010As]2.8-0.9+0.4ASZ1.60-0.76+0.77Ωm0.265-0.021+0.026σ80.807-0.030+0.034zre10.3-2.5+2.6H071.4-2.4+2.2Figure 11
Marginalized posterior likelihood contours for the parametersAIR and Δ (the feature amplitude and spacing, resp.) of the multiple-bump model. Black and grey regions correspond to values allowed at 95.4% and 99.7% confidence levels, respectively.From the likelihood plot, Figure11, there is evidently a preferred value of the feature spacing, roughly Δ~0.75, for which the constraints are the weakest. This can be understood as follows. For very high density Δ→0 the localized features from IR cascading smear together, and the total power spectrum (106) is P(k)~As(k/k0)ns-1+C where the size of the constant deformation scales linearly with the density of features: C∝Δ-1. Therefore, the upper bound on the amplitude AIR should scale linearly with Δ. Indeed, this linear trend is very evident from Figure 11 in the small-Δ regime. This linear behaviour must break down at some point since as the features become infinitely widely spaced the constraint on AIR must go to zero. This explains the bump in the likelihood plot, Figure 11, near Δ~0.75.In passing, notice that the behaviourP(k)~As(k/k0)ns-1+C for Δ≪1 also explains why the best-fit As in Table 3 is somewhat lower than the standard value and why the spectral tilt ns-1 is somewhat more red.To get some intuition for the kinds of multibump distortions that are allowed by the data, we consider an example withAIR=1×10-9, Δ=0.75 and fix the vanilla parameters to As=1.04×10-9, ns=0.93. This choice of parameters is consistent with the data at 2σ and corresponds to a reasonable coupling g2~0.02. In Figure 12 we plot the primordial power spectrum P(k) and also the CMB TT angular power spectrum for this example.(a) shows a sample multiple-bump distortion with amplitudeAIR=1×10-9 which corresponds to g2~0.02. The feature spacing is Δ=0.75. We take the vanilla parameters to be As=1.04×10-9, ns=0.93 so that the scale of inflation is slightly lower than in the standard scenario and the spectral tilt is slightly redder. (b) shows the CMB angular TT power spectrum for this example.
(a)(b)


## 7. Particle Physics Models

In previous sections of this paper we have studied in some detail the observational signatures of inflationary particle production in model (4), and we have also derived observational constraints on the model parameters g2 and ϕ0. We now aim to provide some explicit microscopic derivations of model (5) from popular theories of particle physics including string theory and SUSY. Our examples are meant to be illustrative rather than exhaustive. From the low energy perspective interactions of type (4) are completely generic, and hence we expect that many microscopic derivations may be possible.


### 7.1. Open String Inflation Models

String theory inflation models may be divided into two classes depending on the origin of the inflaton: closed string models and open string models. In the former case the inflaton is typically a geometrical modulus associated with the compactification manifold (examples include racetrack inflation [169], Kähler modulus inflation [170], Roulette inflation [171], and fibre inflation [172]). In the latter case the inflaton is typically the position modulus of some mobile D-brane (One notable exception is inflation driven by the open string tachyon, for example nonlocal string field theory models [42–44, 47–50].) living in the compactification manifold (examples include brane inflation [173, 174] such as the warped KKLMMT model [175], D3/D7 inflation [176], and DBI inflation [51, 52]). In open string inflation models there may be, in addition to the mobile inflationary brane, some additional “spectator” branes. If the mobile brane collides with any spectator brane during inflation, then some of the stretched string states between these branes will become massless at the moment when the branes are coincident [54, 80], precisely mimicking interaction (4). Thus, we expect particle production, IR cascading, and the bump-like features described above to be a reasonably generic prediction of open string inflation.


### 7.2. String Monodromy Models

A concrete example of the heuristic scenario discussed in the last subsection is provided by the brane monodromy and axion monodromy string theory inflation models proposed in [77–79]. In the original brane monodromy model [77] one considers type IIA string theory compactified on a nil manifold that is the product of two twisted tori. The metric on each of these twisted tori has the form(108)ds2α′=Lu12du12+Lu22du22+Lx2(dx′+Mu1du2)2,
where x′=x-(M/2)u1u2, and M is an integer flux number. The dimensionless constants Lu1, Lu2, and Lx determine the size of the compactification.Inflation is realized by the motion of a D4-brane along the directionu1 of the internal manifold. The D4 spans our large 3 dimensions and wraps 1 cycle along the direction u2 of the internal space. The size of this 1 cycle, in string units, is given by(109)L=Lu22+Lx2M2u12.
Hence, the brane prefers to minimize its world-volume by moving to the location u1=0 where this 1 cycle has minimal size. This preference gives a potential to the D4-brane position which goes like V~u1 in the large u1 regime that is relevant for large field inflation.In [54] it was shown that this scenario allows for the inclusion of a number of spectator branes stabilized at positions u1=i/M (with i integer) along the inflationary trajectory. As the mobile inflationary D4 rolls through these points, particles (strings) will be produced, and the resulting distribution of features will look nearly identical to the simple multi-bump scenario studied in Section 6.4. To see this, let us now determine the distribution of features that is predicted from brane monodromy inflation. The canonical inflaton ϕ can be related to the position of the mobile D4 as(110)ϕ=Bu11/p,
where B is a constant with dimensions of mass that depends on model parameters. Hence, the effective potential during inflation has the power-law form(111)V(ϕ)=μ4-pϕp.
For the simplest scenario described above one has p=2/3. However, formulas (110) and (111) still hold for the variant considered in [77] with p=2/5 as long as one replaces u1 by a more complicated linear combination of coordinates. These relations also hold for axion monodromy models [78] with p=1 and u1 replaced by the axion, c, arising from a 2-form RR potential C(2) integrated over a 2-cycle Σ2. For all models of the form of (111) the number of e-foldings N from ϕ=ϕ(N) to the end of inflation (which occurs at ϕ=p/2 when the slow-roll parameter ϵ(ϕ)=1) is given by(112)N=12pϕ2(N)Mp2-p4=12pB2Mp2u12/p-p4.
Since the spectator branes are located at u1=i/M, the bursts of particle production must occur at times N=Ni where(113)Ni=12pB2Mp2(iM)2/p-p4.
The location k=ki of the ith feature is defined, roughly, by the scale leaving the horizon at the moment N=Ni. Hence, the distribution of features for brane/axion monodromy models is given by(114)ln[kiH]≅B̃2i2/p-p4
with p=2/3 or p=2/5 for brane monodromy and p=1 for axion monodromy. In (114) the dimensionless number B̃ depends on model parameters.Although the distribution of features (114) is not exactly the same as the evenly space distribution considered in Section 6.4, the two are essentially indistinguishable over the range of scales which are probed observationally (corresponding to roughly 10 e-foldings of inflation). The reason for this is simple: the inflaton is nearly constant during the first 10 e-foldings of inflation and hence δN~δϕ~δu1 within the observable region. It follows that ki≅e(i-1)Δk1 to very good approximation for a huge class of models. To see this more concretely in the case at hand, let us compute dN/du1 from (112). It is straightforward to show that(115)dNdu1=1pp1[2ϵ(ϕ)]1-p/2(BMp)p,
where (116)ϵ(ϕ)≡Mp22(V′V)2=p22(Mpϕ)2
is the usual slow-roll parameter. Observational constraints on the running of the spectral index imply that ϵ(ϕ) cannot change much over the observable 10 e-foldings of inflation. Since dN/du1≅const to very high accuracy it follows trivially that N=N(u1) is very close to linear and ki≅e(i-1)Δk1 as desired.In the context of axion monodromy inflation models [78] the multiple bump features discussed here will be complimentary to the oscillatory features described in [79] which result from the sinusoidal modulation of the inflaton potential by instanton effects. If the bursts of particle production are sufficiently densely spaced, then signal from IR cascading may appear oscillatory; however, it differs from the effect discussed in [79] in both physical origin and also functional form.Let us now estimate the effective value of the couplingsgi2 appearing in the prototype interaction (105) that are predicted from the simplest brane monodromy model. A complete calculation would involve dimensionally reducing the DBI action describing the brane motion and requires knowledge of the full 10-dimensional geometry with the various embedded branes. For our purposes, however, a simple heuristic estimate for the collision of two D4-branes will suffice. When N D-branes become coincident the symmetry is enhanced from U(1)N to a U(N) Yang Mills gauge theory. The gauge coupling for this Yang Mills theory is given by(117)gYM2=gs(2π)2L,
where L is the volume of the 1 cycle that the D4 branes wrap and is given by (5). If the inflationary brane is at position u1 and the ith spectator brane is at u1,i, then the distance between the two branes is given by(118)d2=α′Lu12(u1-u1,i)2.
The mass of the gauge bosons corresponding to the enhanced symmetry is(119)Mi2=gYM2d2(2π)2(α′)2=gsLu12(u1-u1,i)2α′Lu22+Lx2M2u12.
To put this in the prototype form Mi2=gi2(ϕ-ϕi)2 we must first convert to the canonical variable ϕ using formula (110) with p=2/3 and(120)B=M1/26π2Lu1Lx1/2gsα′.
Next, we must Taylor expand the resulting equation about the minimum ϕ=ϕi. We find(121)Mi2≅gi2(ϕ-ϕi)2+⋯,(122)gi2=16gs2π4MLxu1,i1Lu22+Lx2M2u1,i2=16gs2π4Lxi1Lu22+Lx2i2,
where on the second line of (122) we have used the fact that u1,i=i/M (with i integer) in the simplest models. We see that the effective couplings gi2 become larger as the D4 unwinds during inflation. (The apparent divergence for u1,i=0 in formula (122) is an artifact of the fact that relation (110) is not valid at small values of u1. This will not concern us here since inflation has already terminated at the point that our formulas break down.)To compute the amplitude of the bump-like feature produced by brane monodromy inflation we should take into account also combinatorial factors. When two branes become coincident the symmetry is enhanced fromU(1)2 to U(2) so there are 22-2=2 additional massless spin-1 fields appearing at the brane collision. Thus, using (107), the amplitude of the feature that will be imprinted in the CMB is(123)Ai,eff=2×(22-2)×[10-6·gi15/4],
where the extra factor of 2 counts the polarizations of the massless spin-1 fields. This combinatorial enhancement can be much larger if the inflationary brane collides with a stack of spectators.The above discussion is predicated on the assumption that the original brane monodromy setup 1 is supplemented by additional spectator branes. This may seem like an unnecessary contrivance; however, in order for this model to reheat successfully it may benecessary to include spectator branes. For example, with the reheating mechanism proposed in [177], semirealistic particle phenomenology can be obtained by confining the standard model (SM) to a D6 brane which wraps the compact space. In order to reheat into SM degrees of freedom, we orient this brane so that its world-volume is parallel to the mobile (inflationary) D4. In this case the end of inflation involves multiple oscillations of the D4 about the minimum of its potential. At each oscillation the D4 collides with the D6 and SM particles are produced by parametric resonance preheating [92]. However, due to the periodic structure of the compactification, D4/D6 collisions will necessarily occur also during inflation, leading to IR cascading features in the CMB.The timing of these D4/D6 collisions was computed in [177] for the minimal p=2/3 brane monodromy model, assuming the same choices of parameters used in [77]. For this particular case there is only one collision (and hence one feature) during the first 10 e-foldings of inflation, and the phenomenology is essentially the same as that considered in Section 6.3. What is the amplitude of this feature? Assuming, again, the parameters employed in [77] and noting that the first collision takes place at i=13 [177], (122) gives g12≅0.001. From (123) we find the effective amplitude of the feature to be A1,eff/As≅0.01. This value is well within the observational bounds derived in Section 6.3.We stress that the conclusions in the previous paragraph applyonly for the particular choice of model parameters employed in [77]. There exist other consistent parameter choices for which the simplest brane monodromy model predicts a much higher density of features with much larger amplitude.Note that both brane and axion monodromy models may be used to realize trapped inflation [54]. Here we are restricting ourselves to the large-field regime where the potential V=μ4-pϕp is flat enough to drive inflation without the need for trapping effects. For a given choice of parameters one should verify that this classical potential dominates over the quantum corrections from particle production.


### 7.3. A Supersymmetric Model

Another microscopic realization of multiple bursts of particle production and IR cascading during inflation which does not rely on string theory can be obtained from the so-called “distributed mass” model derived in [76] with warm inflation [82] in mind; however, the theory works equally well for our scenario. This model is based on 𝒩=1 global SUSY and allows for the inclusion of multiple points along the inflationary trajectory where both scalar degrees of freedom and also their associated fermion super-partners become massless. The distribution of features in this set-up is essentially arbitrary.


## 7.1. Open String Inflation Models

String theory inflation models may be divided into two classes depending on the origin of the inflaton: closed string models and open string models. In the former case the inflaton is typically a geometrical modulus associated with the compactification manifold (examples include racetrack inflation [169], Kähler modulus inflation [170], Roulette inflation [171], and fibre inflation [172]). In the latter case the inflaton is typically the position modulus of some mobile D-brane (One notable exception is inflation driven by the open string tachyon, for example nonlocal string field theory models [42–44, 47–50].) living in the compactification manifold (examples include brane inflation [173, 174] such as the warped KKLMMT model [175], D3/D7 inflation [176], and DBI inflation [51, 52]). In open string inflation models there may be, in addition to the mobile inflationary brane, some additional “spectator” branes. If the mobile brane collides with any spectator brane during inflation, then some of the stretched string states between these branes will become massless at the moment when the branes are coincident [54, 80], precisely mimicking interaction (4). Thus, we expect particle production, IR cascading, and the bump-like features described above to be a reasonably generic prediction of open string inflation.


## 7.2. String Monodromy Models

A concrete example of the heuristic scenario discussed in the last subsection is provided by the brane monodromy and axion monodromy string theory inflation models proposed in [77–79]. In the original brane monodromy model [77] one considers type IIA string theory compactified on a nil manifold that is the product of two twisted tori. The metric on each of these twisted tori has the form(108)ds2α′=Lu12du12+Lu22du22+Lx2(dx′+Mu1du2)2,
where x′=x-(M/2)u1u2, and M is an integer flux number. The dimensionless constants Lu1, Lu2, and Lx determine the size of the compactification.Inflation is realized by the motion of a D4-brane along the directionu1 of the internal manifold. The D4 spans our large 3 dimensions and wraps 1 cycle along the direction u2 of the internal space. The size of this 1 cycle, in string units, is given by(109)L=Lu22+Lx2M2u12.
Hence, the brane prefers to minimize its world-volume by moving to the location u1=0 where this 1 cycle has minimal size. This preference gives a potential to the D4-brane position which goes like V~u1 in the large u1 regime that is relevant for large field inflation.In [54] it was shown that this scenario allows for the inclusion of a number of spectator branes stabilized at positions u1=i/M (with i integer) along the inflationary trajectory. As the mobile inflationary D4 rolls through these points, particles (strings) will be produced, and the resulting distribution of features will look nearly identical to the simple multi-bump scenario studied in Section 6.4. To see this, let us now determine the distribution of features that is predicted from brane monodromy inflation. The canonical inflaton ϕ can be related to the position of the mobile D4 as(110)ϕ=Bu11/p,
where B is a constant with dimensions of mass that depends on model parameters. Hence, the effective potential during inflation has the power-law form(111)V(ϕ)=μ4-pϕp.
For the simplest scenario described above one has p=2/3. However, formulas (110) and (111) still hold for the variant considered in [77] with p=2/5 as long as one replaces u1 by a more complicated linear combination of coordinates. These relations also hold for axion monodromy models [78] with p=1 and u1 replaced by the axion, c, arising from a 2-form RR potential C(2) integrated over a 2-cycle Σ2. For all models of the form of (111) the number of e-foldings N from ϕ=ϕ(N) to the end of inflation (which occurs at ϕ=p/2 when the slow-roll parameter ϵ(ϕ)=1) is given by(112)N=12pϕ2(N)Mp2-p4=12pB2Mp2u12/p-p4.
Since the spectator branes are located at u1=i/M, the bursts of particle production must occur at times N=Ni where(113)Ni=12pB2Mp2(iM)2/p-p4.
The location k=ki of the ith feature is defined, roughly, by the scale leaving the horizon at the moment N=Ni. Hence, the distribution of features for brane/axion monodromy models is given by(114)ln[kiH]≅B̃2i2/p-p4
with p=2/3 or p=2/5 for brane monodromy and p=1 for axion monodromy. In (114) the dimensionless number B̃ depends on model parameters.Although the distribution of features (114) is not exactly the same as the evenly space distribution considered in Section 6.4, the two are essentially indistinguishable over the range of scales which are probed observationally (corresponding to roughly 10 e-foldings of inflation). The reason for this is simple: the inflaton is nearly constant during the first 10 e-foldings of inflation and hence δN~δϕ~δu1 within the observable region. It follows that ki≅e(i-1)Δk1 to very good approximation for a huge class of models. To see this more concretely in the case at hand, let us compute dN/du1 from (112). It is straightforward to show that(115)dNdu1=1pp1[2ϵ(ϕ)]1-p/2(BMp)p,
where (116)ϵ(ϕ)≡Mp22(V′V)2=p22(Mpϕ)2
is the usual slow-roll parameter. Observational constraints on the running of the spectral index imply that ϵ(ϕ) cannot change much over the observable 10 e-foldings of inflation. Since dN/du1≅const to very high accuracy it follows trivially that N=N(u1) is very close to linear and ki≅e(i-1)Δk1 as desired.In the context of axion monodromy inflation models [78] the multiple bump features discussed here will be complimentary to the oscillatory features described in [79] which result from the sinusoidal modulation of the inflaton potential by instanton effects. If the bursts of particle production are sufficiently densely spaced, then signal from IR cascading may appear oscillatory; however, it differs from the effect discussed in [79] in both physical origin and also functional form.Let us now estimate the effective value of the couplingsgi2 appearing in the prototype interaction (105) that are predicted from the simplest brane monodromy model. A complete calculation would involve dimensionally reducing the DBI action describing the brane motion and requires knowledge of the full 10-dimensional geometry with the various embedded branes. For our purposes, however, a simple heuristic estimate for the collision of two D4-branes will suffice. When N D-branes become coincident the symmetry is enhanced from U(1)N to a U(N) Yang Mills gauge theory. The gauge coupling for this Yang Mills theory is given by(117)gYM2=gs(2π)2L,
where L is the volume of the 1 cycle that the D4 branes wrap and is given by (5). If the inflationary brane is at position u1 and the ith spectator brane is at u1,i, then the distance between the two branes is given by(118)d2=α′Lu12(u1-u1,i)2.
The mass of the gauge bosons corresponding to the enhanced symmetry is(119)Mi2=gYM2d2(2π)2(α′)2=gsLu12(u1-u1,i)2α′Lu22+Lx2M2u12.
To put this in the prototype form Mi2=gi2(ϕ-ϕi)2 we must first convert to the canonical variable ϕ using formula (110) with p=2/3 and(120)B=M1/26π2Lu1Lx1/2gsα′.
Next, we must Taylor expand the resulting equation about the minimum ϕ=ϕi. We find(121)Mi2≅gi2(ϕ-ϕi)2+⋯,(122)gi2=16gs2π4MLxu1,i1Lu22+Lx2M2u1,i2=16gs2π4Lxi1Lu22+Lx2i2,
where on the second line of (122) we have used the fact that u1,i=i/M (with i integer) in the simplest models. We see that the effective couplings gi2 become larger as the D4 unwinds during inflation. (The apparent divergence for u1,i=0 in formula (122) is an artifact of the fact that relation (110) is not valid at small values of u1. This will not concern us here since inflation has already terminated at the point that our formulas break down.)To compute the amplitude of the bump-like feature produced by brane monodromy inflation we should take into account also combinatorial factors. When two branes become coincident the symmetry is enhanced fromU(1)2 to U(2) so there are 22-2=2 additional massless spin-1 fields appearing at the brane collision. Thus, using (107), the amplitude of the feature that will be imprinted in the CMB is(123)Ai,eff=2×(22-2)×[10-6·gi15/4],
where the extra factor of 2 counts the polarizations of the massless spin-1 fields. This combinatorial enhancement can be much larger if the inflationary brane collides with a stack of spectators.The above discussion is predicated on the assumption that the original brane monodromy setup 1 is supplemented by additional spectator branes. This may seem like an unnecessary contrivance; however, in order for this model to reheat successfully it may benecessary to include spectator branes. For example, with the reheating mechanism proposed in [177], semirealistic particle phenomenology can be obtained by confining the standard model (SM) to a D6 brane which wraps the compact space. In order to reheat into SM degrees of freedom, we orient this brane so that its world-volume is parallel to the mobile (inflationary) D4. In this case the end of inflation involves multiple oscillations of the D4 about the minimum of its potential. At each oscillation the D4 collides with the D6 and SM particles are produced by parametric resonance preheating [92]. However, due to the periodic structure of the compactification, D4/D6 collisions will necessarily occur also during inflation, leading to IR cascading features in the CMB.The timing of these D4/D6 collisions was computed in [177] for the minimal p=2/3 brane monodromy model, assuming the same choices of parameters used in [77]. For this particular case there is only one collision (and hence one feature) during the first 10 e-foldings of inflation, and the phenomenology is essentially the same as that considered in Section 6.3. What is the amplitude of this feature? Assuming, again, the parameters employed in [77] and noting that the first collision takes place at i=13 [177], (122) gives g12≅0.001. From (123) we find the effective amplitude of the feature to be A1,eff/As≅0.01. This value is well within the observational bounds derived in Section 6.3.We stress that the conclusions in the previous paragraph applyonly for the particular choice of model parameters employed in [77]. There exist other consistent parameter choices for which the simplest brane monodromy model predicts a much higher density of features with much larger amplitude.Note that both brane and axion monodromy models may be used to realize trapped inflation [54]. Here we are restricting ourselves to the large-field regime where the potential V=μ4-pϕp is flat enough to drive inflation without the need for trapping effects. For a given choice of parameters one should verify that this classical potential dominates over the quantum corrections from particle production.


## 7.3. A Supersymmetric Model

Another microscopic realization of multiple bursts of particle production and IR cascading during inflation which does not rely on string theory can be obtained from the so-called “distributed mass” model derived in [76] with warm inflation [82] in mind; however, the theory works equally well for our scenario. This model is based on 𝒩=1 global SUSY and allows for the inclusion of multiple points along the inflationary trajectory where both scalar degrees of freedom and also their associated fermion super-partners become massless. The distribution of features in this set-up is essentially arbitrary.


## 8. Nongaussianity from Particle Production and IR Cascading

We have seen that the dynamics of particle production, rescattering, and IR cascading during inflation in model (5) leads to a bump-like contribution to the primordial power spectrum of the observable cosmological fluctuations. However, this dynamics must also have a nontrivial impact on nongaussian statistics, such as the bispectrum. Indeed, it is evident from the analytical formalism presented in Sections 4 and 5 that the inflaton fluctuations δϕ generated via the rescattering diagram in Figure 1 may be significantly nongaussian. (This is evident since the particular solution of (49) is bi-linear in the Gaussian field χ.) In this section we characterize this nongaussianity. The results of this section appeared for the first time in [12, 14, 15].


### 8.1. Probability Density Function and Cummulants

As discussed in the introduction, nongaussianity is often characterized by the size, shape, and running of the bispectrumB(ki). A popular measure of the size of the nongaussianity is the magnitude of B(ki) on equilateral triangles. This approach has the advantage of being straightforward; however, it may give misleading results when one wishes to compare bispectra with different shapes [178] or scaling behaviour. This is particularly true in models, such as (5), where the bispectrum is large only for triangles with a characteristic size and where the nongaussian part of ζ is uncorrelated with the gaussian part. A more meaningful single number to compare between models is the skewness (defined below), which contains information about the bispectrum integrated over all scales (up to some UV smoothing scale) and all shape configurations. (See also [179] for a related discussion and alternative methodology.)In what follows it will be useful to define a normalized curvature perturbation(124)ξ≡ζσζ,
where ζ is the usual gauge invariant curvature perturbation, and σζ is the variance(125)σζ2=〈ζ2〉=∫d3k1(2π)3/2d3k2(2π)3/2〈ζk1ζk2〉.
If ζ is generated by the quantum vacuum fluctuations of the inflaton field, then we have σζ≅10-9/2.Let us define the probability density function (PDF),P(ξ), as the probability that the (normalized) curvature perturbation has fluctuations of size ξ. The nth central moment of the PDF is simply(126)〈ξn〉=∫-∞+∞ξnP(ξ)dξ.
In order to quantify the nongaussianity of the PDF it is useful to define the cummulants(127)Ŝn≡〈ξn〉c=〈ζn〉cσζn,
where the subscript c is a reminder that only the connected part of the correlation  function should be included. For a gaussian PDF we have Ŝn=0 for n=3,4,…, hence the cummulants provide a measure of nongaussianity. (Notice that, since ξ is constructed to have unit variance, we do not need to distinguish between the reduced and dimensionless cummulants.)In the case thatζ (and hence ξ) is a free field then the central limit theorem implies that the PDF is gaussian:(128)P(ξ)dξ=dξ2πe-ξ2/2.
This expression admits a simple generalization to the case where P(ξ) is close to (but not exactly) gaussian. This generalization is the Edgeworth expansion(129)P(ξ)dξ≅dξ2πe-ξ2/2×[1+Ŝ36H3(ξ)+(Ŝ424H4(ξ)+Ŝ3272H3(ξ))+⋯]
(See [27] for a review and [180] for an alternative derivation.) Where Hn(ξ) denotes the Hermite polynomials of order n and Ŝn are the dimensionless cummulants, defined by (127). Note that result (129) is an expansion in cummulants. This expression provides an accurate approximation of the true PDF provided the cummulants are well ordered in the sense that(130)1≫|Ŝ3|≫|Ŝ4|≫⋯.
See [181, 182] for more discussion on the structure of the correlation functions and the validity of cosmological perturbation theory.Finally, it remains to relate the cummulantsŜn (n≥3) to the correlation functions (such as the bispectrum) of the primordial curvature perturbation. We can write the dimensionless skewness, Ŝ3, as an integral over the 3-point correlation function [27, 180]:(131)Ŝ3=1σζ3∫d3k1(2π)3/2d3k2(2π)3/2d3k3(2π)3/2〈ζk1ζk2ζk3〉c.
A similar expression may be derived for the dimensionless kurtosis, Ŝ4. From (131) it is clear that the skewness provides a measure of the magnitude of the bispectrum integrated over all triangles. Hence, this provides a reasonable single number to compare the size of nongaussianity between different models. Moreover, the skewnesses (and higher cummulants) are of observational interest since these determine the probability of rare fluctuations [27, 166–168, 183–185]. Observables such as the abundance of collapsed objects may therefore be used to constrain Ŝ3, Ŝ4, and so forth.We have reviewed the derivation of the Edgeworth expansion (129) only to illustrate that the dimensionless skewness (131) encodes information about the bispectrum integrated over all momenta. However, our actual calculation of the PDF in the next subsection will be fully nonperturbative and hence does not rely on the validity of (129) or (130).


### 8.2. Magnitude of the Nongaussianity from IR Cascading

We have argued in the last subsection that the moments of the PDF provide a useful tool for quantifying the size of nongaussianity and drawing comparisons between microscopic models whose bispectra may have very different shapes and scaling properties. Our goal now is to compute the skewness and kurtosis generated by particle production, rescattering and, IR cascading in model (5). We will then compare this to the skewness that would be generated by more familiar constructions, such as the local model ζ=ζg+(3/5)fNLζg2.It is straightforward to extract the PDF numerically from our HLattice simulations by computing the fraction of the simulation box which contains the fieldδϕ∝ζ at a particular value. We again stress that this numerical evaluation of the PDF is non-perturbative; our results do not rely on any prior assumptions about the size or ordering of the cummulants.We present our results for the numerical evaluation of the PDF from IR cascading in Figure13. For illustration we have assumed a chaotic potential V(ϕ)=m2ϕ2/2; however, our results are largely insensitive to the choice of background inflation model.Figure 13
Probability density function ofζ for g2=0.1 at a series of different values of the scale factor, a. The dashed curve shows a Gaussian fit at late time a=6.2, and we have normalized the scale factor so that a=1 at the moment when particle production occurs. For illustration, we have chosen standard chaotic inflation V(ϕ)=(m2/2)ϕ2; however, our qualitative results are independent of the background inflation model.We can understand physically the behaviour of the PDF plotted in Figure13. Shorty after the initial burst of particle production, the curvature perturbation is extremely nongaussian, due to the sudden appearance of the source term J∝χ2 in the equation of motion for the inflaton fluctuations δϕ. Quickly, in less than an e-folding, nonlinear interactions begin to drive the system towards gaussianity. A similar behaviour has been observed in lattice simulations of out-of-equilibrium interacting scalar fields during preheating [186, 187]. In the case of rescattering during preheating, the system will eventually become gaussian when the fields thermalize. However, in our case the universe is still inflating during the process of rescattering and IR cascading. As a result, nongaussian inflaton fluctuations generated by rescattering are stretched out by the quasi-de Sitter expansion and must freeze once their wavelength crosses the Hubble scale. Hence, at late times the PDF does not become completely gaussian but rather freezes-in with some non-trivial skewness. Within a few e-foldings from the moment of particle production the time evolution of the PDF has become completely negligible.Given our numerical results for the PDF, it is straightforward to compute dimensionless cummulants such as the skewness (Ŝ3) and kurtosis (Ŝ4) for various values of the coupling g2. We have summarized our results in Table 4. At small coupling, g2<1, we always have |Ŝ3|>|Ŝ4| which suggests that the cummulants are well ordered. In particular, note that for g2=0.01, the kurtosis Ŝ4 is too small to be measured accurately from our HLattice simulations.Table 4
Dimensionless cummulants for IR cascading.
g2SkewnessKurtosisEquivalent(Ŝ3)(Ŝ4)fNLlocal1-0.510.16-28000.1-0.49-0.08-27000.01-0.033𝒪(10-3)-183In order to give a sense of the magnitude of the nongaussianity from particle production we have also computed an “equivalentfNLlocal” defined by 5Ŝ3/(18σζ). For a given g2 this effective fNLlocal is the magnitude of fNL which would be necessary to reproduce the skewness Ŝ3 of the IR cascading PDF using a local ansatz ζ=ζg+(3/5)fNLζg2. (Our sign conventions for fNL are consistent with WMAP [6]. See [27] for a discussion of various conventions employed in the literature.)From Table4 we see that IR cascading during inflation may generate significant nongaussianity. Even taking g2=0.01 (which is compatible with cosmological data for any choice of ϕ0), we still obtain a skewness Ŝ3=-0.033, which is the same value that would be produced by a local model with fNL=-183. This “equivalent” local nongaussianity is larger than current observational bounds, which suggests that nongaussian features from particle production during inflation might be observable for reasonable values of the coupling g2. Note, however, that observational bounds on fNLlocal cannot be directly applied to our model since the bispectrum in our case is uncorrelated with the vacuum fluctuations and is far from scale invariant.Depending on the value ofϕ0 (and hence the location of the feature in the power spectrum) the nongaussianity from IR cascading may lead to a variety of observational signatures. If the feature is localized on scales relevant for CMB experiments, then the key observable is the primordial bispectrum, B(ki). Efficient computation of the observational constraints on nongaussianity from inflationary particle production will require a separable template for the bispectrum, an issue which we briefly comment on in the next subsection and which will be considered in more detail in [15]. On the other hand, the feature might show up on smaller scales and leave an imprint on large scale structures. In this case there are at least two observables which could be used to constrain our model: higher order correlations of LSS probes (such as the galaxy bispectrum) and the abundance of collapsed objects. In the first case, we require the evolved bispectrum while in the latter case we require the PDF of the evolved density field δρ/ρ. Our model robustly predicts negative skewness for both the curvature perturbation ζ and the density field. Hence we should expect less rare objects (as compared to the gaussian case) on some characteristic scale. (The recent weak lensing measurement of the dark matter mass of the high-redshift galaxy cluster XMMUJ2235.3-2557 [188] has been construed as a possible hint of nongaussian initial conditions [184]. Unfortunately, our model does not produce the correct sign of skewness to explain such observations.)When computing the statistics of density fluctuations on some given length scaleR one should smooth the perturbation field over wave-numbers larger than kmax~R-1. This may be implemented in PDF (129) by the formal substitution ζk→𝒲(k)ζk in expressions (125) and (131) for the variance and skewness. Here 𝒲(k) is some window function which tends to zero for k≫R-1 (e.g., a standard choice of window is a top-hat in real space). This substitution makes any predictions for the statistics of the density fluctuations implicitly dependent on the smoothing scale, R. (It also renders the momentum integrals UV finite.) In our analysis we have computed the PDF numerically by measuring the fraction of the simulation box in which the inflaton field takes a given value. Thus, our cummulants are automatically smoothed on the UV cutoff scale kmax~Δx-1 defined by the lattice spacing Δx.One may wonder how sensitively our results for the dimensionless cummulants depend on this particular choice of smoothing. Recall that the bispectrum from IR cascading (73) is strongly peaked on triangles with a particular size, corresponding to the location of the bump-like distortion in the power spectrum. If we denote this characteristic scale by kbump, then it should be clear that the cummulants in our model are insensitive the smoothing scale, as long as R≲kbump-1.


### 8.3. Shape of the Bispectrum

In Section4 we presented a schematic discussion of how to compute analytically the primordial bispectrum from IR cascading, (73). In [14], this expression is evaluated analytically. We find that the bispectrum (73) from particle production has a unique shape that has not been considered in the previous literature. In order to describe this novel shape we would like to first factor out that strong overall dependence of Bϕ(ki) on the size of the triangle. To this end we define a “shape function” S(ki) in terms of the bispectrum as follows:(132)S(ki)=N-1(k1k2k3)2B(ki),
where B(ki) is related to Bϕ(ki) by (104), and N is a normalization factor to be discussed shortly. The function S(ki) has the advantage that the strong k6 running of the bispectrum is extracted. Hence, any residual scaling behaviour displayed by S(ki) must be a result of nonlinear interactions.Since we expect the 3-point correlation function to be of orderPζ2(k), a natural choice of normalization is N=(2π)4𝒫ζ2 where 𝒫ζ1/2=5×10-5 is the usual amplitude of the scale invariant fluctuations from inflation. With this choice of normalization our function S(ki) coincides with the quantity 𝒢(ki)/(k1k2k3) which was used to study localized nongaussian features from models with steps in the inflaton potential in [59, 60].Symmetry of the bispectrum under permutations of momenta implies that we can focus only on the regionk1≥k2≥k3, to avoid counting the same configuration twice. Moreover, the triangle inequality implies that 1-k2/k1≤k3/k1. Therefore we can completely specify the shape of the bispectrum for a given size of triangle k by plotting S(k,kx2,kx3) in the region x3≤x2≤1 and 1-x2≤x3. Because our bispectrum is very far from being scale invariant, it follows that this shape function is sensitive to the choice of k. Therefore, in Figure 14 we choose several representative choices: ln(k/kbump)=-1,0,1,2.The shape functionS(k,kx2,kx3), defined by (132), as a function of the dimensionless quantities x2,x3 which parametrize the shape of the triangle. (a) corresponds to k=e-1kbump, (b) is k=kbump, (c) is k=e+1kbump, and (d) is k=e+2kbump. In the IR (k≤kbump) the shape of the bispectrum is similar to the equilateral shape; however, there is also some support on flattened triangles near k~e+1kbump. At larger values of k the shape is unlike any other template proposed in the literature.
(a)(b)(c)(d)We see that a rich array of shapes is possible: fork≲kbump, the bispectrum is qualitatively similar to the equilateral model; however, at slightly larger k there is also considerable support on flattened triangles. Note that for k≳7.4kbump the shape of the bispectrum is extremely unusual and is not easily comparable to any shape that has been proposed in the previous literature.We find that the bispectrum from IR cascading is, to a first approximation, factorisable in the sense thatB(k1,k2,k3)≅∏i=13F(ki) where F(k)~C1k3e-C2k2 for some constants C1 and C2. This separable form is important since it permits fast algorithms for forecasting, data analysis, and simulations [189].In this section we have attempted to characterize the nongaussian signature associated with inflationary particle production. This signature is rather unique in the literature. Our model predicts uncorrelated localized nongaussian features with a unique shape of bispectrum. We have quantified the size of this effect by studying the dimensionless cummulants (such as the skewness) and have argued that this nongaussianity can be significant, depending ong2. We leave a detailed study of the observational constraints on such nongaussianities to future studies [15].


### 8.4. Multiple Bursts of Particle Production

Note that our discussion of nongaussianity generalizes easily to the case where there are multiple pointsϕi (i=0,…,n) along the inflaton trajectory where new degrees of freedom χi become massless, such as model (105). Such a construction may be quite natural in the context of brane/axion monodromy inflation [54, 77–79]. In this case the nongaussian features in B(ki) from each burst of particle production may superpose to generate a broad-band signal. Such a construction may (but need not) be associated with trapped inflation [54]. Depending on the spacing of the points ϕi and the couplings gi2, a rich variety of nongaussianities may be possible. In this case that the points ϕi are sufficiently densely spaced, we expect a nongaussian signal which is close to equilateral but also has some support on flattened triangles, consistent with the analysis of [54, 81]. We leave a detailed discussion to future studies.


## 8.1. Probability Density Function and Cummulants

As discussed in the introduction, nongaussianity is often characterized by the size, shape, and running of the bispectrumB(ki). A popular measure of the size of the nongaussianity is the magnitude of B(ki) on equilateral triangles. This approach has the advantage of being straightforward; however, it may give misleading results when one wishes to compare bispectra with different shapes [178] or scaling behaviour. This is particularly true in models, such as (5), where the bispectrum is large only for triangles with a characteristic size and where the nongaussian part of ζ is uncorrelated with the gaussian part. A more meaningful single number to compare between models is the skewness (defined below), which contains information about the bispectrum integrated over all scales (up to some UV smoothing scale) and all shape configurations. (See also [179] for a related discussion and alternative methodology.)In what follows it will be useful to define a normalized curvature perturbation(124)ξ≡ζσζ,
where ζ is the usual gauge invariant curvature perturbation, and σζ is the variance(125)σζ2=〈ζ2〉=∫d3k1(2π)3/2d3k2(2π)3/2〈ζk1ζk2〉.
If ζ is generated by the quantum vacuum fluctuations of the inflaton field, then we have σζ≅10-9/2.Let us define the probability density function (PDF),P(ξ), as the probability that the (normalized) curvature perturbation has fluctuations of size ξ. The nth central moment of the PDF is simply(126)〈ξn〉=∫-∞+∞ξnP(ξ)dξ.
In order to quantify the nongaussianity of the PDF it is useful to define the cummulants(127)Ŝn≡〈ξn〉c=〈ζn〉cσζn,
where the subscript c is a reminder that only the connected part of the correlation  function should be included. For a gaussian PDF we have Ŝn=0 for n=3,4,…, hence the cummulants provide a measure of nongaussianity. (Notice that, since ξ is constructed to have unit variance, we do not need to distinguish between the reduced and dimensionless cummulants.)In the case thatζ (and hence ξ) is a free field then the central limit theorem implies that the PDF is gaussian:(128)P(ξ)dξ=dξ2πe-ξ2/2.
This expression admits a simple generalization to the case where P(ξ) is close to (but not exactly) gaussian. This generalization is the Edgeworth expansion(129)P(ξ)dξ≅dξ2πe-ξ2/2×[1+Ŝ36H3(ξ)+(Ŝ424H4(ξ)+Ŝ3272H3(ξ))+⋯]
(See [27] for a review and [180] for an alternative derivation.) Where Hn(ξ) denotes the Hermite polynomials of order n and Ŝn are the dimensionless cummulants, defined by (127). Note that result (129) is an expansion in cummulants. This expression provides an accurate approximation of the true PDF provided the cummulants are well ordered in the sense that(130)1≫|Ŝ3|≫|Ŝ4|≫⋯.
See [181, 182] for more discussion on the structure of the correlation functions and the validity of cosmological perturbation theory.Finally, it remains to relate the cummulantsŜn (n≥3) to the correlation functions (such as the bispectrum) of the primordial curvature perturbation. We can write the dimensionless skewness, Ŝ3, as an integral over the 3-point correlation function [27, 180]:(131)Ŝ3=1σζ3∫d3k1(2π)3/2d3k2(2π)3/2d3k3(2π)3/2〈ζk1ζk2ζk3〉c.
A similar expression may be derived for the dimensionless kurtosis, Ŝ4. From (131) it is clear that the skewness provides a measure of the magnitude of the bispectrum integrated over all triangles. Hence, this provides a reasonable single number to compare the size of nongaussianity between different models. Moreover, the skewnesses (and higher cummulants) are of observational interest since these determine the probability of rare fluctuations [27, 166–168, 183–185]. Observables such as the abundance of collapsed objects may therefore be used to constrain Ŝ3, Ŝ4, and so forth.We have reviewed the derivation of the Edgeworth expansion (129) only to illustrate that the dimensionless skewness (131) encodes information about the bispectrum integrated over all momenta. However, our actual calculation of the PDF in the next subsection will be fully nonperturbative and hence does not rely on the validity of (129) or (130).


## 8.2. Magnitude of the Nongaussianity from IR Cascading

We have argued in the last subsection that the moments of the PDF provide a useful tool for quantifying the size of nongaussianity and drawing comparisons between microscopic models whose bispectra may have very different shapes and scaling properties. Our goal now is to compute the skewness and kurtosis generated by particle production, rescattering and, IR cascading in model (5). We will then compare this to the skewness that would be generated by more familiar constructions, such as the local model ζ=ζg+(3/5)fNLζg2.It is straightforward to extract the PDF numerically from our HLattice simulations by computing the fraction of the simulation box which contains the fieldδϕ∝ζ at a particular value. We again stress that this numerical evaluation of the PDF is non-perturbative; our results do not rely on any prior assumptions about the size or ordering of the cummulants.We present our results for the numerical evaluation of the PDF from IR cascading in Figure13. For illustration we have assumed a chaotic potential V(ϕ)=m2ϕ2/2; however, our results are largely insensitive to the choice of background inflation model.Figure 13
Probability density function ofζ for g2=0.1 at a series of different values of the scale factor, a. The dashed curve shows a Gaussian fit at late time a=6.2, and we have normalized the scale factor so that a=1 at the moment when particle production occurs. For illustration, we have chosen standard chaotic inflation V(ϕ)=(m2/2)ϕ2; however, our qualitative results are independent of the background inflation model.We can understand physically the behaviour of the PDF plotted in Figure13. Shorty after the initial burst of particle production, the curvature perturbation is extremely nongaussian, due to the sudden appearance of the source term J∝χ2 in the equation of motion for the inflaton fluctuations δϕ. Quickly, in less than an e-folding, nonlinear interactions begin to drive the system towards gaussianity. A similar behaviour has been observed in lattice simulations of out-of-equilibrium interacting scalar fields during preheating [186, 187]. In the case of rescattering during preheating, the system will eventually become gaussian when the fields thermalize. However, in our case the universe is still inflating during the process of rescattering and IR cascading. As a result, nongaussian inflaton fluctuations generated by rescattering are stretched out by the quasi-de Sitter expansion and must freeze once their wavelength crosses the Hubble scale. Hence, at late times the PDF does not become completely gaussian but rather freezes-in with some non-trivial skewness. Within a few e-foldings from the moment of particle production the time evolution of the PDF has become completely negligible.Given our numerical results for the PDF, it is straightforward to compute dimensionless cummulants such as the skewness (Ŝ3) and kurtosis (Ŝ4) for various values of the coupling g2. We have summarized our results in Table 4. At small coupling, g2<1, we always have |Ŝ3|>|Ŝ4| which suggests that the cummulants are well ordered. In particular, note that for g2=0.01, the kurtosis Ŝ4 is too small to be measured accurately from our HLattice simulations.Table 4
Dimensionless cummulants for IR cascading.
g2SkewnessKurtosisEquivalent(Ŝ3)(Ŝ4)fNLlocal1-0.510.16-28000.1-0.49-0.08-27000.01-0.033𝒪(10-3)-183In order to give a sense of the magnitude of the nongaussianity from particle production we have also computed an “equivalentfNLlocal” defined by 5Ŝ3/(18σζ). For a given g2 this effective fNLlocal is the magnitude of fNL which would be necessary to reproduce the skewness Ŝ3 of the IR cascading PDF using a local ansatz ζ=ζg+(3/5)fNLζg2. (Our sign conventions for fNL are consistent with WMAP [6]. See [27] for a discussion of various conventions employed in the literature.)From Table4 we see that IR cascading during inflation may generate significant nongaussianity. Even taking g2=0.01 (which is compatible with cosmological data for any choice of ϕ0), we still obtain a skewness Ŝ3=-0.033, which is the same value that would be produced by a local model with fNL=-183. This “equivalent” local nongaussianity is larger than current observational bounds, which suggests that nongaussian features from particle production during inflation might be observable for reasonable values of the coupling g2. Note, however, that observational bounds on fNLlocal cannot be directly applied to our model since the bispectrum in our case is uncorrelated with the vacuum fluctuations and is far from scale invariant.Depending on the value ofϕ0 (and hence the location of the feature in the power spectrum) the nongaussianity from IR cascading may lead to a variety of observational signatures. If the feature is localized on scales relevant for CMB experiments, then the key observable is the primordial bispectrum, B(ki). Efficient computation of the observational constraints on nongaussianity from inflationary particle production will require a separable template for the bispectrum, an issue which we briefly comment on in the next subsection and which will be considered in more detail in [15]. On the other hand, the feature might show up on smaller scales and leave an imprint on large scale structures. In this case there are at least two observables which could be used to constrain our model: higher order correlations of LSS probes (such as the galaxy bispectrum) and the abundance of collapsed objects. In the first case, we require the evolved bispectrum while in the latter case we require the PDF of the evolved density field δρ/ρ. Our model robustly predicts negative skewness for both the curvature perturbation ζ and the density field. Hence we should expect less rare objects (as compared to the gaussian case) on some characteristic scale. (The recent weak lensing measurement of the dark matter mass of the high-redshift galaxy cluster XMMUJ2235.3-2557 [188] has been construed as a possible hint of nongaussian initial conditions [184]. Unfortunately, our model does not produce the correct sign of skewness to explain such observations.)When computing the statistics of density fluctuations on some given length scaleR one should smooth the perturbation field over wave-numbers larger than kmax~R-1. This may be implemented in PDF (129) by the formal substitution ζk→𝒲(k)ζk in expressions (125) and (131) for the variance and skewness. Here 𝒲(k) is some window function which tends to zero for k≫R-1 (e.g., a standard choice of window is a top-hat in real space). This substitution makes any predictions for the statistics of the density fluctuations implicitly dependent on the smoothing scale, R. (It also renders the momentum integrals UV finite.) In our analysis we have computed the PDF numerically by measuring the fraction of the simulation box in which the inflaton field takes a given value. Thus, our cummulants are automatically smoothed on the UV cutoff scale kmax~Δx-1 defined by the lattice spacing Δx.One may wonder how sensitively our results for the dimensionless cummulants depend on this particular choice of smoothing. Recall that the bispectrum from IR cascading (73) is strongly peaked on triangles with a particular size, corresponding to the location of the bump-like distortion in the power spectrum. If we denote this characteristic scale by kbump, then it should be clear that the cummulants in our model are insensitive the smoothing scale, as long as R≲kbump-1.


## 8.3. Shape of the Bispectrum

In Section4 we presented a schematic discussion of how to compute analytically the primordial bispectrum from IR cascading, (73). In [14], this expression is evaluated analytically. We find that the bispectrum (73) from particle production has a unique shape that has not been considered in the previous literature. In order to describe this novel shape we would like to first factor out that strong overall dependence of Bϕ(ki) on the size of the triangle. To this end we define a “shape function” S(ki) in terms of the bispectrum as follows:(132)S(ki)=N-1(k1k2k3)2B(ki),
where B(ki) is related to Bϕ(ki) by (104), and N is a normalization factor to be discussed shortly. The function S(ki) has the advantage that the strong k6 running of the bispectrum is extracted. Hence, any residual scaling behaviour displayed by S(ki) must be a result of nonlinear interactions.Since we expect the 3-point correlation function to be of orderPζ2(k), a natural choice of normalization is N=(2π)4𝒫ζ2 where 𝒫ζ1/2=5×10-5 is the usual amplitude of the scale invariant fluctuations from inflation. With this choice of normalization our function S(ki) coincides with the quantity 𝒢(ki)/(k1k2k3) which was used to study localized nongaussian features from models with steps in the inflaton potential in [59, 60].Symmetry of the bispectrum under permutations of momenta implies that we can focus only on the regionk1≥k2≥k3, to avoid counting the same configuration twice. Moreover, the triangle inequality implies that 1-k2/k1≤k3/k1. Therefore we can completely specify the shape of the bispectrum for a given size of triangle k by plotting S(k,kx2,kx3) in the region x3≤x2≤1 and 1-x2≤x3. Because our bispectrum is very far from being scale invariant, it follows that this shape function is sensitive to the choice of k. Therefore, in Figure 14 we choose several representative choices: ln(k/kbump)=-1,0,1,2.The shape functionS(k,kx2,kx3), defined by (132), as a function of the dimensionless quantities x2,x3 which parametrize the shape of the triangle. (a) corresponds to k=e-1kbump, (b) is k=kbump, (c) is k=e+1kbump, and (d) is k=e+2kbump. In the IR (k≤kbump) the shape of the bispectrum is similar to the equilateral shape; however, there is also some support on flattened triangles near k~e+1kbump. At larger values of k the shape is unlike any other template proposed in the literature.
(a)(b)(c)(d)We see that a rich array of shapes is possible: fork≲kbump, the bispectrum is qualitatively similar to the equilateral model; however, at slightly larger k there is also considerable support on flattened triangles. Note that for k≳7.4kbump the shape of the bispectrum is extremely unusual and is not easily comparable to any shape that has been proposed in the previous literature.We find that the bispectrum from IR cascading is, to a first approximation, factorisable in the sense thatB(k1,k2,k3)≅∏i=13F(ki) where F(k)~C1k3e-C2k2 for some constants C1 and C2. This separable form is important since it permits fast algorithms for forecasting, data analysis, and simulations [189].In this section we have attempted to characterize the nongaussian signature associated with inflationary particle production. This signature is rather unique in the literature. Our model predicts uncorrelated localized nongaussian features with a unique shape of bispectrum. We have quantified the size of this effect by studying the dimensionless cummulants (such as the skewness) and have argued that this nongaussianity can be significant, depending ong2. We leave a detailed study of the observational constraints on such nongaussianities to future studies [15].


## 8.4. Multiple Bursts of Particle Production

Note that our discussion of nongaussianity generalizes easily to the case where there are multiple pointsϕi (i=0,…,n) along the inflaton trajectory where new degrees of freedom χi become massless, such as model (105). Such a construction may be quite natural in the context of brane/axion monodromy inflation [54, 77–79]. In this case the nongaussian features in B(ki) from each burst of particle production may superpose to generate a broad-band signal. Such a construction may (but need not) be associated with trapped inflation [54]. Depending on the spacing of the points ϕi and the couplings gi2, a rich variety of nongaussianities may be possible. In this case that the points ϕi are sufficiently densely spaced, we expect a nongaussian signal which is close to equilateral but also has some support on flattened triangles, consistent with the analysis of [54, 81]. We leave a detailed discussion to future studies.


## 9. Conclusions

In this paper we have considered the possibility that some non-inflation (iso-curvature) particles were produced during the observable range ofe-foldings of inflation. Inflationary particle production might occur as a result of a phase transition, parametric resonance, or other non-adiabatic processes. In order to illustrate the basic physics of inflationary particle production we restricted our analysis to a simple prototype model with coupling (g2/2)(ϕ-ϕ0)2χ2 between the inflaton ϕ and iso-inflaton χ. However, we expect that our qualitative results will apply also to SUSY models, gauged interactions, higher spin iso-inflaton fields, and phase transitions.Models of the type we study have attracted considerable interest recently in connections with trapped inflation, trans-Planckian effects, and observable features/nongaussianity in the primordial curvature fluctuations from inflation. Moreover, such models are quite natural from the microscopic perspective and may be obtained in popular models of open string inflation, such as brane/axion monodromy.We have shown that inflationary particle production in the model(g2/2)(ϕ-ϕ0)2χ2 leads to a new mechanism for generating cosmological fluctuations. This mechanism is qualitatively different from previous proposals in that this approach does not rely on the quantum vacuum fluctuations of light fields during inflation. Rather, the scenario involves the production of massive χ particles during inflation, which subsequently rescatter off the slow-roll condensate ϕ(t) to emit bremsstrahlung radiation of light inflaton fluctuations δϕ. We have studied this dynamics using classical lattice field theory simulations, analytical QFT computations, and also second-order cosmological perturbation theory. All of these approaches yield consistent results. We have found that rescattering proceeds with a time scale short compared to the expansion time. Moreover, the emission of long-wavelength inflaton fluctuations is very energetically inexpensive. The combination of these two effects leads to a rapid build-up of power in IR inflation fluctuations shortly after the moment of particle production. This dynamical process is called IR cascading.Our numerical and analytical studies of rescattering and IR cascading during inflation may have relevance for trapped inflation, preheating, moduli trapping, and also non-equilibrium QFT more generally. For instance, we have seen that, even with a small number of out-of-equilibriumχ particles, multiple rescatterings can nevertheless generate long-wavelength δϕ fluctuations with huge occupation numbers. We have also observed, for the first time, the dynamical approach to the scaling regime discovered in [96, 97].IR cascading during inflation leads to observable features in the primordial cosmological fluctuations. In particular, this process generates a bump-like contribution to the primordial scalar power spectrum. This signature is very different from what would be obtained in a model with transient violation of slow roll during inflation (such as a step-like feature inV(ϕ)), contrary to some claims in the literature.We have studied the observational constraints on bump-like features from inflationary particle production during inflation. We found that relatively large distortions, of order10% of the usual scale invariant vacuum fluctuations, are compatible with current data. We have derived observational bounds on the coupling g2 for a given ϕ0, which play a crucial role in determining the detectability of nongaussianity from particle production. Our observational bounds on particle production during inflation have implications for brane/axion monodromy inflation models and other microscopic constructions.IR cascading also has a nontrivial impact on nongaussian statistics, such as the bispectrum. The model(g2/2)(ϕ-ϕ0)2χ2 leads to a very novel nongaussian signature: uncorrelated, localized nongaussian features with a unique shape of bispectrum. For reasonable values of the coupling, g2≲0.01, this new kind of nongaussianity may be detectable in future missions.The nongaussian signature predicted by inflationary particle production is rather unusual, as compared to other models of inflation which are frequently studied in the literature. However, the underlying field theory description is extremely simple and rather generic from the low-energy perspective. In order to obtain observable nongaussianity it was not necessary to fine-tune the inflaton trajectory or appeal to re-summation of an infinite series of high-dimension operators. Indeed, the only “tuning” which is required for our signal to be observable is the requirement thatϕ=ϕ0 during the observable range of e-foldings.There are a variety of directions for future studies. From the theoretical perspective, it would be interesting to explicitly generalize our results to more complicated models with particle production during inflation (such as SUSY models, fermionχ fields, and phase transitions). There are also a wide range of interesting phenomenological possibilities. Varying the location of the feature we can have a variety of possible signatures for the CMB and LSS. We expect that IR cascading will also have implications for the spectrum of gravity waves from inflation and also primordial black holes. We can also imagine superposing multiple bursts of particle production to obtain an even richer variety of signatures. We leave these possibilities for future investigation.



---

*Source: 156180-2010-08-30.xml*