Systematic treatment of non-linear effects in Baryon Acoustic Oscillations

In this contribution we will discuss the non-linear effects in the baryon acoustic oscillations and present a systematic and controllable way to account for them within time-sliced perturbation theory.


Introduction
Baryon acoustic oscillations (BAO) are widely used to establish the distance-redshift relation enabling to constrain the expansion history and the composition of the Universe [1,2]. BAO measurements aim at providing (sub-) percent precision in the near future and thus it is imperative to provide the best possible theoretical control over the BAO in order to fully exploit the potential of future surveys.
The BAO are most prominent in the 2-point correlation function of matter density in position space where they form a peak at r BAO ∼110 Mpc/h in comoving coordinates. It has been observed long ago that the shape of the BAO peak retrieved from N-body simulations differs significantly from the prediction of linear theory (see Fig. 1, left panel). Even though r BAO is significantly larger than the characteristic scale of non-linear clustering 2πk −1 NL ∼ 20 Mpc/h, the leading correction computed in Eulerian standard perturbation theory (SPT) failed to capture the behavior seen in N-body data (see Fig. 1, left panel). This disagreement is caused by large-scale bulk flows whose interaction with short modes is amplified if the distribution of matter has a feature. At leading order the effect of bulk motions is to worsen the correlation between galaxies at separations of order r BAO , which results in the suppression of the BAO peak.
Many approaches have been put forward in order to account for bulk flows. From the numerical side one can mention BAO reconstruction [3], which undoes bulk motions directly in the data and yields a sharper BAO peak with the signal-to-noise ratio improved by a factor of 2. From the theoretical side among the most successful approaches we would like to mention renormalized perturbation theory [4,5], Lagrangian perturbation theory [6] and IR -resummed effective field theory of large scale structures [7,8]. In this short contribution we will present a new way to systematically describe the non-linear evolution of BAO in the framework of time-sliced perturbation theory (TSPT) [9,10].  Figure 1. Left panel: Two-point matter correlation functions in position space at redshift zero, ξ(r) ≡ δ(r)δ(0) : N-body data from the Horizon Run 3 simulation [14] (red points), the prediction of linear theory (black dotdashed line), and the 1-loop SPT result (blue solid line). Right panel: Ratio of oscillatory (wiggly) part P w of the linear power spectrum to the smooth part P s . The ΛCDM cosmological parameters have been chosen as in [14].

Interaction with bulk flows: a simple physical picture
In this section we present a simple intuitive picture behind the non-linear evolution of the BAO (see Refs [5][6][7][8][9][10][11][12] for more details). We start with the density contrast field in Fourier space, δ(k, η) ≡ (ρ(η, k) −ρ(η))/ρ(η), where η is conformal time, ρ(η, k) is the local density andρ(η) is the average matter density in the Universe. δ(η, k) is a random stochastic field and one is typically interested in its 2-point correlation function (power spectrum) defined as P(η; k)δ (3) Dirac (k + k ) ≡ δ(η, k)δ(η, k ) . The linear power spectrum in our Universe can be decomposed into the smooth and wiggly parts (see Fig. 1, right panel), P lin (η; k) = P s (η; k) + P w (η; k) , and these two parts receive non-linear corrections from long-wavelength perturbations (IR modes) in very different ways. Let us first focus on the situation in which the linear power spectrum has only the smooth part, i.e. is featureless. Consider two galaxies on top of a large scale bulk flow. If the separation between these galaxies is shorter than the long wavelength associated to the flow, then at zeroth order its effect is a uniform acceleration, a = −∇Φ L (x) const (Φ L (x) is a gravitational potential induced by a longwavelength perturbation). Then, according to the equivalence principle, at a given time the effect of this acceleration can be totally removed by a proper coordinate transformation, so that it cannot affect equal-time physical observables such as the power spectrum. At next-to-leading order one can take into account inhomogeneity in the acceleration, in other words, the fact that the flow is not uniform but rather "converging" or "diverging". In that case two galaxies are likely to be found at separations shorter ("converging" flows) or larger ("diverging" flows) than in linear theory. However, since the smooth power spectrum in the real Universe is close to scale-invariant, the effect of the motions due to long modes should be small. In SPT the latter fact is manifest if one considers the 1-loop correction to the power spectrum in the IR limit (q k, where q is the loop momentum), where σ 2 S = q max k 0 d 3 qP lin (q) 1 is the large scale variance of the density field. From Eq.(2) it is clear that the interaction between the mode with momentum k and the IR modes is suppressed provided that the derivatives of the linear power spectrum are small, k n P (n) lin ∼ P lin . Thus, expanding over the coupling to the IR modes (∼ σ 2 S ) yields a converging perturbative expansion in this case. As the reader might have already guessed, this is not the case for a power spectrum with a feature. Now consider what happens in the presence of the wiggly component. In that case the effect of non-uniform accelerations is not small anymore. Indeed, the presence of a feature at r BAO ≡ k −1 osc means, at the linear level, that two galaxies have bigger probability to be found at this separation. At the non-linear level bulk flows make these galaxies move to shorter ("converging" flows) or larger ("diverging" flows) separations, which results in the degradation of the correlation between them at r BAO once these motions are averaged over long modes. Since matter is collapsing, "converging" flows are more common and two galaxies should be typically found at separations systematically smaller than r BAO . This effect is known as the physical shift of the BAO peak [12,13] and this is a next-to-leading order effect as compared to the broadening due to bulk flows. Assuming P w (k) = f s (k) sin(k/k osc ) ( f s (k) is a monotonic function) the SPT 1 loop result yields (here we assume q ≤ q max k osc ), From (3b) one clearly sees the mentioned effects on acoustic oscillations: the suppression of their amplitude and the shift in their phase w.r.t linear theory. The effective coupling constant to the IR modes ∼ σ 2 S (k 2 /k 2 osc ) is parametrically enhanced compared to the featureless case and in fact is of order unity for realistic cosmology. Thus, the perturbative expansion as (3b) does not make sense unless one can resumm relevant contributions to all orders in perturbation theory. We will show now how this procedure, called IR-resummation can be systematically performed in TSPT.

Overview of Time-Sliced Perturbation Theory
Time -sliced perturbation theory is a novel approach to large-scale structure (LSS) formation aimed at studying equal-time correlators in perturbation theory [9]. It is instructive to focus first on the SPT-like approach, in which one studies the following partition function 1 , where the ensemble average is taken over the initial field that has a Gaussian distribution. The timedependence of the field δ(η, k) is governed by the system of the continuity, Euler and Poisson equations for a perfect pressurless fluid (see, for example, Eq (2.1) of [9]). In SPT one finds the solution to these equations as a series over the initial field δ 0 , where g(η) is the time propagator (growth factor), F 2 is the SPT integral kernel. Then one can plug the expansion (5) into Eq.(4) and find equal-time correlation functions by taking functional derivatives w.r.t the source J. It is well known that the expansion (5) breaks the invariance under the transformations related to the equivalence principle, which eventually produces spurious IR enhancements in individual SPT diagrams [15,16]. The latter complicates the analysis of physical IR effects (bulk flows) within SPT.
Since only equal-time correlators typically have practical value, it seems natural to study the partition function for cosmological fields at particular times. To do this, one can transform the measure in (4) δ 0 → δ η (δ 0 ), which will define a time-dependent distribution function of δ η , The equation that governs the time evolution of P[δ η , η] is a classical Liouville continuity equation, One can write P[δ η , η] ≡ exp −Γ[δ η ] and assume, in the spirit of perturbation theory, Then, using the fluid equations for δ η , one can eliminate ∂δ η /∂η from Eq. (7) and arrive at the the chain of equations defining the time -evolution of the TSPT vertices Γ tot n . Remarkably, these equations can be solved exactly in the case of Einstein-de Sitter Universe without requiring any truncation of hierarchy. It is convenient to split Γ tot n into a regular part Γ n and a singular 'counterterm' C n , The terms C n are related to the singular Jacobian that appears as a result of the measure transformation in Eq. (6). The vertices Γ n represent physical non-linear n-point couplings and are non-local in space. Having imposed Gaussian initial conditions for these vertices one finds that the time evolution factorizes out, andΓ n are time-independent. Then one can find cosmological equal-time correlators by performing a perturbative expansion of (6) completely analogous to that of quantum field theory. In this expansion the growth factor g 2 (η) is treated as a coupling constant, the initial power spectrum P(k) =Γ −1 2 (k, −k) plays the role of the free propagator, and higher verticesΓ n (n > 2) represent connected parts of the amputated tree-level n-point correlation functions.
Remarkably, all TSPT vertices are finite in the soft limit, i.e. when some of their arguments are soft, lim →0Γn (η, k 1 , ..., k m , q 1 , ..., q n−m ) = O( 0 ). This property is related to the fact that TSPT operates only equal-time quantities, which are protected by the equivalence principle. Thus, each individual TSPT diagram is IR safe in contrast to SPT.
Given that the vertex¯ 0 2 generates all the higher vertices by the recursion relation (A.5b), one can introduce a similar decomposition for all the¯ 0 n vertices 0 n =¯ 0s n + ¯ 0w n + O( 2 ) .
(2.15) The¯ 0w n (¯ 0s n ) vertices are computed using¯ 0w 2 (¯ 0s 2 ) as an input in (A.5b). The decomposition (2.15) can be introduced back in the partition function (2.7). Since the C n and K n kernels are not functionals of the linear power spectrum [24] (see also Appendix A), they are not subject to the wiggly-smooth decomposition. The leading order O( 0 ) corresponds to the smooth correlation functions. The O( ) results include the wiggly contribution. In terms of diagrams, they can be summarized in the form of a new wiggly propagator g 2 P w (represented by a wiggly line) and by wiggly vertices 0w n (represented by a dashed circle), see Fig. 4. We use small dots to depict the smooth vertices and straight lines to represent the smooth power spectrum. The terms O( 2 ) are quadratic in P w and will be neglected.
The graphs with the wiggly elements are loosely referred to as wiggly graphs. For Hard loop

Soft loop
The loop integration in each diagram can be either hard, q are counted at order 2 h , or soft, q < kS, and are of order can be IR enhanced, so we focus on them for the momen of (3.15) are never IR-enhanced, i.e. they are at most o not contain a wiggly vertex. On the other hand, the di do receive an IR enhancement. The first diagram in th according to (3.12), (3.7) of order 1/" ⇥ 2 S . The last (3.12), (3.7) we find that it is of order 1/" 2 ⇥ 2 S and thus diagram. At leading order in " it is given bȳ [dq]Ps(q)DqD qPw(⌘ where in the last step we defined the operator S which ca In terms of our power counting As discussed previously, the product 1/" 2 ⇥ 2 S can be of O this one-loop contribution can be comparable to the linea In order to identify and eventually resum such terms the order of an arbitrary L-loop diagram in our power n-point correlation function, with external wavenumbers ki kS, kosc. Given a TSPT diagram with L loops (i.e.   soft can be divided into the hard and soft parts by introducing the hard and k L , such that k H k L . We will call the loop variable of the definition it is of the same order with the external momentum k, i.e. and |q 0 |, |q 0 k| k H . The loop variable of the soft domain, as in t will be called q.
In general, one expects the following scaling of diagrams with l s a hard loop, Let us identify the contributions of this type at two loop order, i.e. l The first set of enhanced diagrams are 1-particle reducible (1 obtained by dressing the external wiggly propagators of the first five a soft loop, The loop integration in each diagram can be either hard, q > k S , in whic are counted at order 2 h , or soft, q < k S , and are of order 2 S . Only the can be IR enhanced, so we focus on them for the moment. The diagra of (3.15) are never IR-enhanced, i.e. they are at most of order " 0 ⇥ not contain a wiggly vertex. On the other hand, the diagrams with do receive an IR enhancement. The first diagram in the second line according to (3.12), (3.7) of order 1/" ⇥ 2 S . The last diagram conta (3.12), (3.7) we find that it is of order 1/" 2 ⇥ 2 S and thus is the most IR diagram. At leading order in " it is given bȳ where in the last step we defined the operator S which can be written In terms of our power counting As discussed previously, the product 1/" 2 ⇥ 2 S can be of O(1) at low red this one-loop contribution can be comparable to the linear wiggly spec In order to identify and eventually resum such terms we now discus the order of an arbitrary L-loop diagram in our power counting. Th n-point correlation function, with external wavenumbers k i around the k i k S , k osc . Given a TSPT diagram with L loops (i.e. scaling as g 2L 1. choose for each propagator and each vertex whether it is smooth are interested in diagrams that contain one power of P w , at most propagator or vertex) can be wiggly. To obtain the full answer, o to sum over all possibilities to choose a vertex or propagator to b 2. assign each loop to be either hard (q > k S ) or soft (q < k S ). be done by splitting the linear input spectrum into two parts k S )P lin (q) + ✓(k S q)P lin (q), and calling a loop hard if all propa along the loop contain only power spectra of the former type 4 . T 4 Strictly speaking, the splitting into hard and soft contributions is only necess contain a wiggly vertex. For other loops, since 2 S + 2 h ⇡ 2 h in a realistic case, it whether we make this split or not.

IR -resummation in TSPT
The advantageous IR properties of TSPT can be readily used to account for bulk flows. As we already pointed out, the linear power spectrum in our Universe can be decomposed into the smooth and wiggly parts. Similar decomposition can also be done for the TSPT vertices, Thus, the TSPT loop expansion is split into the smooth and wiggly elements separately (see left panel of Fig. 2), which allows to isolate and resumm the relevant for the BAO contributions. One can find that the wiggly vertices are enhanced in the soft limit, e.g.
while the smooth vertices do not receive any enhancement in agreement with physical expectations. In order to identify and resumm the enhanced contributions at the diagrammatic level one can establish power-counting rules that determine the order of an arbitrary L -loop diagram. For that one has to set up a separation scale k S which splits the loop integrals into the soft and hard parts. This scale is a priory arbitrary and any residual dependence on it should be considered a theoretical error. Then, for each diagram containing a wiggly vertex (other diagrams are not IR-enhanced) one has to: (1) Assign each loop to be hard (q > k S ) or soft (q < k S ). Its value then can be associated with the variance of the density field taken over the corresponding domain, (2) Count the number of soft legs which are attached to a wiggly vertex. According to Eq.(12) this number is equal to the level of enhancement due to IR modes. Examples of our power counting rules are shown in the right panel of Fig. 2. It is easy so see that the leading contributions (∼ (σ 2 S /ε 2 ) L , L-loop order) to the wiggly power spectrum are given by the so-called 'daisy' graphs, whose resummation has a simple diagrammatic representation shown in Fig.3. One can show that the daisy diagrams exponentiate, However, as emphasized in [8], the Taylor expansion of P w (|k + q|) becomes unreliable for k osc . q ⌧ k. This means that non-linear corrections to the correlation functions at scale k receive large corrections from IR modes q within this range. In this work we identify these contributions for all n vertices, and establish a power counting scheme to compute corrections to the most enhanced terms. The leading contributions to the oscillatory part of the power spectrum are given by a set of 'daisy' diagrams, and their resummation is represented diagrammatically in the following form (see Sec. 4 for details), and j 0 (x), j 2 (x) are spherical Bessel functions. Remarkably, at leading order the result is the same for exact dynamics and the Zel'dovich approximation.

Comparison with N-body data
In this section we compare our theoretical predictions with the N-body simulations Horizon Run 2 and 3 [14]. For that we evaluate the matter correlation function in position space, using the LO and NLO IR -resummed power spectra Eq.(13) and Eq. (14) as an input. The result is shown in Fig. 4. On the left panel we display the LO IR -resummed result. The dependence on the separation scale, which we assume to be a theoretical error, is quite strong in this case. In the right panel of Fig. 4 we show the NLO IR -resummed correlation function. The agreement with the data is significantly improved as compared to LO. Also the dependence on the separation scale is reduced. This is an important observation because any dependence on k S must vanish, in principle, in the exact result. Thus, the reduction of the theoretical error when going to next-to-leading order ensures the consistency of our resummation scheme. It is useful to compare the TSPT NLO predictions with the  Zel'dovich approximation (ZA), which is known to be a quite accurate model for the BAO feature. In Fig. 5 we show the NLO correlation function divided by ZA, along with the N-body data Horizon Run 2 and 3. The difference between our results and ZA is of order 5% at the BAO peak, and it seems that the NLO result performs slightly better there. We notice, however, that the difference is comparable with the uncertainties in the N-body data in this region 2 , thus an apparent improvement over ZA seen in the Horizon Run data should be taken with a grain of salt. Nevertheless, the difference between TSPT NLO and ZA is larger than the forecasted precision of future surveys. The TSPT framework can also be used to systematically incorporate the NNLO corrections and thus provides a reliable tool for going beyond ZA.
Another non-linear effect on the BAO is the shift of the BAO peak [12,13]. This effect is quite small ∼ −0.3% for the ΛCDM, but can be measured with future observations [17]. The shift of the BAO peak may be big in modified gravity models [18], which makes it an interesting observable for constraining the deviations from ΛCDM. The shift is induced by terms which are off-phase with the linear BAO oscillations (see (3b)). Those correspond to NLO IR contributions in our power counting (scale as ∼ (σ 2 S /ε 2 ) L ε), and, in principle, can be explicitly extracted from Eq. (14) to yield δr BAO /r BAO ≈ −(0.3 ÷ 0.5)% which is consistent with values measured in simulations, e.g. [19].

Conclusions
In these notes we discussed the non-linear effects in the BAO and sketched the way how one can systematically take them into account within time-sliced perturbation theory. We outlined the physical picture behind the interactions with large scale bulk flows and argued the need for IR resummation if one works within the Eulerian framework. Then we made a short introduction into the TSPT formalism and discussed its key virtues relevant for IR resummation: a clear way to separate the perturbative expansion into the smooth and wiggly components, and manifest IR safety of the TSPT loop integrands. We introduced the power counting rules which were used to identify and resumm relevant sets of diagrams at leading and next-to-leading orders. Finally, we compared our results with N-body data and found good agreement within data errors.
We point out that the TSPT framework can be easily extended to incorporate higher-order corrections due to non-linear clustering, as well as new physics, e.g. the effects of neutrino masses or primordial non-gaussianity.