Impact of shell crossing and scope of perturbative approaches in real and redshift space.

Impact of shell crossing and scope of perturbative approaches in real and redshift space

Key Words.:
gravitation; cosmology: theory – large-scale structure of Universe

Abstract

Context:

Aims:We study the effect of nonperturbative corrections associated with the behavior of particles after shell crossing on the matter power spectrum. We compare their amplitude with the perturbative terms that can be obtained within the fluid description of the system, to estimate the range of scales where such perturbative approaches are relevant.

Methods:We use the simple Zeldovich dynamics as a benchmark, as it allows the exact computation of the full nonlinear power spectrum and of perturbative terms at all orders. Then, we introduce a “sticky model” that coincides with the Zeldovich dynamics before shell crossing but shows a different behavior afterwards. Thus, their power spectra only differ in their nonperturbative terms. We consider both the real-space and redshift-space power spectra.

Results:We find that the potential of perturbative schemes is greater at higher redshift for a CDM cosmology. For the real-space power spectrum, one can go up to order of perturbation theory at , and to order at , before the nonperturbative correction surpasses the perturbative correction of that order. This allows us to increase the upper bound on where systematic theoretical predictions may be obtained by perturbative schemes, beyond the linear regime, by a factor at and at . This provides a strong motivation to study perturbative resummation schemes, especially at high redshifts .

In the context of cosmological reconstruction methods, the Monge-Ampère-Kantorovich scheme appears to be close to optimal at . There also seems to be little room for improvement over current reconstruction methods of the baryon acoustic oscillations at . This can be understood from the small number of perturbative terms that are relevant at , before nonperturbative corrections dominate.

We also point out that the rise of the power spectrum on the transition scale to the nonlinear regime strongly depends on the behavior of the system after shell crossing.

We find similar results for the redshift-space power spectrum, with characteristic wavenumbers that are shifted to lower values as redshift-space distortions amplify higher order terms of the perturbative expansions while decreasing the resummed nonlinear power at high .

Conclusions:

1 Introduction

The growth of large-scale structures in the Universe through the amplification of small primordial fluctuations by gravitational instability is a key ingredient of modern cosmology (Peebles 1980), and it can be used to constrain cosmological parameters through the dependence of the matter power spectrum on scale and redshift. On very large scales or at high redshifts, where the amplitude of the density fluctuations is small, it is sufficient to use linear theory, whereas on small scales, in the highly nonlinear regime, one must use numerical simulations or phenomenological models, such as the halo model (Cooray & Sheth 2002), which are also calibrated on simulations. In the weakly nonlinear regime one expects perturbative approaches to provide a useful tool, as they allow going beyond linear theory in a systematic and controlled fashion. Several observational probes, such as weak lensing surveys (Massey et al. 2007; Munshi et al. 2008) or measures of acoustic baryonic oscillations (Eisenstein et al. 1998, 2005), are mostly sensitive to these intermediate scales, and to meet the accuracy of future observations we need a theoretical accuracy of about . Phenomenological models typically have an accuracy of in this range, while numerical simulations may suffer from finite resolution and finite size effects and require a long computational time if we need to obtain the power spectra over a fine grid of cosmological parameters.

This has led to a renewed interest in perturbative approaches, as it may be possible to improve over the standard perturbation theory (Goroff et al. 1986; Bernardeau et al. 2002) by using resummation schemes that allow systematic partial resummations of higher order terms . Thus, Crocce & Scoccimarro (2006b, a) present a partial resummation of the diagrammatic series associated with the response function (propagator), within a high- limit, which provides improved predictions for the density power spectrum (Crocce & Scoccimarro 2008). On the other hand, Valageas (2007a) describes a path-integral formalism that allows applying the tools of field theory, such as large- expansions, to compute the power spectrum and higher order statistics like the bispectrum (Valageas 2008). One of these large- expansions was recovered by Taruya & Hiramatsu (2008), as a “closure theory” where one closes the hierarchy of equations satisfied by the many body correlations at the third order, following the “direct interaction approximation” introduced in hydrodynamics (Kraichnan 1959). This also improves the predictions for the power spectrum on the scales probed by the baryonic acoustic oscillations (Taruya et al. 2009; Valageas & Nishimichi 2010). Other approaches have been proposed by Matarrese & Pietroni (2007), using the dependence on a running high- cutoff, by Pietroni (2008), using a truncation of the hierarchy satisfied by the many-body correlations, by Matsubara (2008), within a Lagrangian framework, and by McDonald (2007), using a renormalization group technique. Most1 of these approaches start from the fluid description of the system, where the density and velocity fields obey hydrodynamical equations of motion. This corresponds to a single-stream approximation that neglects shell crossing. Then, the domain of validity of most such perturbative schemes is limited to wavenumbers where shell-crossing effects are negligible, even if we could sum all perturbative terms.

This problem with the impact of shell crossing also appears in the context of cosmological reconstruction, where one attempts to follow the matter distribution observed in a given galaxy survey back in time (Peebles 1989). Thus, an efficient algorithm for building such a reconstruction (which can then be used to estimate the velocity field) is provided by the Monge-Ampère-Kantorovich method, which neglects multistreaming (Brenier et al. 2003; Mohayaee et al. 2006). Similar methods are used to reconstruct the baryon acoustic oscillations (BAO) of the density power spectrum, in order to improve the accuracy of cosmological distance measurements and tighten the constraints on cosmological parameters such as the amount and evolution of dark energy (Eisenstein et al. 2007; Seo et al. 2010). Then, by comparing the amplitude of nonperturbative corrections with the perturbative terms obtained within the fluid description, one can estimate the scale down to which these reconstruction schemes can be used.

Thus, to estimate the potential of such approaches, it is necessary to evaluate the effect of shell crossing on the matter power spectrum. This question has already been investigated by Afshordi (2007) by comparing the phenomenological halo model with a modified variant where halos are collapsed to pointlike masses. In this paper we revisit this problem in a more systematic fashion, within the framework of the Zeldovich dynamics. Thus, we compare the Zeldovich dynamics with a second model (named the “sticky model” in the following), which only differs after shell crossing. Therefore, both power spectra have the same perturbative expansions and only differ by nonperturbative terms. Taking advantage of our being able to explicitly compute the full nonlinear spectra, as well as perturbative terms at all orders and these nonperturbative corrections, we can compare their respective amplitudes in detail. This allows a more detailed discussion of the importance of shell-crossing effects and of the scope of perturbative approaches. This also enables us to distinguish the dependence on the behavior of the dynamics after shell crossing of the rise of the density power spectrum on mildly nonlinear scales, in the intermediate regime where the logarithmic power goes from to . In addition, we can perform the same analysis for the redshift-space power spectrum, which is actually the quantity most directly observed in galaxy surveys.

This article is organized as follows. We first recall in Sect. 2.1 the nonlinear real-space power spectrum associated with the Zeldovich dynamics, and its perturbative expansions in Sects. 2.2 and 2.3. We present our “sticky model”, which only differs from the Zeldovich dynamics after shell crossing, in Sect. 2.4, and we give the associated nonperturbative correction. Then, we describe the numerical results obtained for a CDM cosmology in Sect. 3, comparing the various perturbative and nonperturbative terms. Finally, we extend our analysis to the redshift-space power spectrum in Sect. 4, focusing on wavenumbers that are aligned with the line of sight. We conclude in Sect. 5.

2 Computation of the density power spectrum in real space

2.1 Nonlinear Zeldovich power spectrum

As is well known, the Zeldovich approximation (Zeldovich 1970) sets the Eulerian position, , of the particle of Lagrangian coordinate , equal to the position given by the linear displacement field, ,

 x(q,t)=q+ΨL(q,t)with∇q⋅ΨL=−δL(q,t), (1)

where is the linear growing mode of the density contrast, which is defined by

 δ(x,t)=ρ(x,t)−¯¯¯ρ¯¯¯ρ. (2)

Here is the mean matter density of the Universe, we work in comoving coordinates and as usual we only consider the linear growing mode. It is well known (Schneider & Bartelmann 1995; Taylor & Hamilton 1996) that the explicit expression of the matter power spectrum can be derived from Eq.(1) by using the conservation of matter, which reads as

 ρ(x)dx=¯¯¯ρdq% whence1+δ(x)=∣∣∣det(∂x∂q)∣∣∣−1. (3)

For an arbitrary displacement field , this also reads as

 1+δ(x)=∫dqδD[x−q−Ψ(q)], (4)

where is the Dirac distribution, and this yields in Fourier space (for , that is, disregarding a term ):

 ~δ(k)=∫dx(2π)3e−ik⋅xδ(x)=∫dq(2π)3e−ik⋅(q+Ψ). (5)

Defining the density power spectrum as

 ⟨~δ(k1)~δ(k2)⟩=δD(k1+k2)P(k1), (6)

we obtain from Eq.(5), using statistical homogeneity,

 P(k)=∫dq(2π)3⟨eik⋅[x(q)−x(0)]⟩. (7)

Equation (7) is quite general since we have not used the Zeldovich approximation (1) yet. It shows how the density power spectrum is related to the statistical properties of the displacement field, for any mapping (which can include some shell-crossing, as in the Zeldovich case). The great simplification provided by the Zeldovich approximation (1) is that in this case the quantity is a Gaussian random variable, so that the mean can be computed at once as

 P(k)=∫dq(2π)3eik⋅qe−12⟨(k⋅[ΨL(q)−ΨL(0)])2⟩. (8)

Next, using the second relation (1) to compute the average , we obtain the explicit expression

 P(k)=∫dq(2π)3eik⋅qe−∫dw[1−cos(w⋅q)](k⋅w)2w4PL(w). (9)

It is convenient to perform the integration over the angles of the wavenumber by expanding over spherical harmonics (Schneider & Bartelmann 1995),

 ∫dweiw⋅q(k⋅w)2w4PL(w)=k2I0(q)+k2(1−3μ2)I2(q), (10)

where , and we introduced

 Iℓ(q)=4π3∫∞0dwPL(w)jℓ(qw), (11)

where is the spherical Bessel function of order . In particular, the variance of the one-dimensional displacement field (or of the linear velocity field, up to a time-dependent multiplicative factor), reads as

 σ2v=13⟨|ΨL|2⟩=I0(0). (12)

Therefore, the expression (9) also writes as

 P(k)=∫dq(2π)3cos(kqμ)e−k2[σ2v−I0(q)−(1−3μ2)I2(q)]. (13)

Following Schneider & Bartelmann (1995), we can perform the integration over the angles of by expanding part of the exponential, and using the property

 ∫10dμcos(kqμ)(1−μ2)ℓ=ℓ!(2kq)ℓjℓ(kq), (14)

which gives

 P(k) = Missing or unrecognized delimiter for \right

2.2 Standard perturbative expansion

In the standard perturbative approach (Goroff et al. 1986; Bernardeau et al. 2002), one writes the nonlinear density contrast as a series over powers of the linear growing mode,

 ~δ(k)=∞∑n=1~δ(n)(k)with~δ(n)(k)∝(~δL)n; (16)

that is,

 ~δ(n)(k) = ∫dw1..dwnδD(w1+..+wn−k)Fn(w1,..,wn) (17) ×~δL(w1)..~δL(wn).

Substituting this expansion (and the one associated with the velocity field) into the equations of motion one obtains a recursion relation for the kernels , which allows terms of increasing order to be computed in a sequential manner. As is well known, for the Zeldovich dynamics, where , we do need to follow this route, since by expanding the exponential (5) over , and using , we obtain at once all terms,

 Fn(w1,..,wn)=1n!k⋅w1w21...k⋅wnw2n. (18)

Then, substituting the expansion (16) into the definition (6) of the power spectrum and taking the Gaussian average, one obtains the standard perturbative series

 P(k)=∞∑n=1P(n)(k)withP(n)(k)∝(PL)n. (19)

In particular, the two lowest order terms are

 P(1)(k)=PL(k),P(2)(k)=P22(k)+P13(k), (20)

where and arise from terms of the form and , with

 P22(k) = ∫dw1dw2δD(w1+w2−k)(k⋅w1)2(k⋅w2)22w41w42 (21) ×PL(w1)PL(w2),

and

 P13(k)=−k2σ2vPL(k). (22)

In fact, for the Zeldovich dynamics it is not necessary to use this procedure, since by expanding the exponential in the exact result (9) or (LABEL:Pkjn) we obtain all terms at once,

 P(n)(k) = ∫∞0dqq22π2n∑p=01p![k2(I0(q)−2I2(q)−σ2v)]p (23) ×(6kI2(q)q)n−pjn−p(kq).

2.3 Renormalized perturbative expansion

As pointed out by Crocce & Scoccimarro (2006b), it is useful to keep the -independent exponential term in Eq.(13) and to define a “renormalized” perturbative expansion

 P(k)=e−k2σ2v∞∑n=1P(n)σv(k)withP(n)σv(k)∝(PL)n. (24)

Indeed, while the standard perturbative terms grow increasingly fast at high and show large cancellations, the “renormalized” terms are positive and show a Gaussian decay at high . Thus, each term peaks on a well-defined range of and one can clearly see the contribution of diagrams of a given order to the full nonlinear power spectrum. In fact, as noticed in Crocce & Scoccimarro (2006b), these terms also write as

 P(n)σv(k) = n!∫dw1..dwnδD(w1+..+wn−k) (25) ×Fn(w1,..,wn)2PL(w1)..PL(wn).

This can be seen at once by expanding the expression (9) and recognizing the square of the kernel given in Eq.(18). Thus, each new term is associated with the kernel that couples linear modes, so that one can follow the contribution of higher order mode couplings. Again, from Eq.(LABEL:Pkjn) we obtain all terms at once,

 P(n)σv(k) = ∫∞0dqq22π2n∑p=01p![k2(I0(q)−2I2(q))]p (26) ×(6kI2(q)q)n−pjn−p(kq),

while the two lowest order terms are

 P(1)σv(k)=PL(k),P(2)σv(k)=P22(k), (27)

where was given in Eq.(21).

For numerical purposes, it is convenient to obtain the standard perturbative terms (19) from the renormalized ones (24). Expanding the prefactor in Eq.(24) gives

 P(n)(k)=n−1∑p=01p!(−k2σ2v)pP(n−p)σv(k). (28)

On the other hand, to compute the full nonlinear power spectrum (LABEL:Pkjn) it is convenient to subtract the two terms , so as to avoid integrating over badly behaving terms for .

2.4 Nonperturbative correction

The expressions recalled in the previous sections give the nonlinear power spectrum and its perturbative expansions for the usual Zeldovich dynamics, where particles follow the linear trajectories (1). This includes shell crossing and leads to a decay of the nonlinear power spectrum at high , as particles freely escape to infinity, and these random trajectories erase small-scale features (Schneider & Bartelmann 1995; Taylor & Hamilton 1996; Valageas 2007b; Bernardeau & Valageas 2010a). However, because of the simple nature of the Zeldovich dynamics, the full nonlinear result (13) can be obtained by resumming2 the perturbative expansions (19) or (24). The latter can be obtained from the usual hydrodynamical equations of motion in Eulerian space, which actually break down at shell crossing.

Our goal in this work is to compare the perturbative terms (19) and (24) with the nonperturbative terms that are associated with the physics that takes place after shell crossing. Within the framework of the Zeldovich dynamics considered in this article, this means that we wish to compare the previous results with those that would be obtained for a second dynamics, which coincides with the Zeldovich dynamics until shell crossing. Then, the difference between both power spectra would give us the amplitude of the correction due to shell-crossing effects, which are generically nonperturbative.

To introduce this second model, let us first recall that the linear growing mode of the displacement field, , and the associated peculiar velocity field, , are curlfree (Peebles 1980), and are derived from a velocity potential, . Using the Poisson equation we can see that this potential is equal to the linear gravitational potential , up to a time-dependent factor (Peebles 1980; Vergassola et al. 1994). In particular, the Lagrangian-space to Eulerian-space mapping defined by the Zeldovich dynamics (1) derives from a Lagrangian potential ,

 x(q,t)=∂φL∂q, (29)

with

 φL(q,t)=|q|22−a(t)4πG¯¯¯ρΦL(q,t), (30)

where is the linear gravitational potential, which obeys the Poisson equation,

 ΔqΦL=4πG¯¯¯ρa(t)δL. (31)

Here is Newton’s constant and the scale factor. Then, from Eq.(3) the nonlinear density contrast is given by the Hessian determinant of ,

 1+δ(x)=∣∣∣det(∂x∂q)∣∣∣−1=∣∣ ∣∣det(∂2φL∂qi∂qj)∣∣ ∣∣−1. (32)

At early times, , (and if there is not too much power on small scales), the Lagrangian potential is dominated by the first term in Eq.(30) and . Then, the Hessian matrix is definite positive and goes to the identity matrix, whence . As time increases and structures form, the Lagrangian potential becomes increasingly sensitive to the fluctuations in the linear gravitational potential, and the Hessian determinant deviates from unity. However, the eigenvalues remain strictly positive until shell crossing, which means that the Lagrangian potential remains a strictly convex function. Then, at shell crossing one eigenvalue goes through zero and becomes negative (generically the collapse proceeds at different rates along the three axes). Thus, as is known (Vergassola et al. 1994; Brenier et al. 2003; Bernardeau & Valageas 2010b), the onset of shell crossing is associated with the change in sign of the Hessian determinant of and with the loss of convexity of the Lagrangian potential (the Hessian matrix is no longer positive-definite).

A well-studied dynamics that agrees with the Zeldovich dynamics until shell crossing is provided by the “adhesion model” (Gurbatov et al. 1989, 1991; Vergassola et al. 1994). More precisely, the “geometrical adhesion model” (Bernardeau & Valageas 2010b; Valageas & Bernardeau 2010) leads to replacing the linear Lagrangian potential , which defines the mapping through (29), by a nonlinear Lagrangian potential given by the convex hull of , that is, . Then, particles no longer cross each other but form shocks. Moreover, this second dynamics still coincides with the Zeldovich dynamics outside of shocks, where coincides with its convex hull.

In this article, we do not compute the density power spectrum associated with this “adhesion model”, or another explicit dynamical system, which is a difficult task. Since we are merely interested in the amplitude of the effects associated with shell crossing, we simply make use of the property that shell crossing is associated with the loss of convexity of the Lagrangian potential . Then, we note that as long as is strictly convex its restriction along a line that goes through two arbitrary points and is also strictly convex, which implies (see also Noullez & Vergassola (1994); Brenier et al. (2003)),

 before shell crossing, for any qA≠qB: (33) [x(qB)−x(qA)]⋅(qB−qA)>0.

Indeed, we can choose a coordinate system so that and , and from strict convexity we obtain , whence along the first axis and , which reads as (33) in its general form3. This means that the projection of the Eulerian separation vector, , onto the Lagrangian vector, , is positive, so that lies in the forward half-space delimited by the plane orthogonal to .

Going back to the general expression (7), the density power spectrum is fully determined by the mean , for any , where we note . To compute this quantity we can take along the first axis, that is, with . Then, before shell crossing we have from the constraint (33) the property . Therefore, we consider the “sticky model” defined by:

 sticky model'', forq=|q|e1:Δx1=max(ΔxL1,0), (34) Δx2=ΔxL2,Δx3=ΔxL3,

where is the linear Eulerian separation, as given by Eq.(1). Thus, this second model only differs from the Zeldovich dynamics when the parallel linear Eulerian separation, , is negative, in which case we set it equal to zero. This is thus a simplified version of the “adhesion model”, as once reaches zero, it remains equal to zero forever. However, the model (34) cannot be explicitly derived from the “adhesion model”, since we take neither transverse directions nor larger scales into account. Therefore, we use the more generic name “sticky model”, to refer to this sticking along one direction for the pair separation.

It is clear that the condition , where the two models differ, is only a sufficient condition for shell crossing, but it is not necessary. Thus, it is a local condition that does not take the “cloud-in-cloud” problem into account : even though no shell crossing seems to have appeared on scale yet, it may happen that this region is enclosed within a larger domain of size that has already collapsed, so that particles in the smaller domain have actually experienced shell crossing (Bond et al. 1991). In terms of the Lagrangian potential , which defines the Lagrangian mapping, , through Eq.(29), the absence of shell crossing on a small domain of size means that is equal to its convex hull in this domain (Vergassola et al. 1994; Bec & Khanin 2007; Bernardeau & Valageas 2010b). However, the construction of the convex hull is a global problem, as one must consider the behavior of over all the space, thereby taking into account the “cloud-in-cloud” problem, while in the definition of the model (34) we only check a weaker condition, since we only consider the two points and . This means that we somewhat underestimate the effects of shell crossing, but we can expect to obtain a reasonable estimate of their amplitude because the probability of collapse decreases on larger scales and we perform a statistical integration over the angles of in Eq.(7).

From the previous discussions, the “sticky model” (34) and the Zeldovich dynamics (1) coincide before shell crossing, since then we have . This implies that both theories coincide at all orders of the perturbation theory; that is, they show the same expansions (23) and (26) over powers of . However, they differ through nonperturbative terms, which arise from their different behaviors after shell crossing.

From Eq.(34) we obtain , for , and

 ⟨eik⋅Δx⟩qe1=eik⋅qe−∫dw[1−cos(w⋅q)](k⋅w)2w4PL(w) (35) ×e12k21σ2∥∫∞−∞dΔΨL1√2πσ∥e−(ΔΨL1)2/(2σ2∥)eik1ΔΨ1,

where we factorized the result associated with the usual Zeldovich dynamics (8) in the first two terms and we introduced the variance of the linear longitudinal displacement,

 σ2∥(q) = ⟨(ΔΨL1)2⟩=2∫dw[1−cos(w1q)]w21w4PL(w) (36) = 2σ2v−2I0(q)+4I2(q). (37)

Separating the contribution from , we obtain for the last two terms of Eq.(35),

 e12k21σ2∥⟨eik1ΔΨ1⟩ = 1+e12k21σ2∥∫−q−∞dΔΨL1√2πσ∥e−(ΔΨL1)2/(2σ2∥) (39) ×(e−ik1q−eik1ΔΨL1) =1+12e12k21σ2∥−ik1qerfc(q√2σ∥)−12erfc⎛⎝q+ik1σ2∥√2σ∥⎞⎠,

where is the complementary error function (extended to the complex plane),

 erfc(z)=2√π∫∞zdte−t2. (40)

Then, substituting into Eqs.(35) and (7), we can see that the density power spectrum of the “sticky model”, , is equal to the usual Zeldovich power spectrum obtained in Sect. 2.1 plus a correction term ,

 Psticky(k)=PZel(k)+Ps.c.(k), (41)

with

 Ps.c.(k) = 12∫dq(2π)3e−k2(1−μ2)[σ2v−I0(q)−I2(q)]e−q2/(2σ2∥) (42) ×⎧⎨⎩w(iq√2σ∥)−w⎛⎝iq−kμσ2∥√2σ∥⎞⎠⎫⎬⎭.

 w(z)=e−z2erfc(−iz), (43)

which satisfies the asymptotic expansion (Abramowitz & Stegun 1970)

 |arg(z)|<3π4,z→∞:w(iz)∼1√πz(1−12z2+..). (44)

 Ps.c.(k) = ∫∞0dqq2(2π)2e−q2/(2σ2∥)∫10dμe−k2(1−μ2)[σ2v−I0−I2] (45) ×Re⎧⎨⎩w(iq√2σ∥)−w⎛⎝iq−kμσ2∥√2σ∥⎞⎠⎫⎬⎭.

As expected, we can check on Eqs.(42) and (45), using the behavior (44), that the correction due to shell-crossing effects is nonperturbative, in the sense that, because of the term , its expansion over powers of the amplitude of the linear power spectrum is identically zero4.

At low the shell-crossing contribution (45) behaves as . This agrees with the generic behavior5 associated with small-scale redistributions of matter (Peebles 1974). Thus, the “sticky model” provides an explicit example of a nonperturbative power spectrum, which satisfies these generic behaviors.

3 Numerical results in real space

We now describe the numerical results we obtain for the Zeldovich and “sticky model” power spectra, Eqs.(13) and (41), as well as their common perturbative expansions. We consider a CDM cosmology, with , , , and .

3.1 Logarithmic power

We show in Fig. 1 the power per logarithmic interval of , defined as

 Δ2(k)=4πk3P(k), (46)

for redshifts and . We compare the Zeldovich nonlinear power spectrum (13) with its perturbative expansions (19) and (24) and the nonperturbative correction (45). As is well known, higher order terms of the standard perturbative expansion (19) grow increasingly fast at high with changes in sign and large cancellations between different orders. The “renormalized” perturbative expansion (24) gives positive terms (see Eq.(25)) that peak on a well-defined range (Crocce & Scoccimarro 2006b) and are much easier to distinguish. Thus, while we only plot the first three orders of the standard expansion, and the absolute values and , we plot the first five orders of the “renormalized” expansion, to , as well as , , , , and (always multiplied by the factor of Eq.(46)). As we go to higher orders, contributions become narrower and more densely packed, which implies that as we go deeper in the nonlinear regime, we need increasingly more perturbative terms per logarithmic interval of wavenumber.

We also plot a few partial sums of expansion (24), that is, , with , and . We can check that they agree with the full nonlinear power spectrum (13) until the wavenumber associated with the peak of , after which they follow the Gaussian decay associated with the prefactor . These partial sums are also slightly more efficient than those obtained from the standard expansion (19) at the same order (not shown in the figure).

As expected, the nonperturbative correction (45) is very small on quasi-linear scales, so that there is indeed a range where higher order terms of the perturbative expansions (19) and (24) (i.e. beyond the first order associated with the linear regime) are relevant. At higher , the nonperturbative correction becomes dominant and the perturbative expansions become irrelevant, since one can no longer neglect the physics beyond shell crossing.

In the highly nonlinear regime, we recover the well-known decay of the Zeldovich logarithmic nonlinear power spectrum (Schneider & Bartelmann 1995; Taylor & Hamilton 1996; Valageas 2007b). This is because particles escape to infinity after shell crossing (i.e. keep moving on their straight trajectories), and these random trajectories (in the sense of random initial conditions) erase small-scale features in the density field. This is expressed by the Gaussian decaying factor in Eq.(13), where the quantity within brackets is always positive as can be seen from expression (9). Even though this leads to a decay at high , the falloff is not Gaussian because of the integration over the Lagrangian distance (Taylor & Hamilton 1996). In particular, for power-law initial power spectra, with , one finds that (Valageas 2007b) (the Zeldovich dynamics is not well-defined for because of ultraviolet divergences). Contrary to the gravitational case, the nonlinear power spectrum decreases faster at high for higher values of . This is due to the greater smearing out of small-scale features by the larger amplitude of the random linear displacements at small wavelengths.

It is interesting to note that, thanks to its shell-crossing correction (45), the nonlinear power spectrum (41) of the “sticky model” does not show this fast decay, and its logarithmic power still increases in the nonlinear regime. This is due to the prescription (34), which in a sense prevents particles from escaping to infinity in one direction, as they stick together. This approximately models the formation of Zeldovich pancakes, associated with collapse along one axis (although there is no precise relationship, since Eq.(34) is only a statistical model and does not consider the “cloud-in-cloud” problem). Thus, small-scale structures are no longer erased, since we keep a trace of planar features. This is expressed by the factor in the exponential argument in expression (45), which suppresses the Gaussian decaying term of the form for . This gives a width , hence at high , as would be the case for a density field where planar objects are the relevant nonlinear structures (i.e. bi-dimensional structures, as opposed to pointlike masses or lines for instance). Contrary to the Zeldovich power spectrum there is no dependence on for the high- slope, for power-law linear power spectra with . This gives the universal asymptote for the nonlinear logarithmic power of the “sticky model” at high .

3.2 Perturbative and nonperturbative contributions

To clearly see the range of scales where perturbative schemes are relevant, we plot in Fig. 2 the ratios of the successive terms of the “renormalized” perturbative expansion (24) with respect to the nonlinear power spectrum of the “sticky model”, , for to , as well as , and . We also plot the ratio , to compare with the amplitude of the nonperturbative correction associated with shell-crossing effects. We consider the four redshifts , and . Here we focus on the “renormalized” perturbative expansion (24), rather than on the standard expansion (19), to avoid the interferences brought by the changes in sign and the cancellations between various terms.

We can see that the potential of perturbative expansions grows at higher redshift for a CDM power spectrum, because of the change in slope of . Thus, at we can go up to the order before the perturbative term becomes subdominant with respect to the nonperturbative correction, while at the crossover takes place at , as the term is the lowest order one that is fully below the nonperturbative correction . This agrees with the observation that in the gravitational case perturbative schemes (and resummation approaches) seem to fare better at higher (Carlson et al. 2009). We can see that perturbative schemes are relevant over roughly one decade over . Thus, if we require an accuracy of at , linear theory is sufficient up to Mpc, while higher order perturbative terms allow reaching Mpc. At higher one must take the nonperturbative correction associated with shell-crossing effects into account, which implies going beyond the fluid approximation and requires new approaches.

To help the reader, we give in Table 1 the wavenumber , below which the linear term is enough to reach a or accuracy, for the four redshifts shown in Fig. 2. We also give the wavenumber beyond which the nonperturbative correction is required to reach an accuracy of , or (in units of ). Thus, the interval gives the range where perturbation theories based on the fluid description are relevant. Of course, this range shifts to higher at higher redshift. It is interesting to note that this range is also broader at higher redshift, as the slope of the CDM linear power spectrum on the relevant scales changes slowly. The last column gives the last order, , of the “renormalized” perturbative expansion that is not fully below the nonperturbative term. As noticed above in Fig. 2, this expansion order is significantly higher at higher redshift. This corresponds to a greater potential for perturbative schemes. However, grows faster than the logarithmic width of the perturbative range, . Indeed, as shown in Fig. 2, peaks associated with higher order perturbative terms are increasingly narrow on the axis. This implies that to multiply the upper wavenumber , defined by a fixed accuracy, by a given amount, one needs to add an increasingly greater number of new perturbative terms.

In practice, we do not expect that perturbative terms will be computed up to such high orders, since in the case of the gravitational dynamics this would involve multidimensional integrals that are beyond the reach of current numerical possibilities. However, resummation schemes allow one to consider parts of such higher order terms (actually, an infinite number of diagrams that contribute to terms of all orders). Then, the hope is that such methods can efficiently sum most of the contributions of higher order terms and accelerate the convergence on weakly nonlinear scales. The comparison displayed in Fig. 2 shows the potential of such methods (i.e. the best result one can obtain, for the “sticky model” considered here), which appears to be quite significant at higher redshifts, .

Since the Zeldovich dynamics provides a reasonable approximation of the gravitational dynamics down to weakly nonlinear scales (Coles et al. 1993; Pauls & Melott 1995), and its accuracy can be improved by implementing the “adhesion model” that only differs after shell crossing (Weinberg & Gunn 1990; Melott et al. 1994; Sathyaprakash et al. 1995), we can expect that to a large extent these results still apply to the gravitational case. In particular, while we find that at the nonperturbative correction to the density power spectrum is around at Mpc and at Mpc, Afshordi (2007) finds the wavenumbers Mpc and Mpc with a phenomenological “sticky halo model”. It is comforting that these two very different approaches give similar results. Then, the property that the range where perturbative schemes are relevant is greater at than at , with a higher order , should remain valid.

This is confirmed by the analysis of unified models that combine perturbation theories with halo models, for the gravitational case. Thus, as shown in Valageas & Nishimichi (2010), at higher redshift the intermediate range, where the power spectrum departs from one-loop perturbation theory but is not yet well described by the “one-halo” contribution, becomes wider. There, the “one-halo” contribution plays the role of the nonperturbative contribution (45), as it is also fully nonperturbative and decays at low as (if we only consider matter conservation, see the discussion in Valageas & Nishimichi (2010)). Then, the broadening of this intermediate range also suggests that higher orders of perturbation theory become relevant, as explicitly found in Figs. 1 and 2.

In the context of cosmological reconstruction6, a comparison with -body simulations (Mohayaee et al. 2006) shows that the Monge-Ampère-Kantorovich method is able to recover the nonlinear displacement field down to Mpc at , which corresponds roughly to Mpc. We can see from Table 1 that this is a very good result, as going to smaller scales requires taking shell crossing into account (in fact, at Mpc nonperturbative corrections have already started to dominate). Thus, the Monge-Ampère-Kantorovich method appears to be close to optimal at , because it goes as far as any scheme that disregards shell-crossing effects can be expected to go. This can be understood partly from the fact that relatively few orders of perturbations theory are relevant at (since ), so that it may not help much to explicitly include the effects of higher order terms.

Reconstruction techniques are also used to sharpen the acoustic peak of the real-space correlation function or to restore the harmonics of the oscillations of the power spectrum, in order to improve cosmological distance measurements and constraints on dark energy (Eisenstein et al. 2007; Seo et al. 2010). Then, one can read in Table 1 the wavenumber up to which one can hope to recover these baryon acoustic oscillations. At we can see that present schemes, which are based on the linear displacement field and manage to reach Mpc (not necessarily for the amplitude but at least for the shape and position of the oscillations) are not far from the upper bound, as could be expected from only a few orders of perturbation theory being relevant (). At it seems that one could go much beyond present schemes (which do not go much farther than Mpc), in agreement with , which means that higher orders of perturbation theory are relevant. However, for the specific purpose of measuring cosmological distances from the baryon acoustic oscillations, the potential is limited by the relative amplitude of the oscillations of the linear power spectrum itself decreasing at higher , so that even a very good reconstruction would not greatly enhance the signal-to-noise ratio. Nevertheless, pushing to higher orders (e.g., through resummation schemes) remains useful for other purposes, such as weak-lensing studies.

3.3 Rise of power at the transition

The increase of power on the transition scale to nonlinearity shown in Fig. 1 (especially in the left panel at ) is reminiscent of a similar feature observed for the gravitational dynamics (Hamilton et al. 1991; Peacock & Dodds 1996). This is usually interpreted from a Lagrangian point of view inspired by the spherical collapse dynamics. Thus, Padmanabhan (1996) argues that, on these intermediate scales, one has in real space with . More generally, in Fourier space one writes for the nonlinear power per logarithmic interval of wavenumber, , the parametric system (Peacock & Dodds 1996)

 kL = [1+Δ2(k)]−1/3k, (47) Δ2(k) = f[Δ2L(kL)], (48)

with a function to be determined. These relations express the conservation of matter, since the Lagrangian scale collapses down to the Eulerian scale . The linear regime implies that for , whereas in the highly nonlinear regime the stable-clustering ansatz (Peebles 1982) gives the scaling for . At the transition one observes a sharper growth, which is consistent with (Padmanabhan 1996). In practice, one builds a fitting formula for to match numerical simulations and to account for the dependence on the shape of the linear power spectrum.

In any case, such models usually estimate the shape of the nonlinear two-point correlation function or of the nonlinear power spectrum by considering the collapse of a “typical” overdensity (Padmanabhan 1996; Valageas & Schaeffer 1997) (or merely obtaining from simulations without further interpretation). It is interesting to note that this collapse also takes place within the Zeldovich dynamics studied here. In particular, the function that describes the spherical collapse is no longer given by cycloids (Peebles 1980) but by the simple expression (Bernardeau & Kofman 1995; Valageas 2009). (Collapse to a point is delayed from to , since the motion does not accelerate as the gravitational potential well becomes deeper.)

We show in Fig. 3 the functions defined by the system (47)-(48) for the Zeldovich dynamics and the “sticky model”. Indeed, from the knowledge of the nonlinear power we obtain the Lagrangian wavenumber from Eq.(47) and next . For a given Eulerian wavenumber , the Lagrangian wavenumbers obtained for both models, Zeldovich dynamics and “sticky model”, are different. Then, Fig. 3 shows that, within the Zeldovich dynamics, this growth of the nonlinear density contrast through the spherical collapse is not sufficient to build up the increase in power on these mildly nonlinear scales, although one can see a modest rise at . Indeed, the shape of the nonlinear power spectrum is quickly governed by the decay that takes place at higher , due to the escape of particles beyond shell crossing that erases small-scale features, as recalled in Sect. 3.1. In contrast, such a growth is clearly seen in the “sticky model”.

This means that such arguments, based on the evolution of a “typical” overdensity through the spherical collapse dynamics, are not good enough to explain the shape of the power spectrum on these scales since they clearly fail for the Zeldovich dynamics. Thus, one needs to consider the behavior of a whole range of typical density fluctuations (i.e., truly perform the Gaussian average over the initial conditions), which includes a significant number of configurations where shell crossing has already taken place. Then, one cannot neglect the dependence on their behavior after shell crossing, and whereas the free escape associated with the Zeldovich dynamics is sufficient to erase any growth of power on the transition scale, a trapping of particles allows such a growth. In the “sticky model” studied here, this trapping after shell crossing is a simple sticking along one direction, as defined in (34), which eventually leads to a universal tail at very high for , as noticed in Sect. 3.1. In terms of the function , this corresponds to a large- scaling from the system (47)-(48).

In the gravitational case, this regime corresponds to the virialization within 3D bound structures, with an asymptotic tail at high that is poorly known. In particular, although overdensities seem to form halos with an almost universal profile (Navarro et al. 1997), there is no first-principle derivation of the high- exponent of the density power spectrum, and its degree of universality is not well known. However, especially in the case of CDM-like power spectra that decay as at high , halos are rather smooth with a low amount of substructures. Then, the nonlinear logarithmic power spectrum does not seem to show a strong high- tail, such as the universal forms , or , associated with pointlike masses, lines, or sheets, but appears to be consistent with a logarithmic or shallow slope asymptote (Smith et al. 2003). This means that the steep rise in the function