Average and dispersion of the luminosity-redshift relation in the concordance model
I. Ben-Dayan, M. Gasperini, G. Marozzi, F. Nugier, G. Veneziano
BBA-TH/666-12, CERN-PH-TH/2012-362, LPTENS-13/01, DESY 13-011
Average and dispersion of the luminosity-redshift relationin the concordance model
I. Ben-Dayan , M. Gasperini , , G. Marozzi , , F. Nugier and G. Veneziano , , Deutches Elektronen-Synchrotron DESY, Theory Group, D-22603 Hamburg, Germany Dipartimento di Fisica, Universit`a di Bari, Via G. Amendola 173, 70126 Bari, Italy Istituto Nazionale di Fisica Nucleare, Sezione di Bari, Bari, Italy Coll`ege de France, 11 Place M. Berthelot, 75005 Paris, France Universit´e de Gen`eve, D´epartement de Physique Th´eorique and CAP,24 quai Ernest-Ansermet, CH-1211 Gen`eve 4, Switzerland Laboratoire de Physique Th´eorique de l’ ´Ecole Normale Sup´erieure,CNRS UMR 8549, 24 Rue Lhomond, 75005 Paris, France CERN, Theory Unit, Physics Department, CH-1211 Geneva 23, Switzerland Center for Cosmology and Particle Physics, Department of Physics, New York University4 Washington Place, New York, NY 10003, USA
Starting from the luminosity-redshift relation recently given up to second order in the Poissongauge, we calculate the effects of the realistic stochastic background of perturbations of the so-called concordance model on the combined light-cone and ensemble average of various functions ofthe luminosity distance, and on their variance, as functions of redshift. We apply a gauge-invariantlight-cone averaging prescription which is free from infrared and ultraviolet divergences, makingour results robust with respect to changes of the corresponding cutoffs. Our main conclusions, inpart already anticipated in a recent letter for the case of a perturbation spectrum computed inthe linear regime, are that such inhomogeneities not only cannot avoid the need for dark energy,but also cannot prevent, in principle, the determination of its parameters down to an accuracy oforder 10 − − − , depending on the averaged observable and on the regime considered for thepower spectrum. However, taking into account the appropriate corrections arising in the non-linearregime, we predict an irreducible scatter of the data approaching the 10% level which, for limitedstatistics, will necessarily limit the attainable precision. The predicted dispersion appears to bein good agreement with current observational estimates of the distance-modulus variance due toDoppler and lensing effects (at low and high redshifts, respectively), and represents a challenge forfuture precision measurements. PACS numbers: 98.80-k, 95.36.+x, 98.80.Es
I. INTRODUCTION
In a recent letter [1] we have presented the main ideas and most significant results of a preliminary study of theeffects of a stochastic background of inhomogeneities on the determination of the dark-energy parameters in thecontext of modern precision cosmology. The main conclusions of that analysis, based on the use of a perturbationspectrum valid in the linear regime, were as follows. On the one hand, such kind of perturbations cannot simulatea substantial fraction of dark energy: their contribution to the averaged flux-redshift relation is both too small (atlarge values of the redshift z ) and has the wrong z -dependence. On the other hand, stochastic fluctuations add anew and relatively important dispersion with respect to the prediction of the homogeneous and isotropic Friedmann-Lemaˆıtre-Robertson-Walker (FLRW) cosmology. This dispersion is independent of the experimental apparatus, ofthe observational procedure, of the intrinsic fluctuations in absolute luminosity, and may prevent a determination ofthe dark-energy parameter Ω Λ ( z ) down to the percent level – at least if we are using the luminosity-redshift relationalone. Another important conclusion was that (light-cone averages of) different functions of the same observableget biased in different ways, with the energy flux sticking out as the observable which gets minimally affected byinhomogeneities, irrespectively of the redshift binning utilized . We should recall here that other possible sources ofuncertainty, bias and scatter in the Hubble diagram have been studied in many previous papers (see e.g. [3]-[7]).The method we have followed, in order to arrive at the above-mentioned conclusions, consists of two differentsteps. In the first step we start from an exact expression for the luminosity-redshift relation in the special “geodesiclight-cone” gauge introduced in [8]. We then transform this expression, up to second order in perturbation theory, to Also redshift binning reduces biases and selects the flux as the preferred variable [2]. a r X i v : . [ a s t r o - ph . C O ] S e p another gauge in which perturbations are known up to that order, the so-called Poisson gauge (PG) (see e.g. [9]). Thesecond step consists of performing the relevant light-cone and ensemble averages, and in inserting a realistic powerspectrum of stochastic perturbations. The light-cone average procedure appropriate to this context was formulatedand discussed in [8], by extending to null hypersurfaces the gauge-invariant procedure for space-like domains previouslydefined in [10, 11] and also applied in [12].Details on the first stage of this two-step process have been presented in a recent paper [13], while in this work weprovide a detailed implementation of the second step. The computation method is basically the same as the one alreadyused in [14], but it will involve the full second-order results obtained in [13]. We will consider in detail both a ColdDark Matter (CDM) model and a ΛCDM one (the so-called concordance model). We will reproduce, in particular,the results reported in [1] based on the use of the power spectrum of [15], valid in the linear perturbative regime.However, we will also extend our treatment by adding the effects of baryons and by considering two parametrizationsof the HaloFit model [16, 17], describing the density power spectrum in the non-linear regime.The paper is organized as follows. In Sect. II we recall, for the sake of completeness, the results of [13] for thecontribution of scalar perturbations to the light-cone average of the flux-redshift relation, to second order, in thePoisson gauge. We also reorganize the many different contributions in a convenient form for the actual estimatesto be carried out. In Sect. III we present a few important aspects and consequences of the process of combininglight-cone and ensemble averages, considering in particular the luminosity distance and its phenomenologically mostrelevant functions. We also introduce a convenient spectral parametrization of the inhomogeneous averaged terms. InSect. IV we discuss some relevant results about the dynamical evolution of scalar perturbations, up to second order,required for the computation of their averaged contribution. In Sect. V, as a warm-up exercise, we apply our methodsto a simple perturbed CDM model, where most calculations – except for the explicit mode integration over the givenpower spectrum – can be done analytically. In the relevant range of z , the full leading result can be written in termsof an explicit (and simple) function of z times a particular moment of the spectrum, and this shows that such a modelbadly fails in explaining the data, both in magnitude (in particular, at large redshift) and in z -dependence.We then turn to the case of a perturbed ΛCDM model, where calculations are more involved. They are simplifiedby restricting our attention to the so-called “enhanced terms” (i.e the dominant ones in the relevant range of z ),already identified in the CDM case. In order to discuss a realistic perturbation background, we will also considera power spectrum which includes the contribution of baryonic matter and takes into account the “Silk-damping”effect. In Sect. VI we restrict our computations to the linear power-spectrum proposed in [15], and we evaluate theimpact of stochastic inhomogeneities not only on the averaged flux-redshift relation, but also on other functions ofthe luminosity distance (in particular, on the distance modulus used in the analyses of the Supernovae data). Wealso discuss the dispersion induced by the presence of the perturbations. In Sect. VII we take into account the effectsof the non-linear regime by using the HaloFit parametrizations of [16, 17]. We also compute the variance/dispersionexpected in the distance modulus and attempt a first comparison with available data and phenomenological fits, inparticular for what concerns lensing at large redshifts. Finally, in Sect. VIII we summarize our result and draw someconclusions. In addition, in Appendix A we discuss why vector and tensor perturbations, although interesting ontheir own, do not contribute to the averages discussed in this paper. In Appendix B we report the explicit results forvarious spectral coefficients used in the analysis of the CDM model. II. AVERAGING THE LUMINOSITY FLUX AT SECOND ORDER
In this section we recall previous results on the light-cone average of the luminosity flux, (cid:104) Φ (cid:105) ∼ (cid:104) d − L (cid:105) (where d L is the luminosity distance), first computed in a generally inhomogeneous metric background, and then specialized tothe case of a spatially flat FLRW metric perturbed, to second order, by the presence of small fluctuations of scalar,vector and tensor type. A. Exact expression for (cid:104) d − L (cid:105) in the geodesic light-cone gauge When describing the propagation of light emitted by sources lying on the past light-cone of a given observer, itis convenient to identify the null hypersurfaces along which the photons reach the observer with those on which anull coordinate takes constant values. For this reason we have introduced in [8] an adapted system of coordinates –defining what we have called the “geodesic light-cone” (GLC) gauge – in which several quantities greatly simplify, whilekeeping all the required degrees of freedom for applications to general geometries. The coordinates x µ = ( τ, w, (cid:101) θ a )(with a = 1 , (cid:101) θ = (cid:101) θ , (cid:101) θ = (cid:101) φ ), specifying the metric in the GLC gauge, correspond to a complete gauge fixing of theso-called observational coordinates, defined e.g. in [18, 19] (see also [20]).The GLC metric depends indeed on six arbitrary functions (a function Υ, a two-dimensional “vector” U a and asymmetric matrix γ ab ), and its line-element takes the form : ds GLC = Υ dw − dwdτ + γ ab ( d (cid:101) θ a − U a dw )( d (cid:101) θ b − U b dw ) . (2.1)In matrix form, the metric and its inverse are then given by: g GLCµν = − Υ (cid:126) − Υ Υ + U − U b (cid:126) T − U Ta γ ab , g µνGLC = − − Υ − − U b / Υ − Υ − (cid:126) − ( U a ) T / Υ (cid:126) T γ ab , (2.2)where (cid:126) , U b = ( U , U ), while the 2 × γ ab and γ ab lower and raise the two-dimensional indices.Clearly w is a null coordinate (i.e. ∂ µ w∂ µ w = 0), and a past light-cone hypersurface is specified by the condition w = const. We also note that u µ ∼ ∂ µ τ defines a geodesic flow, i.e. that ( ∂ ν τ ) ∇ ν ( ∂ µ τ ) = 0 (as a consequence of g ττ = − k µ = g µν ∂ ν w = g µw = − δ µτ Υ − , meaning that photons travel at constant w and (cid:101) θ a . Thismakes the calculation of the redshift and of the area distance particularly easy in this gauge.In fact, let us denote by the subscripts o and s , respectively, a quantity evaluated at the observer and sourcespace-time position, and consider a light ray emitted by a geodetic source (with four-velocity u µ = − ∂ µ τ ) lying atthe intersection between the past light-cone of a given geodetic observer (defined by the equation w = w o ) and thespatial hypersurface τ = τ s , where τ s for the moment is a constant parameter. The light ray will be received by ourstatic geodetic observer at τ = τ o > τ s . The redshift z s associated with this light ray is then given by [8]:(1 + z s ) = ( k µ u µ ) s ( k µ u µ ) o = ( ∂ µ w∂ µ τ ) s ( ∂ µ w∂ µ τ ) o = Υ( w o , τ o , (cid:101) θ a )Υ( w o , τ s , (cid:101) θ a ) . (2.3)We shall be interested in averaging the luminosity at fixed redshift, hence on the two-dimensional surface Σ( w o , z s )(topologically a sphere) which lies on our past light-cone ( w = w o ) and is associated with a fixed redshift ( z = z s ). Interms of the τ coordinate such a surface corresponds to the equation τ = τ s ( w o , z s , ˜ θ a ) enforcing Eq. (2.3). Hereafter τ s will denote this (in general angle-dependent) quantity.Also the area distance d A , related to the luminosity distance d L of a source at redshift z s by one of Etherington’srelations [21]: d A = (1 + z s ) − d L , (2.4)takes a particularly simple form in the GLC gauge. A direct derivation [8] starts from its general definition [22]: d A = dS s d Ω o , (2.5)where d Ω o is the infinitesimal solid angle subtended by the source at the observer position, and dS s is the area elementon the surface orthogonal to both the photon momentum and to the source 4-velocity at the source’s position. It iseasy to check that the surface dS s is characterized by having constant w and τ , and that the induced 2-metric on itis nothing but γ ab [13]. Therefore : d A = d ˜ θ √ γd ˜ θ sin ˜ θ = √ γ sin ˜ θ , (2.6)where we have used the fact that photons travel at constant ˜ θ a . Our averaging surface Σ( w o , z s ), being one of constant z , differs from the one of constant w and τ but – amusingly – the same formula holds, locally, for the area elementon it, so that Eq. (2.6) can also be used for our light-cone averages [13]. We have put tildas on the GLC gauge θ a coordinates in order to be consistent with our previous notations in [1, 13, 14]. This result has been checked meanwhile [23] by explicitly constructing the relevant Jacobi map in the GLC gauge and by using its knownrelation to the area distance [24]. This method confirms the results discussed below (and in our previous papers), modulo a Lorentztransformation of d Ω o connected to the peculiar velocity of the observer as measured in the longitudinal gauge. A similar correction isalso needed, of course, in order to take into account the peculiar motion of our galaxy hence, strictly speaking, our unintegrated resultsonly hold without these peculiar-velocity-related effects. In practice, the correction exactly vanishes for the averaged flux (Eq. (2.7)below) and is numerically negligible for the other averages discussed in our papers. The above result singles out the received luminosity flux, Φ ∼ d − L = (1 + z s ) − d − A , as an important – andextremely simple – observable to average over the 2-sphere Σ( w o , z s ) embedded in the light-cone. In fact (see [13] formore details): (cid:104) d − L (cid:105) ( w o , z s ) = (1 + z s ) − (cid:82) dS d Ω o dS (cid:82) dS = (1 + z s ) − (cid:82) d Ω o (cid:82) dS = (1 + z s ) − π A ( w o , z s ) , (2.7)where A ( w o , z s ) = (cid:90) Σ( w o ,z s ) d ˜ θ a √ γ (2.8)is the proper area of Σ( w o , z s ) computed with the metric γ ab , and expressed in terms of internal coordinates ( w o , z s )parametrizing the deformed 2-sphere Σ( w o , z s ).Eq. (2.7) holds non-perturbatively for any space-time geometry, and is the starting point for the computation ofthe average flux summarized in [1] and presented in details here. By using the notations introduced in [1] we canwrite, in particular, (cid:104) d − L (cid:105) ( w o , z s ) = (1 + z s ) − (cid:34)(cid:90) d ˜ θ a π √ γ ( w o , τ s ( w o , z s , ˜ θ a ) , ˜ θ a ) (cid:35) − ≡ (cid:0) d F LRWL (cid:1) − I − φ ( w o , z s ) , (2.9)where we have defined I φ ( w o , z s ) ≡ A ( w o , z s )4 π (cid:104) a ( η (0) s )∆ η (cid:105) , (2.10)and where d F LRWL = (1 + z s ) a ( η (0) s )∆ η is the luminosity distance for the unperturbed FLRW geometry, with scalefactor a ( η ). Here η is the conformal time coordinate, ∆ η = η o − η (0) s , and we have denoted with η (0) s the backgroundsolution of the equation for the source’s conformal time η s = η s ( z s , (cid:101) θ a ) (see [13, 14]). Note that, according to theabove equation, the interpretation of I φ ( w o , z s ) is straightforward: it is simply the ratio of the area of the 2-sphereat redshift z s on the past light-cone (deformed by inhomogeneities), over the area of the corresponding homogeneous2-sphere. B. Second-order expression for (cid:104) d − L (cid:105) in the Poisson gauge Let us consider a space-time geometry that can be approximated by a spatially flat FLRW metric distorted by thepresence of scalar, vector and tensor perturbations. In the so-called Poisson gauge (PG) [9] (a generalization of thestandard Newtonian gauge beyond first order), the corresponding metric (in Cartesian coordinates) takes the form: ds P G = a ( η ) (cid:8) − (1 + 2Φ) dη + 2 ω i dηdx i + [(1 − δ ij + h ij ] dx i dx j (cid:9) . (2.11)Here Φ and Ψ are scalar perturbations, ω i is a transverse vector perturbation ( ∂ i ω i = 0) and h ij is a transverse andtraceless tensor perturbation ( ∂ i h ij = 0 = h ii ). This metric depends on six arbitrary functions, hence it is completelygauge fixed. By including first-order and second-order contributions, the (generalized) Bardeen potentials Φ and Ψcan be defined as follows: Φ ≡ ψ + 12 φ (2) , Ψ ≡ ψ + 12 ψ (2) , (2.12)where we have assumed the absence of anisotropic stress in order to set Ψ = Φ = ψ at first order.It is important to stress, at this point, that for the purpose of this paper we can safely restrict our subsequentdiscussion to the case of pure scalar perturbations. In fact, it is true that at second order different perturbations getmixed: vector and tensor perturbations are automatically generated from scalar perturbations (see e.g. [25, 26]), whilesecond-order scalar perturbations are generated from first-order vector and tensor perturbations. However, a singlevector or tensor perturbation does not contribute to our angular averages on Σ( w o , z s ) (see Appendix A for furtherdetails). Furthermore, we will treat ω i and h ij as second order quantities. In other words, we assume the first-orderperturbed metric to be dominated by scalar contributions (which is indeed the case if perturbations are generated bya phase of standard slow-roll inflation, see e.g. [27, 28]). As a result, we shall also neglect the contributions induced,at second order, by first-order vector and tensor perturbations.We have already established in [13] a connection between the second-order perturbative expression of the luminositydistance d L and of the integrand of I φ (controlling (cid:104) d − L (cid:105) ), both written in terms of the PG perturbations and of theobserver’s angles ˜ θ a (remember that photons reach the observer traveling at constant ˜ θ a ). In particular, by defining d L ( z s , ˜ θ a )(1 + z s ) a o ∆ η = d L ( z s , ˜ θ a ) d F LRWL ( z s ) = 1 + δ (1) S ( z s , ˜ θ a ) + δ (2) S ( z s , ˜ θ a ) , (2.13)and I φ ( z s ) = (cid:90) d ˜ θ a π sin ˜ θ (1 + I + I , + I ) , (2.14)we have found that [13] I = 2 δ (1) S + (t . d . ) (1) , I , + I = 2 δ (2)S + ( δ (1)S ) + (t . d . ) (2) , (2.15)where the (t . d . ) (1 , denote total derivatives terms w.r.t. the ˜ θ a angles, giving vanishing contributions either byperiodicity in (cid:101) φ or by the vanishing of the integrand at (cid:101) θ = 0 , π . We only recall here, for later use, the first-ordertotal derivative: (t . d . ) (1) = 2 J (1)2 , J (1)2 = 1∆ η (cid:90) η o η (0) s dη η − η (0) s η o − η ∆ ψ ( η, η o − η, (cid:101) θ a ) , (2.16)where ∆ = ∂ (cid:101) θ + cot (cid:101) θ∂ (cid:101) θ + sin − (cid:101) θ∂ (cid:101) φ . It is also convenient to rewrite the PG metric using spherical coordinates (butstill considering that photons travel at constant ˜ θ a ), and define the following quantities [13]: P ( η, r, (cid:101) θ a ) = (cid:90) ηη in dη (cid:48) a ( η (cid:48) ) a ( η ) ψ ( η (cid:48) , r, (cid:101) θ a ) , Q ( η + , η − , (cid:101) θ a ) = (cid:90) η − η + dx ˆ ψ ( η + , x, (cid:101) θ a ) , Ξ s = 1 − H s ∆ η , J = ([ ∂ + Q ] s − [ ∂ + Q ] o ) − ([ ∂ r P ] s − [ ∂ r P ] o ) . (2.17)Here H = d ln a/dη , and the lower limit η in represents an early enough time when the perturbation (or better theintegrand) was negligible. When ambiguities may occur the superscript (0) denotes the background solution of agiven quantity (similarly, the superscripts (1), (2) will denote, respectively, the first- and second-order perturbedvalues of that quantity). In the above equations we have also introduced the useful (zeroth-order) light-cone variables η ± = η ± r , with corresponding partial derivatives: ∂ η = ∂ + + ∂ − , ∂ r = ∂ + − ∂ − , ∂ ± = ∂∂η ± = 12 ( ∂ η ± ∂ r ) . (2.18)We shall use hereafter a hat to denote a quantity expressed in terms of the ( η + , η − , (cid:101) θ a ) variables, so that, for instance,ˆ ψ ( η + , η − , (cid:101) θ a ) ≡ ψ ( η, r, (cid:101) θ a ). Finally, in order to understand the physical meaning of the above quantities, it may behelpful to recall that the radial gradient of P is related to the Doppler effect (due to peculiar velocities of sourceand observer), while the gradient of Q with respect to ∂ + represents the Sachs-Wolfe and the integrated Sachs-Wolfeeffect. The last term J corresponds to a combination of the three mentioned effects (see [13], Sect. 4, for a moredetailed discussion).The results obtained in [13] can then be reported in the following form: I = (cid:88) i =1 T (1) i ; I , = (cid:88) i =1 T (1 , i ; I = (cid:88) i =1 T (2) i ; (2.19) Note the simplified notation with respect to [13]. We have also omitted the indication that w = w o . As already mentioned in a previous footnote, the terms appearing on the r.h.s. of Eq.(2.13) should be corrected, in principle, for thechange in d Ω o stemming from the peculiar velocity of the observer. However, such a correction has no effect on the integral I φ itself,since it can be compensated by a Lorentz transformation of the angular variables. where I , I , , I are, respectively, the first-order, quadratic first-order, and genuine second-order contributions ofour stochastic fluctuations, and where: T (1)1 = − ψ ( η (0) s , r (0) s , (cid:101) θ a ) ; T (1)2 = 2Ξ s J ; T (1)3 = − η Q s ; (2.20) T (1 , = Ξ s (cid:2) ψ s − ψ o (cid:3) ; T (1 , = Ξ s (cid:0) ([ ∂ r P ] s ) − ([ ∂ r P ] o ) (cid:1) ; T (1 , = − s ( ψ s + [ ∂ + Q ] s ) [ ∂ r P ] s ; T (1 , = 12 Ξ s ( γ ab ) s (2 ∂ a P s ∂ b P s + ∂ a Q s ∂ b Q s − ∂ a Q s ∂ b P s ) ; T (1 , = − Ξ s lim r → (cid:2) γ ab ∂ a P ∂ b P (cid:3) ; T (1 , = 2Ξ s Q s (cid:18) ∂ r ψ o + 2 ∂ η ψ o − ∂ r ψ s + 2 (cid:90) η η (0) s dη (cid:48) ∂ r ψ (cid:16) η (cid:48) , η o − η (cid:48) , (cid:101) θ a (cid:17) + [ ∂ r P ] s (cid:19) ; T (1 , = 2Ξ s J H s ([ ∂ η ψ ] s + H s [ ∂ r P ] s ) ; T (1 , = 2Ξ s J H s [ ∂ r P ] s ; T (1 , = − Ξ s (cid:90) η (0) s η in dη (cid:48) a ( η (cid:48) ) a ( η (0) s ) ∂ r (cid:2) − ψ + ( ∂ r P ) + γ ab ∂ a P ∂ b P (cid:3) ( η (cid:48) , ∆ η, ˜ θ a ) ; T (1 , = Ξ s (cid:90) η o η in dη (cid:48) a ( η (cid:48) ) a ( η o ) ∂ r (cid:2) − ψ + ( ∂ r P ) + γ ab ∂ a P ∂ b P (cid:3) ( η (cid:48) , , ˜ θ a ) ; T (1 , = 2Ξ s (cid:90) η (0) − s η (0)+ s dx ∂ + (cid:20) ˆ ψ ∂ + Q + 14 ˆ γ ab ∂ a Q∂ b Q (cid:21) ( η (0)+ s , x, ˜ θ a ) ; T (1 , = (cid:20) Ξ s − H s ∆ η (cid:18) − H (cid:48) s H s (cid:19)(cid:21) J ; T (1 , = − ψ s J ; T (1 , = 2Ξ s (cid:26) ψ o + [ ∂ r P ] o − Q s ∆ η (cid:27) J ; T (1 , = − J − ψ s ) Q s ∆ η ; T (1 , = − ψ s − ∂ + Q s ) Q s ∆ η ; T (1 , = (cid:18) Q s ∆ η (cid:19) ; T (1 , = 1 H s ( γ ab ) s ∂ a Q s ∂ b J ; T (1 , = 12 ( γ ab ) s ∂ a Q s ∂ b Q s ; T (1 , = 2 J H s ( − [ ∂ η ψ ] s + [ ∂ r ψ ] s ) ; T (1 , = 2 Q s [ ∂ r ψ ] s ; T (1 , = − η (cid:90) η (0) − s η (0)+ s dx (cid:20) ˆ ψ ∂ + Q + 14 ˆ γ ab ∂ a Q ∂ b Q (cid:21) ( η (0)+ s , x, ˜ θ a ) ; T (1 , = 18 sin ˜ θ ∂∂ ˜ θ cos ˜ θ (cid:32)(cid:90) η (0) − s η (0)+ s dx [ˆ γ b ∂ b Q ]( η (0)+ s , x, ˜ θ a ) (cid:33) ; (2.21)and T (2)1 = −
12 Ξ s (cid:16) φ (2) s − φ (2) o (cid:17) ; T (2)2 = 12 Ξ s (cid:16) ψ (2) s − ψ (2) o (cid:17) ; T (2)3 = − Ξ s (cid:90) η (0) s η in dη (cid:48) a ( η (cid:48) ) a ( η (0) s ) [ ∂ r φ (2) ]( η (cid:48) , ∆ η, ˜ θ a ) ; T (2)4 = Ξ s (cid:90) η o η in dη (cid:48) a ( η (cid:48) ) a ( η o ) [ ∂ r φ (2) ]( η (cid:48) , , ˜ θ a ) ; T (2)5 = 12 Ξ s (cid:90) η (0) − s η (0)+ s dx ∂ + (cid:104) ˆ φ (2) + ˆ ψ (2) (cid:105) ( η (0)+ s , x, ˜ θ a ) ; T (2)6 = − ψ (2) s ; T (2)7 = − η (cid:90) η (0) − s η (0)+ s dx (cid:34) ˆ φ (2) + ˆ ψ (2) (cid:35) ( η (0)+ s , x, ˜ θ a ) . (2.22)In the above equations γ ab = diag( r − , r − sin − (cid:101) θ ), and for T (1 , we have used the following identity: − [ ∂ Q ] s + [ ∂ + ˆ ψ ] s = 2 ∂ r ψ o + 2 ∂ η ψ o − ∂ r ψ s + 2 (cid:90) η η (0) s dη (cid:48) ∂ r ψ (cid:16) η (cid:48) , η o − η (cid:48) , (cid:101) θ a (cid:17) . (2.23)Let us point out, finally, that the last term in Eq.(2.21) corresponds to a total derivative, and thus to a boundarycontribution, that superficially looks non vanishing. We believe that this is the result of a naive treatment of theangular coordinate transformation, which becomes singular near the poles of the 2-sphere. This contribution, indeed,has the same form as the irrelevant one coming from an overall SO (3) rotation. III. COMBINING THE LIGHT-CONE AND ENSEMBLE AVERAGE OF d L AND OF ITS FUNCTIONS
In the cosmological model we are considering, the deviations from the standard FLRW geometry are sourced bya stochastic background of primordial perturbations satisfying ψ = 0, ψ (cid:54) = 0, where the bar denotes statistical(or ensemble ) average. Hence, non-trivial effects on the ensemble average of d L , or of a generic function of it, canonly originate either from quadratic and higher-order perturbative corrections, or from the spectrum of correlationfunctions such as d L ( z, (cid:101) θ a ) d L ( z (cid:48) , (cid:101) θ (cid:48) a ) (see [29]). In this paper, rather than considering the ensemble average of d L , weshall deal with that of (cid:104) d L (cid:105) , where the angular brackets refer to the light-cone average defined in [8] and presentedin Sect. II (see [30–32] for previous attempts of combining ensemble and space-time averages in the case of spacelike hypersurfaces).As already stressed in [14], given the covariant (light-cone) average of a perturbed (inhomogeneous) observable S , the average of a generic function of this observable differs, in general, from the function of its average, i.e. (cid:104) F ( S ) (cid:105) (cid:54) = F ( (cid:104) S (cid:105) ) (as a consequence of the nonlinearity of the averaging process). Expanding the observable tosecond order as S = S + S + S , one finds [1] (cid:104) F ( S ) (cid:105) = F ( S ) + F (cid:48) ( S ) (cid:104) S + S (cid:105) + F (cid:48)(cid:48) ( S ) (cid:104) S / (cid:105) , (3.1)where in general (cid:104) S (cid:105) (cid:54) = 0 as a consequence of the so-called “induced backreaction” terms, arising from the couplingbetween the inhomogeneity fluctuations of S and those of the integration measure (see [14]). The overall correction to (cid:104) F ( S ) (cid:105) thus depends not only on the intrinsic inhomogeneity of the observable S , but also on the covariance propertiesof the adopted averaging procedure. Eq. (3.1) implies, in our case, that different functions of the luminosity distance(or of the flux) may be differently affected by the process of averaging out the inhomogeneities, and may requiredifferent “subtraction” procedures for an unbiased determination of the relevant observable quantities.Let us consider, in particular, the luminosity flux Φ ∼ d − L (not to be confused with the Bardeen potential!),computed in Sect. II. Performing the stochastic average of Eq. (2.9) (and using Eq. (2.14)) we obtain (cid:104) d − L (cid:105) ( z ) = ( d F LRWL ) − ( I Φ ( z )) − ≡ ( d F LRWL ) − [1 + f Φ ( z )] , (3.2)where: f Φ ( z ) ≡ (cid:104)I (cid:105) − (cid:104)I , + I (cid:105) . (3.3)We can now apply the general result (3.1) to the flux variable, by setting S = Φ and considering two importantfunctions of the flux: F (Φ) = Φ − / ∼ d L , and F (Φ) = − . Φ + const ∼ µ (the distance modulus). They will beconsidered in the following sections, together with the flux. For the luminosity distance we can introduce a fractionalcorrection f d , in analogy with Eq. (3.2), such that: (cid:104) d L (cid:105) ( z ) = d F LRWL [1 + f d ( z )] . (3.4)Then, by using the general expression (3.1), we find: f d = − f Φ + 38 (cid:104) (Φ / Φ ) (cid:105) , (3.5)where, in terms of the quantities defined in (2.13), we have (cid:104) (Φ / Φ ) (cid:105) = 4 (cid:104) ( δ (1) S ) (cid:105) , and where f Φ is defined by Eq.(3.3). For the distance modulus we obtain, instead, (cid:104) µ (cid:105) − µ F LRW = − . e ) (cid:104) f Φ − (cid:104) (Φ / Φ ) (cid:105) (cid:105) . (3.6)We can also consider, for any given averaged variable (cid:104) S (cid:105) , the associated dispersion σ S controlling how broad isthe distribution of a perturbed observable S around its mean value (cid:104) S (cid:105) . This dispersion is due to both the geometricfluctuations of the averaging surface and to the statistical ensemble fluctuations, and is defined, in general, by [14]: σ S ≡ (cid:115)(cid:28)(cid:16) S − (cid:104) S (cid:105) (cid:17) (cid:29) = (cid:114) (cid:104) S (cid:105) − (cid:16) (cid:104) S (cid:105) (cid:17) . (3.7)The dispersion associated with the flux is thus given by: σ Φ = (cid:114) (cid:104) (Φ / Φ ) (cid:105) − (cid:16) (cid:104) Φ / Φ (cid:105) (cid:17) = (cid:113) (cid:104) (Φ / Φ ) (cid:105) , (3.8)while for the distance modulus we find: σ µ = (cid:114) (cid:104) µ (cid:105) − (cid:16) (cid:104) µ (cid:105) (cid:17) = 2 . e ) (cid:113) (cid:104) (Φ / Φ ) (cid:105) . (3.9)The above results will be applied to the case of a realistic background of cosmological perturbations of inflationaryorigin in the following sections.Let us conclude this section by introducing a convenient spectral parametrization to be used for the various termscontributing to the fractional corrections of our observables, and to the corresponding dispersions. We start byrecalling that the simplest way to implement the ensemble average of a stochastic background of scalar perturbations ψ is to consider its Fourier decomposition in the form: ψ ( η, (cid:126)x ) = 1(2 π ) / (cid:90) d k e i(cid:126)k · (cid:126)x ψ k ( η ) E ( (cid:126)k ) , (3.10)where – assuming that the fluctuations are statistically homogeneous and isotropic – E is a unit random variablesatisfying E ∗ ( (cid:126)k ) = E ( − (cid:126)k ), as well as the ensemble -average conditions E ( (cid:126)k ) = 0 and E ( (cid:126)k ) E ( (cid:126)k ) = δ ( (cid:126)k + (cid:126)k ). As asimple illustrative example one has: (cid:104) ψ s ψ s (cid:105) = (cid:90) d k d k (cid:48) (2 π ) E ( (cid:126)k ) E ( (cid:126)k (cid:48) ) (cid:90) d Ω4 π (cid:104) ψ k ( η (0) s ) e ir(cid:126)k · ˆ x (cid:105) r = η − η (0) s (cid:104) ψ k (cid:48) ( η (0) s ) e ir (cid:126)k (cid:48) · ˆ x (cid:105) r = η − η (0) s = (cid:90) d k (2 π ) | ψ k ( η (0) s ) | (cid:90) − d (cos θ )2 (cid:2) e ik ∆ η cos θ (cid:3) (cid:2) e − ik ∆ η cos θ (cid:3) = (cid:90) d k (2 π ) | ψ k ( η (0) s ) | = (cid:90) ∞ dkk P ψ ( k, η (0) s ) , (3.11) (cid:104) ψ s (cid:105) (cid:104) ψ s (cid:105) = (cid:90) d k d k (cid:48) (2 π ) E ( (cid:126)k ) E ( (cid:126)k (cid:48) ) (cid:20)(cid:90) d Ω4 π ψ k ( η (0) s ) e ir(cid:126)k · ˆ x (cid:21) r = η − η (0) s (cid:20)(cid:90) d Ω (cid:48) π ψ k (cid:48) ( η (0) s ) e ir (cid:126)k (cid:48) · ˆ x (cid:48) (cid:21) r = η − η (0) s = (cid:90) d k (2 π ) | ψ k ( η (0) s ) | (cid:20)(cid:90) − d (cos θ )2 e ik ∆ η cos θ (cid:21) (cid:20)(cid:90) − d (cos θ (cid:48) )2 e − ik ∆ η cos θ (cid:48) (cid:21) = (cid:90) d k (2 π ) | ψ k ( η (0) s ) | (cid:18) sin( k ∆ η ) k ∆ η (cid:19) = (cid:90) ∞ dkk P ψ ( k, η (0) s ) (cid:18) sin( k ∆ η ) k ∆ η (cid:19) , (3.12)where in the second line of both terms we made use of isotropy (i.e. ψ k only dependent on k = | (cid:126)k | ), and defined θ and θ (cid:48) as the angles between (cid:126)k and (cid:126)x ≡ r ˆ x and between (cid:126)k (cid:48) and (cid:126)x (cid:48) ≡ r ˆ x (cid:48) . We recall that ∆ η = η o − η (0) s . We have alsointroduced the (so-called dimensionless) power spectrum of ψ : P ψ ( k, η ) ≡ k π | ψ k ( η ) | . (3.13)To give a slightly more involved example, we can consider the average of a term like ψ s Q s / ∆ η : (cid:28) ψ s Q s ∆ η (cid:29) = − η (cid:90) d kd k (cid:48) (2 π ) E ( (cid:126)k ) E ( (cid:126)k (cid:48) ) (cid:90) d Ω4 π (cid:20) ψ k ( η (0) s ) e i(cid:126)k · ˆ x ∆ η (cid:90) η o η (0) s dη ψ k (cid:48) ( η ) e i(cid:126)k (cid:48) · ˆ x ( η o − η ) (cid:21) ⇒ − (cid:90) ∞ dkk P ψ ( k, η o ) SinInt( k ∆ η ) k ∆ η , (3.14) (cid:104) ψ s (cid:105) (cid:28) Q s ∆ η (cid:29) = − η (cid:90) d kd k (cid:48) (2 π ) E ( (cid:126)k ) E ( (cid:126)k (cid:48) ) (cid:20)(cid:90) d Ω4 π ψ k ( η (0) s ) e i(cid:126)k · ˆ x ∆ η (cid:21) (cid:20)(cid:90) d Ω (cid:48) π (cid:90) η o η (0) s dη ψ k (cid:48) ( η ) e i(cid:126)k (cid:48) · ˆ x (cid:48) ( η o − η ) (cid:21) ⇒ − (cid:90) ∞ dkk P ψ ( k, η o ) sin( k ∆ η ) k ∆ η SinInt( k ∆ η ) k ∆ η , (3.15)where the arrows refer to the case of a time-independent fluctuation mode andSinInt( x ) ≡ (cid:90) x dyy sin y . (3.16)As one can see, the angular average is making the results completely different in the two cases. We remark that thepresence of the SinInt function is a direct consequence of the integration over time in Q s . Consequently, the non-localnature of the backreaction terms is reflected in the form of the corresponding spectral coefficients.Following this approach, all the relevant contributions to the averaged functions of the luminosity redshift relation,at second order, can be parameterized in the form: (cid:104) X (cid:105) = (cid:90) ∞ dkk P ψ ( k, η o ) C X ( k, η o , η (0) s ) , (3.17) (cid:104) X (cid:48) (cid:105)(cid:104) Y (cid:48) (cid:105) = (cid:90) ∞ dkk P ψ ( k, η o ) C X (cid:48) ( k, η o , η (0) s ) C Y (cid:48) ( k, η o , η (0) s ) , (3.18)valid for any given model of perturbation spectrum. Here, X (cid:48) and Y (cid:48) ( X ) are first (second) order generic functionsof ( η, r, θ a ), and the C are the associated spectral coefficients. In the particularly simple case of a CDM-dominatedgeometry the spectral distribution of sub-horizon scalar perturbations is time-independent ( ∂ η ψ k = 0) and, as we shallsee later, all the spectral coefficients C can be calculated analytically. We should stress, however, that when performingnumerical calculations the above integration limits will be replaced by appropriate cut-off values determined by thephysical range of validity of the considered spectrum. IV. DYNAMICAL EVOLUTION OF SCALAR PERTURBATIONS UP TO SECOND-ORDER
For a full computation of the fractional correction f Φ what we need, at this point, is the combined angular and ensemble averages of the three basic quantities, I , I , and I . We must evaluate, in particular, the spectralcoefficients {C T (1) i } . {C T (1 , i } and {C T (2) i } , related to the terms defined in Eqs. (2.20-2.22) in terms of the first andsecond-order Bardeen potential. For this purpose we need to know the dynamical evolution of the scalar fluctuations,at first order for the computation of I , I , and at second order for I .Let us consider, first of all, a general model with cosmological constant plus dust sources. For the evolution of thescalar degrees of freedom in the Poisson gauge we will follow the analysis performed, up to second order, in [25]. Ina general ΛCDM model the linear scalar perturbation obeys the evolution equation ψ (cid:48)(cid:48) + 3 H ψ (cid:48) + (cid:0) H (cid:48) + H (cid:1) ψ = 0 . (4.1)Considering only the growing mode solution we can set ψ ( η, (cid:126)x ) = g ( η ) g ( η o ) ψ o ( (cid:126)x ) , (4.2)where g ( η ) is the so-called “growth-suppression factor”, or – more precisely – the least decaying mode solution, and ψ o is the present value of the gravitational potential. This growth factor can be expressed analytically in terms ofelliptic functions [33] (see also [34]), and it is well approximated by a simple function of the critical-density parametersof non-relativistic matter (Ω m ) and cosmological constant (Ω Λ ) as follows: g = 52 g ∞ Ω m Ω / m − Ω Λ + (1 + Ω m / Λ / . (4.3)Here g ∞ represents the value of g ( η ) at early enough times when the cosmological constant was negligible, and is fixedby the condition g ( η o ) = 1.The second-order potentials obey a similar evolution equation, containing, however, an appropriate source term.Their final expression in terms of ψ o has been given in [25] and reads: ψ (2) ( η ) = (cid:18) B ( η ) − g ( η ) g ∞ −
103 ( a nl − g ( η ) g ∞ (cid:19) ψ o + (cid:18) B ( η ) − g ( η ) g ∞ (cid:19) O ij ∂ j ψ o ∂ i ψ o + B ( η ) O ij ∂ j ψ o ∂ i ψ o + B ( η ) O ij ∂ j ψ o ∂ i ψ o , (4.4) φ (2) ( η ) = (cid:18) B ( η ) + 4 g ( η ) − g ( η ) g ∞ −
103 ( a nl − g ( η ) g ∞ (cid:19) ψ o + (cid:34) B ( η ) + 43 g ( η ) (cid:18) e ( η ) + 32 (cid:19) − g ( η ) g ∞ (cid:35) ×O ij ∂ j ψ o ∂ i ψ o + B ( η ) O ij ∂ j ψ o ∂ i ψ o + B ( η ) O ij ∂ j ψ o ∂ i ψ o , (4.5)0where O ij = ∇ − (cid:18) δ ij − ∂ i ∂ j ∇ (cid:19) , O ij = ∂ i ∂ j ∇ , O ij = δ ij , (4.6)and where we have introduced the functions B A ( η ) = H − o [ l ( η o ) + 3Ω m ( η o ) / − ˜ B A ( η ), with A = 1 , , ,
4, and withthe following definitions:˜ B ( η ) = (cid:90) ηη m d ˜ η H (˜ η )( l (˜ η ) − C ( η, ˜ η ) , ˜ B ( η ) = 2 (cid:90) ηη m d ˜ η H (˜ η ) (cid:104) l (˜ η ) − − m (˜ η ) (cid:105) C ( η, ˜ η ) , (4.7)˜ B ( η ) = 43 (cid:90) ηη m d ˜ η (cid:18) e (˜ η ) + 32 (cid:19) C ( η, ˜ η ) , ˜ B ( η ) = − (cid:90) ηη m d ˜ η C ( η, ˜ η ) , (4.8)and with C ( η, ˜ η ) = g (˜ η ) a (˜ η ) a ( η o ) (cid:104) g ( η ) H (˜ η ) − g (˜ η ) a (˜ η ) a ( η ) H ( η ) (cid:105) , e ( η ) = l ( η ) / Ω m ( η ) , l ( η ) = 1 + g (cid:48) / ( H g ) . (4.9)Here η m denotes the time when full matter domination starts [25]. Its precise value is irrelevant since the region ofintegration around η m is strongly suppressed. Finally, a nl is the so-called non-gaussianity parameter (see [25]), whichapproaches unity in the standard inflationary scenario.For further use let us now evaluate the ensemble (and angular/light-cone) average of the different operators definedin Eq.(4.6), when applied to ∂ i ψ o ∂ j ψ o . Considering first the ensemble average of O ij ∂ i ψ o ∂ j ψ o , and Fourier-expanding ψ o , we get (see [35]): O ij ∂ i ψ o ∂ j ψ o = (cid:90) d q d k (2 π ) δ (3) ( (cid:126)q ) e i(cid:126)q · ˆ xr ψ | (cid:126)k | ψ ∗| (cid:126)k − (cid:126)q | (cid:34) − (cid:126)k · (cid:126)q ) + | (cid:126)k | | (cid:126)q | + 3 ( (cid:126)k · (cid:126)q ) | (cid:126)q | (cid:35) . (4.10)(from this point, and up to the end of this section, we will neglect all suffixes “ o ” present in terms inside the integrals).By using the Taylor expansion of ψ ∗| (cid:126)k − (cid:126)q | around (cid:126)q = 0 we have: ψ ∗| (cid:126)k − (cid:126)q | (cid:39) ψ ∗ k − (cid:126)k · (cid:126)qk ∂ k ψ ∗ k + 12 (cid:40)(cid:32) q k − ( (cid:126)k · (cid:126)q ) k (cid:33) ∂ k ψ ∗ k + ( (cid:126)k · (cid:126)q ) k ∂ k ψ ∗ k (cid:41) + O ( q ) , (4.11)where k ≡ | (cid:126)k | , q ≡ | (cid:126)q | , and where the latter terms have been obtained by using the Hessian matrix H ij = ∂ k i ∂ k j ψ ∗ k = ∂ k i [( k j /k ) ∂ k ψ ∗ k ]. Combining the last two results, writing the integral over k as (cid:82) ∞ πk dk (cid:82) +1 − d cos α ( ... ), where (cid:126)k · (cid:126)q = kq cos α , and integrating over cos α , we obtain: O ij ∂ i ψ o ∂ j ψ o = (cid:90) d q (2 π ) δ (3) ( (cid:126)q ) e i(cid:126)q · ˆ xr (cid:90) ∞ πk dk (cid:20) k ψ k ∂ k ψ ∗ k + 4 k ψ k ∂ k ψ ∗ k (cid:21) . (4.12)Note that the integrand’s dependence on the angle θ between (cid:126)x and (cid:126)q only arises from the exponential term exp( i(cid:126)q · ˆ xr ),which disappears in the presence of δ (3) ( (cid:126)q ). As a consequence, the angular average has no impact on this particularterm and we get: (cid:104)O ij ∂ i ψ o ∂ j ψ o (cid:105) = O ij ∂ i ψ o ∂ j ψ o = (cid:90) ∞ dk k π (cid:20) k ψ k ∂ k ψ ∗ k + 2 k ψ k ∂ k ψ ∗ k (cid:21) . (4.13)Let us note also that, quite generally, ∂ k ψ k ∼ ψ k /k , and thus the above term is of the same order as (cid:104) ψ o ψ o (cid:105) .By repeating exactly the same procedure for the O ij ∂ i ψ o ∂ j ψ o term, we find: O ij ∂ i ψ o ∂ j ψ o = (cid:90) d q (2 π ) δ (3) ( (cid:126)q ) e i(cid:126)q · ˆ xr (cid:90) ∞ πk dk (cid:20) k | ψ k | + kq ψ k ∂ k ψ ∗ k + k q ψ k ∂ k ψ ∗ k (cid:21) , (4.14)where the last two contributions are vanishing because of the q -integration. We are thus left with the following simpleresult (insensitive to the angular average, as before) (cid:68) O ij ∂ i ψ o ∂ j ψ o (cid:69) = O ij ∂ i ψ o ∂ j ψ o = (cid:90) ∞ dkk (cid:20) k P ψ ( k, η o ) (cid:21) . (4.15)Finally, the last term O ij ∂ i ψ o ∂ j ψ o = | (cid:126) ∇ ψ o | is trivial, and gives (cid:68) O ij ∂ i ψ o ∂ j ψ (cid:69) = O ij ∂ i ψ o ∂ j ψ o = (cid:90) ∞ dkk (cid:2) k P ψ ( k, η o ) (cid:3) . (4.16)1 TABLE I: The spectral coefficients C T (1) i for the (cid:68) T (1) i (cid:69) (cid:68) T (1) j (cid:69) terms. T (1) i C T (1) i ( k, η , η s ) T (1)1 − sin k ∆ ηk ∆ η T (1)2 − (cid:16) − H s ∆ η (cid:17) (cid:16) − sin k ∆ ηk ∆ η (cid:17) + 2 (cid:16) − H s ∆ η (cid:17) f s ∆ η (cid:16) cos k ∆ η − sin k ∆ ηk ∆ η (cid:17) T (1)3 4 k ∆ η SinInt( k ∆ η ) V. AN ILLUSTRATIVE EXAMPLE: THE CDM MODEL
In the previous sections we have computed, up to the second perturbative order, the general form of the correctionsinduced by the stochastic fluctuations of the geometry on the averaged luminosity flux, for a spatially flat FLRWmetric and for a generic spectrum of metric perturbations. We have found that such corrections are controlled by thecombined angular and ensemble averages of three basic quantities, I , I , and I . In this section we will evaluatesuch averages for the particular case of the standard CDM model, as a simple illustrative example. The results wewill obtain, will give us useful information about the dominant terms to be selected also for the discussion of thephenomenologically more relevant case, the ΛCDM model, presented in the forthcoming sections. A. The quadratic first-order contributions (cid:104)I (cid:105) and (cid:104)I , (cid:105) Let us now explicitly calculate the spectral coefficients for the first two backreaction terms (induced by the averaging)contributing to Eq. (3.3). The genuine second-order term (cid:104)I (cid:105) will be discussed in the next subsection. We startconsidering the first term, (cid:104)I (cid:105) . From Eq.(2.20) we obtain: (cid:104)I (cid:105) = (cid:88) i =1 3 (cid:88) j =1 (cid:68) T (1) i (cid:69) (cid:68) T (1) j (cid:69) = (cid:90) ∞ dkk P ψ ( k, η o ) (cid:88) i =1 3 (cid:88) j =1 C T (1) i ( k, η o , η s ) C T (1) j ( k, η o , η s ) (5.1)(notice that, from now on, the background solution η (0) s will be simply denoted by η s ). For a dust-dominated phasethe spectral distribution of sub-horizon scalar perturbations is time independent ∂ η ψ k = 0, and the scale factor a ( η )can be written as a ( η ) = a ( η o )( η/η o ) . We can also define f o,s ≡ (cid:90) η o,s η in dη a ( η ) a ( η o,s ) = η o,s − η in η o,s (cid:39) η o,s (5.2)(recall that η in satisfies, by definition, η in (cid:28) η o,s ). All the spectral coefficients of Eq. (5.1) can then be easilycalculated, and the result is reported in Table I . In a similar way, the second contribution to Eq. (3.3) can beexpressed as: (cid:104)I , (cid:105) = (cid:90) ∞ dkk P ψ ( k, η o ) (cid:88) i =1 C T (1 , i , (5.3)and the explicit form of the spectral coefficients C T (1 , i , for the CDM case, is presented in Appendix B.Considering all contributions generated by I and I , , we find that the dominant contributions are all contained in (cid:104)I , (cid:105) , and are characterized by spectral coefficients proportional to k (such dominant terms have been emphasized,in the Appendix B, by enclosing them in a rectangular box). They correspond, in particular, to the terms T (1 , i with We take the opportunity to point out two misprints appearing in Table 2 of [14] where the spectral coefficients of both A and A shouldhave the opposite sign. { i = 2 , , , , , , , , } . Including only such dominant contributions we find that, to leading order, (cid:34) (cid:88) i =1 C T (1 , i (cid:35) Lead = ℵ s f s + f o k + 2(Ξ s − H s f s k + Ξ s f s k − Ξ s f o k = f o k (cid:18) ℵ s − Ξ s (cid:19) + f s k (cid:18) ℵ s s − (cid:19) = − k H (cid:101) f , ( z ) , (5.4)where we have defined ℵ s = Ξ s − H s ∆ η (cid:18) − H (cid:48) s H s (cid:19) , (5.5)and we have used the relation H s (cid:39) / (3 f s ). Also, we have included into the function (cid:101) f , ( z ) all the z -dependence ofthese leading contributions. After some simple algebra we find: (cid:101) f , ( z ) = 10 − √ z + 5 z (cid:0) √ z (cid:1)
27 (1 + z ) (cid:0) − √ z (cid:1) . (5.6)This function (and thus the corresponding backreaction) flips sign around z ∗ = 0 .
205 (as illustrated also in Fig. 1).It should be noted, finally, that some of the genuine second-order terms contained into I are also associated tospectral coefficients proportional to k . However, as we shall see in the next subsection, the corresponding contributionsto Eq. (3.3) turn out to be roughly an order of magnitude smaller than the above ones, because of approximatecancellations. B. The genuine second-order contribution (cid:104)I (cid:105) In order to complete the calculation of f Φ ( z ) we still have to consider the genuine second-order backreaction term (cid:104)I (cid:105) . In particular, we must evaluate the spectral coefficients {C T (2) i } , corresponding to the various contributionsdefined in Eq. (2.22), in terms of the first-order Bardeen potential. By using the results of Sect. IV, and starting fromEqs. (4.4) and (4.5), we can easily see that the only possible k -enhanced contributions arising from the coefficients {C T (2) i } (which only contain functions of φ (2) and ψ (2) ) should correspond to the term B ( η ) ∇ − ∂ i ∂ j ( ∂ i ψ o ∂ j ψ o ) + B ( η ) ∂ i ψ o ∂ i ψ o . Therefore, we should obtain φ (2) (cid:39) ψ (2) at leading order.Let us discuss and estimate all possible contributions for the CDM model we are considering in this section. In thissimple case l ( η ) = 1 and B ( η ) = B ( η ) = 0 (see Eq.(4.7)). Furthermore, restricting our attention to the standardinflationary scenario, we can set a nl = 1. Eqs.(4.4) and (4.5) thus reduce to: ψ (2) = − ψ o − O ij ∂ i ψ o ∂ j ψ o + B ( η ) O ij ∂ i ψ o ∂ j ψ o + B ( η ) O ij ∂ i ψ o ∂ j ψ o , (5.7) φ (2) = 2 ψ o + 2 O ij ∂ i ψ o ∂ j ψ o + B ( η ) O ij ∂ i ψ o ∂ j ψ o + B ( η ) O ij ∂ i ψ o ∂ j ψ o . (5.8)Finally, we can use (as before) the scale factor a ( η ) = a ( η o ) ( η/η o ) , and obtain (according to Eq. (4.8)): B ( η ) = 2021 1 H , B ( η ) = −
27 1 H . (5.9)We are now in the position of evaluating the genuine second-order terms, by exploiting the results given in Eqs.(4.13),(4.15) and (4.16). Following the classification of Eq.(2.22) we can see that the first two terms C ( T (2)1 ) and C ( T (2)2 )exactly cancel for the CDM case (while they cancel only at leading order for the case of a ΛCDM model). For CDMwe have, in particular, C ( T (2)1 ) = −C ( T (2)2 ) = − Ξ s
252 ( η s − η o ) k . (5.10)Another interesting simplification concerns the terms T (2)3 , T (2)4 and T (2)5 , namely those terms for which the in-tegrand contains ∂ r ψ (2) or ∂ r φ (2) . From our previous results, in particular from Eqs. (4.13), (4.15) and (4.16), it3is easy to see that the ψ (2) and φ (2) contributions are unchanged when one averages over the 2-sphere (i.e. over θ by isotropy), since the θ -dependence is removed by the presence of δ (3) ( (cid:126)q ). On the other hand, the presence ofthe r -derivative brings a further factor | (cid:126)q | cos θ , and the q -integration gives zero, even before performing the angularaverage. It follows that, for a general model, C ( T (2)3 ) = 0 ; C ( T (2)4 ) = 0 . (5.11)The contribution of T (2)5 , on the contrary, is nonvanishing because of the presence of the partial derivative ∂ η , actingon the B A ( η ) coefficients. For the CDM model we find, in particular, C ( T (2)5 ) = Ξ s
126 ( η s − η o ) k . (5.12)Finally, for the last two terms T (2)6 and T (2)7 we obtain: C ( T (2)6 ) = − − k η s
126 + 32 k ψ k ∂ k ψ ∗ k | ψ k | + 8 k ψ k ∂ k ψ ∗ k | ψ k | , C ( T (2)7 ) = η o − η s η k + 16 k ψ k ∂ k ψ ∗ k | ψ k | + 4 k ψ k ∂ k ψ ∗ k | ψ k | . (5.13)The sum of all contributions then leads to: (cid:88) i =1 C ( T (2) i ) = − (cid:20) Ξ s f s − f o ) − f s − f s − f o ∆ η (cid:21) k + 16 k ψ k ∂ k ψ ∗ k | ψ k | + 4 k ψ k ∂ k ψ ∗ k | ψ k | . (5.14)All the above spectral coefficients are now to be numerically evaluated by using the power spectrum of the CDMmodel. We can easily check, however, that the leading k -contributions of these coefficients are given by: (cid:34) (cid:88) i =1 C ( T (2) i ) (cid:35) Lead = 17 (cid:20) Ξ s f s − f o ) − f s − f s − f o ∆ η (cid:21) k = − k H o (cid:101) f ( z ) , (5.15)where: (cid:101) f ( z ) = − − √ z + z (cid:0) − √ z (cid:1) (1 + z )( √ z − . (5.16)Such second-order contributions turn out to be about one order of magnitude smaller than the leading contributionsof the squared first-order terms of Sect. V A (as can be easily checked, for instance, by comparing the plots of ˜ f , and ˜ f ). C. Full numerical results for the CDM model
At this point, in order to perform the numerical computations, we need to insert the explicit form of the powerspectrum. Limiting ourselves to sub-horizon perturbations we can simply obtain ψ k , for the CDM model, by applyingan appropriate, time-independent transfer function to the primordial (inflationary) spectral distribution (see e.g. [27]).The power spectrum of the Bardeen potential is then given by: P ψ ( k ) = (cid:18) (cid:19) ∆ R T ( k ) , ∆ R = A (cid:18) kk (cid:19) n s − , (5.17)where T ( k ) is the transfer function which takes into account the sub-horizon evolution of the modes re-enteringduring the radiation-dominated era, and ∆ R is the primordial power spectrum of curvature perturbations outside thehorizon. The typical parameters of such a spectrum, namely the amplitude A , the spectral index n s and the scale k ,are determined by the results of recent WMAP observations [36]. In our computations we will use, in particular, thefollowing approximate values: A = 2 . × − , n s = 0 . , k = 0 .
002 Mpc − . (5.18)4Finally, since our main purpose here is to present an illustrative example, it will be enough for our needs to approximate T ( k ) by the effective shape of the transfer function for density perturbations without baryons, namely T ( k ) = T ( k ),where [15]: T ( q ) = L L + q C ( q ) , L ( q ) = ln(2 e + 1 . q ) , C ( q ) = 14 . . q , q = k . k eq , (5.19)and where k eq (cid:39) .
07 Ω m h M pc − is the scale corresponding to matter-radiation equality, with h ≡ H / (100 km s − Mpc − ).We can easily check that the above transfer function goes to 1 for k (cid:28) k eq , while it falls like k − log k for k (cid:29) k eq .For the numerical estimates we will use h = 0 . a o = 1, Ω m = 1. In that case we obtain k eq (cid:39) .
036 Mpc − (see [15]), and we can more precisely define the asymptotic regimes of our transfer function as T (cid:39) k < ∼ − Mpc − , and T ∼ k − log k for k > ∼ . − .Following the results of the previous subsections we can now set f Φ ( z ) = (cid:90) ∞ dkk P ψ ( k, z = 0) (cid:104) f , ( k, z ) + f ( k, z ) (cid:105) , (5.20)where f , and f are complicated –but known for the CDM case– analytic functions of their arguments. However, asalready stressed, the leading contributions in the range of z relevant to dark-energy phenomenology are sourced byterms of the type f ( k, z ) ∼ ( k/ H o ) (cid:101) f ( z ). In that range of z we can thus write, to a very good accuracy: f Φ ( z ) (cid:39) (cid:104) (cid:101) f , ( z ) + (cid:101) f ( z ) (cid:105) (cid:90) ∞ dkk (cid:18) k H o (cid:19) P ψ ( k, z = 0) , (5.21)where (cid:101) f , ( z ) and (cid:101) f ( z ) are given, respectively, by Eqs. (5.6) and (5.16).To proceed, we need to insert a power spectrum as well as infrared (IR) and ultraviolet (UV) cutoffs in (5.20). Theformer can be identified with the present horizon H − o ; however, considering the used spectra, larger scales give acompletely negligible contribution. On the other hand, in spite of the fact that our expressions converge in the UVfor any reasonable power spectrum, some mild sensitivity to the actual UV cutoff will be shown to occur in certainobservables. The absolute value (and sign) of f Φ ( z ) for the CDM model, obtained from both Eq. (5.20) and (5.21),are illustrated in Fig. 1, where we can explicitly check the accuracy of the leading order terms (5.21). The figure alsoconfirms that the backreaction of a realistic spectrum of stochastic perturbations induces negligible corrections to theaveraged flux at large z (the larger corrections at small z , due to “Doppler terms”, have been already discussed alsoin [14]). In addition, it shows that such corrections have the wrong z -dependence (in particular, they change sign atsome z ) for simulating even a tiny dark-energy component. VI. THE Λ CDM MODEL: POWER SPECTRUM IN THE LINEAR REGIME
We will now extend the procedure of the previous section to the case of the so-called “concordance” cosmologicalmodel, using first a power spectrum computed in the linear regime, and then adding the effects of non-linearitiesfollowing the parametrizations proposed in [16] and [17]. In both cases, we will restrict our attention only to the k -enhanced terms already identified in the CDM model for the light-cone average of the flux variable, as well as tothe k -enhanced terms which, as we will see, will appear in the variance, or in the averages of other functions of d L ( z ).Hereafter in all numerical computations we will use, in particular, the following numerical values: Ω Λ0 = 0 . m = 0 .
27, Ω b = 0 . h = 0 . A. Second-order corrections to the averaged luminosity flux
The power spectrum of the ΛCDM model is, in general, time-dependent. Considering for the moment only thelinear regime, the scalar power spectrum can be written, starting from Eq.(4.2), as P ψ ( k, η ) = (cid:20) g ( η ) g ( η o ) (cid:21) P ψ ( k, η = η o ) , (6.1)5 CDM (cid:45) f (cid:70) f (cid:70) (cid:45) (cid:45) (cid:45) z F r ac ti on a l C o rr ec ti on s FIG. 1: The fractional correction f Φ of Eq. (5.20) (solid curve), compared with the same quantity given to leading order byEq. (5.21) (dashed curve), in the context of an inhomogeneous CDM model. We have used for the spectrum the one definedin Eq. (5.17). The plotted curves refer, as an illustrative example, to an UV cutoff k UV = 1Mpc − . and in this case we can easily extend the results previously obtained for the CDM model, concerning the leading( k -enhanced) contributions to the averaged-flux integral I φ . Using the general definitions of the T (1 , i terms, wefirst notice that such enhanced contributions arise from the following particular (sub)-terms appearing in Eq. (2.21): T (1 , ,L = Ξ s (cid:0) ([ ∂ r P ] s ) − ([ ∂ r P ] o ) (cid:1) ; T (1 , ,L = 12 Ξ s ( γ ab ) s (2 ∂ a P s ∂ b P s ) ; T (1 , ,L = − Ξ s lim r → (cid:2) γ ab ∂ a P ∂ b P (cid:3) ; T (1 , ,L = − s ([ ∂ r P ] s ) ; T (1 , ,L = 2Ξ s H s (cid:18) ψ s − (cid:90) η o η s dη (cid:48) ∂ r ψ ( η (cid:48) , η o − η (cid:48) , ˜ θ a ) (cid:19) [ ∂ r P ] s ; T (1 , ,L = (cid:20) Ξ s − H s ∆ η (cid:18) − H (cid:48) s H s (cid:19)(cid:21) (cid:0) ([ ∂ r P ] s ) + ([ ∂ r P ] o ) (cid:1) ; T (1 , ,L = 2Ξ s ([ ∂ r P ] o ) ; T (1 , ,L = − H s ( γ ab ) s ∂ a Q s ∂ b [ ∂ r P ] s ; T (1 , ,L = − H s [ ∂ r P ] s [ ∂ r ψ ] s , (6.2)where we have added the suffix “ L ” to stress that we are reporting here only the leading terms.In order to calculate the corresponding spectral coefficients, we note that their time dependence can be factorizedwith respect to the k -dependence whenever the time variable does not appear in the exponential factor exp( i(cid:126)k · (cid:126)x )present in our integrals (see for instance Eqs. (3.11), (3.12), (3.14), (3.15)). This is indeed the case for the terms T (1 , ,L , T (1 , ,L , T (1 , ,L , T (1 , ,L , T (1 , ,L , T (1 , ,L and T (1 , ,L . In that case the previous CDM results for the leading spectralcoefficients can be simply generalized to the ΛCDM case through the following procedure: ( i ) by inserting a factor g ( η s ) /g ( η o ) whenever ψ s is present in the initial term; and ( ii ) by replacing the f o,s factors (see Eq.(5.2)), arisingfrom the presence of P o,s terms, by: (cid:101) f o,s ≡ (cid:90) η o,s η in dη a ( η ) a ( η o,s ) g ( η ) g ( η o ) . (6.3)6For the remaining two terms T (1 , ,L and T (1 , ,L the integrals are performed along the path r = η o − η , and thetime dependence cannot be fully factorized. Therefore, the evaluation of the double integrals over η and k is muchmore involved than in the previous cases. However, a good approximation of the exact result can be obtained byreplacing ψ ( η (cid:48) , η o − η (cid:48) , ˜ θ a ), appearing in the integrands of T (1 , ,L and T (1 , ,L , with ψ ( η s , η o − η (cid:48) , ˜ θ a ) (this is so sincethe leading contributions to the time integral arise from a range of values of η approaching η s ). By adopting such anapproximation we can follow the same procedure as before, and we obtain: C ( T (1 , ,L ) = 23 Ξ s (cid:101) f s H s k , C ( T (1 , ,L ) = − (cid:101) f s H s k . (6.4)This is formally the same result as in the CDM case (see Appendix B for the leading terms of T (1 , ,L and T (1 , ,L ), withthe only difference that f s is replaced by (cid:101) f s .Let us now move to the evaluation of the leading contributions present in the genuine second-order part I . Thefinal results of Sect. V, for the particular case of a CDM model, can be easily generalized to the ΛCDM casestarting from the observation that the leading contributions can only arise from terms containing the operators O ij and O ij in Eqs.(4.4) and (4.5). As a consequence, the first two terms T (2)1 and T (2)2 will give a subleading overallcontribution, while the general result that the terms T (2)3 and T (2)4 give identically zero still holds. The remainingleading contributions can be easily obtained, using the results in Eqs. (4.15) and (4.16), as follows: C ( T (2)5 ,L ) = − Ξ s (cid:20)
13 ( B ( η o ) − B ( η s )) + ( B ( η o ) − B ( η s )) (cid:21) k , (6.5) C ( T (2)6 ,L ) = − (cid:18) B ( η s ) + B ( η s ) (cid:19) k , (6.6) C ( T (2)7 ,L ) = 2∆ η (cid:90) η o η s dη (cid:48) (cid:18) B ( η (cid:48) ) + B ( η (cid:48) ) (cid:19) k . (6.7)These leading contributions can be now evaluated using Eqs. (4.8) and (4.9) and moving to redshift space, where wecan write: H ( z ) = H o z (cid:2) Ω m (1 + z ) + Ω Λ0 (cid:3) / , Ω m = Ω m (1 + z ) Ω m (1 + z ) + Ω Λ0 , Ω Λ = Ω Λ0 Ω m (1 + z ) + Ω Λ0 , (6.8)where the suffix “0” appended to Ω m and Ω Λ denotes the present value of those fractions of critical density. We notethat, as in the CDM case, the leading genuine second-order contributions are about one to two orders of magnitudesmaller than the leading squared first-order contributions, evaluated above.We need now to insert the explicit form of the power spectrum. Considering the general solution of Eq. (4.2) wecan also re-express the z -dependence of the power spectrum, in the linear regime, as follows: P ψ ( k, z ) = (cid:18) (cid:19) ∆ R T ( k ) (cid:18) g ( z ) g ∞ (cid:19) , (6.9)where the previous CDM result (5.17), based on the transfer function T ( k ) given in [15], is modified by the presenceof the factor g ( z ) /g ∞ originating from the time dependence of the gravitational perturbations. Another modificationwith respect to the CDM result, implicitly contained into the transfer function T ( k ), concerns the different numericalvalue of the equilibrium scale k eq (which now turns out to be lower because of the lower value of Ω m ). The effectsof such modifications are illustrated in Fig. 2 by comparing the ΛCDM spectrum of Eq. (6.9), at different values of z , with two z -independent spectra: the primordial spectrum of scalar perturbations (3 / ∆ R , and the “transferred”spectrum of the CDM model, introduced in Sect. V (notice that, as expected, the ΛCDM and CDM spectra tendto coincide at large enough values of z ). In the ΛCDM case the solid curves are obtained with a transfer functionwhich takes into account the presence of baryonic matter, while the dashed curves correspond to the transfer function T ( k ) = T ( k ) of Sect.V C (without baryons). Here, however, we are always using a spectrum evaluated in the linearregime.As illustrated in Fig. 2, the effect of neglecting the baryonic fraction of Ω m (and thus the associated Silk-dampingeffect) may lead to an overestimation up to 40% of the corresponding transfer function for scalar perturbations, inthe range k > ∼ . h Mpc − (see [15]). In order to take into account this baryonic contribution we have to replace the7 (cid:45) (cid:215) (cid:45) (cid:215) (cid:45) Λ CDM CDMz=0z=0.5z=1.5 P r i m o r d i a l Sp e c t r u m P Ψ ( k ) k [ h Mpc − ] FIG. 2: A comparison of the primordial inflationary spectrum (long-dashed curve) with the spectrum of the CDM modelneglecting baryons (thick solid curve) and of a ΛCDM model (thin solid curves), at various values of z . The dotted curves forthe ΛCDM case describe the spectrum obtained by neglecting the baryon contribution (hence without taking into account theSilk-damping effect). value of q used in the previous section, i.e. q = k/ (13 . k eq ), with the more accurate value given by [15]: q = k × M pch Γ , (6.10)Γ = Ω m h (cid:18) α Γ + 1 − α Γ . ks ) (cid:19) , (6.11) α Γ = 1 − .
328 ln(431Ω m h ) Ω b Ω m + 0 .
38 ln(22 . m h ) (cid:18) Ω b Ω m (cid:19) , (6.12) s = 44 . . / Ω m h ) (cid:112) b h ) / Mpc . (6.13)Here Ω b is the baryon density parameter, s is the sound horizon and Γ is the k -dependent effective shape parameter.We have compared the above transfer function to the one which includes baryon acoustic oscillations (BAO) [15],and to a transfer function calculated numerically by using the so-called “code for anisotropies in the microwavebackground” (CAMB) [37]. We have checked, in particular, that the above simple form of transfer function isaccurate to within a few percent compared to the one calculated numerically by CAMB, for all scales of interest. Inaddition, the effect of including BAO only produces oscillations of the spectrum around the above value. Since weare considering here integrals over a large range of k , the presence of BAO has a negligible effect on our final results,and will be neglected in the rest of the paper.We are now in the position of computing the fractional corrections to the averaged flux variable in a perturbedΛCDM geometry. As discussed before, there are complicated average integrals which can be performed only by usingsome approximations. Once this is done, the remaining integration over k can be done numerically, exactly as in thecase of the CDM model.In a ΛCDM context we may generally expect smaller corrections to the averaged flux, due to the fact that theperturbation spectrum P ψ is suppressed by the presence of g ( z ). In addition (and as already stressed) the transferfunction [15] turns out to be suppressed, at large k , because of a smaller value of the parameter k eq (see Eq. (5.19)).These expectations are fully confirmed by an explicit numerical computation of | f Φ | , which we have performed withand without the inclusion of the baryon contributions into the transfer function T ( k ).The results of such a computation are illustrated in Fig. 3, and a comparison with Fig. 1 clearly shows that | f Φ | (cid:76) CDM f d (cid:45) f (cid:70) f (cid:70) (cid:45) (cid:45) (cid:45) z F r ac ti on a l C o rr ec ti on s (cid:76) CDM f d (cid:45) f (cid:70) f (cid:70) (cid:45) (cid:45) (cid:45) z F r ac ti on a l C o rr ec ti on s FIG. 3: The fractional correction to the flux f Φ of Eq. (3.3) (thin curves) is plotted together with the fractional correction tothe luminosity distance f d of Eq. (3.5) (thick curves), for a ΛCDM model with Ω Λ0 = 0 .
73. We have used two different cutoffvalues: k UV = 0 . − (dashed curves) and k UV = 1Mpc − (solid curves). The left panel shows the results obtained witha linear spectrum without baryon contributions. The right panel illustrates the effects of including baryons (we have used, inparticular, Ω b = 0 . is smaller in the ΛCDM case than in the CDM case, and further (slightly) depressed when we take into accountthe presence of a small fraction of baryon matter (the curves presented in the right panel). In any case, the smallvalues of | f Φ | at relatively large z , for a realistic ΛCDM scenario, lead us to conclude that the averaged flux is aparticularly appropriate quantity for extracting from the observational data the “true” cosmological parameters. Aswe will discuss now, the situation is somewhat different for other functions of d L . B. Second-order corrections to other observables and dispersions
Let us now consider other observables, beyond the flux, to see how the impact of the inhomogeneities may change.We will treat, in particular, the two important examples introduced in Sect. III, namely the luminosity distance d L and the distance modulus µ . The fractional corrections to their averages (see Eqs. (3.4–3.6)) are qualitatively differentfrom those of the averaged flux (represented by f Φ ), because of the presence of extra contributions, unavoidable forany non-linear function of the flux and proportional to the square of the first order fluctuation (Φ / Φ ) .For a better understanding of such contributions let us start with the results obtained in [13], concerning thesecond-order perturbative expansion of the luminosity distance, d L = d (0) L + d (1) L + d (2) L , and summarized in Sect. II B(see in particular Eqs .(2.13) and (2.15)). Using those results we obtain:Φ Φ = − d (1) L d (0) L = −I + (t . d . ) (1) = 2 (cid:18) − Ξ s J + Q s ∆ η + ψ s + J (1)2 (cid:19) . (6.14)In order to determine the leading corrections we first notice that, by applying the procedure of Sect. V, and computingthe averages for the CDM case, we obtain, to leading order, (cid:42)(cid:18) Φ Φ (cid:19) (cid:43) L = 4 (cid:40)(cid:28)(cid:16) J (1)2 (cid:17) (cid:29) + Ξ s (cid:104) (cid:104) ([ ∂ r P ] s ) (cid:105) + (cid:104) ([ ∂ r P ] o ) (cid:105) (cid:105)(cid:41) . (6.15)Working in the context of a ΛCDM model, considering only these leading terms, and limiting ourselves to the linearregime, we find that the terms multiplying Ξ s on the right hand side of the above equation can be calculated withoutapproximations (as seen in the previous subsection). Also, we find that their contribution is controlled by the spectralfactor k P ψ ( k, η o ). The first term, on the contrary, is due to to the so-called “lensing effect”, dominates at large z (as already shown in [14] for a CDM model), and has leading spectral contributions of the type k P ψ ( k, η o ).For such term, however, the integrals over time cannot be factorized with respect to the k integrals, and we mustuse the approximation already introduced in the previous subsection (namely, we have to replace in the integrands ψ ( η (cid:48) , η o − η (cid:48) , ˜ θ a ) with ψ ( η s , η o − η (cid:48) , ˜ θ a )). The full result for the spectral coefficient can then be finally written as9 (cid:45) (cid:45) z Μ (cid:45) Μ M (cid:45) (cid:45) z Μ (cid:45) Μ M FIG. 4: The averaged distance modulus (cid:104) µ (cid:105) − µ M of Eq. (3.6) (thick solid curve), and its dispersion of Eq. (3.9) (shadedregion) are computed for Ω Λ0 = 0 .
73 and compared with the homogeneous value for the unperturbed ΛCDM models with,from bottom to top, Ω Λ0 = 0 .
69, 0 .
71, 0 .
73, 0 .
75, 0 .
77 (dashed curves). We have used k UV = 1 Mpc − . The left panel showsthe results obtained with a linear spectrum without baryon contributions. The right panel illustrates the effects of includingbaryons, with Ω b = 0 . follows: C ((Φ / Φ ) L ) (cid:39)
43 Ξ s (cid:16) (cid:101) f s + (cid:101) f o (cid:17) k + 415 (cid:20) g ( η s ) g ( η o ) (cid:21) ∆ η k SinInt( k ∆ η ) . (6.16)Because of the new term (due to lensing) affecting the average of µ and d L – but not of the flux Φ – we may expectlarger fractional corrections for these variables, as well as for other functions of Φ, at higher redshifts. This is indeedconfirmed by the plots presented Fig. 3 reporting the results of an explicit numerical integration and comparing, inparticular, the value of f d with the absolute value of f Φ . We obtain | f Φ | (cid:28) f d , at large values of z where the lensingterm dominates, both in the presence and in the absence of the baryon contribution to the total energy density. Itshould be stressed, however, that also the new k -enhanced contributions are free from IR and UV divergences, atleast for the class of models we are considering.Let us now discuss to what extent the enhanced corrections due to the square of the first-order flux fluctuationcan affect the determination of the dark-energy parameters, if quantities other than the flux are used to fit theobservational data. To this purpose we may consider the much used (average of the) distance modulus given in Eq.(3.6), referring it, as usual, to a homogeneous Milne model with µ M = 5 log [(2 + z ) z/ (2 H )]. Considering Eqs. (3.6)and (3.9), where the averaged value of the distance modulus and its dispersion are given in function of f Φ and of (cid:68) (Φ / Φ ) (cid:69) , we have investigated the magnitude of the effect in this case. In particular, in Fig. 4 we have comparedthe averaged value (cid:104) µ (cid:105) − µ M with the corresponding expression for homogeneous ΛCDM models with different valuesof Ω Λ0 . We have also illustrated the expected dispersion around the averaged result, represented by the dispersionpreviously reported in Eq. (3.9) (and already computed in [14] for the CDM case).We have found that the given inhomogeneities, on the average, may affect the determination of Ω Λ0 obtained fromthe measure of the distance modulus, at large z , only at the third decimal figure (at least if the spectral contributionsare computed in the linear regime). As we can see from Fig. 4, the curves for (cid:104) µ (cid:105) and for the correspondingunperturbed value µ FLRW (with the same Ω Λ0 ) practically coincide at large enough z . It should be stressed, also,that the dispersion on the distance modulus computed from Eq. (3.9) reaches, at large redshift, a value which iscomparable with a change of about 2% in the dark energy parameter Ω Λ0 (cf. [1], considering however that baryoniceffects have now been added in the transfer function). We shall see in the next section that this effect is enhanced bythe use of non-linear power spectra. VII. THE Λ CDM MODEL: POWER SPECTRUM IN THE NON-LINEAR REGIME
The linear spectra considered so far are sufficiently accurate only up to scales of order 0 . h Mpc − . If we want tobetter study the effect of shorter-scale inhomogeneities on our light-cone averages we need to go beyond such linearapproximation, taking into consideration the non-linear evolution of the gravitational perturbations. This can be doneby using the so-called “HaloFit” models, which are known to reproduce quite accurately the results of cosmological0 N -body simulations. In particular, the HaloFit model of [16] and its recent upgrade of [17] provide an accurate fittingformula for the power spectrum up to wavenumber k (cid:39) h Mpc − .In our previous paper [1] the analysis was limited to k < − , but with the use of the non-linear power spectrawe can now extend our analysis up to the maximum scale of the mentioned HaloFit models. The main difficulty withsuch an extension is that the time (i.e. z ) dependence of the spectrum becomes more involved than in the linear case,since different scales no longer evolve independently. Hence, the need for introducing approximations in performingthe integrals becomes even more essential in this case.Let us start by recalling a few details of the HaloFit model. The fractional density variance per unit ln k isrepresented by the variable ∆ ( k ), defined by [16]: σ ≡ δ ( x ) δ ( x ) = (cid:90) d k (2 π ) | δ (cid:126)k | = (cid:90) ∆ ( k ) d ln k , (7.1)which implies ∆ ( k ) = k π | δ k | . (7.2)where δ ( x ) = δρ ( x ) /ρ is the fractional density perturbation of the gravitational sources, and δ k is the associatedFourier component. On the other hand, the power spectrum of scalar perturbations, P Ψ ( k, z ), is related to ∆ ( k, z )by the Poisson equation, holding at both the linear ( L ) and non-linear ( N L ) level [34]: P L , NL ψ ( k, z ) = 94 Ω m H k (1 + z ) ∆ , NL ( k, z ) . (7.3)The linear part of the spectrum is used to introduce the normalization equation, defining the non-linearity lengthscale k − σ ( z ), as follows: σ ( k − σ ) ≡ (cid:90) ∆ ( k, z ) exp( − ( k/k σ ) ) d ln k ≡ . (7.4)It is then obvious that the scale k σ is redshift-dependent. Since the non-linear power spectra obtained from the models[16] and [17] are using such a scale, they will be characterized by an implicit z -dependence which is more complicatedthan the one following from the usual growth factor of Eq. (4.3), and which cannot be factorized.Besides k σ , the linear spectrum also determines two additional parameters, important for the construction of theHaloFit model: the effective spectral index n eff and the parameter C , controlling the curvature of the spectral indexat the scale k σ . They are defined by:3 + n eff ≡ − d ln σ ( R ) d ln R (cid:12)(cid:12)(cid:12)(cid:12) R = k − σ , C ≡ − d ln σ ( R ) d ln R (cid:12)(cid:12)(cid:12)(cid:12) R = k − σ . (7.5)Once the values of k σ , n eff and C are determined, they can be inserted into a given HaloFit model [16, 17] to producethe non-linear power spectrum ∆ ( k, z ). One finally goes back to P NL ψ ( k, z ) using again the Poisson equation (7.3).The non-linear spectra obtained in this way, and based on the previous linear spectrum with baryon contributions,are illustrated in Fig. 5 for the two HaloFit models [16] and [17]. The results are compared with the spectrum of thelinear regime, for different values of the redshift ( z = 0 and z = 1 . k > ∼ . h Mpc − ), as a result of the intricate redshift dependence. We can also observe thatthe two HaloFit models lead to the same result at low values of k (including baryons, but without BAO). A. Numerical results and comparison with the linear regime
In order to evaluate our integrals using the non-linear power spectra we need to face the additional (alreadymentioned) problem due to the fact that the non-linear spectra – unlike the linear ones – cannot be factorized as afunction of k times a function of z . In that case the full two-dimensional integration is highly non-trivial, and we havethus exploited a further approximation. In the presence of time-integrals of the mode functions (like those appearing,for instance, in the leading terms of (cid:104)I , (cid:105) ), we have parametrized the non-linear power spectrum in a factorized formas follows: P NLΨ ( k, z ) = g ( z ) g ( z ∗ ) P NLΨ ( k, z ∗ ) , (7.6)1 (cid:45) (cid:45) (cid:45) z=1.5 z=0 k P NLΨ ( k ) k P LΨ ( k ) k P L , N L Ψ ( k ) [ h M p c − ] k [ h Mpc − ] FIG. 5: The linear spectrum P L ψ (dotted curves) and the non-linear spectrum P NL ψ for the HaloFit model of [16] (dashed curves)and of [17] (solid curves). In all three cases the spectrum is multiplied by k (for graphical convenience) and is given for z = 0(thin curves) and z = 1 . b = 0 . and we have chosen z ∗ = z s / . ≤ z ≤
2, the above approximation is most inadequate only when we consider z s = 2 (i.e. z ∗ = 1), and our mode functions are evaluated inside the integrals at the values of z most distant from z ∗ (namely, z = 0 .
015 and z = 2). In that case we are lead to underestimate the spectrum by about a 40% factor for z = 0 . z = 2. However, this only occurs in two narrow bands of k centeredaround the values k = 1 h Mpc − and k = 2 h Mpc − , while, outside these bands, our approximation is good. Sincethese errors are limited to only a part of the region of integration, both in z and in k , we can estimate an overallaccuracy at the 10% level, at least, for the results given by the adopted approximation. We have also checked thesensitivity of the numerical results to changes in z ∗ , such as z ∗ = z s , and checked that our final results are only weaklydependent (typically at the 1% level) on the choice of z ∗ .Once we have established the range of validity of the parametrization (7.6), we proceed in evaluating the z (or η )integrals as we did in Sect. VI. The final results of this procedure are illustrated by the curves plotted in Figs. 6and 7, computed with the non-linear power spectrum following from the HaloFit model of [17] and including baryoncontributions . As illustrated for instance in Fig. 6, the fractional correction to d L turns out to be of order of a fewparts in 10 − around z = 2, and smaller in the rest of the intermediate redshift range relevant for cosmic acceleration.On the other hand, in the same redshift range, the fractional correction to the flux is about two orders of magnitudesmaller.Comparing with the results obtained with the linear spectrum (see Figs. 3, 4), we can see that taking into accountthe non-linearity distortions (and using higher cut-off values) enhances the backreaction effects on the consideredfunctions of the luminosity distance, but not enough to reach a (currently) observable level. On the other hand, thedispersion, already large in the linear case, is further enhanced when non-linearities are included. In particular, thedispersion on µ due to inhomogeneities is of order 10% around z = 2, which implies that the predictions of homogeneousmodels with Ω Λ0 ranging from 0 .
68 to 0 .
78 lie inside one standard deviation with respect to the averaged predictions The two HaloFit models quoted above give similar results, and we have chosen to present here, for simplicity, only those obtained withone of them. Λ0 = 0 . (cid:76) CDM f d (cid:45) f (cid:70) f (cid:70) (cid:45) (cid:45) (cid:45) z F r ac ti on a l C o rr ec ti on s FIG. 6: The fractional correction to the flux ( f Φ , thin curves) and to the luminosity distance ( f d , thick curves), for a perturbedΛCDM model with Ω Λ0 = 0 .
73. Unlike in Fig. 3, we have taken into account the non-linear contributions to the power spectrumgiven by the HaloFit model of [17] (including baryons), and we have used the following cutoff values: k UV = 10 h Mpc − (dashedcurves) and k UV = 30 h Mpc − (solid curves). (cid:45) (cid:45) z Μ (cid:45) Μ M (cid:45) (cid:45) z Μ (cid:45) Μ M FIG. 7: The averaged distance modulus (cid:104) µ (cid:105) − µ M of Eq. (3.6) (thick solid curve), and its dispersion of Eq. (3.9) (shadedregion), for a perturbed ΛCDM model with Ω Λ0 = 0 .
73. Unlike Fig. 4, we have taken into account the non-linear contributionsto the power spectrum given by the HaloFit model of [17] (including baryons), and used the cut-off k UV = 30 h Mpc − . Theaveraged results are compared with the homogeneous values of µ predicted by unperturbed ΛCDM models with (from bottomto top) Ω Λ0 = 0 .
68, 0 .
69 0 .
71, 0 .
73, 0 .
75, 0 .
77, 0 .
78 (dashed curves). The right panel simply provides a zoom of the samecurves, plotted in the smaller redshift range 0 . ≤ z ≤ B. Comparing theory and observations via the intrinsic dispersion of the data
Let us now consider in more detail our prediction for the dispersion σ µ induced by the presence of the inflationaryperturbation background, and compare it with the intrinsic dispersion of the distance modulus that can be inferredfrom SNe Ia data. Our results for the dispersion are already implicitly contained in Fig. 7 but, for the sake ofclarity, we have separately plotted our value of σ µ in Fig. 8, where the thick solid curve represents the value of σ µ z . We can see from the figure that σ µ has a characteristic z -dependence, with a minimal value of about 0 .
016 reached around z = 0 . z ), and the leading lensing contribution obtained fromthe first term of Eq. (6.15) (represented by the dashed curve approaching zero at small z ).The total variance σ obs µ associated with the observational data, on the other hand, can be decomposed in generalas follows (see e.g. [38–40]): ( σ obs µ ) = ( σ fit µ ) + ( σ zµ ) + ( σ int µ ) . (7.7)Here σ fit µ is the statistical uncertainty due, for instance, to the method adopted for fitting the light curve (e.g. theso-called SALT-II method [41]), but also to the uncertainty in the modeling of the supernova process. The term σ zµ represents instead the uncertainty in redshift due to the peculiar velocity of the supernova as well as to the precisionof spectroscopic measurements. Finally, σ int µ is an unknown phenomenological quantity, needed to account for theremaining dispersion of the data with respect to the chosen homogeneous model. This part of the dispersion canbe subsequently redefined whenever we are able to estimate some of the possible contributions it contains. Thecontribution we are mainly interested in here is the one originating from the lensing effect, which is dominant at largeredshift. We can thus write, at large z : ( σ int µ ) = ( (cid:100) σ int µ ) + ( σ lens µ ) , (7.8)where (cid:100) σ int µ is the remaining source of intrinsic dispersion.Given the typical precision of current data [41, 42], a reasonable fit of the Hubble diagram does not seem to requirea strong z -dependence of the parameter σ int µ : for instance, a nearby sample gives σ int = 0 . ± .
02, to be comparedwith the value σ int = 0 . ± .
02 obtained for distant supernovae. On the other hand, as illustrated in Fig. 8, theresults of our computations at z > ∼ . σ lens µ ( z ) = 0 . z . We should also note that this contribution stays below 0 .
12 up to z ∼
2, which makes it perfectlycompatible with observations so far performed.It is remarkable that the (above mentioned) simple linear fit of our curve for σ µ ( z ) at z > ∼ . σ lens µ = (0 . ± . z . Also, such a fit is wellcompatible with the results of [44], namely σ lens µ = (0 . +0 . − . ) z . These two observational estimates of the lensingdispersion, with the relative error bands, are illustrated by the shaded areas of Fig. 8. Our result is also in relativelygood agreement with the simulations carried out in [45] which predicted an effect of 0 . z . Other authors ([46–48]),however, have found no indication of this z -dependence of σ µ , and it is true that such a signal waits for a betterobservational evidence.It is likely that future improvements in the accuracy of SNe data will detect (or disprove) this effect at a higherconfidence level. In this respect, our results for σ µ ( z ) stand out as a challenging prediction , which, we believe, couldrepresent a further significative test of the concordance model. VIII. CONCLUSIONS
Starting from the result of a previous paper [13], where the luminosity-redshift relation has been computed to secondorder in the Poisson gauge, we have proceeded here to the evaluation of the effects of a realistic stochastic backgroundof perturbations on the determination of dark energy parameters. The basic tool we have used is the gauge-invariantlight-cone averaging procedure proposed in [8], applied to different functions of the luminosity distance d L averagedover a constant-redshift surface lying on our past light-cone.As already explained in [14], different functions of d L differ by the sensitivity of their light-cone averages to fluctu-ations. Remarkably, a directly observable variable like the luminosity flux, Φ ∼ d − L turns out to be the least sensitiveto perturbations. Its averaged expression is also the simplest and, fortunately enough, other averages are readilycomputed once the one of the flux is given. Similarly, calculation of the dispersion is straightforward. We remark, incidentally, that our prediction for the Doppler-related dispersion at small z is also consistent with previous findings [49]on the so-called Poissonian peculiar-velocity contribution to σ µ . We have checked that our result for the Doppler contribution (see Fig.8) is well fitted by an inverse power law: σ µ ( z ) ∼ . z − in very good agreement, shape-wise, with the corresponding result in [49].The fact that our prediction is about a factor 1.7 larger than the one of [49] is due, presumably, to the use of somewhat different powerspectra (linear or non-linear, with or without baryons). z Σ Μ FIG. 8: The z -dependence of the total dispersion σ µ is illustrated by the thick solid curve, and it is separated into its “Doppler”part (dashed curve dominant at low z ) and “lensing” part (dashed curve dominant at large z ). The slope of the dispersion inthe lensing-dominated regime is compared with the experimental estimates of Kronborg et al. [43] (dark shaded area), and ofJ¨onsson et al. [44] (light shaded area). As far as modeling inhomogeneities goes, we have used a concordance ΛCDM model with an arbitrary Ω Λ0 , uponwhich we have added a realistic spectrum of stochastic scalar perturbations up to second order (and we have alsoexplained why, to this order, we have no contribution from vector and tensor perturbations). The perturbations havebeen taken as those originating from a quasi-scale-invariant primordial spectrum with a realistic transfer function(including baryons and the relative Silk-damping) after it undergoes a non-linear evolution according to the so-calledHaloFit model for structure formation. This model appears to agree well with numerical N -body simulations as wellas with large-scale structure data.Our main conclusions (already succinctly presented in [1] for the case of a perturbation spectrum computed inthe linear regime) are that the effect of perturbations on the averaged flux are extremely small, typically of order10 − at z ∼ O (1). Thus the average flux stands out as an extremely safe observable for determining dark-energyparameters using the simplest FLRW geometry. Such observable is also practically insensitive (see Fig. 6) to theshort-distance behaviour of the power spectrum. Other variables (like d L and the commonly used distance modulus µ ) receive corrections that are typically from two to three orders of magnitude larger (at large values of the redshift z ), but still small-enough for allowing dark-energy measurements at the percent level without invoking theoreticalcorrections. On the other hand, they are more sensitive to the chosen value of the UV cutoff, and to the correctionsarising from the power spectrum in the non-linear regime. The enhanced bias in these flux-related variables is simplydue to the scatter in the flux at fixed z (compare Eqs. (3.5) and (3.6) with (3.8)). It is an effect that can (and needsto) be taken into account in the analyses of SNe data. The absence of this enhanced bias for the flux itself confirmsit as the best variable for all observational purposes.We find, however, that the predicted intrinsic dispersion (or scatter) of the data due to just stochastic inhomo-geneities –and not to other well-known sources of dispersion– is considerably larger than their effect on averages. Theyimply that data should fluctuate in a band which, at large redshifts, covers the FLRW luminosity curves correspondingto a spread in Ω Λ0 of nearly 10% (see Fig. 7). For limited statistics this irreducible dispersion will limit the precisionwith which dark-energy parameters can be extracted from the data. Particularly interesting, at large redshift, is thescatter due to lensing. Such an effect has been observed and a linear phenomenological fit to σ µ has been proposed[43, 44] with a slope dσ µ /dz ∼ .
05 but with large errors. Our theoretical prediction is well described (for the consid-ered range of z ) by the linear behaviour σ µ ( z ) ∼ . z which is not only consistent with the above phenomenologicalfits but also provides an interesting test of the concordance model if and when a more precise determination of σ lens µ will become available. Also at small redshifts our (Doppler-induced) scatter, obeying an approximate inverse powerlaw σ µ ( z ) ∼ . z − , looks compatible with observations and with previous theoretical estimates [49]. As a resultof both effects we find the intrinsic dispersion σ µ ( z ) to have a minimum of about 0 .
016 at z ∼ . ACKNOWLEDGMENTS
IBD would like to thank Eric Switzer and Pascal Vaudrevange for fruitful discussions. GM would like to thank R.Durrer, E. Di Dio and V. Marra for useful discussions. FN wishes to thank J. Guy and D. Hardin for interestingdiscussions about dispersion in the SNe Ia data. FN and GV would like to thank R. Scoccimarro for interestingdiscussions about non-linear spectra and the HaloFit model.IBD would like to acknowledge the hospitality of the Hebrew University where part of this work was carried out.GM has enjoyed the hospitality of the Department of Mathematics at Rhodes University and of the Astrophysics,Cosmology and Gravity Centre of the University of Cape Town during the completion of this work. GV and MGacknowledge the hospitality of the University of Geneva during the last stages of this work.The research of IBD is supported by the German Science Foundation (DFG) within the Collaborative ResearchCenter 676 ”Particles, Strings and the Early Universe”. GM is supported by the Marie Curie IEF, Project NeBRiC -“Non-linear effects and backreaction in classical and quantum cosmology”.
Appendix A. Second-order vector and tensor perturbations
We have already stressed in Sect. II B that vector and tensor perturbations automatically appear, at second order,sourced by the squared first-order perturbation terms. Hence, vector and tensor perturbations must be includedin a consistent second-order computation of the luminosity distance, even if their contributions is negligible at firstorder (as expected, in particular, for a background of super-horizon perturbations generated by a phase of slow-rollinflation).Working in the Poisson gauge, and moving to spherical coordinates x i = ( r, θ, φ ), we can rewrite the relevant partof the PG metric (2.11) as follows: ds P G = a (cid:2) − dη + 2 v i dηdx i (cid:3) + a [( γ ) ij + χ ij ] dx i dx j , (A.1)so that g µνP G ( η, r, θ a ) = a − (cid:32) − v i v j γ ij − χ ij (cid:33) , (A.2)where γ ij = diag(1 , r − , r − sin − θ ), and where we have called v i and χ ij the vector and tensor perturbations writtenin spherical polar coordinates. They satisfy the conditions ∇ i v i = 0 = ∇ i χ ij and γ ij χ ij = 0, where ∇ i is the covariantgradient of three-dimensional Euclidean space in spherical coordinates.Following the same procedure as in the scalar case we can now evaluate the vector and tensor contributions to thecoordinate transformation connecting Poisson and GLC gauge, and then express the perturbed GLC metric, up tosecond order, including the vector and tensor variables v i and χ ij . Such a detailed computation has already beenperformed, and its results presented in [13]. For the purpose of this paper it will be enough to recall here the vectorand tensor contributions to the coordinate transformation between θ and (cid:101) θ :˜ θ a = ˜ θ a (0) + ˜ θ a (2) = θ a + 12 (cid:90) η − η + dx (cid:32) ˆ v a ( η + , x, θ a ) − ˆ χ ra ( η + , x, θ a ) + ˆ γ ab ( η + , x, θ a ) (cid:90) xη + dy ∂ b ˆ α r ( η + , y, θ a ) (cid:33) , (A.3)6and to the 2 × γ ab appearing in the GLC metric: a ( η ) γ ab = γ ab − χ ab ++ (cid:34) γ ac (cid:90) η − η + dx ∂ c (cid:32) ˆ v b ( η + , x, θ a ) − ˆ χ rb ( η + , x, θ a ) + ˆ γ bd ( η + , x, θ a ) (cid:90) xη + dy ∂ d ˆ α r ( η + , y, θ a ) (cid:33) + ( a ↔ b ) (cid:35) , (A.4)where α r ≡ ( v r / − ( χ rr / v i and χ ij , and compute d L according to Eqs. (2.4) and (2.5), we find that the right-handside of Eq. (2.13) has to be modified by the addition of a new term, δ (2) V,T ( z s , (cid:101) θ a ), representing the effect of the vectorand tensor part of the perturbed geometry (see [13] for its explicit expression). No modification is induced, however,on the corresponding equation for I φ ( z s ) controlling the light-cone average of d − L , so that Eq. (2.14) holds even inthe presence of vector and tensor perturbations.The sought average, in fact, is proportional to the proper area of the deformed two-sphere Σ( w o , z s ), and is given by I φ ∼ (cid:82) d (cid:101) θ √ γ (see Eq. (2.9)). Considering the vector and tensor contributions to γ − = det γ ab (obtained from Eq.(A.4)), and computing from Eq. (A.3) the Jacobian determinant | ∂ (cid:101) θ/∂θ | , we can express I φ as an angular integralover the two-sphere with unperturbed measure d Ω = sin θdθdφ . In that case many terms cancel among each other,and we end up with the result: I φ ( w o , z s ) − π (cid:90) π sin θdθ (cid:90) π dφ f ( η, r, θ, φ ) , (A.5)where the integrand f ( η, r, θ, φ ) is a simple expression proportional to the components of the vector and tensorperturbations.The above angular integrals are all identically vanishing, as we can check by expanding the perturbations in Fouriermodes v ik and χ ijk . For each mode (cid:126)k we can choose, without loss of generality, the x axis of our coordinate systemaligned along the direction of (cid:126)k . Considering, for instance, tensor perturbations, we can then write the most generalperturbed line-element, in Cartesian coordinates (omitting, for simplicity, the Fourier index), as follows: h ij dx i dx j = h + ( η, x )( dx dx − dx dx ) + 2 h × ( η, x ) dx dx (A.6)(we have called h + and h × , as usual, the two independent polarization modes). After transforming to sphericalcoordinates, using the standard definitions x = r sin θ cos φ , x = r sin θ sin φ , x = r cos θ , we easily obtain: χ rr = sin θ ( h + cos 2 φ + h × sin 2 φ ) ; χ rθ = sin 2 θ r ( h + cos 2 φ + h × sin 2 φ ) ; χ rφ = 1 r ( − h + sin 2 φ + h × cos 2 φ ) ; χ θθ = cos θr ( h + cos 2 φ + h × sin 2 φ ) ; χ θφ = cos θr sin θ ( − h + sin 2 φ + h × cos 2 φ ) ; χ φφ = − r sin θ ( h + cos 2 φ + h × sin 2 φ ) . (A.7)Since h + = h + ( η, r cos θ ), h × = h × ( η, r cos θ ), all perturbation components depend on φ only through cos 2 φ or sin 2 φ ,so that their contribution averages to zero when inserted into Eq. (A.5). The same is true for the case of vectorperturbations, with the only difference that the φ dependence of v i , in spherical coordinates, is through cos φ or sin φ (corresponding to waves of helicity one instead of helicity two as in the tensor case). Also the vector contributionthus averages to zero when inserted into Eq. (A.5). Appendix B. Computation of the C ( T (1 , i ) spectral coefficients of (cid:104)I , (cid:105) We give here the result for the C ( T (1 , i ) spectral coefficients of (cid:104)I , (cid:105) computed in the CDM case. We haveintroduced the convenient notation l = k ∆ η , and we have enclosed in a box the leading contributions. Finally, we7have defined Sinc( l ) = sin( l ) /l . C ( T (1 , ) = 0 , C ( T (1 , ) = Ξ s f s − f o ∆ η l , C ( T (1 , ) = 0 , C ( T (1 , ) = 2Ξ s f s ∆ η l s l (cid:8) − l + ( l −
2) cos( l ) + l sin( l ) + l SinInt( l ) (cid:9) +4Ξ s f s ∆ η (cos l − l ) + l SinInt( l )) , C ( T (1 , ) = − s f o ∆ η l , C ( T (1 , ) = 4Ξ s [1 − Sinc( l )] − s f s ∆ η (cos( l ) − Sinc( l )) , C ( T (1 , ) = − s (cid:40) f s ∆ η l − f s ∆ η (cos( l ) − Sinc( l )) − f s ∆ η H o ∆ η [2 cos( l ) + ( l − l )] (cid:41) , C ( T (1 , ) = 2Ξ s H s f s H o ∆ η (cid:110) H o ∆ η l + 2[ − H o ∆ η + ( l − l ) − H o ∆ η )( − l )Sinc( l ) (cid:111) , C ( T (1 , ) = 0 , C ( T (1 , ) = 0 , C ( T (1 , ) = 8Ξ s (cid:26) − l + (cid:20) l − (cid:21) cos( l ) −
16 Sinc( l ) − l l ) (cid:27) , C ( T (1 , ) = (cid:20) Ξ s − H s ∆ η (cid:18) − H (cid:48) s H s (cid:19)(cid:21) (cid:26) − Sinc( l )] + 2 f o − f s ∆ η [ cos ( l ) − Sinc( l )]+ f o + f s ∆ η l − f o f s ∆ η [2 cos( l ) + ( l − l )] (cid:41) , C ( T (1 , ) = 4 (cid:26) − Sinc( l ) + f o ∆ η [ cos ( l ) − Sinc( l )] (cid:27) , C ( T (1 , ) = 2Ξ s (cid:104) f o − f s ∆ η (cos l − Sinc( l )) + f o η l − f o f s ∆ η (cid:0) l ) + ( l − l ) (cid:1) + (cid:18) f o + f s ∆ η (cid:19) (1 − Sinc( l )) (cid:105) , C ( T (1 , ) = 4 Ξ s f o + f s ∆ η { − Sinc( l ) } − l SinInt( l ) , C ( T (1 , ) = 8 l SinInt( l ) , C ( T (1 , ) = 8 l {− l ) + l SinInt( l ) } , C ( T (1 , ) = − H s ∆ η (cid:40) − cos( l ) − Sinc( l ) − l SinInt( l ) + 2 f s ∆ η (cid:32) l − cos( l ) − l SinInt( l ) (cid:33) + 2 f o ∆ η [1 − Sinc( l )] (cid:41) , C ( T (1 , ) = 43 l (cid:8) − l + ( l −
2) cos( l ) + l sin( l ) + l SinInt( l ) (cid:9) , C ( T (1 , ) = − H s ∆ η (cid:40) − cos( l ) + Sinc( l ) + f s η l − f o ∆ η [2 cos( l ) + ( l − l )] (cid:41) , C ( T (1 , ) = 4 [1 − Sinc( l )] , C ( T (1 , ) = − l (cid:2) − l ) + (4 + l ) cos( l ) + l sin( l ) + l (6 + l )SinInt( l ) (cid:3) , C ( T (1 , ) = 0 . [1] I. Ben-Dayan, M. Gasperini, G. Marozzi, F. Nugier and G. Veneziano, Phys. Rev. Lett. , 021301 (2013).[2] Y. Wang, Astrophys. J. , 531 (2000).[3] J. A. Frieman, astro-ph/9608068.[4] D. E. Holz, astro-ph/9806124.[5] R. B. Metcalf, Mon. Not. Roy. Astron. Soc. , 746 (1999).[6] C. Porciani and P. Madau, Astrophys. J. , 522 (2001).[7] K. Bolejko, C. Clarkson, R. Maartens, D. Bacon, N. Meures and E. Beynon, Phys. Rev. Lett. , 021302 (2013).[8] M. Gasperini, G. Marozzi, F. Nugier and G. Veneziano, JCAP , 008 (2011).[9] E. Bertschinger, Cosmological Dynamics in Cosmology and Large Scale Structure, proc. Les Houches Summer School,Session LX, ed. R. Schaeffer, J. Silk, M. Spiro and J. Zinn-Justin (Amsterdam:Elsevier Science), 273-347, 1996;L. Bombelli, W. E. Couch and R. J. Torrence, Class. Quant. Grav. , 139 (1994). [10] M. Gasperini, G. Marozzi and G. Veneziano, JCAP , 011 (2009).[11] M. Gasperini, G. Marozzi and G. Veneziano, JCAP , 009 (2010).[12] G. Marozzi, JCAP (2011) 012.[13] I. Ben-Dayan, G. Marozzi, F. Nugier and G. Veneziano, JCAP , 045 (2012).[14] I. Ben-Dayan, M. Gasperini, G. Marozzi, F. Nugier and G. Veneziano, JCAP , 036 (2012).[15] D. J. Eisenstein, W. Hu, Astrophys. J. , 605 (1998).[16] R. E. Smith et al. [Virgo Consortium Collaboration], Mon. Not. Roy. Astron. Soc. (2003) 1311.[17] R. Takahashi, M. Sato, T. Nishimichi, A. Taruya and M. Oguri, Astrophys. J. , 152 (2012).[18] R. Maartens, PhD Thesis, University of Cape Town, South Africa (1980).[19] G. F. R. Ellis, S. D. Nel, R. Maartens, W. R. Stoeger and A. P. Whitman, Phys. Rep. , 315 (1985).[20] C. Clarkson and R. Maartens, Class. Quantum Grav. , 124008 (2010).[21] I. M. H. Etherington, Phil. Mag. , 761 (1933).[22] J. Kristian and R. K. Sachs, Ap. J. , 379 (1966).[23] G. Fanizza, M. Gasperini, G. Marozzi and G. Veneziano, arXiv:1308.4935 [astro-ph.CO].[24] P. Schneider, J. Ehlers and E. E. Falco, Gravitational Lenses (Springer-Verlag editor, 1992).[25] N. Bartolo, S. Matarrese and A. Riotto, JCAP , 010 (2006).[26] S. Matarrese, S. Mollerach and M. Bruni, Phys. Rev. D , 043504 (1998).[27] R. Durrer, The Cosmic Microwave Background (Cambridge University Press, Cambridge, UK 2008).[28] M. Gasperini,
Elements of String Cosmology (Cambridge University Press, Cambridge, UK 2007).[29] C. Bonvin, R. Durrer and M. A. Gasparini, Phys. Rev. D , 023523 (2006) [Erratum-ibid. D , 029901 (2012)].[30] N. Li and D. J. Schwarz, Phys. Rev. D , 083531 (2008).[31] C. Clarkson, K. Ananda and J. Larena, Phys. Rev. D 80 , 083525 (2009).[32] C. Clarkson and O. Umeh, Class. Quantum Grav. , 164010 (2011).[33] D. J. Eisenstein, [astro-ph/9709054].[34] P. Peter and J.-P. Uzan, Primordial cosmology (Oxford Univ. Press, England, 2009).[35] G. Marozzi and G. P. Vacca, Class. Quant. Grav. (2012) 115007.[36] E. Komatsu et al. , Astrophys. J. Suppl. , 18 (2011); G. Hinshaw, et al. , arXiv:1212.5226 [astro-ph.CO].[37] A. Lewis, A. Challinor and A. Lasenby, Astrophys. J. , 473 (2000).[38] M. C. March, R. Trotta, P. Berkes, G. D. Starkman and P. M. Vaudrevange, Mon. Not. Roy. Astron. Soc. , 2308(2011).[39] A. Conley, J. Guy, M. Sullivan, N. Regnault, P. Astier, C. Balland, S. Basa and R. G. Carlberg et al. , Astrophys. J. Suppl. (2011) 1.[40] R. G. Vishwakarma and J. V. Narlikar, Res. Astron. Astrophys. , 1195 (2010).[41] J. Guy, P. Astier, S. Baumont, D. Hardin, R. Pain, N. Regnault, S. Basa and R. G. Carlberg et al. , Astron. Astrophys. (2007) 11.[42] M. Kowalski et al. [Supernova Cosmology Project Collaboration], Astrophys. J. , 749 (2008).[43] T. Kronborg, D. Hardin, J. Guy, P. Astier, C. Balland, S. Basa, R. G. Carlberg and A. Conley et al. , arXiv:1002.1249[astro-ph.CO].[44] J. Jonsson, M. Sullivan, I. Hook, S. Basa, R. Carlberg, A. Conley, D. Fouchez and D. A. Howell et al. , arXiv:1002.1374[astro-ph.CO].[45] D. E. Holz and E. V. Linder, Astrophys. J. (2005) 678.[46] N. V. Karpenka, M. C. March, F. Feroz and M. P. Hobson, arXiv:1207.3708 [astro-ph.CO].[47] L. L. R. Williams and J. Song, Mon. Not. Roy. Astron. Soc. (2004) 1387.[48] B. Menard and N. Dalal, Mon. Not. Roy. Astron. Soc. (2005) 101.[49] L. Hui and P. B. Greene, Phys. Rev. D 73 , 123526 (2006).[50] E. W. Kolb, Class. Quant. Grav.28