A fitting formula for the non-Gaussian contribution to the lensing power spectrum covariance

A fitting formula for the non-Gaussian contribution to the lensing power spectrum covariance

J. Pielorz Argelander-Institut für Astronomie (AIfA), Universität Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
11email: pielorz@astro.uni-bonn.de
   J. Rödiger Argelander-Institut für Astronomie (AIfA), Universität Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
11email: pielorz@astro.uni-bonn.de
   I.Tereno Argelander-Institut für Astronomie (AIfA), Universität Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
11email: pielorz@astro.uni-bonn.de
   P. Schneider Argelander-Institut für Astronomie (AIfA), Universität Bonn, Auf dem Hügel 71, 53121 Bonn, Germany
11email: pielorz@astro.uni-bonn.de
Received; accepted
Key Words.:
gravitational lensing – Methods: -body simulations – Cosmology: theory – large-scale structure of the Universe

Context:Weak gravitational lensing is one of the most promising tools to investigate the equation-of-state of dark energy. In order to obtain reliable parameter estimations for current and future experiments, a good theoretical understanding of dark matter clustering is essential. Of particular interest is the statistical precision to which weak lensing observables, such as cosmic shear correlation functions, can be determined.

Aims:We construct a fitting formula for the non-Gaussian part of the covariance of the lensing power spectrum. The Gaussian contribution to the covariance, which is proportional to the lensing power spectrum squared, and optionally shape noise can be included easily by adding their contributions.

Methods:Starting from a canonical estimator for the dimensionless lensing power spectrum, we model first the covariance in the halo model approach including all four halo terms for one fiducial cosmology and then fit two polynomials to the expression found. On large scales, we use a first-order polynomial in the wave-numbers and dimensionless power spectra that goes asymptotically towards for , i.e., the result for the non-Gaussian part of the covariance using tree-level perturbation theory. On the other hand, for small scales we employ a second-order polynomial in the dimensionless power spectra for the fit.

Results:We obtain a fitting formula for the non-Gaussian contribution of the convergence power spectrum covariance that is accurate to for the off-diagonal elements, and to for the diagonal elements, in the range and can be used for single source redshifts in WMAP5-like cosmologies.


1 Introduction

Weak gravitational lensing by the large-scale structure, or cosmic shear, is an important tool to probe the mass distribution in the Universe and to estimate cosmological parameters. The constraints it provides are independent and complementary to those found by other cosmological probes such as cosmic microwave background (CMB) anisotropies, supernovae (SN) type Ia, baryon acoustic oscillations (BAO) or galaxy redshift surveys. The cosmic shear field quantifies the distortion of faint galaxy images that is induced by continuous light deflections caused by the large-scale structure in the Universe (e.g., Bartelmann & Schneider 2001; Schneider 2006). Since this effect is too small to be measured for a single galaxy, large surveys with millions of galaxies are required to detect it in a statistical way. The cosmic shear signal has been successfully measured in various surveys, since the first detections of Bacon et al. (2000); Kaiser et al. (2000); Van Waerbeke et al. (2000); Wittman et al. (2000). Most recently, shear two-point correlation functions were measured in the Canada-France-Hawaii Telescope Legacy Survey (CFHTLS) and were used to constrain the amplitude of dark matter clustering, , with uncertainty (Fu et al. 2008).

The next generation of galaxy surveys will greatly improve the precision with which weak lensing effects can be measured (Albrecht et al. 2006) enabling us to obtain, with accurate redshift information and tomographic measurements, precise constraints on the evolution of dark energy. However, the expected improvement in future data only leads to a significant improvement of the precision and accuracy of the cosmological interpretation if the systematic errors, the underlying physics, and the statistical precision of cosmic shear estimators are well understood. Systematics currently identified arise mainly from non-cosmological sources of shear correlations, i.e., intrinsic alignments of galaxies (e.g., Schaefer 2008, for a review), and biases on the shear measurement (Massey et al. 2007; Semboloni et al. 2008). This paper addresses the issue of the statistical precision of cosmic shear estimators, determined by the covariance of the estimator. Since much of the scales probed by cosmic shear lie in the non-linear regime, being affected by non-linear clustering, the covariance depends on non-Gaussian effects and has a non-Gaussian, as well as a Gaussian, contribution. Indeed, even though the non-Gaussianity of the shear field is weaker than that of the matter field due to the projection along the line-of-sight, various studies indicate that the non-Gaussian contribution to the covariance cannot be neglected when constraining cosmological parameters with weak lensing (Scoccimarro et al. 1999; White & Hu 2000; Cooray & Hu 2001; Kilbinger & Schneider 2005; Semboloni et al. 2007; Takada & Jain 2009).

Most cosmic shear results are based on the measurement of two-point correlation functions of the shear field. Since, in general, the number of independent measurements is insufficient to infer the complete covariance directly from observations, one may derive it from ray-tracing maps of numerical -body simulations. This, however, requires a large number of realizations and, in addition, is very time-consuming if an exploration of the covariance in the parameter space is needed. An alternative is to compute the covariance with an analytic approach. For shear two-point correlation functions, Schneider et al. (2002) derived an expression for the Gaussian contribution to the covariance. Semboloni et al. (2007) fitted the ratio between that expression and a covariance computed with -body simulations, containing both Gaussian and non-Gaussian contributions, providing thus a formula to compute the total covariance from the Gaussian term. In Fourier space, large-scale modes are independent and, differently from real space, the Gaussian contribution to the covariance of the convergence power spectrum (i.e., of the Fourier transform of the two-point shear correlation function) is diagonal and can be computed from the convergence power spectrum alone (Kaiser 1992; Joachimi et al. 2008), whereas the non-Gaussian contribution can be computed from the trispectrum of the convergence (Scoccimarro et al. 1999). The trispectrum, on large scales, can be accurately derived in tree-level perturbation theory111We refer by tree-level perturbation theory to the lowest, non-vanishing order of the considered quantity in perturbation theory., and, on small scales, is well represented by the one-halo term of a halo model approach. A non-Gaussian part of the covariance consisting of a perturbation theory term and a one-halo term was used, e.g., in Takada & Jain (2009).

This paper aims at producing an accurate expression for the non-Gaussian contribution of the covariance of the convergence power spectrum that is fast to compute, contributing thus to accurate estimates of cosmological parameters. Following Scoccimarro et al. (1999) and Cooray & Hu (2001), we start from a canonical estimator of the dimensionless convergence power spectrum and use it to derive an analytic expression for the corresponding covariance. The various spectra involved are evaluated using the halo model approach of dark matter clustering (Seljak 2000; Ma & Fry 2000; Scoccimarro et al. 2001; Cooray & Sheth 2002). The halo model approach assumes that all dark matter in the Universe is bound in spherical halos, and uses results from numerical -body simulations to characterize halo properties such as their profile, abundance and clustering behavior.

The evaluation of the covariance of the convergence power spectrum in the halo model approach is time-consuming. In addition, it may be needed to repeat it for different cosmological models for the purpose of parameter estimation. To allow for a faster computation, we construct a fitting formula for the non-Gaussian part of the convergence power spectrum covariance. On small scales, we fit the halo model result with a polynomial in the non-linear dimensionless convergence power spectrum. On large scales, we fit the ratio between the halo model covariance and the perturbation theory covariance. We stress that it is a fit to the halo model covariance, not involving a covariance computed from -body simulations. The result is, however, calibrated by -body simulations, since they determine the halo model parameters.

The paper is organized as follows. We define in Sect. 2 the reference cosmology, considering the growth of matter perturbations. We introduce the convergence spectra, construct an estimator for the dimensionless convergence power spectrum, and derive an expression for its covariance in Sect. 3. In Sect. 4, we describe the halo model approach, and compute the covariance of the power spectrum. The covariance depends on the values of halo model parameters, which are also defined here. It also depends on the power spectrum, bispectrum and trispectrum of the correlations of halo centers. Expressions for these spectra, in the framework of perturbation theory, are given in the Appendix. Section 5 tests the accuracy of the halo model predictions, for both the power spectrum and its covariance, against two sets of ray-tracing simulations. Section 6 presents the fitting formula for the non-Gaussian contribution to the covariance where its coefficients are given as function of source redshift. We conclude in Sect. 7.

2 Structure formation in a CDM cosmology

Throughout this work we assume a spatially flat cold dark matter model with a cosmological constant (), as supported by the latest 5-year data release of WMAP results (Komatsu et al. 2009). The expansion rate of the Universe, , in such models is described by the Friedmann equation , where is the Hubble constant, denotes the combined contributions from dark matter and baryons today in terms of the critical density , and is the density parameter of the cosmological constant. The comoving distance to a source at is then


where the scale factor is related to the redshift via the relation using the convention today.

In structure formation, the central quantity is the Fourier transform of the density contrast , which describes the relative deviation of the local matter density to the comoving average density of the Universe at time . We suppress the time dependence of in the following. In this way, the mean density contrast is by definition zero, and we can describe matter perturbations in the early Universe as zero-mean Gaussian random fields. In this case, the statistical properties of the Fourier transformed density field,


are completely characterized by the power spectrum


where is the ensemble average and denotes the Dirac delta distribution. Note that throughout this paper the tilde symbol is used to denote the Fourier transform of the corresponding quantity.

In linear perturbation theory, which is valid on large scales, the power spectrum at a scale factor is characterized by


where the amplitude is normalized in terms of , denotes the spectral index of the primordial power spectrum, and is the transfer function. Note that all Fourier modes of the matter density grow at the same rate, i.e., , where


is the growth factor which we normalize as . In the non-linear regime, i.e., on small scales, different Fourier modes couple and the Gaussian assumption cannot be maintained. Thus we have to consider higher-order moments of the density field to describe its statistical properties. In perturbation theory, it is possible to find analytic expressions for these moments, which hold up to the quasi-linear regime. In Appendix A, we derive the expressions for the bispectrum and trispectrum in tree-level perturbation theory, which are the Fourier transforms of the three- and four-point-correlation functions, respectively.

3 Covariance of the convergence power spectrum

A central quantity in weak lensing applications is the two-dimensional projection of the density contrast on the sky, which is known as effective convergence . It is obtained by projecting the density contrast along the backward-directed light-cone of the observer according to


where denotes the redshift-dependent comoving distance, is the comoving distance to the horizon and the weight function takes into account the distribution of source galaxies along the line-of-sight. We assume for simplicity that all background sources are situated at a single comoving distance , such that the weight function has the form


where denotes the Heaviside step function. To take advantage of the Fourier properties, we analyze the statistical properties of the Fourier counterpart of . For the theoretical consideration of the convergence power spectrum covariance we need the second- and the fourth-order moments, as will become apparent later. They are defined by


where the subscript ‘c’ refers to the connected part of the corresponding moment and is a sum of Fourier wave-vectors. The convergence power spectrum and trispectrum are calculated using the flat-sky and Limber’s approximation (Kaiser 1998; Scoccimarro et al. 1999; Bernardeau et al. 2002):


where and are the corresponding three-dimensional matter power spectrum and trispectrum (Fourier transform of the four-point correlation function).

We are interested in estimating the dimensionless convergence power spectrum


and the corresponding covariance for wave-vectors of different length . A natural choice for the estimator of the dimensionless convergence power spectrum is (Scoccimarro et al. 1999; Cooray & Hu 2001; Takada & Bridle 2007)


which is unbiased in the limit of infinitesimal small bin sizes, since . Here, denotes the solid angle of a survey with a fractional sky coverage of , and the integration is performed over the Fourier modes lying in the annulus defined by , where is the width of the -th bin. We denote the integration area formed by the annulus as .

The evaluation of the covariance of the estimator, Eq. (13), results in an expression of the form


where is given by Eq. (12), and denotes the Kronecker delta. The second term in square brackets is the bin-averaged convergence trispectrum,


where is the non-linear convergence trispectrum as defined in Eq. (11). To derive Eq. (14), we made use of the definitions of the convergence power spectrum and trispectrum in Eqs. (8) and (9), and of the discrete limit of the delta distribution . The derived expression consists of two terms: a Gaussian part , which scales as the convergence power spectrum squared and only contributes to the diagonal of the covariance matrix (see, e.g. Joachimi et al. 2008), and a non-Gaussian part , which scales as the dimensionless bin-averaged convergence trispectrum and introduces correlations between the wave-vectors of different bins (see, e.g. Scoccimarro et al. 1999; Cooray & Hu 2001). Both terms are inversely proportional to the survey area , but have a different behavior with respect to the bin width . While the Gaussian term decreases with increasing bin size, the non-Gaussian term is independent of the binning, since the bin area cancels out after the integration.

We can analytically perform one of the integrations of the bin-averaged trispectrum in Eq. (15). First, note that it only depends on the parallelogram configuration of the convergence trispectrum, i.e., setting , and in Eq. (11). Also, if we choose an appropriate coordinate system for the integration over the wave-vectors, the problem becomes symmetric under rotations and we can parametrize the convergence trispectrum by the length of the two sides of the parallelogram and and the angle between them, . Hence, we define


Making use of the symmetry properties of this problem, one angular integration becomes trivial and the integration in Eq. (15) simplifies to


If the bin-width is sufficiently small (), the integration area is , and we can make use of the mean value theorem. In this way, we approximate the integral in Eq. (17) by an angular average


Note that if is independent of the angle between and , an approximation of the covariance can be calculated without having to perform an integration at all. In particular, this is the case for the 1-halo term of the three-dimensional matter trispectrum, as we will see later in Eq. (42).

4 Halo Model

We have seen that the covariance of the dimensionless convergence power spectrum estimator consists of two terms: a Gaussian part, which is proportional to the dimensionless convergence power spectrum squared and a non-Gaussian part, which is the bin-averaged dimensionless convergence trispectrum (see Eq. 14). We will compute these terms using the halo model approach (Seljak 2000; Ma & Fry 2000; Scoccimarro et al. 2001; see also the comprehensive review by Cooray & Sheth 2002).

4.1 Overview

With the assumption that all dark matter is bound in spherically-symmetric, virialized halos, the halo model provides a way to calculate the three-dimensional polyspectra of dark matter in the non-linear regime. In Sect. 4.5 below, we summarize the equations one obtains for the dark matter power spectrum and trispectrum.

In the halo model description, the density field at an arbitrary position in space is given as a superposition of all halo density profiles such that


where denotes the density profile of the -th halo with center of mass at and is the normalized profile. By parametrizing the halo profile in this way, we assume that the shape of the -th halo depends only on the halo mass and the halo concentration parameter , which we define below. The dark matter polyspectra of the density field follow then from taking the ensemble averages, , of products of the density at different points in space. Assuming that the number of halos is , where is the average number density of halos and the considered volume, we compute the ensemble averages by integrating over the joint probability density function (PDF) for the halos that form the field (Smith & Watts 2005), i.e.,


where denotes the PDF. If one considers that position and mass of a single halo are independent random variables, the PDF factorizes as


where is the halo mass function and is the concentration probability distribution for halos given a mass .

4.2 Ingredients

The halo model approach provides a scale-dependent description of the statistical properties of the large-scale structure. On small scales, the correlation of dark matter is governed by the mass profiles of the halos, whereas on large scales the clustering between different halos determines the nature of the correlation. As there are a multitude of models to describe the behavior on different scales, and an even larger number of parameters one has to set judiciously, there exists no such thing as a unique halo model. In order to have reproducible results, it is therefore necessary to specify ones choice of parameters. For this work, we will adopt the following parameters for the halo model:

  1. The average mass of a halo is defined as the mass within a sphere of virial radius as , where denotes the overdensity of the virialized halo with respect to the average comoving mass density in the Universe. Typically, values for are derived in the framework of the non-linear spherical collapse model (e.g. Gunn & Gott 1972). Expressions valid for different cosmologies are summarized in Nakamura & Suto (1997). In our implementation, we use the results which are valid for a flat CDM-Universe, i.e.,


    where . We find for our fiducial WMAP5-like cosmology .

  2. -body simulations suggest that the density profile of a halo follows a universal function. We choose to use the NFW profile (Navarro et al. 1997), which is in good agreement with numerical results and has an analytical Fourier transform. It is given by


    where is the amplitude of the density profile and characterizes the scale at which the slope of the density profile changes. For small scales () the profile scales with , whereas for large scales it behaves as . The Fourier transform of the NFW profile is


    where , we truncated the integration at in the second step, and introduced the concentration parameter in the third step. Additionally, we use for the sine- and cosine integrals the definitions

  3. The abundance of halos of mass at a redshift is given by


    where we introduced the dimensionless variable ,


    where denotes the redshift-dependent growth factor, is the smoothed variance of the density contrast, and denotes the value of a spherical overdensity that collapses at a redshift as calculated from linear perturbation theory. In our work, we use the expression from Nakamura & Suto (1997), which is valid for a CDM Universe,


    where . The quantity has only a weak dependence on redshift and we find for our fiducial model.

    The advantage of introducing is that part of the mass function can be expressed by the multiplicity function , which has a universal shape, i.e., is independent of cosmological parameters and redshift. In this work, we employ the Sheth and Tormen mass function (Sheth & Tormen 1999)


    which is an improvement over the original Press-Schechter formulation (Press & Schechter 1974). We use the parameter values , , and amplitude , which follows from mass conservation.

  4. The concentration parameter characterizes the form of the halo profile. From -body simulations one finds that the average, , depends on the halo mass (Bullock et al. 2001) like


    where is the characteristic mass defined within the Press-Schechter formalism as . In the following, we will use the values and as proposed by Takada & Jain (2003). This implies that more massive halos are less centrally concentrated than less massive ones. However, results from numerical -body simulations (Jing 2000; Bullock et al. 2001) indicate that there is a significant scatter in the concentration parameter for halos of the same mass. Furthermore, Jing (2000) proposes that such a concentration distribution can be described by a log-normal distribution


    Typical values for the concentration dispersion range from to (Jing 2000; Wechsler et al. 2002). Note that the width of the distribution is independent of the halo mass. The variation of the halo concentration can be attributed to the different merger histories of the halos (Wechsler et al. 2002). We will analyze the impact of this effect on different spectra in Sect. 4.7. When we use only the mean concentration parameter, we have to replace the probability distribution of the concentration, needed for example in Eq. (21), by a Dirac delta distribution

  5. On large scales, the correlation of the dark matter density field is governed by the spatial distribution of halos. Since the clustering behavior of halos and matter density differ, one introduces the bias factors such that


    In this way, the halo density contrast, , is expressed as a Taylor expansion of the matter density contrast, . The bias parameters are in general derived based on the Sheth-Tormen mass function introduced above. For the linear halo bias one obtains then


    where and match the values used in the mass function. Expressions for higher-order bias factors can be found, e.g., in Scoccimarro et al. (2001). Since they only have a small impact on the quantities employed here, we take into account only the first-order bias. In Fourier space we may then write

  6. To obtain the final correlation function, one has to perform integrations along the halo mass and optionally along the halo concentration, with limits formally extending from to . In practice, we use the mass limits and . Masses smaller than give no significant contribution to the considered quantities, while, due to the exponential cut-off in mass, masses larger than are rare. For the concentration, we employ the integration limits and .

  7. Due to the cut-off in mass, the consistency relation (Scoccimarro et al. 2001)


    does not hold. To cure this problem we consider a rescaled linear bias such that , where is the result of the integral in Eq. (36). In this way, one ensures that the halo term with the largest contribution to the correlation equals the perturbation theory expression on large scales (see Fig. 1).

4.3 Building Blocks

Using the ingredients described in the previous section, it is possible to define building blocks, which simplify significantly the notation for expressing the polyspectra (Cooray & Hu 2001):


In the case , we additionally define , for consistency.

4.4 Power spectrum

We can now compute the power spectrum from Eq. (20). The result consists of two terms, the 1-halo and the 2-halo terms, (Seljak 2000). They are given by


where denotes the linear perturbation theory power spectrum, defined in Eq. (4), and we use the ingredients summarized earlier-on. Note that, for convenience, we omit the redshift-dependence in the notation. Using the building blocks from Eq. (37), these terms can be written in the following compact form :


The 1-halo term, , denotes correlations in space between two points in the same halo, whereas the 2-halo term, , takes into account correlations between two different halos. Hence, the 1-halo term is dominant on small scales and the 2-halo term is dominant on large scales. Note that the 2-halo term converges to the linear power spectrum on large scales because of the consistency relation of the first-order halo bias factor (see Eq. 36) and of the limit for .

4.5 Trispectrum

We compute now the dark matter trispectrum in the halo model approach. As discussed in Sect. 3, only parallelogram configurations of the trispectrum wave-vectors contribute to the covariance of the convergence power spectrum. Restricting our calculations to these configurations, we obtain four different halo term contributions,


which are simpler than in the general case. In addition, we neglect terms involving higher-order halo bias factors since, on most scales, they provide only a small correction (Ma & Fry 2000; Takada & Jain 2003). We further note that the perturbative expansion of halo centers, used in the calculations, was shown to become inaccurate on non-linear scales (Smith et al. 2007). The contributions to the trispectrum take the following forms, using the compact notation of the building blocks (see Cooray & Sheth 2002 for the expression of the halo model trispectrum including higher-order bias factors): The 1-halo term, dominant on the smallest scales, is


The 2-halo term has two contributions, , consisting of a term , which corresponds to correlations of three points within one halo and a fourth point in a second halo, and a term , which describes correlations involving two points in a first halo and the other two points in the second halo. They read,


The 3-halo term is given by


Finally, the 4-halo term, dominant on large scales, is


and describes correlations of points distributed in four different halos. Note that, like for the power spectrum, the 2-halo term is computed from the linear power spectrum. On the other hand, the 3- and 4-halo terms depend on and , respectively, which are the lowest-order, non-vanishing, perturbation theory contributions to the bispectrum and trispectrum. Both spectra are derived in Appendix A.2.

4.6 Convergence spectra

Figure 1: Square configuration of the dimensionless convergence trispectrum against wave-number for the range for . The upper panel displays the four individual halo terms, the sum of all four trispectrum halo terms (solid line) and the tree-level perturbation theory trispectrum , as indicated in the key. The lower panel shows the ratio between the indicated contributions and the complete trispectrum. The double dashed line illustrates the corresponding ratio of the approximation . Note that we consider the 4-halo term only in the upper panel since it resembles the term on large scales.

The convergence power spectrum and trispectrum, needed to evaluate the covariance in Eq. (14), are computed by projecting and , according to Eqs. (10) and (11). Fig. 1 shows the dimensionless convergence trispectrum, defined as for a square configuration where all wave-numbers have a length , where denotes the contribution from the corresponding halo term. The plot shows the individual contributions to the dimensionless (the projections of Eqs. 42-46), illustrating on which scales the individual terms are important, as well as the total contribution of all halo terms . We show, in addition, the dimensionless projected , which closely follows the 4-halo term. We see that the commonly used approximation is accurate for large wave-numbers () but has a deviation of about from the complete trispectrum for small wave-numbers ().

4.7 Stochastic halo concentration

The previous results were computed using the deterministic concentration-mass relation of Eq. (30). We now analyze the impact of scatter in the halo concentration parameter on the covariance of the convergence power spectrum, using the stochastic concentration relation given by the the log-normal concentration distribution of Eq. (31).

Cooray & Hu (2001), analyzed the effect of a stochastic concentration on the three-dimensional power spectrum and trispectrum and found that the behavior of the corresponding 1-halo terms were increasingly sensitive to the width of the concentration distribution for smaller wave-numbers . Furthermore, the effect of a stochastic concentration relation was more pronounced for the trispectrum than for the power spectrum, since the tail of the concentration distribution is weighted more strongly in higher-order statistics.

Performing the same analysis for the projected power spectrum and trispectrum, we find a similar trend as in the three-dimensional case, but with a smaller sensitivity to the concentration width of the distribution on small scales. For , we find, in the case of the 1-halo term of the power spectrum, a deviation from a deterministic concentration relation of about for wave-numbers larger than , whereas, for the 1-halo term of the trispectrum, the deviation is of the order of in the same -range. Thus, when considering the covariance of the convergence power spectrum, one should take into account the concentration dispersion in the 1-halo term of the trispectrum but can safely neglect it for the power spectrum. Additionally, we find that a stochastic concentration has only a small impact on the 2-halo terms of the power spectrum and trispectrum (Pielorz 2008).

From this analysis, we expect the effect of a concentration distribution to be the strongest on the non-Gaussian part of the covariance, which depends on the trispectrum. To directly infer the impact of a concentration distribution on the covariance, we calculate the 1-halo contribution to the non-Gaussian covariance, i.e., we perform the bin averaging of Eq. (15),


for two different concentration dispersions , using our halo model implementation. Fig. 2 shows the ratio


where denotes the deterministic concentration relation, for two values in two contour plots. In agreement with the previous results, the largest impact of the concentration dispersion occurs for wave-numbers larger than . For , the deviation of the covariance from the original deterministic concentration relation is small, of about . The effect becomes non-negligible for , where we find a deviation of (for ) to (for ).

Figure 2: Contour plots of the ratio between the non-Gaussian covariance 1-halo term contribution (see Eq. 48) computed with a concentration dispersion and with a deterministic concentration, against wave-numbers ranging from to . The left panel is for concentration dispersion , whereas in the right panel .

5 Comparison with -body simulations

In Sect. 4, we computed the power spectrum and trispectrum of the dark matter fluctuations. Projecting them according to Eqs. (10) and (11), and inserting the result in Eqs. (12), (14), and (15), we obtained the covariance of the convergence power spectrum estimator.

To test the accuracy of the halo model predictions for the statistics of the dark matter density field, we compare the dimensionless convergence power spectrum and the corresponding covariance, calculated in the halo model approach, with results from two different ray-tracing simulations.

Virgo 0.3 0.7 0.7 0.0 0.9 1.0 0.21 1 (2) BBKS
Gems 0.25 0.75 0.7 0.04 0.78 1.0 0.14 1 (2) EH
Table 1: Cosmological parameters used to set up the initial power spectrum, which determines how the simulation particles are distributed initially. The simulations employ either the BBKS (Bardeen et al. 1986) or the EH (Eisenstein & Hu 1998) transfer function. Convergence maps with source galaxies situated at and were produced for both sets of simulations.
Virgo 141.3 200 0.25
Gems 150.0 220 16.00
Table 2: Parameters used for generating the two -body simulations considered in this paper and for producing the resulting convergence maps with the multiple-lens-plane ray-tracing algorithm (see e.g., Jain et al. 2000; Hilbert et al. 2008). From left to right these are: the side length, , of the cubic simulation box, the number of particles, , used for the simulation, their mass, , and the number of available convergence maps, , with area .

5.1 Virgo and Gems simulation

For our comparison, we chose one simulation from Jenkins et al. (1998) and ten simulations from Hartlap et al. (2009), which we denote in the following as Virgo and Gems simulation, respectively. The Virgo simulation was carried out in 1997 by the Virgo-Consortium for a CDM cosmology (see Tab. 1) with particles in a periodic box of comoving side length (see Tab. 2). It uses the PP-/PM-code HYDRA, which places subgrids of higher resolution in strongly clustered regions. Structures on scales larger than can be considered as well resolved. The Gems simulations were set up in cubic volumes of comoving side length with particles (see Tab. 2). The cosmology chosen reflects the WMAP5 results (Komatsu et al. 2009) and thus has a lower value for than the Virgo simulation. It uses the GADGET-2 code to simulate the evolution of dark matter particles (Springel 2005) and has a softening length of .

5.2 Ray-tracing

The output of numerical simulations are three-dimensional distributions of particles in cubic boxes over a range of redshift values. In order to compare the results with the predicted convergence power spectrum from the halo model, we make use of the multiple-lens-plane ray-tracing algorithm (see e.g., Jain et al. 2000; Hilbert et al. 2008) to construct effective convergence maps. The basic idea is to introduce a series of lens planes perpendicular to the central line-of-sight of the observer’s backward light cone. In this way, the matter distribution within the light cone is sliced and can be projected onto the corresponding lens plane. By computing the deflection of light rays and its derivatives at each lens plane, one simulates the photon trajectory from the observer to the source by keeping track of the distortions of ray bundles. In this way, the continuous deflection of light rays is approximated by a finite number of deflections at the lens planes. As a result, one obtains the Jacobian matrix for the lens mapping from source to observer and can construct convergence maps.

For both simulations, a similar number of around 200 effective convergence maps were produced, with source galaxies situated at a single redshift of either or (see Tab. 2). The maps produced with the Gems simulation have an area of 16 , while the ones from Virgo are much smaller, with 0.25 .

Figure 3: Dimensionless convergence power spectrum against wave-number for galaxy sources at redshift (upper panels) and (lower panels). Points with error bars show measurements from the two sets of numerical simulations: Gems (left plots) and Virgo (right plots). Both simulations use a similar particle setup but differ in the area of the -maps (see Tabs. 1 and 2). They are compared with the corresponding results obtained with fitting formulae from Smith et al. (long-dashed line), Peacock-Dodds (short-dashed line), and the halo model predictions for a deterministic halo concentration (solid line). See text for a discussion.

5.3 Convergence power spectrum

To test the accuracy of the halo model approach in describing the non-linear evolution of dark matter, we compare the dimensionless projected power spectrum, computed in the halo model approach, to the ones estimated from the numerical -body simulations.

The dimensionless convergence power spectra of the simulations are measured from the real-space two-dimensional convergence maps of length and grid-size . For this, we first apply a Fast Fourier Transform222For the FFT we use an algorithm from the GNU scientific library (see http://www.gnu.org/software/gsl/ for more details). to obtain on each grid-point. Then, we estimate the power spectrum at a wave-number for the -th convergence map by averaging over all Fourier modes in the band }, i.e.,


where denotes the number of modes in each band of width . Finally, the ensemble average of the dimensionless power spectrum is obtained by averaging over the results of the various convergence maps, i.e., . The error bars for the power spectrum estimate are computed from the dispersion over the convergence maps used, and are due to sample variance. Since modes with a length similar to the side length of the convergence map are only poorly represented, the sampling variance is largest for small . The effect is stronger for the Virgo simulation than for the Gems simulation, since the Gems convergence maps have a much larger area .

Fig. 3 shows a good agreement between the halo model and the -body simulation convergence power spectra. We also include, in Fig. 3, the convergence power spectra computed using the fitting formulae of Peacock & Dodds (1996) and Smith et al. (2003). Both, the halo model and the fitting formulae, show a better agreement with simulations for the lower source redshift case. On small scales, the halo model and the Smith et al. fitting formula agree well with the simulations, whereas the Peacock-Dodds fitting formula has too little power, in particular in the case of the Gems simulation. On intermediate scales (), the halo model is less accurate than the Smith et al. fitting formula, suggesting that the halo model suffers from the halo exclusion problem on these scales, as described, e.g., in Tinker et al. (2005). This means that, while in simulations halos are never separated by distances smaller than the sum of their virial radii, this is not accounted for in the framework of our halo model and is probably the cause for the deviation. The good agreement of the Smith et al. formula is not surprising, since it is based on simulations of similar resolution than the ones we consider here. In contrast, a similar comparison using the convergence power spectrum estimated with the Millennium Run (Hilbert et al. 2008) clearly favors the halo model prediction over of the two fitting formulae, with both fitting formulae strongly underestimating the power on intermediate and small scales.

The good overall accuracy of the halo model results was expected since its ingredients, such as the mass function and the halo profile, were obtained from -body simulations. We note that, in this analysis, we used a deterministic concentration parameter, since the use of a stochastic one has only a small effect on the small scales of the convergence power spectrum.

Figure 4: Relative difference (see Eq. 51) between the halo model prediction for the covariance of the convergence power spectrum and the results from the Virgo simulation () against wave-numbers . The binning scheme is linear with a constant bin width of ranging from to . The left panel displays the full covariance, including all four halo terms for the trispectrum, and uses a deterministic concentration-mass relation denoted by . The right panel illustrates the same covariance but with a stochastic concentration distribution of width for the 1-halo term of the trispectrum.
Figure 5: Relative difference (see Eq. 51) between the halo model prediction for the covariance of the convergence power spectrum and the results from the Gems simulation () against wave-numbers . The binning scheme is linear with a constant bin width of ranging from to . The left panel displays the full covariance, including all four halo terms for the trispectrum, and uses a deterministic concentration-mass relation denoted by . The right panel illustrates the same covariance but with a stochastic concentration distribution of width for the 1-halo term of the trispectrum.

5.4 Covariance of the convergence power spectrum

The similarity between simulation and halo model power spectra was to some extent expected. The ability to make an accurate description of higher-order correlations provides a stronger test of the halo model. Due to its important role in calculating the error of the power spectrum and for parameter estimates, we focus in this section on the accuracy of the covariance of the dimensionless convergence power spectrum.

We need an appropriate estimator for the power spectrum covariance of the simulations. As each simulation provides different -maps, we have realizations of the power spectrum. From these we estimate the covariance by applying the unbiased sample covariance estimator which has for our purpose the form:


where is the projected power spectrum estimate of the -th effective convergence map at a wave-number (see Eq. 49). We evaluate for both, Virgo and Gems simulations, for the case .

The halo model results, , are calculated as described in Sect. 4 and include all terms of the non-Gaussian covariance. In their computation, we use our fiducial halo model with the ingredients summarized in Sect. 4.2, and use the cosmological parameters values corresponding to each simulation, given in Tab. 1, for the case . We consider both, a deterministic concentration-mass relation and a stochastic one, with dispersion for the 1-halo term of the trispectrum.

We compare the halo model with the simulations covariance matrices considering their relative deviation


We note that although the number of available convergence maps for our simulations () is too small to obtain a percentage level accuracy for the covariance estimate (Takahashi et al. 2009), we assume that the covariance from the simulations is the reference one, and put it in the denominator of Eq. (51). The comparisons are displayed in Figs. 4 and 5 for the Virgo and Gems simulation, respectively. In the case of the Virgo simulation with , the best agreement is a relative deviation of found on intermediate scales between and . For small scales the halo model underestimates the simulation by or more. On large scales, i.e., for wave-numbers , the halo model overestimates the simulation by around . However, since the simulation covariances suffer from a large sampling variance due to the small size of the convergence maps, the comparison is not meaningful on these scales. The agreement improves for , in this case there is a larger range of scales where the deviation is small. For the Gems simulation the agreement is much better. There is now an interval of scales () where the deviation is in the range . Throughout all the scales probed, the deviation on the diagonal of the covariances is around . However, for off-diagonal components, at small scales, the agreement is still poor. Finally, the improvement of considering is less significant than in the Virgo case.

Although the halo model predictions strongly deviate from the simulation estimates of the covariance on small scales, this analysis does not necessarily imply a poor accuracy of those predictions. Indeed, the simulation covariances estimated in this analysis have a strong scatter. In particular, in the case of the Gems simulation, we found from bootstrap subsamples of 50 convergence maps that the resulting covariances can deviate up to from the average covariance of the complete sample (Pielorz 2008). This is in agreement with the results by Takahashi et al. (2009) who need to use 5000 simulations to obtain an estimate of the matter power spectrum covariance at a sub-percent level accuracy.

There are however very recent indications that Eq. (14) indeed underestimates the covariance of the convergence power spectrum on small scales (Sato et al. 2009). In particular, sample variance in the number of halos in a finite field is not accounted for by Eq. (14). Indeed, the mass function yields a mean number density, but there are fluctuations on the number of halos due to the large-scale mass fluctuations. Sample variance in the number of clusters in a volume-limited survey (Hu & Kravtsov 2003) has been considered in cluster abundance studies (e.g., Vikhlinin et al. 2009). In the halo model framework, the sample variance in the number of halos was derived in Takada & Bridle (2007) and its contribution to the covariance of the convergence power spectrum was found, in Sato et al. (2009), to boost the non-Gaussian errors of a 25 survey by one order of magnitude on scales . The increase is reduced for larger survey areas.

6 Fitting formula for the covariance of the convergence power spectrum

Figure 6: Ratio