Perturbed redshifts from N-body simulations

In order to keep pace with the increasing data quality of astronomical surveys the observed source redshift has to be modeled beyond the well-known Doppler contribution. In this letter I want to examine the gauge issue that is often glossed over when one assigns a perturbed redshift to simulated data generated with a Newtonian N-body code. A careful analysis reveals the presence of a correction term that has so far been neglected. It is roughly proportional to the observed length scale divided by the Hubble scale and therefore suppressed inside the horizon. However, on gigaparsec scales it can be comparable to the gravitational redshift and hence amounts to an important relativistic effect.


I. INTRODUCTION
While the standard analysis of redshift space perturbations [1] only contains the leading Doppler term, recent and future galaxy surveys have or will have sufficient statistical power to detect subleading terms such as gravitational redshift, weak lensing, transverse Doppler shift, time delays etc. These include general relativistic effects that are interesting targets for testing gravity at cosmological scales. For instance, a first detection of gravitational redshift in galaxy clusters has been claimed in [2]. Perturbation theory can be used to understand these effects analytically (e.g. [3]), and many results have been derived in longitudinal gauge. However, their application to N-body simulations requires due care with respect to how the data are mapped to the relativistic spacetime.
One possibility is to run relativistic simulations directly in longitudinal gauge [4], but the use of Newtonian N-body codes is by far the more common practice. Coming from the Newtonian world it is not immediately clear how their results should be interpreted in a relativistic context (e.g. [5][6][7]). In [8] it was finally realized that one can specify a gauge, the so-called "N-body gauge," 1 in which the relativistic equations are formally identical to the Newtonian ones, therefore providing a framework for a fully relativistic interpretation of Newtonian simulations. The same authors went on to analyze initial conditions in this framework, showing that the usual recipes are consistent with general relativity as well [9].
The issue has so far mostly been discussed in relation to the matter power spectrum, and the implications for the perturbed redshift and some other observables have therefore not yet been fully appreciated. In the recent literature the redshift formula of longitudinal gauge is of- * Electronic address: julian.adamek@obspm.fr 1 In fact, [7] had already discovered "half of" the necessary coordinate transformation, leading to what they called "Newtonian matter gauge." The spatial part of the transformation is however different from the N-body gauge, leading to a non-vanishing volume perturbation.
ten directly applied to N-body simulation data, thereby disregarding the fact that these are not provided in the appropriate coordinate system. For instance, [10,11] do not account for this aspect, but I will show that the correction is fortunately very small on the scales they are interested in. While the gauge issue is noted in [12] where a correction term is applied to the density perturbation (see their equation 17), the redshift of individual sources is still computed without such a correction (see their equations 9 and 10). In the following I wish to clarify this issue by studying the coordinate transformations involved, and finally by deriving the correct formulas for the perturbed redshift in N-body gauge. I show that a small correction term due to the spatial coordinate transformation appears and should in principle be included in the analysis. However, this term only becomes relevant on extremely large scales.

II. WEAK-FIELD METRIC
Astrophysical objects with high compactness exist on scales 0.01 parsec (the largest known supermassive black holes), while on extremely large scales 100 megaparsec the Universe can be described entirely in terms of linear equations. In between these two extremes there lies a vast range of scales where the distribution of matter can be very inhomogeneous but the gravitational fields are weak, and the geometry is therefore only weakly perturbed. This empirical fact is confirmed every time we point a telescope at the sky and observe that most rays seem to propagate along almost straight paths.
It is therefore possible to describe the dynamics on those scales in terms of the nonlinear evolution of matter on a geometry with linearly perturbed metric. Of course, as for any geometry, there also exist many coordinate systems in which the metric takes a wildly nonlinear form. A particularly well-known coordinate system in which the smallness of the geometric perturbation is carried into effect is the one of the longitudinal gauge. For the purpose of a more general discussion I adopt the notation introduced in [13] and write the generic line ele-ment for a metric with linear scalar perturbations on top of a Friedmann-Lemaître model as where a(τ ) is the scale factor, τ and x i are conformal time and comoving coordinates on the spacelike hypersurface, respectively, and A, B, H L , H T are scalar functions describing perturbations. Vector and tensor perturbations shall not be discussed here. A setting is said to be "weakfield" if there exists a coordinate system for which all the above perturbation variables are small, i.e. |A|, |H L |, One can then find other such coordinate systems by making a small change of coordinates, generated by two scalar fields T and L, so that τ → τ +T and x i → x i +∇ i L. This shows that the four scalar perturbations only contain two physical modes. The longitudinal gauge is obtained by choosing T and L such that B = H T = 0, leaving only the lapse perturbation A = Ψ and the volume perturbation H L = Φ, where the two potentials Ψ and Φ denote precisely the two physical modes when written as firstorder gauge-invariant expressions. The typical amplitude of these perturbations is ∼ 10 −5 in our Universe, except for the vicinity of black holes or neutron stars where the weak-field description breaks down.
In longitudinal gauge the Hamiltonian constraint reads where H is the conformal Hubble rate, a prime denotes partial derivative with respect to τ , and the equation is linearized in Φ and Ψ but not in δρ. In fact, the matter perturbation δρ can be very large, but it is computed and evolved on a linearly perturbed geometry. This equation is however not the one used in a Newtonian N-body code. First, such a code uses counting densities that are not corrected for perturbations of the volume element -Newtonian theory assumes Euclidean geometry. Second, the Poisson equation lacks some of the terms featured on the left-hand side. From now on I restrict the discussion to the case where gravitational fields are sourced exclusively by nonrelativistic matter. In this case anisotropic stress can be neglected, implying Φ = −Ψ.

III. PERTURBED REDSHIFT IN LONGITUDINAL GAUGE
The effect of geometry (and perturbations thereof) on observables can be understood by studying the geodesics of photons that reach the observation event. Let me denote the tangent vector of such a geodesic as k µ . The condition g µν k µ k ν = 0 implies that where n i is the unit vector (δ ij n i n j = 1) pointing in the direction the photon is traveling. The geodesic equation can then be expressed as an evolution equation for the energy, and an equation describing the deflection of the ray, A reference clock is specified through a unit timelike vec- is the peculiar (coordinate) velocity of the clock's rest frame. In such a frame, the measured photon energy is −g µν u µ k ν . A first integral of eq. (3) yields the following expression for the observed redshift between source (src) and observer (obs): Here dχ is a conformal distance element along the photon path. The redshift perturbations are easily identified as the Doppler shift due to peculiar motion, the gravitational redshift due to time dilation, and the Rees-Sciama effect (also known as integrated Sachs-Wolfe effect in the context of linear theory). The boundary terms in this expression have to be evaluated at the coordinate time at which the photon geodesic actually intersects the world lines of source and observer. The time of flight of the photon is affected by the Shapiro delay, which for a coordinate distance χ between source and observer is given by the following relation:

IV. FROM LONGITUDINAL TO N-BODY GAUGE
Given a linearly perturbed geometry specified in longitudinal gauge I now make a small change of coordinates to set H L = 0. Considering the Lie derivative of the metric tensor the required transformation is generated by T , L that satisfy I furthermore choose L ′ = 0 such that velocities are not transformed. In this new coordinate system one finds B = T and H T = −L. It is important to verify that these new perturbations do not violate the weak-field conditions 2 . As explained below, L shall be chosen such that ∆L is of order Ψ, and hence the same is true for the perturbations generated by H T . One can then easily convince oneself that ∇ i T is of the order of a peculiar velocity -in fact, it coincides with the Zel'dovich approximation thereof -and therefore √ ∇ i B∇ i B ∼ v 10 −3 at low redshift. The shift perturbation in the new coordinate system is therefore substantially larger than Ψ, but still comfortably within the weak-field regime. The lapse perturbation becomes and as explained shortly I arrange that it vanishes at leading order. The Hamiltonian constraint becomes and as a consequence of H L = 0 the density perturbation δρ can be obtained simply by counting the mass elements per coordinate volume, in accordance with the procedure relevant for Newtonian codes. For nonrelativistic particles the geodesic equation reduces to Since I assume that matter is nonrelativistic and hence the pressure perturbation can be neglected, the spatial trace of Einstein's equations yields This shows that an appropriate choice of boundary conditions will set A = 0. Furthermore, with such a choice eqs. (9), (10) are formally identical to the ones of Newtonian gravity if eqs. (7), (8) are used, Keeping in mind that the above equations remain valid even if δρ/ρ becomes large I now want to set the boundary conditions at early times when matter perturbations are still linear. According to eq. (8) the condition A = 0 is satisfied when 2 For instance, had I made a change of coordinates that sets A = B = 0 instead, like in a synchronous gauge, I would have found L ′ = T = −a −1´a (Ψ/H)dã and that H L receives a contribution ∼ ∆Ψ/H 2 ∼ δρ/ρ. The volume perturbation therefore does not remain small everywhere in such a coordinate system, rendering the weak-field treatment inconsistent.
The linear solution of Ψ is constant in matter domination, and with eq. (7) the corresponding choice of L is given by ∆L = 5Ψ in . Here I introduce Ψ in to denote the linear initial condition for Ψ in matter domination. With this choice one can see that ∇ i B = v i in the linear regime, which (together with H L = 0) is the original gauge condition used in [8] for the N-body gauge. So even though my gauge condition H ′ T = −L ′ = 0 is different, the resulting coordinate system is the same whenever only nonrelativistic matter is present. The advantage of my condition is that it does not explicitly refer to a matter perturbation and can hence be easily extended into the nonlinear regime of matter as long as gravitational fields remain weak. In a slight abuse of terminology I shall therefore always call this system of coordinates the one of N-body gauge, as the state of a Newtonian N-body simulation is given in precisely these coordinates.
I shall now discuss the repercussions of this change of coordinates from longitudinal to N-body gauge. The null condition is solved in N-body gauge by

and the photon geodesic equation can be written as
and Considering how the coordinate transformation acts on u µ one sees that u µ = a −1 (1, v i ) in the new coordinates. Thus, a first integral of eq. (14) gives the following new expression for the observed redshift: In order to recover all the terms of eq. (5) I used (HB) ′ = Ψ ′ and HB + B ′ = −Ψ, but evidently a new boundary term HB appears. Noting that B = T this boundary term can be understood as the result of the change of coordinates acting on the background term a obs /a src . In other words, the term has to appear because the equaltime hypersurfaces in N-body gauge do not coincide with the ones of longitudinal gauge. The coordinate time of a Newtonian N-body simulation is the one of N-body gauge, and hence this boundary term needs to be taken into account.
Let me now inspect the time of flight for the photon in N-body gauge, where I again use the gauge conditions to recover the terms known from longitudinal gauge. Compared to eq. (6) there are two new boundary terms. These are expected from the gauge transformation, since the coordinate time transforms as τ → τ + T = τ + B, and the coordinate distance transforms as For a photon trajectory with fixed endpoints, the coordinate time of emission and observation therefore transforms such that the boundary terms due to the shift perturbation in eqs. (16) and (17) cancel exactly. The other boundary term in eq. (17) gives precisely the change in the coordinate distance due to the spatial transformation between longitudinal and N-body gauge. Therefore, the perturbed redshift for the trajectory remains invariant.

V. DISCUSSION
The precedent analysis clarifies that in order to use the longitudinal gauge for computing the perturbed redshifts with N-body simulation data one should transform the coordinates appropriately. In particular, the coordinate distance between observer and sources changes according to a spatial transformation that is independent of time. This has already been pointed out in [5], and I explicitly show how to recover this result in the relativistic framework provided by the N-body gauge. Alternatively the computation of the perturbed redshift can also be carried out directly in N-body gauge. In this case the coordinates of sources are directly taken from the simulation, and the effect appears as a modification of the Shapiro delay.
In order to estimate the amplitude of the correction, let me compute the typical change δχ = n i ∇ i L| obs − n i ∇ i L| src of the coordinate distance. Using the relation ∆L = 5Ψ in the variance of δχ is given by where ∆ Ψ k is the dimensionless power spectrum of Ψ in . Unfortunately the integral has an infrared divergence for nearly scale invariant spectra, but in practice this divergence is regulated by the finite size of a simulation. Imposing a cutoff close to the Hubble scale one finds that δχ/χ ∼ 10 −4 almost independent of scale χ and precise value of the cutoff. Considering how a change in the coordinate distance affects the time of flight one sees that the typical correction to the redshift due to this coordinate effect is δz/(1 + z) ∼ Hδχ. Therefore the effect becomes of the order of the gravitational redshift for trajectories at or above the gigaparsec. At these extreme scales the gauge correction is of the same order as all other relevant terms and should be taken into account.
As suggested in [2] the gravitational redshift may be measured statistically by looking for "excess redshift" of the brightest galaxies at the center of clusters when compared to the fainter galaxies in the outskirts. For simulating such a measurement the relevant correction is given by the difference in time of flight between the different sources, and therefore the effect is suppressed by the small ratio between the scale of the cluster and the Hubble scale. One expects that in this case the correction is typically less than 1% of the signal and can therefore safely be neglected.