Time-sliced perturbation theory with primordial non-Gaussianity and effects of large bulk flows on inflationary oscillating features

TL;DR

This paper extends time-sliced perturbation theory to include primordial non-Gaussianity, revealing how large bulk flows damp primordial oscillating features through systematic IR resummation.

Contribution

It introduces a new formulation of TSPT that incorporates non-Gaussian initial conditions and clarifies the IR structure of non-Gaussian vertices.

Findings

New TSPT formulation with non-Gaussian initial conditions

Damping of primordial oscillations due to large bulk flows

Derivation of damping factors via IR resummation

Abstract

We extend the formalism of time-sliced perturbation theory (TSPT) for cosmological large-scale structure to include non-Gaussian initial conditions. We show that in such a case the TSPT interaction vertices acquire new contributions whose time-dependence factorizes for the Einstein-de Sitter cosmology. The new formulation is free from spurious infrared (IR) enhancements and reveals a clear IR structure of non-Gaussian vertices. We use the new technique to study the evolution of oscillating features in primordial statistics and show that they are damped due to non-linear effects of large bulk flows. We derive the damping factors for the oscillating primordial power spectrum and bispectrum by means of a systematic IR resummation of relevant Feynman diagrams.

Equations

$$⟨\zeta (\mathbf{k})\zeta ({\mathbf{k}}^{\prime })⟩=(2\pi {)}^3{\delta }_D^{(3)}(\mathbf{k}+{\mathbf{k}}^{\prime })P_{\zeta }(k),$$

$$⟨\zeta ({\mathbf{k}}_1)\zeta ({\mathbf{k}}_2)\zeta ({\mathbf{k}}_3)⟩=(2\pi {)}^3{\delta }_D^{(3)}({\mathbf{k}}_1+{\mathbf{k}}_2+{\mathbf{k}}_3)B_{\zeta }({\mathbf{k}}_1,{\mathbf{k}}_2,{\mathbf{k}}_3),$$

$$P_{\zeta }(k)\equiv \frac{{\mathcal{A}}_{\zeta }^2}{k^3}{\left(\frac{k}{k_{\ast }}\right)}^{n_s-1},$$

$${\mathcal{A}}_{\zeta }^2=4.14\cdot 10^{-8},n_s=0.965.$$

$$\displaystyle B^{\text{loc}.}_{\zeta}(k_{1},k_{2},k_{3})=\frac{6}{5}f_{NL}^{\text{loc}.}\Big{(}P_{\zeta}(k_{1})P_{\zeta}(k_{2})+P_{\zeta}(k_{2})P_{\zeta}(k_{3})+P_{\zeta}(k_{3})P_{\zeta}(k_{1})\Big{)}\,,$$

$$B_{\zeta }^{\text{eq}.}(k_1,k_2,k_3)=\frac{18}{5}{f}_{NL}^{\text{eq}.}{\mathcal{A}}_{\zeta }^4\left(-\frac{1}{k_1^3{k}_2^3}-\frac{1}{k_2^3{k}_3^3}-\frac{1}{k_3^3{k}_1^3}-\frac{2}{k_1^2{k}_2^2{k}_3^2}+\sum _{i\ne j\ne l}\frac{1}{k_i{k}_j^2{k}_l^3}\right),$$

$$B_{\zeta }^{\text{orth}.}(k_1,k_2,k_3)=\frac{18}{5}{f}_{NL}^{\text{orth}.}{\mathcal{A}}_{\zeta }^4\left(-\frac{3}{k_1^3{k}_2^3}-\frac{3}{k_2^3{k}_3^3}-\frac{3}{k_3^3{k}_1^3}-\frac{8}{k_1^2{k}_2^2{k}_3^2}+3\sum _{i\ne j\ne l}\frac{1}{k_i{k}_j^2{k}_l^3}\right),$$

$$f_{NL}^{\text{loc}.}=-0.9\pm 5.1,f_{NL}^{\text{eq}.}=-26\pm 47,f_{NL}^{\text{orth}.}=-38\pm 24.$$

$${\Theta }_L(\mathbf{k})={\delta }_L(\mathbf{k})\equiv T(t_0;\mathbf{k})\zeta (\mathbf{k}).$$

$$⟨{\delta }_L(\mathbf{k}){\delta }_L({\mathbf{k}}^{\prime })⟩=(2\pi {)}^3{\delta }_D^{(3)}(\mathbf{k}+{\mathbf{k}}^{\prime })P_L(k),$$

$$P_L(k)=[T(t_0;k){]}^2{P}_{\zeta }(k).$$

$$\begin{array}{rl}{B}_L(k_1,k_2,k_3)=& \sqrt{\frac{P_L(k_1)}{P_{\zeta }(k_1)}\frac{P_L(k_2)}{P_{\zeta }(k_2)}\frac{P_L(k_3)}{P_{\zeta }(k_3)}}{B}_{\zeta }(k_1,k_2,k_3).\end{array}$$

$$\frac{⟨{\delta }_L^3⟩}{(⟨{\delta }_L^2⟩{)}^{3/2}}\sim \frac{B_L}{(P_L{)}^{3/2}}\sim f_{NL}{\mathcal{A}}_{\zeta }\ll 1.$$

$$\frac{\partial \delta }{\partial t}+\nabla \cdot [(1+\delta )\mathbf{u}]=0,$$

$$\frac{\partial \mathbf{u}}{\partial t}+\mathcal{H}\mathbf{u}+(\mathbf{u}\cdot \nabla )\mathbf{u}=-\nabla \varphi ,$$

$${\nabla }^2\varphi =\frac{3}{2}{\mathcal{H}}^2{\Omega }_m\delta .$$

$$\eta =\ln D(t),$$

$$\Theta =-\frac{\nabla \cdot \mathbf{u}}{\mathcal{H}f},\text{with}f=\frac{d\ln D}{d\ln a}.$$

$${\delta }_{\mathbf{k}}\equiv \delta (\mathbf{k})=\int d^{3}x\delta (\mathbf{x})e^{-i\mathbf{k}\cdot \mathbf{x}},\delta (\mathbf{x})=\int \frac{d^{3}k}{(2\pi {)}^3}{\delta }_{\mathbf{k}}{e}^{i\mathbf{k}\cdot \mathbf{x}}.$$

$${\partial }_{\eta }{\delta }_{\mathbf{k}}-{\Theta }_{\mathbf{k}}=\int [dq{]}^2(2\pi {)}^3{\delta }_D^{(3)}(\mathbf{k}-{\mathbf{q}}_{12})\alpha ({\mathbf{q}}_1,{\mathbf{q}}_2){\Theta }_{{\mathbf{q}}_1}{\delta }_{{\mathbf{q}}_2},$$

$${\partial }_{\eta }{\Theta }_{\mathbf{k}}+\frac{1}{2}{\Theta }_{\mathbf{k}}-\frac{3}{2}{\delta }_{\mathbf{k}}=\int [dq{]}^2(2\pi {)}^3{\delta }_D^{(3)}(\mathbf{k}-{\mathbf{q}}_{12})\beta ({\mathbf{q}}_1,{\mathbf{q}}_2){\Theta }_{{\mathbf{q}}_1}{\Theta }_{{\mathbf{q}}_2},$$

$$[dq{]}^n\equiv \frac{d^3{q}_1...d^3{q}_n}{(2\pi {)}^{3n}},{\mathbf{q}}_{\mathrm{1...}n}={\mathbf{q}}_1+...+{\mathbf{q}}_n,$$

$$\alpha ({\mathbf{k}}_1,{\mathbf{k}}_2)\equiv \frac{({\mathbf{k}}_1+{\mathbf{k}}_2)\cdot {\mathbf{k}}_1}{k_1^2},\beta ({\mathbf{k}}_1,{\mathbf{k}}_2)\equiv \frac{({\mathbf{k}}_1+{\mathbf{k}}_2{)}^2({\mathbf{k}}_1\cdot {\mathbf{k}}_2)}{2k_1^2{k}_2^2}.$$

$$\displaystyle\delta_{\bf k}=\sum_{n=1}^{\infty}\big{(}D(\eta)\big{)}^{n}\int[dq]^{n}\,(2\pi)^{3}\delta^{(3)}_{D}({\bf k}-{\bf q}_{1...n})F_{n}({\bf q}_{1},...,{\bf q}_{n})\delta_{L}({\bf q}_{1})...\delta_{L}({\bf q}_{n})\,,$$

$$\displaystyle\Theta_{\bf k}=\sum_{n=1}^{\infty}\big{(}D(\eta)\big{)}^{n}\int[dq]^{n}\,(2\pi)^{3}\delta^{(3)}_{D}({\bf k}-{\bf q}_{1...n})G_{n}({\bf q}_{1},...,{\bf q}_{n})\delta_{L}({\bf q}_{1})...\delta_{L}({\bf q}_{n})\,,$$

$$\begin{array}{rl}& ⟨{\delta }_L({\mathbf{k}}_1){\delta }_L({\mathbf{k}}_2){\delta }_L({\mathbf{k}}_3)⟩=(2\pi {)}^3{\delta }_D^{(3)}({\mathbf{k}}_1+{\mathbf{k}}_2+{\mathbf{k}}_3)B_L(k_1,k_2,k_3).\end{array}$$

$$\displaystyle B_{\rm SPT}^{\rm tree}(\eta;k_{1},k_{2},k_{3})=2\big{(}D(\eta)\big{)}^{4}\sum_{1\leq i<j\leq 3}F_{2}({\bf k}_{i},{\bf k}_{j})P_{L}(k_{i})P_{L}(k_{j})+\big{(}D(\eta)\big{)}^{3}B_{L}(k_{1},k_{2},k_{3})\,,$$

$$\displaystyle P_{\text{SPT}}^{\text{1-loop}}(\eta;k)=\big{(}D(\eta)\big{)}^{4}\bigg{(}2\int[dq]P_{L}(q)P_{L}(|{\bf k}-{\bf q}|)\big{(}F_{2}({\bf q},{\bf k}-{\bf q})\big{)}^{2}$$

$$\displaystyle\qquad+6P_{L}(k)\int[dq]P_{L}(q)F_{3}({\bf k},{\bf q},-{\bf q})\bigg{)}+2\big{(}D(\eta)\big{)}^{3}\int[dq]B_{L}(q,|{\bf k}-{\bf q}|,k)F_{2}({\bf q},{\bf k}-{\bf q})\,.$$

$$F_{2}({\bf k},{\bf q})\Big{|}_{|{\bf q}|\to 0}\simeq\frac{{\bf k}\cdot{\bf q}}{2q^{2}}\;.$$

$$\delta({\bf k})\equiv\delta[\Theta;\eta,{\bf k}]=\sum_{n=1}^{\infty}\frac{1}{n!}\int[dq]^{n}(2\pi)^{3}\delta^{(3)}_{D}\big{(}{\bf k}-{\bf q}_{1\ldots n}\big{)}\,K_{n}(\eta;{\bf q}_{1},...,{\bf q}_{n})\prod_{j=1}^{n}\Theta({\bf q}_{j})\,.$$

$$\dot{\Theta}({\bf k})=\mathcal{I}[\Theta;\eta]\equiv\sum_{n=1}^{\infty}\frac{1}{n!}\int[dq]^{n}(2\pi)^{3}\delta^{(3)}_{D}\big{(}{\bf k}-{\bf q}_{1\ldots n}\big{)}\,I_{n}(\eta;{\bf q}_{1},...,{\bf q}_{n})\prod_{j=1}^{n}\Theta({\bf q}_{j})\,,$$

$$Z[J,J_{\delta};\eta]=\int[\mathcal{D}\Theta]\;{\mathcal{P}}[\Theta;\eta]\;\exp\bigg{\{}\int[dk]\Theta({\bf k})J(-{\bf k})+\int[dk]\delta[\Theta;\eta,{\bf k}]J_{\delta}(-{\bf k})\bigg{\}}\,.$$

$$(2\pi)^{3}\delta^{(3)}_{D}({\bf k}+{\bf k}^{\prime})\,P(\eta;k)=\frac{\delta^{2}Z}{\delta J_{\delta}(-{\bf k})\delta J_{\delta}(-{\bf k}^{\prime})}\bigg{|}_{J=J_{\delta}=0}\,.$$

Full text

CERN-TH-2019-091

INR-TH-2019-012

TTK-19-22

1 Introduction

Large scale structure (LSS) provides us with a powerful tool to study the dynamical and statistical properties of our universe. This structure has been formed during the epoch of matter domination and its dynamics on large scales is governed by dark matter and baryons acting as single pressureless perfect fluids. The three-dimensional distribution of LSS can potentially supersede the cosmic microwave background (CMB) measurements in the number of Fourier modes available for observations. On the other hand, the relevant cosmological information encoded in LSS is concealed by the non-linear effects related to gravitational clustering. Disentangling the cosmological information from these effects presents a highly non-trivial challenge that must be addressed in order to efficiently use the data of the ongoing and future surveys.

The standard tool to model dark matter clustering in the non-linear regime is N-body simulations, which are progressively becoming more robust and accessible. On the other hand, reaching the accuracy level required by future surveys still remains computationally expensive [1]. Besides, implementing non-minimal cosmological scenarios or initial conditions (e.g. primordial non-Gaussianity) still requires non-trivial extensions of the standard N-body codes.

An alternative way is the development of analytic methods based on perturbation theory [2]. Although limited to a range of large scales, this approach has a number of advantages. First, analytic methods allow for a shorter computational time. Second, they give valuable insights into the physics of non-linear clustering. Third, they are flexible and can be easily extended to models beyond the standard $\Lambda$CDM cosmology.

One of the goals of future LSS surveys is probing primordial non-Gaussianity (PNG), which can shed light on the dynamics of inflation or its alternatives. This task is complicated by a number of factors. First, generic non-Gaussian phenomena are encoded in higher-order statistics, which are hard to measure and interpret. Second, the primordial non-Gaussian signal must be disentangled from that generated by non-linear gravitational instability even from initially Gaussian fields. Perturbation theory provides a natural and computationally efficient way to do this [3]. Perturbative calculations including primordial non-Gaussianity were performed for the power spectrum [4] and bispectrum [5] of matter and biased tracers. The analysis of galaxy bias was later refined in [6], see [7] for a review. The Effective Field Theory (EFT) of LSS [8, 9] has been applied to compute the power spectrum and bispectrum with non-Gaussian initial conditions [10]. Based on these developments, an analysis of the ability to extract the PNG signal from LSS was carried out in [12, 11, 15, 13, 14].

A particularly interesting question is existence of oscillatory features in the statistics of primordial density fluctuations. Such resonance features naturally appear in models inspired by the axion monodromy inflation [16] where the inflaton potential is modulated by small sinusoidal oscillations [17, 18]. Alternatively, they can be generated by the interaction between the inflaton and massive particles, thereby encoding information about the particle spectrum at inflation [19, 20, 21]. The latter type of non-Gaussianity has got the name of cosmological collider signatures.

Oscillating features in the primordial power spectrum have been constrained from the CMB data by the Planck collaboration [22, 23]. The prospects to improve these constraints using future LSS and 21 cm surveys were studied in [24, 25, 26, 27, 28], see also the recent white paper [29] and references therein. The Planck collaboration has also performed searches for oscillating primordial non-Gaussianity in CMB with the null results [30]. The imprint of resonant non-Gaussianity on the scale dependent halo bias was discussed in [31], whereas Ref. [28] analyzed the sensitivity of 21 cm intensity mapping to this type of PNG. Refs. [32, 33] presented forecasts for the reach of LSS and 21 cm surveys in probing cosmological collider signatures. As pointed out in Ref. [34], one expects primordial oscillations in LSS statistics to get damped at late times due to the non-linear effects of large bulk flows, similarly to the damping of baryon acoustic oscillations (BAO). Indeed, the displacement of modes with Fourier momentum $k$ caused by non-linear coupling to long modes with momentum $q\ll k$ washes out the features with a characteristic period of oscillations in Fourier space $k_{osc}\ll q$ smaller that the long-mode momentum. However, this effect was not discussed beyond the power spectrum.

In this paper we study the non-linear damping of primordial oscillations in the inflationary 3-point function (bispectrum). We work out a systematic procedure of infrared (IR) resummation that captures the effect of large bulk flows. Our method is inspired by techniques developed in the context of the BAO and is based on Time-Sliced Perturbation Theory (TSPT) [35]. The latter has been proven to be a powerful tool for IR resummation of the power spectrum and higher-order statistics [36, 37].

TSPT is a novel way of performing cosmological perturbation theory utilizing the machinery of quantum field theory and statistical physics. The key idea of this method is to study the time-dependent probability distribution function of density and velocity fields, instead of these fields themselves as done in the standard cosmological perturbations theory (SPT). The central object is the generating functional for cosmological correlators whose time evolution is governed by a Liouville equation of statistical mechanics. A perturbative expansion of the generating functional leads to Feynman rules similar to that of a 3-dimensional Euclidean quantum field theory, in which time plays a role of an external parameter. This expansion produces equal-time correlation functions of density and velocity at fixed time (redshift) slices, which explains the name to the method.

TSPT deals directly with the statistical quantities, such as equal time n-point correlation functions, and thus provides a natural framework for studying non-Gaussian initial conditions. We will show that the incorporation of such conditions is straightforward and results in new contributions to the TSPT vertices. The new contributions take particularly simple form for the Einstein-de Sitter (EdS) cosmology, where time dependence of the TSPT vertices factors out. We introduce a new convenient diagrammatic technique for the primordial non-Gaussian contributions. As suggested by observations, the latter should be treated as small perturbations on top of the standard correlation function corresponding to Gaussian initial conditions. The eventual expressions for TSPT correlation functions at fixed order of perturbation theory coincide with the SPT result, although the individual diagrams are quite different. In particular, we will show that the non-Gaussian TSPT vertices do not contain any singularities originating from dynamical non-linear mode-coupling. In this sense they are IR-safe. Of course, they may still contain physical IR singularities which may be contained in the primordial statistics themselves.

Finally, we show how the new reformulation allows for a systematic IR resummation of the large bulk flow effects on the evolution of oscillatory primordial statistics. We show that IR resummation in this case closely resembles that of the BAO for the Gaussian initial statistics and leads to a smearing of oscillating features. In particular, the relevant contributions to be resummed have a graphical representation of daisy diagrams. We derive explicit expressions for the non-linear damping factors of the power spectrum and bispectrum. In the latter case the damping depends on the full configuration of the wavevectors, including their relative angles.

The paper is organized as follows. In Sec. 2 we set up conventions and discuss the standard treatment of non-Gaussian initial conditions in cosmological perturbation theory. In Sec. 3 we recapitulate the Gaussian TSPT framework. In Sec. 4 we extend it by including the non-Gaussian initial conditions. Section 5 contains the proof of the IR safety of the non-Gaussian TSPT vertices. In Sec. 6 we present a framework for IR resummation of non-Gaussian vertices and apply it to the non-linear evolution of the resonant power spectrum and bispectrum from axion monodromy inflation. We conclude in Sec. 7. Some additional material is given in Appendices.

In the numerical calculations throughout the paper we use the Planck 2018 cosmology [38] with massless neutrinos. The linear power spectrum is computed with the Boltzmann code CLASS [39].

2 Preliminaries

In this section we set up the conventions and discuss the standard treatment of non-Gaussian initial conditions in LSS. Throughout the paper we will mostly focus on the non-trivial initial 3-point function (bispectrum). However, our results can be straightforwardly generalized to the case of higher non-Gaussian statistics. We will comment on this point at the relevant places in the text.

2.1 Primordial and linear statistics

We will distinguish between primordial and linear bispectra. We adopt the inflationary paradigm and call primordial the quantities referring to the curvature perturbation $\zeta$ generated during inflation. Its statistical properties are encoded in the power spectrum and the bispectrum,

[TABLE]

where $\delta^{(3)}_{D}$ is the 3-dimensional Dirac delta-function. The power spectrum of the primordial curvature perturbations can be written as,

[TABLE]

where $k_{*}$ is a pivot scale, $\mathcal{A}^{2}_{\zeta}$ is the amplitude and $(n_{s}-1)$ is the tilt of the spectrum. Their mean values from the latest Planck CMB temperature and polarization measurements [38] for $k_{*}=0.05\,$Mpc*-1* are555Note that $\mathcal{A}_{\zeta}^{2}$ is related to the scalar amplitude $A_{s}$ used in the Planck analysis by $\mathcal{A}^{2}_{\zeta}=2\pi^{2}A_{s}$.

[TABLE]

Let us outline some properties of $B_{\zeta}({\bf k}_{1},{\bf k}_{2},{\bf k}_{3})$ in typical inflationary models, see Ref. [40] for more details. Due to momentum conservation and rotational invariance the bispectrum is a function of three independent variables only. These can be chosen to be either the three momenta norms $k_{1},k_{2},k_{3}$, or the norms of a pair of momenta, say $k_{1},k_{2}$, and the cosine $\mu_{12}$ of the angle between them. Approximate scale invariance implies that $B_{\zeta}$ is a homogeneous function of degree $-6$, symmetric in its arguments. The shapes commonly used in data analysis are666We use the Planck conventions [30] for the bispectrum amplitude $f_{NL}$. Note that the Planck collaboration uses the variable $\Phi\equiv 3\zeta/5$, instead of $\zeta$, in the analysis of non-Gaussianity.

local [41], equilateral [42], and orthogonal [43],

[TABLE]

where in the last two expressions we for simplicity neglected the departures from scale invariance (the tilt). The current observational constraints on $f_{NL}$ (at $68\%$ CL) from the CMB measurements are [30],

[TABLE]

The key objects for the description of LSS are the matter density contrast $\delta=(\rho-\bar{\rho})/\bar{\rho}$ and velocity divergence $\Theta\propto\nabla\cdot{\bm{u}}$ fields. Here $\bar{\rho}$ is the average density of the Universe and ${\bm{u}}$ is the peculiar flow velocity. Initially small, these fields start growing after recombination and become non-linear. This non-linear evolution is captured by the cosmological perturbation theory, which uses as a seed the linear fields $\delta_{L}$, $\Theta_{L}$ evolved up to the present epoch as if the perturbations were always in the linear regime. In the case of adiabatic initial conditions corresponding to a growing mode the two linear fields are identical. They are related to the curvature perturbation $\zeta$ by a transfer function $T(t_{0};{\bf k})$ which encodes the evolution of perturbations from inflation, through recombination, up to the present time $t_{0}$. Thus we write,

[TABLE]

We will understand by linear statistics the properties of the field $\delta_{L}$. Defining the linear matter power spectrum,

[TABLE]

we see that

[TABLE]

From this we immediately deduce the relation between the linear bispectrum of LSS and the primordial one,

[TABLE]

Note that the linear bispectrum is parametrically suppressed by the amplitude of the primordial fluctuations. Indeed, neglecting the shape of the bispectrum one has, schematically,

[TABLE]

The linear power spectrum (2.8) and bispectrum (2.9) serve as the input for the equations describing the non-linear evolution of LSS statistics at times after recombination.

2.2 Standard perturbation theory with non-Gaussian initial

conditions

We are interested in the correlation functions of the overdensity and the velocity divergence fields, whose time-evolution is governed by the following pressureless perfect fluid equations:

[TABLE]

where $\phi$ is the gravitational potential obeying the Poisson equation,

[TABLE]

Here $t$ is the conformal time777Defined as $dt=d\tau/a(\tau)$, where $\tau$ is physical time and $a(\tau)$ is the scale factor of the Friedmann-Lemaitre-Robertson-Walker metric, $ds^{2}=d\tau^{2}-a^{2}(\tau)dx_{i}dx^{i}$., the spatial derivatives are taken with respect to the comoving coordinates, ${\cal H}=aH$ is the rescaled Hubble parameter, $a$ is the scale factor and $\Omega_{m}$ is the matter density fraction.

The single-stream perfect fluid approximation breaks down at short scales due to free-streaming of dark matter particles (also called ‘shell-crossing’). These effects show up when one computes loop contributions generated by hard modes. They are taken into account in the EFT of LSS, where the influence of short scale nonlinearities on physics at large scales is encapsulated by an effective stress tensor added to the r.h.s. of (2.11) [8, 9]. For simplicity, we do not explicitly consider the EFT corrections in this paper, keeping in mind that eventually they will have to be properly included. The main goal of this paper is to study the non-linear effects of very long-wavelength (IR) modes, for which the perfect fluid description is sufficient.

It is well-known [2] that in the case of an EdS universe dominated by non-relativistic matter ($\Omega_{m}=1$) the above equations can be cast in a form free from any explicit time dependence. This is achieved by introducing the time parameter

[TABLE]

where $D$ is the linear growth factor888Following the standard practice, we normalize $D(t)$ to be equal to 1 at the present epoch., and appropriately rescaling the velocity divergence,

[TABLE]

For the realistic $\Lambda$CDM cosmology, the above substitution leaves a small residual time dependence which, however, has little effect on the dynamics [44]. Following conventional practice we will neglect this explicit time dependence in the equations, even though our analysis does not crucially depend on this restriction.

Within the above approximation Eqs. (2.11) can be rewritten in Fourier space as999Our conventions are:

$\delta_{\bf k}\equiv\delta({\bf k})=\int d^{3}x\,\delta({\bf x})e^{-i{\bf k}\cdot{\bf x}}~{},\qquad\delta({\bf x})=\int\frac{d^{3}k}{(2\pi)^{3}}\,\delta_{\bf k}e^{i{\bf k}\cdot{\bf x}}\,.$

Note that they differ from those used in [35, 36] by factors of $(2\pi)^{3}$.

[TABLE]

where we used the notations,

[TABLE]

and introduced the non-linear kernels

[TABLE]

These equations can be solved perturbatively by using the following power series Ansatz:

[TABLE]

where $\delta_{L}$ is the linear density field and the non-linear kernels $F_{n}$ and $G_{n}$ are recursively derived from (2.15). The correlation functions of the fields $\delta_{\bf k}$, $\Theta_{\bf k}$ are computed by averaging the expressions (2.17) using the statistics of the linear field $\delta_{L}$. The resulting expressions are conveniently represented as diagrammatic expansions in the number of loops. This method is known as Eulerian standard perturbation theory (SPT) [2]. Typically, the linear statistics are assumed to be Gaussian and adiabatic, i.e. they are fully characterized by the power spectrum (2.7). Recall that the initial conditions for structure formation are set shortly after recombination when baryons and dark matter started to behave as a single pressureless fluid, but the gravitational instability did not yet have enough time to form non-linear structures. For simplicity, we neglect non-Gaussianity generated by the second order effects at radiation domination and recombination [45, 46, 47].

Suppose now that in addition to the power spectrum the linear distribution of the density contrast is also characterized by a non-trivial 3-point function,

[TABLE]

This generates additional terms in the SPT perturbative expansion, which originate from averages of odd number of fields. For instance, the leading-order (tree-level) matter bispectrum and the 1-loop correction to the matter power spectrum read,

[TABLE]

Note that the primordial non-Gaussian contributions have one less power of the growth factor $D(\eta)$ compared to the terms coming from the non-linear evolution.

Let us discuss an important point. The kernels $\alpha$ and $\beta$ entering the Euler equations have poles if one of their arguments vanishes, see Eqs. (2.16). These IR singularities are inherited by the non-linear kernels $F_{n}$, $G_{n}$; for example,

[TABLE]

This leads, in turn, to IR poles in the integrands of the loop expressions. The origin of the IR singularities can be traced back to large displacements of short-wavelength density perturbations by long-wavelength bulk flows. This effect should, however, cancel in the correlation functions of fields taken at the same moment of time. Indeed, the cancellation of IR singularities in equal-time correlators is a well-known property in Gaussian SPT, see [48] and references therein. For non-Gaussian initial conditions, the cancellation of contributions due to the IR poles can be tracked explicitly at the lowest orders of the perturbation theory; however, we are not aware of a general proof in the literature. Clearly, it is desirable to have a framework where spurious IR singularities are absent altogether. Such framework is provided by the time-sliced perturbation theory [35] described in the subsequent sections.

It is important to stress that the IR poles discussed above should be distinguished from those that may be present in the primordial bispectrum itself and thus are physical. Throughout the paper ‘IR poles’ or ‘IR singularities’ will only refer to the singularities originated from the non-linear mode coupling.

3 Review of Gaussian TSPT

The main idea of the TSPT approach is to substitute the time evolution of the overdensity and velocity divergence fields, $\delta$ and $\Theta$, by that of the their time dependent probability distribution functional (PDF). This idea is natural when one is only interested in equal time correlation functions. For adiabatic initial conditions only one of the two fields is statistically independent. It is convenient to choose the velocity divergence field $\Theta$ as an independent variable and its PDF will be denoted by $\mathcal{P}[\Theta;\eta]$. At any moment in time, the field $\delta$ can be expressed in terms of $\Theta$ as

[TABLE]

Note that, in contrast to the SPT formulae (2.17), we have here on the r.h.s. the full non-linear field $\Theta({\bf q})$. Recursion relations for the kernels $K_{n}$ are given in Appendix A. Equation (3.1) can be used to eliminate the density field from Eq. (2.15b) and obtain the following equation for the velocity divergence,

[TABLE]

with $I_{1}\equiv 1$ corresponding to the adiabatic mode in the perfect fluid approximation. The other kernels $I_{n}$ can be derived from the fluid equations (2.15), see Appendix A. In the EdS approximation both kernels $K_{n}$ and $I_{n}$ are time independent.

One introduces the generating functional,

[TABLE]

Equal-time correlation functions for $\Theta$ and $\delta$ are obtained by taking functional derivatives of $Z$ with respect to the external sources $J$ or $J_{\delta}$, respectively. For example, the matter power spectrum is given by

[TABLE]

The conservation of probability implies the Liouville equation for the probability density functional,

[TABLE]

It is useful to expand $\mathop{\rm ln}\nolimits\mathcal{P}[\Theta;\eta]$ as a power series in $\Theta$,

[TABLE]

where ${\cal N}$ is a normalization factor. Then, using Eq. (3.2) we obtain the following chain of equations for the vertices $\Gamma^{tot}_{n}$,

[TABLE]

where the sum in the second term on the l.h.s. is taken over all permutations $\sigma$ of $n$ indices. It is convenient to decompose the solution of this equation into two pieces,

[TABLE]

where $\Gamma_{n}$ is the solution of Eq. (3.7) with vanishing r.h.s. and with the initial conditions reflecting the initial statistical distribution, whereas $C_{n}$ is the solution of the inhomogeneous equation with vanishing initial conditions. The vertices $\Gamma_{n}$ thus have a physical meaning of 1-particle irreducible contributions to the tree-level correlators with amputated external propagators; the vanishing of the average velocity dispersion, $\langle\Theta({\bf k})\rangle=0$, implies $\Gamma_{1}=0$. On the other hand, $C_{n}$ are counterterms that cancel certain ultraviolet divergences appearing in the loop corrections [35].

The vertices $\Gamma_{n}$ and $C_{n}$ satisfy a hierarchy of equations which replaces the dynamical equations of SPT. From now on we specify to the EdS approximation where the kernels $I_{n}$ are constant in time. In order to find $\Gamma_{n}$, we use an Ansatz which separates time and momentum dependence,

[TABLE]

This leads to the following recursion relations for $\Gamma_{n}^{(l)}$ with $l<n$,

[TABLE]

Note that the r.h.s. contains only vertices $\Gamma^{(l)}_{n^{\prime}}$ with $n^{\prime}$ less than $n$, implying that all $\Gamma_{n}^{(l)}$, $l<n$, are uniquely determined from vertices of lower orders. The remaining function $\Gamma_{n}^{(n)}$ must be fixed from the initial conditions. Not to overload the paper with unnecessary formulae, we do not reproduce here the equations for the counterterms $C_{n}$ which can be found in [35]. The upshot is that $C_{n}$ are completely fixed in terms of the kernels $I_{n}$ and are time independent.

To sum up, the PDF (3.6) is fully specified by providing the ‘diagonal’ vertices $\Gamma_{n}^{(n)}$, $n\geq 2$. These define the early-time asymptotics of the full vertex functions,

[TABLE]

Note that in the last formula we used a convenient trick of sending the initial time $\eta_{i}$ to $-\infty$ assuming the validity of the pressureless fluid description (see Eqs. (2.11)) at all times. Of course, this does not correspond to a realistic Universe dominated by radiation at an early epoch. Physically, the initial data for Eq. (3.7) should be set somewhen after recombination. However, mapping the initial conditions to $\eta_{i}=-\infty$ by formally extrapolating the EdS description towards the past leads to a great simplification of formulas.

For the Gaussian initial conditions the solution to (3.10) is greatly simplified. In that case the PDF ${\cal P}[\Theta;\eta]$ should reduce to a Gaussian weight at early times, upon rescaling the field with the linear growth factor, $\Theta\mapsto\tilde{\Theta}=\Theta/D(\eta)$. In other words, one requires,

[TABLE]

where $P_{L}(k)$ is the linear power spectrum. This implies the following initial conditions for the vertices,

[TABLE]

Then all $\Gamma^{(l)}_{n}$ with $l>2$ vanish, leaving the solution,

[TABLE]

where we have introduced the notations $\bar{\Gamma}_{n}\equiv\Gamma_{n}^{(2)}$ and

[TABLE]

Remarkably, in the case of Gaussian initial conditions and the EdS background, all $\Gamma_{n}$ vertices have universal dependence on time through the factor $g^{-2}(\eta)$. As will be discussed shortly, $g(\eta)$ plays the role of expansion parameter in TSPT. Due to momentum conservation, the vertices are proportional to a $\delta$-function of the sum of their arguments. In what follows we use primes to denote the quantities stripped off such $\delta$-functions,

[TABLE]

Note that the recursion relations (3.14b) and the initial conditions (3.13a) imply that all $\bar{\Gamma}_{n}$ vertices are functionals of the linear power spectrum $P_{L}(k)$.

The TSPT perturbative series is obtained upon expanding the generating functional (3.3) over the Gaussian part of $\mathcal{P}$. Following the standard rules of quantum field theory, this expansion can be represented as a sum of Feynman diagrams. These are built from vertices corresponding to $\Gamma_{n}$, $n\geq 3$, and lines corresponding to the ‘propagator’ $g^{2}P_{L}(k)$, see Fig. 3. One should also include vertices corresponding to counterterms $C_{n}$, $n\geq 1$, in order to subtract certain ultraviolet divergences of loop diagrams. To compute an $n$-point correlation function of the velocity divergence $\Theta$ one needs to draw all diagrams with $n$ external legs. It is straightforward to see that diagrams with higher number of loops are proportional to higher powers of $g(\eta)$. Hence, $g(\eta)$ should be interpreted as the coupling constant of the theory. For the correlators of the density field $\delta$ one uses the expression (3.1) which is reminiscent of expressions for composite operators in quantum field theory. It gives rise to additional vertices proportional to the kernels $K_{n}$; these are denoted by an external arrow.