Testing the conditional mass function of dark matter halos against numerical N-body simulations

TL;DR

This paper evaluates theoretical models of the conditional mass function of dark matter halos against N-body simulations, identifying modifications needed for normalization and providing an analytical bias function for underdense regions.

Contribution

It compares two rescaling methods for the CMF, introduces a normalization modification, and offers an analytical bias function for void studies.

Findings

Neither rescaling method alone provides a normalized CMF.

A modified rescaling improves the match with simulations.

The analytical bias function aids in void statistics computations.

Abstract

We compare the predicted conditional mass function (CMF) of dark matter halos from two theoretical prescriptions against numerical N-body simulations, both in overdense and underdense regions and at different Eulerian scales ranging from $5$ to $30\,h^{-1}\,$Mpc. In particular, we consider in detail a locally-implemented rescaling of the unconditional mass function (UMF) already discussed in the literature, and also a generalization of the standard rescaling method described in the extended Press-Schechter formalism. First, we test the consistency of these two rescalings by verifying the normalization of the CMF at different scales, and showing that none of the proposed cases provides a normalized CMF. In order to satisfy the normalization condition, we include a modification in the rescaling procedure. After this modification, the resulting CMF generally provides a better description of numerical results. We finally present an analytical fit to the ratio between the CMF and the UMF (also known as the matter-to-halo bias function) in underdense regions, which could be of special interest to speed-up the computation of the halo abundance when studying void statistics. In this case, the CMF prescription based on the locally-implemented rescaling provides a slightly better description of the numerical results when compared to the standard rescaling.

Tables

Table 1. Table 1: Values of cosmological parameters for our reference model.

${\Omega }_m$	${\Omega }_b$	${\Omega }_{\Lambda }$	$h$	$n_s$	${\sigma }_8$
0.3	$4.6\times {10}^{-2}$	0.7	0.7	0.96	0.8

Table 2. Table 2: Table summarizing the properties of the two simulations used in this work.

Name	Size	Particle mass	Softening
	$\left(h^{-1}\text{Mpc}\right)$	($h^{-1}{\text{M}}_{\odot }$)	($h^{-1}\text{kpc}$)
LSS	1500	$2.7\times {10}^{11}$	60
SSS	250	$1.2\times {10}^9$	10

Table 3. Table 3: Complete list of fitting terms for computing the parameters A 𝐴 A , b 𝑏 b and c 𝑐 c in the reference cosmology, using eq. ( B.1 ). Index i 𝑖 i runs over lines, index j 𝑗 j runs over columns.

$a_{ij}$	0	1	2	3
0	-0.91	-0.15	0.69	-0.13
1	1.08	0.42	-0.59	0.10
2	-0.12	-0.20	0.17	-0.03
3	-0.00	0.02	-0.02	0.00

Table 4. Table 4: Complete list of fitting terms to transform the A 𝐴 A parameter to the chosen values of z 𝑧 z and σ 8 subscript 𝜎 8 \sigma_{8} , using eqs. ( B.1 ) and ( B.1 ). Index i 𝑖 i selects a line, index j 𝑗 j selects a column.

${\xi }_{ij}^A$	0	1	2	3	4
0	-601.87	515.84	-142.68	11.74	0.28
1	480.91	-349.49	59.17	6.42	-1.73
2	-151.77	89.94	-0.73	-7.22	1.04
3	24.26	-12.77	-1.32	1.58	-0.21
4	-1.67	0.95	0.03	-0.09	0.01

Table 5. Table 5: Complete list of fitting terms to transform the b 𝑏 b parameter to the chosen values of z 𝑧 z and σ 8 subscript 𝜎 8 \sigma_{8} , using eqs. ( B.1 ) and ( B.1 ). Index i 𝑖 i selects a line, index j 𝑗 j selects a column.

${\xi }_{ij}^b$	0	1	2	3	4
0	-5098.88	3004.96	-36.68	-248.98	37.40
1	5948.68	-4065.13	526.96	148.63	-29.76
2	-2635.03	2035.06	-435.59	-6.77	7.41
3	517.88	-440.23	120.16	-8.72	-0.48
4	-37.68	34.46	-10.82	1.24	-0.02

Table 6. Table 6: Complete list of fitting terms to transform the c 𝑐 c parameter to the chosen values of z 𝑧 z and σ 8 subscript 𝜎 8 \sigma_{8} , using eqs. ( B.1 ) and ( B.1 ). Index i 𝑖 i selects a line, index j 𝑗 j selects a column.

${\xi }_{ij}^c$	0	1	2	3	4
0	748.68	-72.63	-320.20	131.82	-14.59
1	-1142.51	618.11	48.90	-73.50	9.99
2	574.04	-434.57	82.90	5.72	-2.02
3	-118.26	104.21	-29.84	2.53	0.06
4	8.61	-8.27	2.77	-0.36	0.01

Table 7. Table 7: Complete list of fitting terms for computing the exponent f ​ ( h , Ω m ) 𝑓 ℎ subscript Ω 𝑚 f(h,\Omega_{m}) using eq. ( 35 ), entering the mass scaling ( 34 ). Index i 𝑖 i selects a line, index j 𝑗 j selects a column.

$f_{ij}$	0	1	2	3	4
0	7824.37	25237.25	3138.79	-323056.91	351919.38
1	55357.30	15986.33	-244637.95	-256695.08	598394.06
2	72768.05	-31325.23	-342372.66	74755.55	391854.50
3	20716.17	6647.99	-158499.00	43559.96	161916.81
4	-1549.50	4181.00	6828.17	-51471.20	55152.71

Equations

$$\frac{\text{d}{n}_{\mathrm{u}}}{\text{d}m}=\frac{{\rho }_{\mathrm{m}}}{m}\frac{\text{d}\ln {\sigma }^{-1}}{\text{d}m}f(\sigma ,z),$$

$${\sigma }^2(r,z)=\frac{D^2(z)}{2{\pi }^2}{\int }_0^{\infty }\text{d}k{k}^{2}P(k)w^2(kr),$$

$$m=\frac{4}{3}\pi {\rho }_{\mathrm{m}}{r}^3,$$

$$w(x)=\frac{3}{x^3}(\sin x-x\cos x).$$

$$f_{\text{TA}}(\sigma )=A\left[{\left(\frac{\sigma }{b}\right)}^{-a}+1\right]e^{-c/{\sigma }^2}.$$

$${\delta }_{\mathrm{m}}=\frac{{\rho }_{\mathrm{c}}-{\rho }_{\mathrm{m}}}{{\rho }_{\mathrm{m}}},$$

$${\delta }_{\mathrm{l}}({\delta }_{\mathrm{m}})=\frac{{\delta }_{\mathrm{c}}}{1.68647}[1.68647-\frac{1.35}{(1+{\delta }_{\mathrm{m}}{)}^{2/3}}$$

$$-\frac{1.12431}{(1+{\delta }_{\mathrm{m}}{)}^{1/2}}+\frac{0.78785}{(1+{\delta }_{\mathrm{m}}{)}^{0.58661}}].$$

$$Q=R(1+{\delta }_{\mathrm{m}}{)}^{1/3}.$$

$${\delta }_{\mathrm{c}}$$

$${\sigma }^2$$

$${\sigma }^{\prime 2}=\left[{\sigma }^2-{\sigma }_1^2\right]{\left[\frac{{\delta }_{\mathrm{c}}}{{\delta }_{\mathrm{c}}-{\delta }_1}\right]}^2,$$

$${\delta }_{\mathrm{c}}^{\prime }$$

$${\sigma }_2^{\prime 2}$$

$$D(q,Q_1,Q_2)=\frac{{\sigma }_{12}(q,Q_1,Q_2)}{{\sigma }_1^2(Q_1)},$$

$${\sigma }_{12}$$

$$=\frac{D(z{)}^2}{2\pi }{\int }_0^{\infty }\text{d}k{k}^{2}P(k)w(k{Q}_1)w(k{Q}_2)\frac{\sin (kq)}{kq},$$

$${\sigma }_2^{\prime 2}=\left[{\sigma }_2^2-D(q,Q_1,Q_2{)}^2{\sigma }_1^2\right]{\left[\frac{{\delta }_{\mathrm{c}}}{{\delta }_{\mathrm{c}}-D(q,Q_1,Q_2){\delta }_1}\right]}^2,$$

$$\frac{\text{d}{n}_{\mathrm{c}}^{\text{avg}}}{\text{d}m}(m,z;Q,{\delta }_1)=\frac{3}{Q^3}{\int }_0^Q\text{d}q{q}^2\frac{\text{d}{n}_{\mathrm{c}}}{\mathrm{d}m}(m,z;q,Q,{\delta }_1).$$

$${\sigma }_2^{\prime 2}=\left[{\sigma }_2^2-D(q,Q_1,Q_2{)}^2{\sigma }_1^2\right]{\left[\frac{\alpha }{\alpha -D(q,Q_1,Q_2){\delta }_1}\right]}^2,$$

$$\alpha =1.25$$

$$n_{\mathrm{c}}^{\text{sim}}(m_j;R,{\delta }_{\mathrm{l}})=\frac{N_j}{N_{\mathrm{s}}{V}_{\mathrm{R}}},$$

$$n_{\mathrm{c}}^{\text{theo}}(m_j;R,{\delta }_{\mathrm{l}})=\frac{1}{N_{\mathrm{s}}{V}_{\mathrm{R}}}\sum _{i=1}^{N_{\mathrm{s}}}{n}_i(m_j;Q,{\delta }_{\mathrm{l}})V_{\text{Q},i}$$

$${\tilde{n}}_{\mathrm{c}}^{\text{theo}}(m_j;R,{\delta }_{\mathrm{l}})=n_{\mathrm{c}}^{\text{theo}}(m_j;R,{\delta }_{\mathrm{l}})\frac{n_{\mathrm{u}}^{\text{sim}}(m_j)}{n_{\mathrm{u}}^{\text{theo}}(m_j)}.$$

$$1+{\delta }_N(m,z;Q,{\delta }_{\mathrm{l}})=[1+{\delta }_{\mathrm{m}}({\delta }_{\mathrm{l}})][1+{\delta }_{\text{ns}}(m,z;Q,{\delta }_{\mathrm{l}})],$$

$$1+{\delta }_{\text{ns}}(m,z;Q,{\delta }_{\mathrm{l}})=\frac{1}{N_{\mathrm{u}}(m,z)}\left[\frac{3}{Q^3}{\int }_0^Q\text{d}q{q}^2{N}_{\mathrm{c}}(m,z;q,Q,{\delta }_{\mathrm{l}})\right].$$

$$N_{\mathrm{u}}(m,z)={\int }_m^{+\infty }\text{d}{m}^{\prime }\frac{\text{d}{n}_{\mathrm{u}}}{\text{d}{m}^{\prime }}(m^{\prime },z),$$

$$N_{\mathrm{c}}(m,z;q,Q,{\delta }_{\mathrm{l}})={\int }_m^M\text{d}{m}^{\prime }\frac{\text{d}{n}_{\mathrm{c}}}{\text{d}{m}^{\prime }}(m^{\prime },z;q,Q,{\delta }_{\mathrm{l}}).$$

$$1+{\delta }_{\text{ns}}(m,z;Q,{\delta }_{\mathrm{l}})=A(m,z,Q)e^{-b(m,z,Q){\delta }_{\mathrm{l}}^2},$$

$$1+{\delta }_{\text{ns}}(m,z;Q,{\delta }_{\mathrm{l}})=A(m,z,Q)e^{-b(m,z,Q){\delta }_{\mathrm{l}}^2+c(m,z,Q){\delta }_{\mathrm{l}}},$$

$$\frac{\text{d}{n}_{\mathrm{u}}}{\text{d}m}(m)={\int }_{-\infty }^X\text{d}{\delta }_{\mathrm{l}}\frac{\text{d}{n}_{\mathrm{c}}}{\text{d}m}(m;Q,{\delta }_{\mathrm{l}})G({\delta }_{\mathrm{l}};{\sigma }_1),m<m(Q),$$

$${\delta }_{\mathrm{c}}\to {\delta }_{\mathrm{c}}^{\prime }={\delta }_{\mathrm{c}}\left[1-D(q,Q_1,Q_2)\frac{{\delta }_1}{\alpha }\right].$$

Full text

††pubyear: 2017††pagerange: LABEL:firstpage–LABEL:lastpage

1 Introduction

The abundance of collapsed dark matter halos as a function of mass, and its dependence on redshift, is one of the basic tools in modern cosmology. The halo mass function was first estimated analytically, as described in Press & Schechter (1974). This formalism constitutes the basis for any later analytic development; it is based on the assumption that initial fluctuations are drawn from a Gaussian distribution, and on the physics of spherical collapse. Bond et al. (1991) developed the so-called ‘excursion set formalism’, in whose framework the Press-Schechter mass function can be derived from the barrier crossing statistics of many independent, uncorrelated random walks. However, this recipe underestimates the abundance of high mass collapsed objects (Efstathiou et al., 1988) due to the fact that in deriving the analytical expression, all mass elements are assumed to be at the centre of the object they belong to (Betancort-Rijo & Montero-Dorta, 2006a). As a consequence of this, the accumulated mass fraction is not a universal function of the linear variance of the density contrast but depends on the shape of the power spectrum (Betancort-Rijo & Montero-Dorta, 2006b).

The treatment of the ellipsoidal collapse was incorporated into the formalism in several works (Bond & Myers, 1996; Sheth et al., 2001; Sheth & Tormen, 2002) by using a moving barrier. Those studies provided an approximated analytical expression for the distribution of first crossings of the moving barrier, which was used in Sheth & Tormen (1999) to provide a new parametrization of the mass function. Comparisons with N-body simulations show that this new formalism significantly improves the Press-Schechter mass function in predicting the halo mass distribution. The Sheth-Tormen parametrization has been tested against N-body simulations in Jenkins et al. (2001), providing a very good agreement. In the same reference an alternative fitting formula was proposed: this functional form for the mass function is capable of achieving agreement below $\sim$ 10–30 % at all masses, and is said to be universal, in the sense that it can be applied to different cosmological models and different redshifts. However, the fit cannot be extrapolated beyond the mass range employed in that work. Warren et al. (2006) improved the fit providing a mass function parametrization which reproduce the halo abundance with an error $\sim 5$ % at redshift $z=0$ and a fixed cosmology; the work in Reed et al. (2007), instead, provides a good fit for the mass function at redshifts 10–30. Notice that in the works cited so far, the fit is done against catalogues obtained using friends-of-friends (FoF) defined halos (Davis et al., 1985), which is computationally efficient and has the advantage that it does not assume any a-priori shape for the halos.

A major issue concernig masses obtained with FoF finders is that they not easily accommodated within the theoretical formalism, which is usually based on spherical overdensity (SO) mass. The corresponding SO algorithm (Lacey & Cole, 1994) identifies objects as spherical regions enclosing a certain overdensity around density peaks. This task has been addressed in Tinker et al. (2008): the proposed fit for the mass function is capable of reproducing the halo mass function at $z=0$ and in the mass range $10^{11}\,h^{-1}\,\text{M}_{\odot}<M<10^{15}\,h^{-1}\,\text{M}_{\odot}$ with an accuracy $\lesssim 5$ %. They claim the halo mass function cannot be universal at this level of accuracy; the dependence on redshift is explored up to $z=2.5$, and also on the value of the overdensity used in the halo definition, in the range $200$ to $3200$ (computed with respect to the mean background density). The fitting parameters depend on both the redshift and the overdensity, in order to account for the non-universality of the mass function, improving the agreement with N-body simulations. This mass function is the standard reference used in most cluster cosmology analyses, and is the one we shall take as reference in the present work as well.

There have been subsequent attempts at improving Tinker’s fit. The impact of the initial conditions for the simulations used to fit the mass function can be found in Crocce et al. (2006). In Crocce et al. (2010), the halo abundance is fitted over five order of magnitude in mass, and its redshift evolution explored up to $z=2$; the fitting formula is the same as in Warren et al. (2006), but with different parameters in order to correct for the underestimation of the most massive halo abundance, resulting in an overall $2\,\%$ fit accuracy. The cosmology and redshift dependence of the mass function is also explored in Courtin et al. (2011), showing the universality is not verified in absolute terms, holding true only for $z=0$ and a limited set of cosmologies; this is due to some cosmology-dependent features in the mass function, like the value of the threshold density for collapse, whose value depends on cosmological parameters. In Angulo et al. (2012) an accurate fit with errors within $5\,\%$ is provided over eight orders of magnitude for the halo mass, extending Tinker’s fitting range; a discussion on how different halo definitions affect the mass function is also provided. Watson et al. (2013) explore in detail the difference between FoF and SO-defined mass functions; they find agreement below $5\,\%$ error with Tinker’s fit, with major differences in the high mass end, and extend the result to higher mass values. Finally, Despali et al. (2016) provide a fit with $5-8\,\%$ error for different redshifts and cosmologies, provided the mass function is defined in terms of the virial overdensity halo mass; in the same work, scaling relations for the fitting parameters are provided in order to compute the mass function for different redshifts and halo mass definitions.

The inclusion of baryon physics in the simulations may also impact the fits of the halo mass function. For instance, the work in Bocquet et al. (2016) shows how the inclusion of baryonic feedback leads to a decrease in the halo masses, which in turn leads to a decrease of the halo abundance at a given mass; the effect is more important at low redshift and low masses, where differences with respect to dark matter-only fits can grow up to $\sim 15\,\%$. The study of baryonic effects on the mass function in other halo mass ranges has been tackled in Velliscig et al. (2014) and Martizzi et al. (2014).

All the works on the mass function described so far consider the abundance of halos originated from density perturbations that grow during their evolution, eventually decoupling from the background expansion and finally ending up into collapsed objects. These perturbations are defined with respect to a background which is considered homogeneous. Such a mass function is usually referred to as ‘unconditional’. We address in the present work the study of the ‘conditional’ mass function, meaning the function describing the abundance of objects inside an environment (the condition) with a definite density contrast with respect to the background. Locally underdense or overdense regions do alter the mass distribution of collapsed halos (e.g., Sheth & Lemson, 1999); this applies, for instance, when studying the halo abundance inside a supercluster or a cosmic void.

Despite all the effort that has been undertaken to provide an accurate parametrization of the unconditional mass function (UMF), the conditional mass function (CMF) has not been so thoroughly studied and not a definite recipe can be found in the literature. Existing work on the subject usually take the UMF as a starting point and try to adapt it to the conditioned environment. In Cole & Kaiser (1989) the abundance of collapsed halos in perturbed regions is computed starting from the Press-Schecter UMF and using the ‘peak-background split’ formalism; their work is particularly focused on the application of the halo bias to astronomical objects like clusters of galaxies, quasars or the study of the diffuse X-ray background. The spatial distribution of dark matter halos and its relation to the underlying mass distribution has also been studied in Mo & White (1996): the authors develop an analytical model, based on the extended Press-Schechter formalism and the spherical collapse assumption, to describe the halo bias function defined as the ratio between the overdensity of halos in a conditioning region and its matter density contrast. Comparison with N-body simulations confirms the accuracy of this model in the case of an Einstein-de-Sitter universe. Other attempts worked in the framework of the spherical collapse, and extended the Press-Schechter mass function, in the excursion set formalism (Bond et al., 1991; Lacey & Cole, 1993). This approach, however, is not adequate to reproduce results from N-body simulations (e.g., Tormen, 1998; Betancort-Rijo & Montero-Dorta, 2006a; Neyrinck et al., 2014). The work in Sheth & Tormen (2002) extends the excursion set formalism by including the physics of ellipsoidal collapse. Their results reasonably agree with N-body simulations; this formalism, however, is used to quantify the difference in halo abundance when considering the mass function evolution from one redshift to another, and cannot be used for computing the CMF at fixed redshift $z$, which is the problem we address in this work. An extension of the work on the CMF can be found in Rubiño-Martín et al. (2008, hereafter RBP08): among the major changes with respect to previous studies, it introduces the dependence of the CMF from the radial position inside the condition (assumed spherical), and explicitly implements the UMF scaling in such a way that the correspondent CMF is properly normalized over the possible values of the density contrast. The resulting CMF shows good agreement with numerical simulation for the conditions explored in that work.

The present work aims at comparing with N-body simulations the prescriptions for the CMF described in RBP08, both at small and large scales. The results enable us to pinpoint the best recipe for computing the CMF starting from the UMF, depending on the condition scale or the value of the density contrast. We also study the normalization of the CMF at different scales, explicitly modifying the UMF scaling in order to satisfy the normalization condition by construction. In order to simplify the computation of the CMF, we also provide a fitting formula for the bias between the CMF and the UMF in the case of large scale underdense regions, which is needed in the study of void statistics.

Whenever cosmological parameters are not specified, it is intended we are using a flat $\Lambda$CDM model with the values for the cosmological parameters listed in table 1.

2 Formalism for the UMF

The unconditional mass function (UMF), which we label $\text{d}n_{\rm u}/\text{d}m$, is the mass distribution of the number density of halos at a given redshift. Its definition is such that $\text{d}n_{\rm u}/\text{d}m(m,z)\,\text{d}m$, is the (comoving) number density of halos with mass in the range $[m,m+\text{d}m]$ at redshift $z$. The UMF can be parametrized as in e.g., Tinker et al. (2008), in the form:

[TABLE]

where $\rho_{\rm m}$ is the comoving matter density of the Universe and $\sigma^{2}(m,z)$ is the mass variance of the linear density field extrapolated to redshift $z$, the mass $m$ also defining the scale over which the variance is computed. As a function of the matter power spectrum $P(k)$, the variance can be computed as:

[TABLE]

in which $w$ is the Fourier transform of the window function and $D(z)$ is the growth factor of linear perturbation normalized to unity at $z=0$. Both the window function and the relation between the scale $r$ and the mass $m$ are determined by the employed smoothing filter. In the case of a real space top-hat filter :

[TABLE]

and

[TABLE]

The term $f(\sigma,z)$ in equation (1) represents the mass fraction contained in collapsed objects per unit $\ln{\sigma^{-1}}$. It depends on the value of $\delta_{\rm c}$, that is the critical density for collapse, obtained by extrapolating the linear overdensity $\delta_{\rm l}$ of a spherical perturbation at the time it collapses. There exist different functional form for $f(\sigma,z)$ in the literature (Jenkins et al., 2001), coming from analytic or semi-analytic models. The original Press & Schechter (1974) UMF provides an analytical expression based on the physics of spherical collapse; however, this simple recipe actually underestimates the abundance of high mass collapsed objects. Improvements to this UMF employ parametric expressions fitted against N-body simulations, like in Sheth & Tormen (1999) or in Warren et al. (2006). In this paper we will use the UMF parametrization by Tinker et al. (2008), hereafter T08, which explicitly separates a power-law regime at small scale and a high mass exponential cut-off:

[TABLE]

The function is fitted against N-body simulations in which halos are identified with a spherical overdensity method (Lacey & Cole, 1994). The improvement in the fit compared to other UMFs comes from the fitting parameters being allowed to depend on redshift and on the value of overdensity $\Delta$ employed by the halo finder (meaning a certain region of the simulation is considered a collapsed object if its density is greater than the mean matter density by a factor $\Delta$). In T08, tables are provided for computing the fitting parameters at different values of $\Delta$ and $z$. In the present work, whenever we talk about the overdensity $\Delta$, it is understood we are considering a density contrast with respect to the matter background density, and not to the critical density, following the convention adopted in T08.

3 Formalism for the CMF

The conditional mass function (CMF), which we label $\text{d}n_{\rm c}/\text{d}m$, is the proper halo mass distribution inside regions with a definite density contrast with respect to the average matter density, where the predictions from the UMF are no longer adequate. We will always consider spherical regions with matter density contrast $\delta_{\rm m}$, defined as:

[TABLE]

where $\rho_{\rm m}$ is the average background matter density and $\rho_{\rm c}$ the matter density inside the condition. The correspondent linear density contrast $\delta_{\rm l}$ can be computed as (Sheth & Tormen, 2002):

[TABLE]

In the framework of the CMF it is useful to employ (comoving) Lagrangian coordinates rather than (comoving) Eulerian coordinates. In fact, an overdense (underdense) region does not follow the Hubble flow and it actually shrinks (dilate) during the evolution of the Universe, meaning its Eulerian size changes (the Lagrangian one is constant by construction). Hereafter, we will use the variables $R$ and $Q$ to denote the condition Eulerian and Lagrangian radius, respectively. The two are related by:

[TABLE]

The use of Lagrangian coordinates ensures we can talk about a unique value of $Q$ which is valid at all redshifts.

Apart from mass and redshift, the conditional mass function also depends on the condition radius and density contrast, $\text{d}n_{\rm c}/\text{d}m(m,z;Q,\delta_{\rm l})$. There exist different approaches in the literature to derive analytically an expression for the CMF, usually based in the rescaling of some quantities entering the definition of the UMF. We will concentrate in particular on two methods, namely the standard rescaling and one that is implemented locally. It is worth mentioning that other ways of computing the CMF have been developed, based on the rescaling of the cosmological parameters rather than on the mass function itself; these methods are based on the assumption that the condition can be treated as an independent universe with a locally-defined cosmology. For instance, Gottlöber et al. (2003) stress the efficiency of computing the halo abundance in voids by rescaling the matter power spectrum via a redefinition of the parameter $\sigma_{8}$; a prescription for computing the local cosmological parameters in voids is derived also in Goldberg & Vogeley (2004), with the aim of improving the efficiency of specialized void numerical simulations. A complete description of these rescalings goes beyond the scope of this work; for more information we redirect to Betancort-Rijo et al. (prep), where the local cosmology is derived not only for cosmic voids but also for overdense regions.

3.1 Standard rescaling

The most straightforward way to obtain the conditional mass function is by taking the same analytical expression for the UMF and express it in terms of the rescaled variables:

[TABLE]

where $\delta_{1}$ is the condition linear density contrast and $\sigma_{1}^{2}=\sigma^{2}(R)$ is the variance of the matter field at the condition scale. We refer hereafter to eqs. (3.1) as the standard rescaling. This recipe is the basis for a number of works on the CMF. For instance, the so called excursion set formalism, also known as extended Press-Schechter (Bond et al., 1991; Lacey & Cole, 1993), which is based on the hypothesis of spherical collapse. An extension to this work, including the ellipsoidal collapse, can be found in (Sheth & Tormen, 2002).

This rescaling can be applied to any UMF parametrizations. However, care must be taken when rescaling the T08 mass function, which does not show any explicit dependence on the value of the critical density $\delta_{\rm c}$. The parameter can be introduced taking into account that in the previous parametrizations the function $f$ depends on $\delta_{\rm c}$ only via the combination $\nu=(\delta_{\rm c}/\sigma)^{2}$ (Sheth & Tormen, 1999). In the case of T08 UMF we can rewrite the variance as $\sigma^{2}=\delta_{\rm c}^{2}/\nu$, and apply the scaling (3.1) to the variable $\nu$. This yields:

[TABLE]

which can be applied to eq. (5).

3.2 Local rescaling

The local rescaling is an attempt at improving the results obtained with the standard recipe. We shall take as a reference for this work the method described in RBP08, where the CMF is obtained from the UMF by a variable rescaling in the style of eqs. (3.1), but with the rescaling being implemented locally. This means the CMF also depends on the radial distance from the centre of the condition $q$: $\text{d}n_{\rm c}/\text{d}m(m,z;q,Q_{1},\delta_{1})$. Let $Q_{2}$ be the Lagrangian radius corresponding to the mass $m$, and $\sigma_{2}=\sigma(Q_{2})$. In formulae, the rescaling reads:

[TABLE]

where

[TABLE]

and $\sigma_{12}$ is the covariance:

[TABLE]

($D(z)$ denotes the normalized growth factor). The full formalism is described in RBP08. Again, in the case of T08 mass function, the dependence on $\delta_{\rm c}$ must be introduced explicitly, yielding the rescaling equation for the variance:

[TABLE]

which can be applied to eq. (5). In RBP08 it is said that the quantity defined in eq. (12), for $Q_{2}\ll Q_{1}$, is practically independent of $Q_{2}$, a property that can be exploited to speed the computation up. In particular, this allowed computing the mass derivative of the variance appearing in eq. (1) without the need of applying each time the rescaling. However, we found that this approximation is not so good in most of the mass range we employ. In this work we therefore always applied the full computation with no approximations.

In RBP08 the local rescaling was compared to the standard one, showing differences in the predicted halo abundance, particularly in the high-mass end and in underdense regions. In order to get rid of the extra dependence on the radial distance $q$, they considered the volume average of the locally-rescaled CMF inside the condition, in the form

[TABLE]

In the following, we shall employ the same definition of averaged CMF when showing the prediction from the local rescaling.

In RBP08 the CMF computed with function (15) was directly tested against the halo abundance extracted from numerical simulations, using two specific regions corresponding to underdense regions, resulting in a reasonable agreement.

3.3 CMF normalization

In RBP08 it is stressed that the conditional mass function should satisfy the normalization condition, which means that by integrating the CMF over the possible values of the linear overdensity inside a condition $Q$, each weighted by its probability distribution, the unconditional mass function must be recovered. In the same reference, it has been shown that the local rescaling as defined in equation (3.2) does not satisfy such a condition. The problem was solved in RBP08 by slightly modifying the rescaling procedure, introducing an additional parameter $\alpha$, which we shall refer to as the normalization parameter. For the case of the T08 UMF, the new proposed rescaling reads:

[TABLE]

which should be used instead of eq. (3.2).

A proper value for $\alpha$ should yield a normalized CMF over all the mass range we are using. In RBP08 the normalization parameter was estimated for the mass function parametrizations by Sheth & Tormen (1999) and Warren et al. (2006), for one particular value of the condition Lagrangian radius. In this work we repeated the same analysis for the T08 mass function and using different values for the condition Lagrangian radius, in the range $[5,30]\,h^{-1}\,\text{Mpc}$, in order to study the dependence of $\alpha$ on the condition size. As shown in RBP08, smaller conditions are much more effective in constraining the value for the normalization parameter, which we found to be:

[TABLE]

for the T08 mass function. This parameter is instead unconstrained by larger conditions; requiring a normalized CMF is thus not enough for constraining the value of $\alpha$ at scales $Q\gtrsim 20\,h^{-1}\,\text{Mpc}$. In order to test whether the same normalization parameter can be used also for large scales, we have to compare the CMF with numerical simulations, as detailed in Section 4. We also repeated the same analysis for the standard rescaling, and found the same results; both the local and the standard rescalings need the introduction of the $\alpha$ parameter for satisfying the normalization condition. More details on the computation of the normalization parameter following the analytical procedure of RBP08 can be found in Appendix A. In addition, we also explore a different normalization condition using simulations in Section 4.3.1, and consistently, we found the same value for the $\alpha$ parameter.

4 Comparison with simulations

In this section we test the validity of the CMF prescriptions by comparing the predicted abundance of halos with numerical simulations. We want to extend the work in RBP08 by employing a considerably larger number of conditions, in order to reduce statistical uncertainties associated to the particular condition choice. We will test the CMFs against counts extracted from simulations at different scales, and both in overdense and underdense regions. We will also explore the effect of employing different values for the normalization parameter $\alpha$ entering the computation of the locally rescaled CMF.

We present in the following the simulation set we considered, the algorithm we employed to identify regions with a given density contrast and build the halo mass distribution, and the results we obtained when comparing such a distribution to the theoretical recipe for the conditional mass function.

4.1 The simulation set

We considered two N-body dark matter-only simulations. Each simulation contains $10^{9}$ particles in a cubic box with periodic conditions at the boundary, in one case of size $250\,h^{-1}\,\text{Mpc}$ , and in the other of size $1500\,h^{-1}\,\text{Mpc}$. We shall refer to these simulations as the Small Scale Simulation (SSS) and the Large Scale Simulation (LSS), respectively. The basic information on these simulations is summarized in table 2. The LSS was selected among the Minerva sample of 100 realisations of the same volume size (Grieb et al., 2016); the SSS, instead, was performed in order to extend the lower halo mass range. Initial conditions were generated starting from an initial glass configuration, and using the 2LPT code implementing the method described in Crocce et al. (2006). We adopted a flat $\Lambda$CDM model with parameters: $h=0.695$, $\Omega_{\rm b}=0$, $\Omega_{\rm m}=0.285$, $n_{\rm s}=0.9632$ and $\sigma_{8}=0.828$ (Sánchez et al. (2013), Table I), and computed the input linear power spectrum with CAMB (Lewis et al., 2000). The initial conditions were evolved from $z=63$ to $z=0$ using GADGET (last described in Springel et al. (2005)). On the simulation snapshot we run a FoF halo finder with linking length $\ell=0.1$. This approximately corresponds to a spherical overdensity $\Delta=1600$, and is used to separate density peaks that would be taken as a single halo by a lower overdensity finder. Then, around each identified peak, halos were searched with a spherical overdensity $\Delta=200$, after processing the FoF output with the SUBFIND code (Springel et al., 2001).

In Figure 1 we show the comparison between the T08 mass function, computed with an overdensity parameter $\Delta=200$, and the halo counts extracted from the simulations. Error bars associated to the simulation data points were computed assuming poissonian statistics inside each mass bin. To avoid boundary effects or resolution effects, the safe mass range we can use for each simulation is $[10^{10.5},10^{14.5}]\,h^{-1}\,\text{Mpc}$ with the SSS, and $[10^{13},10^{15.5}]\,h^{-1}\,\text{Mpc}$ with the LSS. By joining the two simulations we are capable of exploring a range of mass extending over five orders of magnitude. For both simulations, the lower mass limit is a resolution limit: below that scale, the halo finder is no longer able to converge in computing the halo mass with overdensity 200, and all the substructures around the previously identified density peak are lost. When determining the theoretical UMF, we checked that the finite-box effect has a negligible impact on the computation of the variance, remaining below one part in $10^{3}$ also in the SSS.

The first plot in Figure 1 shows that the number density of halos extracted from both simulations is capable of reproducing the overall shape of the T08 mass function, showing the same dependence on mass over the whole range we are exploring. However, the residual plot on the right shows there are still some differences with respect to the theoretical prediction, mostly at the level of a $\sim 10\,\%$ relative error. This is a systematic effect which is due to the use of only one simulation for each box; for instance, the excess of halo abundance from the LSS visible in the high mass end can be removed by averaging many realizations of the same simulation (see Figure 1 of Grieb et al. (2016)); this deviation however regards very high masses where errors become large and are not useful for the following analysis.

Besides, a numerical resolution effect is visible in both simulations: it can be recognised as a change in the mass dependence of the residual halo abundance at masses $\sim 2\times 10^{11}\,h^{-1}\,\text{Mpc}$ in the SSS and $\sim 5\times 10^{13}\,h^{-1}\,\text{Mpc}$ in the LSS. Given the individual particle masses reported in table 2, these values correspond in both cases to halos containing around 200 particles. We actually considered such halos in LSS when performing the analysis described in Section 4.3. In order to test to what extent this numerical issue can affect our results, we explicitly computed the bias as the ratio of the variance of halo number density and the variance of matter, as a function of mass. We found that numerical resolution effects only show up below masses $\sim 6\times 10^{12}\,h^{-1}\,\text{Mpc}$, corresponding to halos made of around 20 particles, which are not considered in our analysis. These resolution effects are therefore not a significant issue.

In order to achieve a better agreement with the theoretical mass function, the average of many simulations in each box should be used. This goes beyond the scope of this work. Here, in order to account for possible systematics in the following analyses, when considering comparisons between simulated and predicted halo abundances in conditioned environments, we will always be comparing the ratios of each CMF to the corresponding UMF, i.e. directly comparing the halo bias function.

4.2 Methodology

Testing the CMF requires comparing the theoretical abundance of halos to the actual counts extracted from regions of the simulation with a definite density contrast. In order to identify such regions we threw a large number of spheres on the simulation snapshot, storing the total number of particles inside each one of them; this allowed us to determine the number density of particles inside each sphere. By knowing the average particle density, the sphere density contrast could be computed using (6). The sphere positions were generated randomly with a uniform distribution inside the simulation box; taking advantage of the periodic boundary conditions, any point of the box was a possible sphere centre. The sphere radii were also generated randomly; in the case of the SSS, which we employed for the smaller conditions, the values for the Eulerian radius $R$ were drawn from a uniform distribution in the interval $[5,15]\,h^{-1}\,\text{Mpc}$. The largest conditions were instead placed over the LSS, with radii drawn from a uniform distribution in the interval $[20,30]\,h^{-1}\,\text{Mpc}$. We employed a total of $4.4\times 10^{5}$ spheres for the LSS, and of $6\times 10^{4}$ spheres for the SSS: the smaller simulation volume ensures the same spatial coverage can be obtained with a lower number of spheres. In both cases the spheres are allowed to overlap; in fact, this is not an issue when building the final halo distributions. In Figure 2 we show a sample slice of the SSS, showing the simulation particles, the identified halos and a few spheres randomly placed over the simulation.

Note that the aforementioned radii are the Eulerian quantities $R$; once the physical overdensity $\delta_{\rm m}$ inside each sphere is computed, its Lagrangian radius can be determined with (8). The correspondent linear overdensity $\delta_{\rm l}$ can be obtained via eq. (3), and this is the last quantity needed to completely characterize each sphere. In Figure 3 we show the final normalized distributions for the quantities $Q$, $\delta_{\rm m}$ and $\delta_{\rm l}$, in both simulations. We see that in the case of the SSS, a consequence of using smaller spheres is that the density contrast is peaked at more negative values compared to the LSS case, and as a result the distribution of Lagrangian radii is rather asymmetric. With the LSS spheres the density contrast distribution peak lies closer to zero, and the Lagrangian radius distribution has a more symmetric shape around the interval allowed by the prior Eulerian radius distribution.

This set of characterized conditions can be used to test the theoretical predictions of the CMF. The comparison between theory and simulation was performed at the level of the halo density as a function of mass, computed as a histogram over a specified set of mass bins. In principle, this comparison could be performed for each one of the spheres we threw over the simulations. Typically, a single sphere contains around $\sim 10^{2}$–$10^{3}$ halos in the SSS and $\sim 10$–$10^{2}$ in the LSS. However, given the high number of spheres we employed, it was convenient to group them into bins of Eulerian radius and linear density contrast, and consider the comparison between the simulated and theoretical histograms averaged inside each $(R,\delta_{\rm l})$ bin. This way we could take advantage of the higher statistics resulting from averaging the halo mass distribution from many conditions with common size and density contrast.

We obtained the simulated histograms as follows. For each sphere, using the information of its position and Eulerian radius, we recorded all collapsed objects retrieved by the halo finder that lie inside the sphere. Let us fix a $(R,\delta_{\rm l})$ bin and let $N_{\rm s}$ be the total number of spheres whose size and density contrast fall in this bin; we joined in a single set all the halos contained inside theses $N_{\rm s}$ spheres, and binned these halos over a specified set of mass bins $(m)_{j}$. We call $N_{j}$ the number of halos whose masses enter the $j$-th mass bin. We defined our halo density histogram extracted from the simulation the quantity:

[TABLE]

where $V_{\rm R}$ is the average eulerian volume of the spheres belonging to the $(R,\delta_{\rm l})$ bin, $V_{\rm R}=4\pi R^{3}/3$. We chose a sufficiently small step $\Delta R=0.5\,h^{-1}\,\text{Mpc}$ between radius bins, in order to ensure the central value $R$ was representative of the size of all spheres in the bin.

We then computed the theoretical prediction for the halo density in the bin $(R,\delta_{\rm l})$. For each sphere in this bin, we computed the theoretical CMF $\text{d}n_{\text{c},i}/\text{d}m(m;Q,\delta_{\rm l})$ using the exact values of $\delta_{\rm l}$ and $Q$ of the sphere; we integrated the CMF over the same mass bins $(m)_{j}$ considered earlier, to transform the distribution into a density as a function of mass, $n_{i}(m_{j};Q,\delta_{\rm l})$. Here the subscript $j$ runs over the mass bins, the subscript $i$ over the different spheres contained inside the $(R,\delta_{\rm l})$ bin, $i=1,...,N_{\rm s}$. The $N_{\rm s}$ histograms obtained this way, one for each sphere, had to be combined to get an average estimate for the $(R,\delta_{\rm l})$ bin. Note that the density of objects extracted from the simulation is Eulerian, while the formalism for the CMF described so far yields the Lagrangian density. So, the theoretical histograms were combined by first multiplying each one by the sphere Lagrangian volume $V_{\text{Q},i}=4\pi Q_{i}^{3}/3$, by summing the correspondent halo numbers, and by translating to final Eulerian density dividing by the mean Eulerian volume of the bin $V_{\rm R}$. In formulae:

[TABLE]

The theoretical histograms were generated using both the CMF from the standard rescaling and from the local rescaling. We also explored the effect of the normalization parameter $\alpha$ in changing the agreement of the CMF with the simulation; we considered in particular the value $\alpha=1.25$, which is the one required for normalizing the T08 mass function, and also an intermediate value $\alpha=1.5$, which lies between 1.25 and $\delta_{\rm c}$. Notice that, although the CMF is naturally expressed in Lagrangian density, we could not bin the spheres in $Q$ and $\delta_{\rm l}$ because these quantities are not independent.

4.3 Results

4.3.1 CMF normalization

The first consistency check to be applied to any analytical prescription to compute the CMF is a normalization condition. The conditional halo mass function must recover the unconditional one when averaged over all possible values of the density contrast. In practice, this is done by evaluating the analytical CMF for all the set of possible conditions $(R,\delta_{\rm l})$ that we have obtained with the spheres described in the previous subsection, and averaging over all of them. Note that this test is independent and complementary to the analysis that we perform in Appendix A. There are two main differences: first, we are using here, by construction, the true probability distribution of $\delta_{\rm l}$ as derived from our simulations (see Figure 3); and second, we are integrating for a fixed Eulerian $R$ value, instead of a fixed Lagrangian $Q$ value.

The normalization condition must hold true for each value of the radius, that is, for all bins in $R$. We show in Figure 4 the results for one radius bin, for both simulations, plotted in terms of the halo density. We plot the analytical predictions from eq. 19 obtained using the two different recipes for the CMFs, namely the local and the standard rescalings of the T08 function, both with no explicit normalization and with normalization parameter $\alpha=1.25$. In order to better visualize the deviations among the different predictions, we show in the lower panel of Figure 4 the corresponding plot of residuals, computed as percent deviations with respect to the unconditional mass function (T08 in this case). Residuals from the analytical determinations of the CMF give us information about the accuracy at which the normalization condition is satisfied. These residuals are clearly reduced when introducing the normalization parameter $\alpha=1.25$ in the UMF rescalings, resulting in a final $\sim 1\,\%$ error with respect to the T08 mass function. This confirms that the modification of the rescaling introduced to accomplish the normalization condition also works when using the proper $\delta_{\rm l}$ distribution from numerical simulation, and for a fixed Eulerian condition radius. These conclusions are the same for both the CMF computed with the standard rescaling and the CMF computed with the local rescaling. Residuals from all theoretical CMFs grow larger than $10\,\%$ above a certain mass threshold, around $\sim 10^{14}\,h^{-1}\,\text{M}_{\odot}$ for the SSS, and $\sim 4\times 10^{14}\,h^{-1}\,\text{M}_{\odot}$ for the LSS. Now, in RBP08 it is underlined that the formalism for computing the CMF is expected to fail when approaching the mass of the condition; in the same reference, a fraction of $1/30$ of the condition mass is chosen as an upper limit setting the mass range where the CMF formalism is valid. According to this prescription, the mass values above which our residuals grow very large correspond to the limiting top masses for conditions of size $Q\sim 20\,h^{-1}\,\text{Mpc}$ for the SSS, and $Q\sim 35\,h^{-1}\,\text{Mpc}$ for the LSS. By looking at the $Q$ probability distributions in Figure 3, we found these radii to fall in the very high tail of the distributions. Errors in normalizations therefore arise in a mass interval dominated by rare massive objects that are found in highly overdense regions, which we are not going to consider in the following analysis.

For completeness, we also added in those panels in Figure 4 the data points corresponding to the average of the actual CMFs obtained from the simulations for those spheres, computed using eq. 18. As expected, the average of all the numerical CMFs recovers the UMF plotted with points in Figure 1 with few percent accuracy, for all $R$ bins, and for both SSS and LSS simulations. This result proves the consistency of our algorithm in placing random spheres and averaging their content. Note that errors grow larger in the high mass end where the box finite-size effect becomes important and the statistics of collapsed objects is poorer.

4.3.2 Halo abundance in overdense and underdense regions

We then checked the predictions of the CMF recipes for individual density contrast bins. Note that results from eq. 19, in this case, are not the quantity that can be compared to the simulation results eq. 18; indeed, according to the discussion in Section 4.1, the abundances extracted from simulations are affected by systematics due to cosmic variance. So, let $n_{\rm u}^{\text{theo}}$ and $n_{\rm u}^{\text{sim}}$ be the halo abundances computed using the unconditional mass functions from theory (T08) and simulations, respectively; both are plotted in Figure 1. In order to take the aforementioned systematics into account, the theoretical prediction can be rescaled as:

[TABLE]

The prediction $\tilde{n}^{\text{theo}}$ is the one we finally compared with the simulation counts $n^{\text{sim}}$ for testing the CMFs recipes. In practice, this is equivalent to directly compare the so called matter-to-halo bias function (Sheth & Lemson, 1999), both from the simulations and the theoretical predictions.

In Figure 5 we show the comparison between the simulations and the CMFs computed using the standard rescaling and the local rescaling with no normalization, in the case of underdense regions; we show results for four bins of $\delta_{\rm l}$. Figure 6 shows two examples for the case of an overdense region. In all plots, a vertical line shows the mass correspondent to $1/30$ of the condition mass. We see that, in overdense regions, both the standard and the local CMF provide a very similar prediction, and in good agreement with the counts from simulations. In underdense regions, however, the agreement with simulations gets poorer, especially when the density contrast is stronger. Note that the oscillations of the theoretical predictions at low masses is simply a consequence of the normalization performed with eq. 20.

In Figures 7 and 8 we consider the same radius and density contrast bins, this time showing the comparison between counts from simulations and the prediction from the local rescaling with explicit normalization. We employed the value $\alpha=1.25$ as required by the T08 mass function normalization; we also show, as a reference, the locally rescaled CMF without normalizing and normalized with an intermediate choice of the normalization parameter, $\alpha=1.5$. In case of underdense regions, the explicit normalization of the CMF implies a clear improvement when reproducing the counts from simulations, particularly when the density contrast is stronger. For the most underdense bin we show, when changing from the not normalized CMF to the one using $\alpha=1.25$, the reduced $\chi^{2}$ changes from $9.7$ to $2.4$ in the SSS, and from $22.6$ to $0.7$ in the LSS. These results suggest that normalizing the local CMF with $\alpha=1.25$ provides a more accurate recipe in underdense regions, both at small and large scales; this is in agreement with the results in RBP08, where the same value for $\alpha$ is employed for normalizing the Sheth & Tormen (1999) mass function. In case of overdense regions, the use of the same value of $\alpha$ does not always produce an improvement in the agreement with simulated counts. When introducing the $1.25$ normalization parameter, the reduced $\chi^{2}$ changes from $1.8$ to $42.7$ in the SSS, actually worsening the fit, while in the LSS we have again improvement with the reduced $\chi^{2}$ changing from $16.5$ to $3.3$.

We repeated the same comparisons using the standard CMF, computed with explicit normalization parameters $\alpha=1.25$ and $\alpha=1.5$. The resulting plots comparing the theoretical prediction with the simulation counts, for the same set of radius and overdensity bins considered so far, are shown in Figures 9 and 10. The results are much alike the ones obtained with the local rescaling, with the explicit normalization substantially improving the agreement with simulation in underdense regions but not always in overdensities. In this case, for the most underdense bin, when introducing the $\alpha=1.25$ parameter the reduced $\chi^{2}$ changes from $10.3$ to $2.4$ in the SSS and from $18.0$ to $2.0$ in the LSS. In the positive $\delta_{\rm l}$ bin, instead, the reduced $\chi^{2}$ changes from $1.7$ to $40.5$ in the SSS, and from $11.4$ to $2.9$ in the LSS.

When comparing these $\chi^{2}$ with the ones quoted above, and for all possible cases, we find that the agreement with simulation is basically the same when using the local CMF or the standard CMF (see e.g., Figures 5 and 6), although the local rescaling provides a better description in underdense and intermediate regions than the standard rescaling. In overdense regions, however, the standard rescaling seems slightly better for this particular simulation.

In conclusion, when trying to reproduce the abundance of halos inside conditioning regions, the best recipe for computing the CMF, in either the standard or the local rescaling, requires the use of the $\alpha$ parameter. When using $\alpha=1.25$, for underdense regions we clearly improve the agreement with numerical simulations, while for the case of overdensities, although there is still an improvement for large scales, we have an overcorrection for small regions. We checked that these conclusions still hold when changing the mass range over which we computed the $\chi^{2}$. The large values for the $\chi^{2}$ obtained with explicit normalisation in the SSS may be due to an underestimation of the errors; what matters here is that in this case the normalised CMF worsen the agreement with simulated abundances. It is not clear whether this issue is due to the particular simulation we are using. The same comparison could be repeated by combining results from many N-body simulations; this goes beyond the scopes of this work. In any case, we conclude that the normalized CMF with $\alpha=1.25$ is the one that works better in the majority of the cases.

5 A fitting formula for the CMF in underdense regions

As discussed in RBP08, one of the main driving goals for developing a proper parametrization of the conditional mass function is the theoretical study of void statistics. A void can be considered, within our framework, as a locally underdense region, and the distribution of collapsed objects within it can be described by the CMF. The formalism for computing the number density of void is discussed in Patiri et al. (2006); voids can be used as a tool to constrain cosmological parameters, as it is shown in Betancort-Rijo et al. (2009).

In those works the probability of finding a void of a given radius requires the computation of the quantity:

[TABLE]

where we made the dependence on the objects mass, on redshift and on the void radius explicit. We can call the first factor a ‘systematic bias’, which is due to the different expansion rate of the void compared to the background. The second factor in eq. (21) is the ‘statistical bias’ (or simply the bias), accounting for the clustering of the protohaloes within the condition before it expands in comoving coordinates. This is usually referred in the literature as the (matter-to-)halo bias function (Sheth & van de Weygaert, 2004; Furlanetto & Piran, 2006; Paranjape et al., 2012; Jennings et al., 2013). Again in Patiri et al. (2006), it is shown that the relation between the bias and the CMF is given by:

[TABLE]

Here $N_{\rm u}(m)$ is simply the cumulative mass function for collapsed objects of mass greater than $m$, in the unconditional case:

[TABLE]

and $N_{\rm c}(m;q,Q,\delta_{\rm l})$ is the same for the conditional mass function:

[TABLE]

In this case the mass range over which the mass function has to be integrated is limited by the mass of the condition; the upper integration extrema should therefore be proper fraction $M$ of the condition mass $m(Q)$. In our computations we always used $M=m(Q)/30$. Notice that the mass function we use in eqs. (22) and (24) is the locally rescaled CMF with normalization parameter $\alpha=1.25$. Indeed, the voids that are best suited for constraining cosmological models are the largest ones, and according to the discussion in Section 4.3 this is the best CMF recipe in the case of underdensities and large conditions.

Eq. (22) is the ratio of the total density of collapsed objects above a given mass in the conditioned environment with respect to the unconditional one. It contains all the information on the conditional mass function, and is the starting point for the study of void statistics. Given all the dependences that show up in eq. (22), the full numerical computation of the bias is very time-consuming. Therefore, efforts have been made to find a numerical fit for the bias, which allows us to avoid heavy numerical computations. In RBP08 a fitting formula in the form:

[TABLE]

is discussed, and a procedure for computing the fitting parameters $A$ and $b$ is provided. The fit (25) is a good approximation for the bias at the scales discussed in RBP08 (around $Q=9\,h^{-1}\,\text{Mpc}$) and for strongly underdense regions. However, as already stated, when dealing with void statistics the largest voids are the ones that best constrain cosmological parameters, on scales around $Q=20\,h^{-1}\,\text{Mpc}$. The fit (25) cannot be applied at those scales, nor when the density contrast approaches zero. In this work, we consider an extension of the fitting formula (25), which can be applied to a broader range of condition sizes and overdensities; the function we propose has the form:

[TABLE]

which differs from (25) in the addition of the extra parameter $c$.

We wanted to check the goodness of this fit in reproducing the bias. To this aim, we implemented the full computation given in eq. (22), using the reference cosmology from table 1, for different values of the condition radius in the range $Q\in[18,28]\,h^{-1}\,\text{Mpc}$. The linear density contrasts we considered span the range $[-1.5,0.0]$. We also considered the dependence on redshift, which we moved in the range $[0,1]$, and on the parameter $\sigma_{8}$, which we allowed to vary in $[0.6,1.0]$. According to the results described in Section 4.3, since we are considering underdense regions, we computed the UMF employing the normalized local rescaling (16), with normalization parameter $\alpha=1.25$. For every case, we fitted both functions (25) and (26) against the bias, as a function of $\delta_{\rm l}$. We show some examples in Figure 11; the use of the additional parameter clearly improves the goodness of the fit, which is capable of reproducing the bias with a relative error in $\sim$ 1–2 %. This result is useful for computing the halo abundance as a function of mass in underdense regions, without implementing the full computation of the local CMF rescaling: by computing the bias with (26) and the integrated UMF with (23), their product already yields the integrated CMF, which can be derived to obtain the proper conditional mass function.

The fitting parameters $A$, $b$ and $c$ inherit the dependence on redshift, on the condition size and on the object mass; in addition, they are dependent on the underlying cosmology. Our goal was to give a fitting formula which is capable of reproducing the bias computed with any chosen cosmology, and any scale and redshift. To this aim, we further fitted the dependence of the parameters on these quantities; the full procedure is described in Appendix B. We first provide fitting formulae for computing the parameters $A$, $b$ and $c$ for a reference cosmology, for any choice of the variables $m$ and $Q$. These parameters can then be scaled to a different cosmology. Since the dependences on redshift and on $\sigma_{8}$ are correlated via the linear growth factor, we first fitted jointly the dependence on these two variables. This way, once the values of $z$, $\sigma_{8}$, $Q$ and $m$ are chosen, it is possible to compute the fit to the bias using the procedure described in Appendix B and using the common values from the tables presented there; the final accuracy at reproducing the bias is always below $\sim 8.5\,\%$, but it most cases is around a few percents. We attempted also at allowing the cosmological parameters $h$ and $\Omega_{\rm m}$ to vary, including their effect in a rescaling of the mass definition; in this case, however, the quality of the fit degrades and the relative error can grow as large as $40\,\%$ in the worst cases.

6 Conclusions

We have compared the predicted CMF of dark matter halos from two theoretical prescriptions against numerical N-body simulations: the standard rescaling of the UMF, and a rescaling locally defined inside the condition (RBP08). Both prescriptions use as a reference the T08 parametrization for the UMF. Drawing from the results of RBP08, we reminded that the locally rescaled mass function requires the introduction of an additional parameter $\alpha$ in order to satisfy the normalization condition, i.e., the integral of the CMF over the weighted density contrast should yield the initial UMF. The problem of normalization thus turns into the problem of determining a proper value for $\alpha$ that satisfies the normalization condition over the whole mass range in which the mass function is to be computed. We repeated the analysis described in RBP08, and found that the standard rescaling also requires the normalization to be implemented explicitly via the parameter $\alpha$. For the first time, we obtained the value of $\alpha$ needed to normalize the T08 mass function, and studied its dependence on the size of the condition. We found that this parameter is well constrained by small regions ($Q\sim 8\,h^{-1}\,\text{Mpc}$), yielding the estimate $\alpha=1.25$ for both the standard and the local rescalings; large conditions ($Q\sim 25\,h^{-1}\,\text{Mpc}$), instead, are much less effective in constraining its value, leaving $\alpha$ undetermined.

In order to test the validity of these CMFs in predicting the abundance of collapsed objects inside conditions, we tested the results from the theoretical rescalings against halo counts extracted from numerical simulations. We employed two different dark matter-only simulations, one at relatively small scales ($250\,h^{-1}\,\text{Mpc}$ box size), the other at larger scales ($1500\,h^{-1}\,\text{Mpc}$ box size). Jointly, the two simulations allowed us to explore the mass range $[10^{10.5},10^{15.5}]\,h^{-1}\,\text{Mpc}$. The correspondent halo catalogues were obtained using a FoF algorithm with $b=0.2$, which corresponds to the $\Delta=200$ spherical overdensity mass function definition used in T08. As a preliminary check, we compared the halo densities as a function of mass, extracted from both simulations, with the prediction computed using the T08 UMF. We found that the simulated counts reproduce the mass dependence of the theoretical halo abundances but are affected by systematic effects at the level of $\sim 10\,\%$ relative error. These offsets are due to the use of a single simulation for each box and have been accounted for in our analyses.

We then placed a large number of spheres over the simulations, with random positions and radii, in order to simulate the conditioning regions. We grouped the spheres in bins of size and density contrast, and obtained the correspondent conditioned halo density as a function of mass. These distributions were then compared with the predictions from the different CMF recipes. For every radius bin, we joined together all spheres with common size and all possible density contrasts, in order to test a new normalization condition. This is an important test since by construction we are using the exact probability distribution of $\delta_{\rm l}$ when averaging over different spheres. We found that the new normalization condition is satisfied within a few percent accuracy by the abundances from simulations, when compared to the UMF from the complete halo catalogue. The theoretical recipes for computing the CMF also achieve that precision in the normalization when compared to the T08 UMF prediction, provided the normalization is explicitly applied with the parameter $\alpha=1.25$. Normalization is no longer achieved in the high mass end, where errors from the simulation grow large and the formalism for computing the CMF begins to fail for getting too close to the condition mass.

The following test was checking the CMF prediction for the halo abundance inside regions with a definite density contrast. We found that the explicitly normalized CMF with $\alpha=1.25$ is the most accurate recipe in reproducing simulated halo counts inside underdense regions, both at small and large scales; the CMF can be computed either with the standard or the local rescalings, both providing very similar results, although the latter provides a slight better fit. In the case of overdense regions, the use of the same normalization parameter results in an improvement of the agreement with simulations at large scales, but in a slight overcorrection at small scales. Apart from this case, the use of the locally rescaled CMF normalised with $\alpha=1.25$ is the recipe that best works in most cases.

Finally, we presented an analytical fit for computing the CMF in underdense regions, in terms of the bias with respect to the UMF. By multiplying the fitting function times the integrated UMF, the integrated CMF can be obtained, avoiding the time-consuming computation required to implement the normalized local rescaling in underdense regions. In RBP08, a possible analytical expression for this fitting function was proposed, which could be applied to extremely underdense regions (voids); we checked that such a fit is actually inadequate at reproducing the bias dependence on the density contrast when $\delta_{\rm l}$ gets closer to zero. We proposed a different fitting function by adding an additional parameter, that clearly improves the agreement with the CMF bias, yielding final errors below $\sim 2\,\%$. The dependence of the fitting parameters on the condition size, the halo mass and the cosmology has also been fitted; it is possible to use the same fitting parameters for computing the bias with different values of the mass, the condition radius, the redshift and $\sigma_{8}$, with a relative error below $\sim 8\,\%$. It is also possible to apply the same fitting parameters moving $h$ and $\Omega_{\rm m}$ from their reference values, computing the bias within a $40\,\%$ error. We present the whole set of fitting parameters in Appendix B.

Acknowledgements

JARM acknowledges financial support from the Spanish Ministry of Economy and Competitiveness (MINECO) under the 2011 Severo Ochoa Program MINECO SEV-2011-0187. This work has been partially funded by the Spanish Ministry of Economy and Competitiveness (MINECO) under the projects ESP2013-48362-C2-1-P and AYA2014-60438-P.

CDV acknowledges financial support from the Spanish Ministry of Economy and Competitiveness (MINECO) under the 2011 and 2015 Severo Ochoa Programs SEV-2011-0187 and SEV-2015-0548, and grants AYA2013-46886 and AYA2014-58308.

Simulations have been performed on the Hydra cluster at the Max Planck Rechenzentrum Garching (RZG) and the COSMA-4 supercomputing facility of the Institute for Computational Cosmology (ICC) of Durham.

This work used the DiRAC Data Centric system at Durham University, operated by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility (www.dirac.ac.uk). This equipment was funded by BIS National E-infrastructure capital grant ST/K00042X/1, STFC capital grant ST/H008519/1, and STFC DiRAC Operations grant ST/K003267/1 and Durham University. DiRAC is part of the National E-Infrastructure.

Appendix A A derivation of the $\alpha$ parameter based on the analytical normalization condition for the CMF

Following RBP08, we can use the following normalization condition

[TABLE]

to derive the optimal value of the $\alpha$ parameter. In this equation, $G(\delta_{\rm l};\sigma_{1})$ stands for a Gaussian with zero mean and variance $\sigma_{1}$. Note that by construction, the expression must be evaluated at masses no larger than the condition mass $m(Q)$. The upper integration limit $X$ is equal to the critical density $\delta_{\rm c}$ for mass functions based on the physics of spherical collapse (Press & Schechter, 1974), and to the value of the moving barrier shape when the ellipsoidal collapse is allowed (Sheth & Tormen, 2002).

We note that this normalization condition has been previously discussed in the literature, mainly in connection to the PS mass function. For example, several authors (Sheth & Lemson, 1999; Musso et al., 2012; Neyrinck et al., 2014) have shown that for the (unphysical) case of an unbound upper limit of the integral (i.e. $X=+\infty$), the normalization is satisfied for the PS mass function. Therefore, when the integration is truncated to the corresponding upper limit, a modification of the normalization is required in order to verify the equation.

In RBP08, this condition was applied to different prescriptions for the unconditional mass function. The method proposed in RBP08 in order to solve the normalization problem is the introduction of a normalization parameter $\alpha$ in the rescaling for the critical density:

[TABLE]

For the T08 parametrization, in which $\delta_{\rm c}$ does not appear explicitly, the new rescaling is to be applied in the form of eq. (16).

Our goal was to find a proper estimate for the normalization parameter and to explore its dependence on $Q$, both for the standard and the local rescaling. We followed the same method described in RBP08. In order to determine the precision with which the normalization is satisfied by the CMF we computed the ratio:

[TABLE]

that quantifies the error in the normalization (how much the integrated conditional mass function differs from the expected unconditional one). In (29) the CMF is to be computed with the new rescaling eq. (28) or (16), depending on the chosen parametrization for the initial UMF. By varying the parameter $\alpha$ we end up with a three-variable function $R(m,Q,\alpha)$. For different values of the condition radius we computed the quantity $(1-R(m,Q,\alpha))\times 100$ and considered the regions in the $(m,\alpha)$ plane where it is at most $10\,\%$, $5\,\%$ and $3\,\%$. The result is plotted in Fig. 12 for the case of local rescaling. We see that it is possible to choose a unique $\alpha$ value for all masses; in the case of $Q=6\,h^{-1}\,\text{Mpc}$, the value $\alpha=1.25$ satisfies normalization with a $3\,\%$ accuracy. This is the value we chose for the comparison of the locally rescaled CMF with simulations. However, the precision contours largely widens when growing the condition, meaning that larger conditions are not effective in constraining the normalization problem. In fact, we see that for $Q\gtrsim 20\,h^{-1}\,\text{Mpc}$ the usual value $\delta_{\rm c}=1.686$ is allowed for the whole mass range. Setting $\alpha=\delta_{\rm c}$ is equivalent to using the original form of the local rescaling, eqs. (3.2) and (14). We repeated the same analysis also for the standard rescaling, and present the normalization contours in Figure 13; we found that the contours evolve in a similar way, constraining the value $\alpha=1.25$ in the smaller conditions.

The conclusion, as anticipated in Section 3.3, is that the normalization condition is effective in constraining the value of $\alpha$ in the case of small conditions, with a radius comparable to the ones used in RBP08, but is not enough when dealing with sufficiently large conditions ($Q\gtrsim 20\,h^{-1}\,\text{Mpc}$).

Appendix B Fitting formulae for the CMF (or halo bias function) in underdense regions

We present in the following the full procedure we employed to parametrize the residual dependence of the bias fitting parameters, and the algorithm for computing them for any choice of cosmology, mass and condition size.

B.1 Description

The fitting parameters $A$, $b$ and $c$ defined in (22) inherit from the bias the dependence on the object mass $m$, the condition size $Q$, the redshift $z$ and the cosmological parameters. It is convenient to split the dependence on the variables $Q$ and $m$ from the cosmological variables (including $z$).

For a given cosmology the dependence on the condition can be fitted by a polynomial function in powers of $Q$ and $\ln{(m)}$. More precisely, the fitting functions we employed are:

[TABLE]

where $m^{\prime}=m/3.51\times 10^{11}\,h^{-1}\,\text{M}_{\odot}$ and $Q^{\prime}=Q/8\,h^{-1}\,\text{Mpc}$. Polynomials are third degree in both radius and mass, resulting in 16 $a_{ij}$ fitting terms for $A$, 16 $b_{ij}$ fitting terms for $b$ and 16 $c_{ij}$ fitting terms for $c$. For our reference cosmology, these terms are listed in table 3. They were obtained allowing the Lagrangian radius to vary in the range $[18.0,34.0]\,h^{-1}\,\text{Mpc}$, and the halo mass in the range $[10^{3},6\times 10^{3}]\,h^{-1}\,\text{Mpc}$; the intervals also fix the limits in which the parameters in table 3 can be applied. Using these tabulated values and eq. (B.1), it is possible to quickly compute the bias with a typical error $\sim$ 1–2 %

The dependence on cosmology is better fitted first in terms of the combination $D(z)\sigma_{8}$, with $D(z)$ the growth factor of linear perturbation at redshift $z$, normalized to unity at $z=0$. Indeed, other combinations of the variables $z$ and $\sigma_{8}$ produce quite a large scatter in the points, making this 1-d fit impossible. This can be understood since the mass function does not depend separately on $z$ and $\sigma_{8}$. In order to fit this dependence, we introduced the variable

[TABLE]

where $\tilde{\sigma}_{8}$ is the value for the reference cosmology. In our case, our reference choice was $\tilde{z}=0$, $\tilde{\sigma}_{8}=0.8$, to which corresponds $\mathcal{C}=1$. We considered redshift values in $[0,1]$ and $\sigma_{8}$ values in $[0.6,1.0]$. Now, given the parameters $A$, $b$ and $c$ for the reference cosmology, we found that a convenient function to fit the parameters change when moving the value of $z$ or $\sigma_{8}$ is:

[TABLE]

The fitting function is the product of a logarithm of $\mathcal{C}$ and a polynomial in $\mathcal{C}$ up to the third power. Note that with this definition of the fitting function, when $\mathcal{C}=1$ the reference values for the parameters are recovered exactly. The fit for each parameter now depends on 4 coefficients, $\xi$, $\eta$, $\zeta$ and $\omega$, which in turn still depend on $m$ and $Q$.

The next step is to fit the residual dependence of these coefficients. We use a fitting function similar to the one in eq. (B.1), but since the dependence of these coefficients on mass and radius is more complex than it was for $A$, $b$ and $c$, we employed a 4-degree polynomial in both $m$ and $Q$. In formulae:

[TABLE]

In this case we have 25 terms $\xi_{i,j}$ for $\xi$, 25 terms $\eta_{i,j}$ for $\eta$ and so on. The whole set must be repeated three times (for $A$, $b$ and $c$), ending up with 300 parameters for this last fit. Tables 4, 5 and 6 list these parameters for $A$, $b$ and $c$ respectively. This set of fitting parameters allows the computation of the bias for any value of $Q$, $m$, $z$ and $\sigma_{8}$ in the specified ranges, with a final error below $\sim 8\,\%$.

The last step is the dependence of fitting parameters $A$, $b$ and $c$ on the matter density parameter $\Omega_{\rm m}$ and the Hubble parameter $h$. We found that the same parameters obtained for the reference cosmology can still be employed, provided the mass is scaled as:

[TABLE]

where $\tilde{h}$ and $\tilde{\Omega}_{\rm m}$ are the values for the reference cosmology, and

[TABLE]

The terms $f_{ij}$ are listed in table 7. We computed them allowing $h$ to vary in $[0.65,0.75]$, and $\Omega_{\rm m}$ in $[0.25,0.35]$. With this recipe it is possible to compute the bias using the same fitting parameters listed in tables 3 to 6. In this case, the final relative error with respect to the bias computed analytically is usually in the range $\sim$ 10–20 %, but in some cases it can grow up to $40\,\%$.

B.2 Implementation

To sum up, in order to determine the value of the bias, given as input a chosen cosmology and the values of $Q$, $m$, $z$, and $\delta_{\rm l}$, one should:

for the chosen values of the Hubble parameter and the matter density parameter, compute the exponent $f(h,\Omega_{\rm m})$ using eq. (35) with the parameters from table 7; 2. 2.

scale the mass according to eq. (34), with the value of $f$ just computed. From now on this is the new mass that has to be used; 3. 3.

compute the parameters $A(m,Q;\tilde{z},\tilde{\sigma_{8}})$, $b(m,Q;\tilde{z},\tilde{\sigma_{8}})$ and $c(m,Q;\tilde{z},\tilde{\sigma_{8}})$ for the reference cosmology, using eq. (B.1) and the parameters listed in table 3; 4. 4.

compute the parameters $(\xi_{\rm A},\eta_{\rm A},\zeta_{\rm A},\omega_{\rm A})$, $(\xi_{\rm b},\eta_{\rm b},\zeta_{\rm b},\omega_{\rm b})$ and $(\xi_{\rm c},\eta_{\rm c},\zeta_{\rm c},\omega_{\rm c})$, using the fit described in eq. (B.1) with parameters from tables 4, 5 and 6; 5. 5.

evolve $A$, $b$ and $c$ to the chosen values of $z$ and $\sigma_{8}$ using eq. (B.1); 6. 6.

compute the final bias with eq. (26).

This computation, though requiring storing a large number of fitting parameters, is much faster than the implementation of the full expression (22). Reading the parameters from files and computing polynomials is a task that can easily be implemented in a numerical code, and is efficiently handled by a processor.