text:
# The mass and temperature functions in a moving barrier model.

## 1 Introduction

In the most promising cosmological scenarios, structure formation is traced back to the hierarchical growth of primordial Gaussian density fluctuations originated from quantum fluctuations (Guth & Pi 1982; Hawking 1982; Starobinsky 1982; Bardeen et al. 1986 - hereafter BBKS). Starting from these fluctuations, collapsed, virialized dark matter haloes condensed out. Within these haloes gas cools and stars form (White & Reese 1977; White & Frenk 1991; Kauffmann et al. 1999). So, the structure of dark matter haloes is of fundamental importance in the study of the formation and evolution of galaxies and clusters of galaxies. From the theoretical point of view, the structure of dark matter haloes can be studied both analytically and numerically. An analytical model that has achieved a wide popularity is the Press-Schecter (1974) (hereafter PS) formula, which allows one to compute good approximations to the mass function. PS and Bond et al. (1991) gave a detailed description of the PS statistics together with an understanding of its dynamical basis. In the PS and Bond et al (1991) papers, the authors described how the statistical properties of the initial density field, assumed to be Gaussian, together with the spherical collapse model, could be used to derive an estimate of the number density of collapsed dark matter haloes at later times, the so called universal “unconditional” mass function. Lacey and Cole (1993), showed how the model could be extended to estimate the merging rate of small objects to form larger ones, thus leading to the possibility of estimating the “conditional” mass function of sub-haloes within parent haloes. Mo & White (1996) applied the model to compute an approximation to the spatial clustering of dark haloes.

Although the analytical framework of the PS model has been greatly refined and extended (as testified by the previous cited papers), it is well known that the PS mass function, while qualitatively correct, disagrees with the results of N-body simulations. In particular, the PS formula overestimates the abundance of haloes near the characteristic mass $`M_{}`$ and underestimates the abundance in the high mass tail (Efstathiou et al. 1988; White, Efstathiou & Frenk 1993; Lacey & Cole 1994; Tozzi & Governato 1998; Gross et al. 1998; Governato et al. 1999). The quoted discrepancy is not surprising since the PS model, as any other analytical model, should make several assumptions to get simple analytical predictions. As previously reported, the main assumptions that the PS model combines are the simple physics of the spherical collapse model with the assumption that the initial fluctuations were Gaussian and small. On average, initially denser regions collapse before less dense ones, which means that, at any given epoch, there is a critical density, $`\delta _c(z)`$, which must be exceeded if collapse is to occur. In the spherical collapse model, this critical density does not depend on the mass of the collapsed object. Taking account of the effects of asphericity and tidal interaction with neighbors, Del Popolo & Gambera (1998) and Sheth, Mo & Tormen (2001) (hereafter SMT), using a parametrization of the ellipsoidal collapse, showed that the threshold is mass dependent, and in particular that of the set of objects that collapse at the same time, the less massive ones must initially have been denser than the more massive, since the less massive ones would have had to hold themselves together against stronger tidal forces. In the second hand, the Gaussian nature of the fluctuation field means that a good approximation to the number density of bound objects that have mass $`m`$ at time $`z`$ is given by considering the barrier crossing statistics of many independent and uncorrelated random walks, where the barrier shape $`B(m,z)`$, is connected to the collapse threshold. Simply changing the barrier shape, SMT showed that it is possible to incorporate the “quoted effects” 111Namely that in the case of objects collapsing at the same time, the less massive regions must initially have been denser than the more massive ones. in the excursion set approach. Moreover, using the shape of the modified barrier in the excursion set approach, it is possible to obtain a good fit to the universal halo mass function. 222Note that at present there is no good numerical test of analytic predictions for the low mass tail of the mass function. As previously reported, the excursion set approach allows one to calculate good approximations to several important quantities, such as the “unconditional” and “conditional” mass functions. Sheth & Tormen (2002) (hereafter ST) provided formulas to calculate these last quantities starting from the shape of the barrier. They also showed that while the “unconditional” and “conditional” mass function is in good agreement with results from numerical simulations, neither the constant nor the moving barrier models (barrier obtained from non-spherical collapse) were able to describe the simulations results at small lookback times, in the case of the rescaled (in terms of $`\nu `$) “conditional” mass function. The reason for this discrepancy is probably due to the excursion set approach’s neglect of correlations between scales (Peacock & Heavens 1990; Bond et al. 1991; ST) or to the too simple parametrization of the ellipsoidal collapse outlined in SMT.

In the present paper, I’ll use the barrier shape obtained in Del Popolo & Gambera (1998), obtained from the parametrization of the nonlinear collapse discussed in that paper, together with the results of ST in order to study the “unconditional” and “conditional” mass function. Finally, I’ll calculate the temperature function for a CDM model by means of the mass functions previously obtained and using an improved version of the M-T relation obtained by Voit (2000), which accounts for the fact that massive clusters accrete matter quasi-continuously, and consequently that the M-T relation evolves, with time, more modestly than what expected in previous models (top-hat model) and taking account of the tidal interaction with neighboring clusters.

The reasons that motivates this study are several:
a) to study the effects of a barrier different from that used by ST on both “unconditional” and “conditional” mass function, as proposed even by ST ;
b) to study how well ST formulas really do work for several barrier shapes;
c) to test if the discrepancies between the temperature function and simulations, observed in several papers (e.g. Governato et al. 1999) are reduced using the mass function and the M-T relation obtained.

The paper is organized as follows: in Sect. 2, I calculate the “unconditional” and “conditional” mass functions. In Sect. 3, I introduce a model for the mass-temperature (M-T) relation and calculate the temperature function. Sect. 4 and 5 are devoted to results and to conclusions, respectively.

## 2 The barrier model and the “unconditional” and “conditional” mass functions

Following Sheth & Tormen (1999) notation, if $`f(m,\delta )dm`$ denotes the fraction of mass that is contained in collapsed haloes that have mass in the range $`m`$-$`m+dm`$, at redshift $`z`$, and $`\delta (z)`$ the redshift dependent overdensity, the associated “unconditional” mass function is:

$$n(m,\delta )dm=\frac{\overline{\rho }}{m}f(m,\delta )dm$$
(1)
where $`\overline{\rho }`$ is the background density. If $`f(m,\delta _1|M,\delta _0)`$ denotes the fraction of the mass of a halo $`M`$ at $`z_0`$ that was in subhaloes of mass $`m`$ at $`z_1`$, ($`z_1>z_0`$), the “conditional” mass function is:

$$N(m,\delta _1|M,\delta _o)dm\frac{M}{m}f(m,\delta _1|M,\delta _o)dm$$
(2)
In the excursion set approach, the average comoving number density of haloes of mass $`m`$ the universal or “unconditional” mass function, $`n(m,z)`$, is given by:

$$n(m,z)=\frac{\overline{\rho }}{m^2}\frac{d\mathrm{log}\nu }{d\mathrm{log}m}\nu f(\nu )$$
(3)
(Bond et al. 1991), where $`\overline{\rho }`$ is the background density, $`\nu =\left(\frac{\delta _\mathrm{c}(z)}{\sigma (m)}\right)^2`$ is the ratio between the critical overdensity required for collapse in the spherical model, $`\delta _\mathrm{c}(z)`$, to the r.m.s. density fluctuation $`\sigma (m)`$, on the scale $`r`$ of the initial size of the object $`m`$. The function $`\nu f(\nu )`$ is obtained by computing the distribution of first crossings, $`f(\nu )d\nu `$, of a barrier $`B(\nu )`$, by independent, uncorrelated Brownian motion random walks. The mass function can be thus calculated once a shape for the barrier is given and the power spectrum is known. In the case of spherical collapse, characterized by a constant barrier (for all $`\nu `$), Bond et al. (1991) obtained:

$$\nu f(\nu )=\left(\frac{\nu }{2\pi }\right)^{\frac{1}{2}}\mathrm{exp}(\frac{\nu }{2})$$
(4)
In the case of a nonspherical collapse, the shape of the barrier is no longer a constant and moreover it depends on mass (Del Popolo & Gambera 1998; SMT). As shown by ST, for a given barrier shape, $`B(S)`$, where $`SS_{}\left(\frac{\sigma }{\sigma _{}}\right)^2=\frac{S_{}}{\nu }`$ and $`\sigma _{}=\sqrt{S_{}}=\delta _{co}`$, the first crossing distribution is well approximated by:

$$f(S)dS=|T(S)|\mathrm{exp}(\frac{B(S)^2}{2S})\frac{dS/S}{\sqrt{2\pi S}}$$
(5)
where $`T(S)`$ is the sum of the first few terms in the Taylor expansion of $`B(S)`$:

$$T(S)=\underset{n=0}{\overset{5}{}}\frac{(S)^n}{n!}\frac{^nB(S)}{S^n}$$
(6)
The previous Eq. (5),(6), reduce to Eq. (4) for constant barriers. In the case of the ellipsoidal barrier shape given in ST:

$$B(\sigma ^2,z)=\sqrt{a}\delta _c(z)\left[1+\frac{\beta }{\left(a\nu \right)^\alpha }\right]$$
(7)
where $`a=0.707`$, $`\delta _\mathrm{c}(z)=1.686(1+z)`$, $`\beta 0.485`$ and $`\alpha 0.615`$, Eqs. (5),(6), give, after truncating the expansion at $`n=5`$ (see ST):

$$f(\nu )d\nu =A\left(1+\frac{0.094}{\left(a\nu \right)^{0.6}}\right)\sqrt{\frac{a\nu }{2\pi }}\mathrm{exp}\{a\nu \left[1+\frac{0.5}{\left(a\nu \right)^{0.6}}\right]^2/2\}$$
(8)
with $`A1`$, which is in good agreement with the fit of the simulated first crossing distribution (ST):

$$f(\nu )d\nu =A_1\left(1+\frac{1}{\left(a\nu \right)^p}\right)\sqrt{\frac{a\nu }{2\pi }}\mathrm{exp}(a\nu /2)$$
(9)
where $`p=0.3`$, $`A_1=0.3222`$ and $`a=0.707`$. 333 Note, that Eq. 9 gives a better fit to Eq. 8 if $`A0.3`$ and $`a0.79`$. Viceversa a smaller value of $`a`$ ($`a0.63`$) and $`A=1.08`$ in Eq. 8 gives a better fit to Eq. 9 (with $`A_1=0.3222`$ and $`a=0.707`$), which was the one ST used to compare model and data.

In the case of the “conditional” mass function, ST showed that an approximation can be obtained making the replacements $`BB(s)B(S)`$ and $`SsS`$ in Eq. (5), (6), that means:

$$f(s|S)ds=|T(s|S)|\mathrm{exp}(\frac{\left[B(s)B(S)\right]^2}{2(sS)})\frac{ds/(sS)}{\sqrt{2\pi (sS)}}$$
(10)
$$T(s|S)=\underset{n=0}{\overset{5}{}}\frac{(Ss)^n}{n!}\frac{^n\left[B(s)B(S)\right]}{s^n}$$
(11)
where as previously reported:

$$N(m,\delta _1|M,\delta _o)dm\frac{M}{m}f(m,\delta _1|M,\delta _o)dm$$
(12)
where $`f(m|M)dm=f(s|S)ds`$, $`\delta _1=\delta (z_1)`$ and $`\delta _0=\delta (z_0)`$ .

Thus, given Eqs. (5)-(6), (10)-(11), it is possible to obtain both the “unconditional” and “conditional” mass function, if the barrier shape and the power spectrum are given. In the following, I’ll use the barrier obtained in Del Popolo & Gambera (1998) to get the mass functions, which shall be compared with those obtained by SMT, ST and with numerical simulations, for several cosmologies. Since the way the barrier is obtained is described in previous papers (see Del Popolo & Gambera 1998, 1999, 2000) the reader is referred to those papers for details. Assuming that the barrier is proportional to the threshold for the collapse, similarly to ST, the barrier can be expressed in the form:

$$B(M)=\delta _\mathrm{c}=\delta _{\mathrm{co}}\left[1+_{r_\mathrm{i}}^{r_{\mathrm{ta}}}\frac{r_{\mathrm{ta}}L^2\mathrm{d}r}{GM^3r^3}\right]\delta _{\mathrm{co}}\left[1+\frac{\beta _1}{\nu ^{\alpha _1}}\right]$$
(13)
where $`\delta _{\mathrm{co}}=1.68`$ is the critical threshold for a spherical model, $`r_\mathrm{i}`$ is the initial radius, $`r_{\mathrm{ta}}`$ is the turn-around radius, $`L`$ the angular momentum, $`\alpha _1=0.585`$ and $`\beta _1=0.46`$. The angular momentum appearing in Eq. (13) is the total angular momentum acquired by the proto-structure during evolution. In order to calculate $`L`$, I’ll use the same model as described in Del Popolo & Gambera (1998, 1999) (more hints on the model and some of the model limits can be found in Del Popolo, Ercan & Gambera (2001)). The CDM spectrum used in this paper is that of Bardeen et al. (1986)(equation (G3)), with transfer function:

$$T(k)=\frac{[\mathrm{ln}\left(1+2.34q\right)]}{2.34q}[1+3.89q+(16.1q)^2+(5.46q)^3+(6.71)^4]^{1/4}$$
(14)
where $`q=\frac{k\theta ^{1/2}}{\mathrm{\Omega }_\mathrm{X}h^2\mathrm{Mpc}^1}`$. Here $`\theta =\rho _{\mathrm{er}}/(1.68\rho _\gamma )`$ represents the ratio of the energy density in relativistic particles to that in photons ($`\theta =1`$ corresponds to photons and three flavors of relativistic neutrinos). The power spectrum was normalized to reproduce the observed abundance of rich cluster of galaxies (e.g., Bahcal & Fan 1998). The barrier given in Eq. (13), differently from that of the spherical collapse is mass dependent. A direct comparison of the threshold given in Eq. (13) and that given in SMT (Eq. 4) is shown in Fig. 1. The dashed line represents $`\delta _\mathrm{c}(\nu )`$ obtained with the present model, while the solid line that of SMT. Both models show that the threshold for collapse decreases with mass, or similarly it increases with $`\sigma `$ since this quantity is a decreasing function of mass.

In other words, this means that, in order to form structure, more massive peaks must cross a lower threshold, $`\delta _c(\nu )`$, with respect to under-dense ones. At the same time, since the probability to find high peaks is larger in more dense regions, this means that, statistically, in order to form structure, peaks in more dense regions may have a lower value of the threshold, $`\delta _c(\nu )`$, with respect to those of under-dense regions. This is due to the fact that less massive objects are more influenced by external tides, and consequently they must be more overdense to collapse by a given time. In fact, the angular momentum acquired by a shell centred on a peak in the CDM density distribution is anti-correlated with density: high-density peaks acquire less angular momentum than low-density peaks (Hoffman 1986; Ryden 1988). A larger amount of angular momentum acquired by low-density peaks (with respect to the high-density ones) implies that these peaks can more easily resist gravitational collapse and consequently it is more difficult for them to form structure. Therefore, on small scales, where the shear is statistically greater, structures need, on average, a higher density contrast to collapse.

Putting Eq. (13) into Eqs. (5)-(6) and truncating the expansion at $`n=5`$, I get the “unconditional” mass function which can be approximated by:

$$f(\nu )d\nu 1.21\left(1+\frac{0.06}{\left(a\nu \right)^{0.585}}\right)\sqrt{\frac{a\nu }{2\pi }}\mathrm{exp}\{a\nu \left[1+\frac{0.57}{\left(a\nu \right)^{0.585}}\right]^2/2\}$$
(15)
where $`a=0.707`$.

In a similar way it is possible to obtain the “conditional” mass function (putting this time Eq. 13 into Eqs. (10)-(11)).

## 3 The temperature function

The haloes mass functions of clusters of galaxies contains information on the structure formation history of the universe: they are a primary input for modeling galaxy formation. The mass function is also a critical ingredient in putting strong constraints on cosmological parameters (principally $`\mathrm{\Omega }_0`$ and $`\mathrm{\Lambda }`$). Observationally the local mass function has been derived from measuring masses of individual clusters from galaxy velocity dispersions or other optical properties by Bahcall and Cen (1993), Biviano et al. (1993), and Girardi et al. (1998). However, the estimated virial masses for individual clusters depend rather strongly on model assumptions. As argued by Evrard et al. (1997) on the basis of hydrodynamical N-body simulations, cluster masses may be presently more accurately determined from a temperature measurement and a mass-temperature relation determined from detailed observations or numerical modeling. Thus alternatively, as a well-defined observational quantity, the X-ray temperature function (XTF) has been measured, which can be converted to the MF by means of the mass-temperature relation. Recent observational improvement on the XTF was made by Markevitch (1998), Henry (2000), Blanchard et al. (2000), and Pierpaoli et al. (2001), using more accurate temperature-measurement results for each cluster with ASCA data (Tanaka et al. 1994).

The evolution of the temperature distribution is very sensitive to the adopted model of structure formation, and its evolution at moderate redshift is considered a crucial test for cosmological models (Kitayama & Suto 1996). The cluster temperature function is defined as:

$$N(T,z)=N(M,z)\frac{dM}{dT}$$
(16)
While the mass function, $`N(M,z)`$, gives the mass and distribution of a population of evolving clusters, the Jacobian $`\frac{dM}{dT}`$ describes the physical properties of the single cluster.

Comparison of the predictions of the PS theory with the SCDM and OCDM cosmologies, performed by Tozzi & Governato (1998) and Governato et al. (1999), have shown discrepancies between PS predictions and N-body simulations, increasing with increasing $`z`$.

In the following, I’ll use the mass function modified as described in the previous section and an improved form of the M-T relation in order to calculate the mass function.

As previously reported, in order to estimate the multiplicity function of real systems one needs to know the temperature-mass (M-T) relation in order to transform the mass distribution into the temperature distribution.

Theoretical uncertainty arises in this transformation because the exact relation between the mass appearing in the PS expression and the temperature of the intra-cluster gas is unknown. Under the standard assumption of the Intra-Cluster (IC) gas in hydrostatic equilibrium with the potential well of a spherically symmetric, virialized cluster, the IC gas temperature-mass relation is easily obtained by applying the virial theorem and for a flat matter-dominated Universe it is given by (Evrard 1990, Evrard et al. 1996, Evrard 1997):

$$T=(6.4h^{2/3}keV)\left(\frac{M}{10^{15}M_{}}\right)^{2/3}(1+z)$$
(17)
The assumptions of perfect hydrostatic equilibrium and virialization are in reality not completely satisfied in the case of clusters. Clusters profile may depart from isothermality, with slight temperature gradients throughout the cluster (Komatsu & Seljak 2001). The X-ray weighted temperature can be slightly different from the mean mass weighted virial temperature. A noteworthy drawback of previous analyses has been stressed by Voit & Donahue (1998) (hereafter V98) and Voit (2000) (hereafter V2000). Using the merging-halo formalism of Lacey & Cole (1993), which accounts for the fact that massive clusters accrete matter quasi-continuously, they showed that the M-T relation evolves, with time, more modestly than what expected in previous models predicting $`T(1+z)`$, and this evolution is even more modest in open universes. Moreover, recent studies have shown that the self-similarity in the M-T relation seems to break at some keV (Nevalainen, Markevitch & Forman (hereafter NMF); Xu, Jin & Wu 2001). By means of ASCA data, using a small sample of 9 clusters (6 at 4 keV and 3 at $`1`$ keV), NMF has shown that $`M_{\mathrm{tot}}T_\mathrm{X}^{1.79\pm 0.14}`$ for the whole sample, and $`M_{\mathrm{tot}}T_\mathrm{X}^{3/2}`$ excluding the low-temperature clusters. Xu, Jin & Wu (2001) has found $`M_{\mathrm{tot}}T_\mathrm{X}^{1.60\pm 0.04}`$ (using the $`\beta `$ model), and $`M_{\mathrm{tot}}T_\mathrm{X}^{1.81\pm 0.14}`$ by means of the Navarro, Frenk & White (1995) profile. Finoguenov, Reiprich & Böeringer (2001), have investigated the T-M relation in the low-mass end finding that $`MT^2`$, and $`MT^{3/2}`$ at the high mass end. This behavior has been attributed to the effect of the formation redshift (Finoguenov, Reiprich & Böeringer 2001) (but see Mathiesen 2001 for a different point of view), or to cooling processes (Muanwong et al. 2001) and heating (Bialek, Evrard & Mohr 2000). Afshordi & Cen (2001) (hereafter AC) have shown that non-sphericity introduces an asymmetric, mass dependent, scatter for the M-T relation altering its slope at the low mass end ($`T3`$ keV).
Clearly this has effects on the final shape of the temperature function. In the following, I’ll use a modified version of the M-T relation obtained improving V98, V2000, to take account of tidal interaction between clusters. This M-T relation is given by:

$$kT8keV\left(\frac{M^{\frac{2}{3}}}{10^{15}h^1M_{}}\right)\frac{\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t}\right)^{\frac{2}{3}}+\frac{K(m,x)}{M^{8/3}}\right]}{\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t_0}\right)^{\frac{2}{3}}+\frac{K_0(m,x)}{M_0^{8/3}}\right]}$$
(18)
(see Appendix for a derivation), where $`t_\mathrm{\Omega }=\frac{\pi \mathrm{\Omega }_0}{H_o\left(1\mathrm{\Omega }_0\mathrm{\Omega }_\mathrm{\Lambda }\right)^{\frac{3}{2}}}`$, $`m=5/(n+3)`$ (being $`n`$ the spectral index), and:

$`K(m,x)`$ $`=`$ $`Fx\left(m1\right)\mathrm{𝐿𝑒𝑟𝑐ℎ𝑃ℎ𝑖}(x,1,3m/5+1)`$ (19)
$`F\left(m1\right)\mathrm{𝐿𝑒𝑟𝑐ℎ𝑃ℎ𝑖}(x,1,3m/5)`$
where $`F`$ is defined in the Appendix and the LerchPhi function is defined as follows:

$$LerchPhi(z,a,v)=\underset{n=0}{\overset{\mathrm{}}{}}\frac{z^n}{(v+n)^a}$$
(20)
and where $`K_0(m,x)`$ indicates that $`K(m,x)`$ must be calculated assuming $`t=t_0`$.

Eq. (18) accounts for the fact that massive clusters accrete matter quasi-continuously, and takes account of tidal interaction between clusters. The obtained M-T relation is no longer self-similar, a break in the low mass end ($`T34\mathrm{k}\mathrm{e}\mathrm{V}`$) of the M-T relation is present. The behavior of the M-T relation is as usual, $`MT^{3/2}`$, at the high mass end, and $`MT^\gamma `$, with a value of $`\gamma >3/2`$ in dependence of the chosen cosmology. Larger values of $`\gamma `$ are related to open cosmologies, while $`\mathrm{\Lambda }`$CDM cosmologies give results of the slope intermediate between the flat case and the open case.

With the previous prescriptions, I am going to calculate the cumulative temperature function, $`N(>kT)`$, and to compare it to simulations.

## 4 Results

In this section, I compare the “unconditional” and “conditional” mass functions with measurements in numerical simulations. The unconditional mass functions in Fig. 2 are taken from SMT, whereas the measurements in Figs. 3 and 4 are from Tozzi & Governato (1998). The “conditional” mass functions in Figs. 5 and 6 are reproduced from ST. The ST and SMT measurements were made in simulations which are a subset of those which were kindly made available to the public by the Virgo collaboration.

In Fig. 2, I plot the “unconditional” mass function, obtained as described in the previous section, for a $`\mathrm{\Lambda }`$CDM model as a function of the scaled variable $`\nu `$ 444Even if not plotted, the “unconditional” mass function was calculted for SCDM and OCDM models. The results in these cases is very similar to the $`\mathrm{\Lambda }`$CDM model plotted and the discussion relative to this case is valid for the other two. The solid line in the figure represents the “unconditional” mass function obtained by SMT, the dashed line the same quantity obtained from the barrier given in Eq. (13). The PS prediction is not plotted since from one hand we are more interested in comparison with the Sheth & Tormen (1999), SMT model, and on the other hand the comparison with the PS formula was performed in Sheth & Tormen (1999) who showed that the PS prediction are not in agreement with simulation at both high and low $`\nu `$.

Filled triangles, open triangles, open squares, filled hexagons and open hexagons show respectively results for $`z=0`$, $`z=0.5`$, $`z=1`$, $`z=2`$ and $`z=4`$ in the so called GIF simulations, a joint effort of astrophysicists from Germany & Israel, (see Kauffmann et al. 1999; SMT).

The simulations refer to a flat model with non-zero cosmological constant, $`\mathrm{\Lambda }`$CDM model ($`\mathrm{\Omega }_0=0.3`$,$`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`h=0.7`$). Fig. 2 show that the Sheth & Tormen (1999) model and that of this paper gives a good description of the mass function at all output times, suggesting, as remarked by Sheth & Tormen (1999) and ST, that the dynamics of the collapse is sensitive to $`\nu `$, and not to the mass scale. So, in the excursion set model the “unconditional” mass function, when expressed as a function of $`\nu `$, is an universal function of $`\nu `$, independent of redshift, cosmology or initial power spectrum.

A comparison between the result of Sheth & Tormen (1999) and that of this paper shows that the model presented here predicts a slightly larger value of the mass function at $`0.1\nu 4`$, and even if not plotted in Fig. 2 555The situation is similar to that of Fig. 10. In Fig. 2, I plotted only values of $`\nu 0.1`$ to compare it directly to Fig. 2 of Sheth & Tormen (1999), a smaller value for small $`\nu `$ ($`\nu <0.1`$) (at $`\nu =0.01`$ there are differences of $`20\%`$). At large $`\nu `$, the mass function of this paper is smaller than that of Sheth & Tormen (1999). The good agreement between Sheth & Tormen (1999) model and that of the present paper is due to the similitude of the barriers of the two papers. In both two, the barrier increases with $`S`$ differently from other models (see Monaco 1997a, b). It is interesting to note that the increasing of the barrier with $`S`$ has several important consequences and these models have a richer structure than the constant barrier model. In the case of non-spherical collapse with increasing barrier, a small fraction of the mass in the universe remains unbound, while for the spherical dynamics, at the given time, all the mass is bound up in collapsed objects. Moreover, incorporating the non-spherical collapse with increasing barrier in the excursion set approach results in a model in which fragmentation and mergers may occur (ST). If the barrier decreases with $`S`$ (Monaco 1997 a,b), this implies that all walks are guaranteed to cross it and so there is no fragmentation associated with this barrier shape.

As expressed in Eq. (3), the comoving number density of dark matter haloes of mass $`m`$ in the interval $`dM`$ and at redshift $`z`$ is connected to $`\nu f(\nu )`$. In the case of a constant barrier, it can be explicitly expressed as (PS):

$$n(m,z)=\sqrt{\frac{2}{\pi }}\frac{\rho }{m^2}\frac{\delta _c}{\sigma D(z)}|\frac{d\mathrm{ln}\sigma }{d\mathrm{ln}m}|\mathrm{exp}(\frac{\delta _c^2}{\sigma ^2D(z)})$$
(21)
where $`D(z)`$ is the linear growth factor normalized to unity at $`z=0`$ (Peebles 1993), and as previously reported $`\delta _c`$ is the linearly evolved density contrast of fluctuations that are virializing at $`z=0`$ ($`\delta _c`$=1.686 for $`\mathrm{\Omega }_0=1`$ and 1.65 for the case $`\mathrm{\Omega }_0=0.3`$). In the case of the barrier of Eq. (13), we have:

$$n(m,z)1.21\frac{\overline{\rho }}{m^2}\frac{d\mathrm{log}(\nu )}{d\mathrm{log}m}\left(1+\frac{0.06}{\left(a\nu \right)^{0.585}}\right)\sqrt{\frac{a\nu }{2\pi }}\mathrm{exp}\{a\nu \left[1+\frac{0.57}{\left(a\nu \right)^{0.585}}\right]^2/2\}$$
(22)
In Fig. 3, I compare the mass function for a SCDM model obtained by means of numerical simulations (Tozzi & Governato 1998) at different redshift (open squares, $`z=0`$; filled squares, $`z=0.5`$; filled triangles, $`z=1`$), with the prediction of the PS formula for the same redshifts (solid line, $`z=0`$; long dashed line, $`z=0.5`$; short dashed line, $`z=1`$). The plot shows a discrepancy between the numerical simulations and the PS predictions. In particular PS underestimates the number density of high mass clusters ($`M>>M_{}`$), discrepancy that grows at higher redshift. For example, at $`z=1`$ the number density is underestimated by a factor $`3`$ for $`M>3\times 10^{14}M_{}`$. The excess of massive clusters at higher redshifts suggests that the cluster mass function for the SCDM model evolves more slowly than the PS mass function. The deficit of low mass haloes ($`M<<M_{}`$) documented in numerical works, e.g. Carlberg & Couchman (1998), has been associated with merger events not accounted for by the PS formalism. As previously reported, incorporating the non-spherical collapse with increasing barrier in the excursion set approach results in a model in which fragmentation and mergers may occur, differently from the PS model. The PS mass function depends sensitively on $`\delta _\mathrm{c}`$ and only changing its value and the $`z`$ dependence it is possible to fit the numerical simulations with the PS formula (Governato et al. 1999).

Fig. 4 gives the same information of Fig. 3, but now the lines represent the predictions of Eq. (22). The plot shows a good agreement with the numerical simulations. So the new mass function is in much better agreement with simulations in agreement with what seen in Fig. 2.

As previously described, an approximation to the conditional mass function can be obtained from Eq. (10), (11). As stressed by ST, the analytical approximation should be more accurate in the case of high-redshift progenitors of massive parents, and less accurate when the parents are not very massive. Fig. 5 plots the “conditional” mass function for a SCDM model ($`\mathrm{\Omega }_0=1`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$, $`h=0.5`$) at a redshift $`z=0.5`$ (left panel) and $`z=4`$ (right panel). The solid line represents the Sheth & Tormen prediction while the dashed line that of the present paper. The symbols show parent haloes with mass in the $`16M/M_{}32`$ range (filled hexagons) at $`z=0.5`$. Fig. 5, shows that the excursion set predictions with the non-spherical collapse are in good agreement with simulations, while the spherical collapse disagrees with simulations both at small and high values of $`m`$ and for different redshifts. A comparison between ST and that of the present paper shows a slight difference between the conditional mass functions, difference which is larger at small $`m`$.

The difference between the mass functions reflects the difference, shown in Fig. 1, between the barrier obtained by SMT and that obtained in Del Popolo & Gambera (1998) and used in this paper (see the following of this section for more insight on this point).

The previous arguments have shown that the excursion set approach with non-spherical collapse gives a good description of the conditional and unconditional mass function. In particular the unconditional mass function is an universal function of $`\nu `$. In the case of the conditional mass function the situation is different, as shown in Figs. 6. In Fig. 6, I plot the rescaled “conditional” mass function for a $`\mathrm{\Lambda }`$CDM model ($`\mathrm{\Omega }_0=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$, $`h=0.7`$) at a redshift $`z=0.5`$ (left panel) and $`z=4`$ (right panel). The solid line represents the Sheth & Tormen prediction while the dashed line that of the present paper. The symbols show the result of rescaling the conditional mass functions in the simulations for parent haloes with mass in the range 1-2 $`M_{}`$ (filled hexagons), and 8-32 $`M_{}`$ (open triangles) at $`z=0`$. As stressed by ST, the symbols in Fig. 6 do not overlap exactly, at fixed $`z`$, the conditional mass functions for different parent haloes do not rescale exactly, and moreover the mass functions at different output times do not rescale either (the band traced out by symbols at $`z=0.5`$ is different from that traced out at $`z=0.5`$). In other words, the conditional mass function is not a universal function of $`\nu `$. A comparison between the simulations (symbols) and the theoretical predictions, (curves), shows that they can provide a good fit only at high redshift (Fig. 6 (right panel)). The situation, even if not represented, is similar for the SCDM model. The reasons of the discrepancy should be the following: one possibility is the neglect of correlations between scales (Peacock & Heavens 1990; Bond et al. 1991; ST), the second is that the parametrization of the collapse described in SMT is too simple and consequently the collapse threshold and the barrier is also too simple. In the model of this paper the parametrization of the collapse is slightly different from that of SMT. In this paper, the spherical collapse is modified to take account of the tidal interaction between protoclusters (see Del Popolo & Gambera 1998 and the appendix of this paper). With the parametrization of this paper the conditional mass function at small $`\nu `$ is lower than SMT (at $`\nu =0.01`$, the difference is $`20\%`$) and the curve gives a better fit of the conditional mass function. At large $`\nu `$ there is very little difference between the two predictions. So the plots show that changing the barrier shape influences the conditional mass function especially at small lookback times. However, if the parametrization of the collapse may have a certain role in final shape of the conditional mass function, the role of correlations between scales is surely important especially at small lookback times (ST). In fact the excursion set approach describes in a reasonably good way the clustering at high redshift, but less accurately at small redshifts. The reason is due to the fact that at large lookback times the largest part of subclumps constitutes a small fraction of the mass of the parent halo. Neglect of correlations between scales introduces a small error since the smoothing scale associated with subclumps is sufficiently different from that of the parent. At smaller lookback times the situation changes, since subclumps and parent are not so well separated and neglecting the correlations introduces a larger error with respect to the previous case (ST). Before going on, it is important to stress that the conditional mass functions in the simulations are not yet sufficiently well determined to make strong statements about whether the simulations prefer the SMT or the barrier predictions of the present paper. Higher resolution simulations are required to see whether one or the other is more accurate. Finally, I have studied the cumulative temperature function. The function I obtained and which I am going to describe differs from predictions based upon the PS formula for two reasons: a) the cluster density is obtained through the excursion set approach but for a non-spherical collapse and a non-constant barrier; b) as described in Sect. (3) the M-T relation is different from the usual self-similar relation $`TM^{2/3}`$. Fig. 7 represents the cumulative temperature function for a SCDM model ($`\mathrm{\Omega }_0=1`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0`$, $`h=0.5`$, $`z=0`$ (left panel) and $`z=0.4`$ (central panel)), and the accretion law for halos of $`5\times 10^{14}h^1M_{}`$ (right panel). The solid line represents the model of this paper for the temperature function, while the open squares Tozzi & Governato (1998) simulations. The dotted lines the PS prediction with $`\delta _\mathrm{c}`$ rescaled as done in Tozzi & Governato (1998), to fit the mass function ($`\delta _\mathrm{c}(z)=1.48\times (1+z)^{0.06}`$). From the plots, it is clear that the PS predicts less massive clusters than simulations. The discrepancy increases with increasing $`z`$ and $`T`$: at $`T10\mathrm{k}\mathrm{e}\mathrm{V}`$ and $`z=0`$ the discrepancy is less than a factor of 3 while at $`T10\mathrm{k}\mathrm{e}\mathrm{V}`$ and $`z=0.4`$, the discrepancy amounts almost to an order of magnitude. This implies that the cumulative temperature function obtained from simulations evolves more slowly than the PS prediction, in the range of $`z`$ and $`T`$ studied. The situation becomes worse at larger $`z`$ and $`T`$. This means that trying to rule out or accept a cosmological model exclusively on the basis of comparisons between data and the PS formula should be taken with caution. As Fig. 7 shows, the corrections introduced by using Eq. (22) for the mass function and Eq. (18) for the M-T relation, noteworthy improves the agreement with numerical simulations. The approach proposed is thus much more realible than the PS model. The right panel of Fig. 7 compares the accretion law calculated according to V98, V2000 model (used in the calculation of the M-T relation) and that obtained by Wechsler at al. 2002. As reported in the Appendix, the accretion law used in the calculation of the M-T relation is that obtained by V98 following Lacey & Cole (1993) prescriptions, namely it is obtained by using the approximation:

$$\frac{dS}{d\omega }\frac{S}{\omega }$$
(23)
(V98), where $`S\sigma ^2(m)`$ and $`\omega (t)`$ is given in V2000 (see also the Appendix). As described in the Appendix, in the case the fluctuation amplitudes can be described by a power-law, we have that:

$$M\omega ^{\frac{3}{n+3}}$$
(24)
Simulations by Wechsler et al. 2002 show that the mass growth rate is better described by an exponential form:

$$M(a)exp(\alpha z),a=1/(1+z)$$
(25)
(Wechsler et al. 2002). I compared the two different growth rates in the right panel of Fig. 7, for the case of halos of mass $`5\times 10^{14}h^1M_{}`$. The solid line plots Eq. (25), while the dashed line represents the prediction of V98, V2000 model. As shown in the plot, Wechsler’s prediction is different from V98, V2000 and then from that of this paper. The difference is due to the fact that I used (as in V98 and V2000) the Extended Press-Schechter (EPS) formalism to evaluate the accretion rate, which can reproduce relative properties of the progenitor halo distribution quite well, but has some difficulties in estimating absolute progenitor masses and the overall conditional mass function (Tormen 1998; Sommerville et al. 2000; Gardner 2001). This is one of the reasons that led many authors to propose alternative expressions for the PS mass function (e.g., SMT, ST). However, note that in the calculation of the cumulative temperature function two factors play a role: the expression for the mass function and the M-T relation. This last introduce a correction of “second order” in the temperature function when compared with the correction introduced by the new form of mass function used (see also Del Popolo & Gambera 1999).

## 5 Conclusions

In this paper, I calculated the unconditional, conditional mass function and temperature function by using the extension of the excursion set model of ST and the barrier shape obtained in Del Popolo & Gambera (1988). I showed that the barrier obtained in Del Popolo & Gambera (1998), which takes account of asphericity and tidal interaction between proto-haloes, is a better description of the mass functions and temperature function than the spherical collapse and is in good agreement with numerical simulations. The results are in good agreement with those obtained by ST, only some differences are observed expecially at the low mass end.

The main results of the paper can be summarized as follows:
1) the non-constant barrier obtained from the non-spherical collapse in Del Popolo & Gambera (1998), taking account of the tidal interaction of proto-clusters with neighboring ones, combined with the ST model gives “unconditional” and “conditional” mass functions, is in reasonably good agreement with results from numerical cosmological simulations.
2) The “unconditional” and “conditional” mass functions obtained with the Del Popolo & Gambera (1998) barrier are slightly different from those obtained by SMT and ST: smaller values at small and large $`\nu `$, with respect to SMT and ST predictions, and larger values for $`0.1\nu 4`$. In the case of the rescaled “conditional” mass function the discrepancy observed by ST for small lookback time is smaller in the model of this paper at small $`\nu `$ ($`0.01\nu 1`$). The discrepancy with simulations is connected to SMT parametrization of the collapse and not only to the neglect of correlations between scales.
3) The mass function in SCDM 666Similar conclusion is valid for the OCDM. is in good agreement with Governato et al. (1999) and Tozzi & Governato (1998) simulations.
4) The temperature function calculated by means of the mass functions obtained in the present paper together with an improved version of the M-T relation, 777As described in the text, the new M-T relation accounts for the fact that massive clusters accrete matter quasi-continuously, and consequently that the M-T relation evolves, with time, more modestly than what expected in previous models (top-hat model) and taking account of the tidal interaction with neighboring clusters is in good agreement with the simulations of SCDM universes of Tozzi & Governato (1998) for different redshifts.
5) The ST formulae really do work for different barrier shapes (at least that used in this paper, that introduced in SMT and that of Monaco (1997a,b)).
6) The behavior of the “unconditional” mass function at small masses is similar to that of Sheth & Tormen (1999), ST, and very different from that proposed by Jenkins et al. (2001) (see ST Fig. 13).

The above considerations show that it is possible to get accurate predictions for a number of statistical quantities associated with the formation and clustering of dark matter haloes by incorporating a non-spherical collapse in the excursion set approach. The improvement is probably connected also to the fact that incorporating the non-spherical collapse with increasing barrier in the excursion set approach results in a model in which fragmentation and mergers may occur, effects important in structure formation.

## 6 Appendix

From the simulations and observations, we know that the mass within a specified density contrast is straightforwardly related to temperature. Two fundamental steps in obtaining the M-T relation are: a) to approximate cluster formation with the evolution of a spherical top-hat perturbation (e.g., Peebles 1993); b) to assume that each cluster we see at a given redshift $`z`$ has just reached the moment of virialization, an assumption known as the late-formation approximation, then $`M_{\mathrm{vir}}T_\mathrm{X}^{3/2}\rho _{\mathrm{cr}}^{1/2}\mathrm{\Delta }_{\mathrm{vir}}^{1/2}`$. Some shortcomings of this approach has been summarized in the previous section (see also Viana & Liddle 1996; Kitayama & Suto 1996; Eke et al. 1996). The late-formation approximation is a good one for many purposes, but a better one can be obtained in the low-$`\mathrm{\Omega }`$ limit. As can be found in the literature, there are two ways of improving the quoted model. One is to define a formation redshift $`z_\mathrm{f}`$ at which a cluster virializes and after the properties of observed clusters at $`z`$ are obtained by integrating over the appropriate distribution of formation redshifts (Kitayama & Suto 1996; Viana & Liddle 1996). The second possibility is that described by V98, V2000. In this approach, the top-hat cluster formation model is substituted by a model of cluster formation from spherically symmetric perturbations with negative radial density gradients. The fact that clusters form gradually, and not instantaneously, is taken into account in the merging-halo formalism of Lacey & Cole (1993). In hierarchical models for structure formation, the growth of the largest clusters is quasi-continuous since these large objects are so rare that they almost never merge with another cluster of similar size (Lacey & Cole 1993). So, Lacey & Cole (1993) approach extends the PS formalism by considering how clusters grow via accretion of smaller virialized objects. To start with, I’ll write the equation governing the collapse of a density perturbation taking account of angular momentum acquisition by proto-structures (see Peebles 1993; Del Popolo & Gambera 1998, 1999; V2000). In the model, the radial acceleration of the particle is:

$$\frac{dv_r}{dt}=\frac{L^2(r)}{M^2r^3}g(r)=\frac{L^2(r)}{M^2r^3}\frac{GM}{r^2}$$
(26)
Assuming a non-zero cosmological constant Eq. (26) becomes:

$$\frac{dv_r}{dt}=\frac{GM}{r^2}+\frac{L^2(r)}{M^2r^3}+\frac{\mathrm{\Lambda }}{3}r$$
(27)
(Peebles 1993; Bartlett & Silk 1993; Lahav 1991; Del Popolo & Gambera 1998, 1999). Integrating Eq. (27) we have:

$$\frac{1}{2}\left(\frac{dr}{dt}\right)^2=\frac{GM}{r}+\frac{L^2}{M^2r^3}𝑑r+\frac{\mathrm{\Lambda }}{6}r^2+ϵ$$
(28)
where the value of the specific binding energy of the shell, $`ϵ`$, can be obtained using the condition for turn-around, $`\frac{dr}{dt}=0`$.

Integrating Eq.(28), I get:

$$t=\frac{dr}{\sqrt{2\left[ϵ+\frac{GM}{r}+_{r_\mathrm{i}}^r\frac{L^2}{M^2r^3}𝑑r+\frac{\mathrm{\Lambda }}{6}r^2\right]}}$$
(29)
A particular shell will collapse if:

$$ϵ+\frac{GM}{r_0}+_{r_\mathrm{i}}^{r_0}\frac{L^2}{M^2r^3}𝑑r+\frac{\mathrm{\Lambda }}{6}r_0^2=0$$
(30)
(see V2000), and the shell reaches its maximum radius (turn-around radius, $`r_{\mathrm{ta}}`$) at a time:

$$t_{\mathrm{ta}}=_0^{r_0}\frac{dr}{\sqrt{2\left[ϵ+\frac{GM}{r}+_{r_\mathrm{i}}^r\frac{L^2}{M^2r^3}𝑑r+\frac{\mathrm{\Lambda }}{6}r^2\right]}}$$
(31)
Eq. (31) can be written in an equivalent form (see Del Popolo & Gambera 1998, 1999; Bartlett & Silk 1993), as:

$`t_{ta}`$ $`=`$ $`{\displaystyle _0^{r_{ta}}}{\displaystyle \frac{dr}{\sqrt{2\left[GM\left(\frac{1}{r}\frac{1}{r_{ta}}\right)+_{r_{ta}}^r\frac{L^2}{M^2r^3}𝑑r+\frac{\mathrm{\Lambda }}{6}(r^2r_{\mathrm{ta}}^2)\right]}}}=`$ (32)
$`2H_0^1\mathrm{\Omega }_0^{1/2}\xi ^{3/2}{\displaystyle _0^1}{\displaystyle \frac{y^2dy}{\sqrt{1+\left[\frac{9\xi ^9}{4\pi ^2\rho ^2r_{ta}^{10}\mathrm{\Omega }H^2}_1^{y^2r_{ta}}\frac{L^2}{y^5}𝑑y\right]y^2+\left(\frac{\mathrm{\Omega }_\mathrm{\Lambda }}{\mathrm{\Omega }_0}\right)\xi ^3y^6\left[1+\left(\frac{\mathrm{\Omega }_\mathrm{\Lambda }}{\mathrm{\Omega }_0}\right)\xi ^3\right]y^2}}}`$
where $`\xi `$ is a parameter given by $`\xi =r_{ta}/x_1`$, $`r_{ta}`$ being the radius of the turn-around epoch, while $`x_1`$ is defined by the relation $`M=4\pi \rho _\mathrm{b}x_1^3/3`$, $`\mathrm{\Omega }_0=\frac{8\pi G\rho _\mathrm{b}}{3H_0^2}`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=\frac{\mathrm{\Lambda }}{3H_0^2}=1\mathrm{\Omega }_0`$, and $`\rho _\mathrm{b}`$ the critical density. The shell collapses at $`t_\mathrm{c}=2t_{\mathrm{ta}}`$. In the case $`L=0`$, Eq. (31) can be analytically integrated and inverted to get $`ϵ(t)`$:

$$ϵ(t)=\frac{1}{2}\left(\frac{2\pi GM}{t}\right)^{2/3}=\frac{1}{2}\left(\frac{2\pi GM}{t_\mathrm{\Omega }}\right)^{2/3}\left(\frac{t_\mathrm{\Omega }}{t}\right)^{2/3}=\frac{1}{2}\left(\frac{2\pi GM}{t_\mathrm{\Omega }}\right)^{2/3}(x1)$$
(33)
where we have defined $`t_\mathrm{\Omega }=\frac{\pi \mathrm{\Omega }_0}{H_o\left(1\mathrm{\Omega }_0\mathrm{\Omega }_\mathrm{\Lambda }\right)^{\frac{3}{2}}}`$ and $`x=1+(\frac{t_\mathrm{\Omega }}{t})^{2/3}`$ which is connected to mass by $`M=M_0x^{3m/5}`$, where $`m=5/(n+3)`$ and $`n`$ is the usual power-law perturbation index (V2000).

Integrating with respect to mass $`ϵ(t)`$, and diving by mass, one gets, as shown by V2000:

$$\frac{E}{M}=\frac{ϵ𝑑M}{M}=\frac{3m}{10(m1)}\left(\frac{2\pi G}{t_\mathrm{\Omega }}\right)^{\frac{2}{3}}M^{\frac{2}{3}}\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t}\right)^{\frac{2}{3}}\right]$$
(34)
In the case $`L0`$, it is numerically possible to integrate and invert Eq. (31) finding the specific energy, which we can indicate with $`ϵ_\mathrm{L}(t)`$ and finally:

$$\frac{E}{M}=\frac{3}{5}\frac{m}{m1}ϵ_\mathrm{L}(t)$$
(35)
An approximate relation for $`ϵ_\mathrm{L}(t)`$, in an Einstein-de Sitter Universe can be obtained as follows. It is possible to obtain the turn-around radius, $`r_{\mathrm{ta}}`$ by solving Eq. (32) (with $`\mathrm{\Lambda }=0`$) for a given mass and a given epoch of interest. This is related to the binding energy of the shell enclosing the mass $`M`$ by Eq. (28) with $`\dot{r}=0`$. In turn, the binding energy is uniquely given by the linear overdensity $`\delta _\mathrm{i}`$ at some arbitrary early time. We may now use the linear theory to “propagate” the overdensity to the “chosen” time to find the linear overdensity at turn-around, $`\delta _\mathrm{c}`$. Using the relation between $`v`$ and $`\delta _\mathrm{i}`$ for the growing mode (Peebles 1980) in Eq. (28), with $`\mathrm{\Lambda }=0`$, at an early time, it is possible to connect $`ϵ`$ and $`\delta _c`$ (see Bartlett & Silk 1993)):

$$ϵ=\frac{5}{6}\mathrm{\Omega }_0H_0^2x_1^2(\delta _\mathrm{i}/a_\mathrm{i})=\frac{1}{2}\left(\frac{2\pi GM}{t}\right)^{2/3}$$
(36)
Using the formula for $`\delta _\mathrm{c}`$ given in Del Popolo & Gambera (1999, Eq. 14), Del Popolo & Gambera (2000, Eq. 6), is possible to write:

$$ϵ_\mathrm{L}=\frac{1}{2}\left(\frac{2\pi GM}{t}\right)^{2/3}\delta _\mathrm{c}=\frac{1}{2}\left(\frac{2\pi GM}{t}\right)^{2/3}\left[1+\frac{r_{\mathrm{ta}}}{GM^3}_0^r\frac{L^2dr}{r^3}\right]$$
(37)
which reduces to Eq. (36) when $`L0`$.

Eq. (37) can also be written as:

$$ϵ_\mathrm{L}=\frac{1}{2}\left(\frac{2\pi GM}{t_\mathrm{\Omega }}\right)^{2/3}(x1)\left[1+\frac{2^{7/3}\pi ^{2/3}\xi \rho _\mathrm{b}^{2/3}}{3^{2/3}H^2\mathrm{\Omega }M^{8/3}}\frac{1}{x1}_0^r\frac{L^2dr}{r^3}\right]$$
(38)
Defining $`M=M_0x^{3m/5}`$ and

$$F=\frac{2^{7/3}\pi ^{2/3}\xi \rho _\mathrm{b}^{2/3}}{3^{2/3}H^2\mathrm{\Omega }}_0^r\frac{L^2dr}{r^3}$$
(39)
Eq. (38) becomes:

$$ϵ_\mathrm{L}=\frac{1}{2}\left(\frac{2\pi G}{t_\mathrm{\Omega }}\right)^{\frac{2}{3}}M^{2/3}\left[\left(\frac{M}{M_o}\right)^{\frac{5}{3m}}1\right]\left[1+\frac{F}{M^{\frac{8}{3}}}\frac{1}{\left(\frac{M}{M_o}\right)^{\frac{5}{3m}}1}\right]$$
(40)
Integrating with respect mass, and dividing again by $`M`$, I get:

$$\frac{E}{M}=\frac{ϵ_\mathrm{L}𝑑M}{M}=\frac{3m}{10(m1)}\left(\frac{2\pi G}{t_\mathrm{\Omega }}\right)^{\frac{2}{3}}M^{\frac{2}{3}}\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t}\right)^{\frac{2}{3}}+\frac{K(m,x)}{M^{8/3}}\right]$$
(41)
where

$`K(m,x)`$ $`=`$ $`Fx\left(m1\right)\mathrm{𝐿𝑒𝑟𝑐ℎ𝑃ℎ𝑖}(x,1,3m/5+1)`$ (42)
$`\left(m1\right)F\mathrm{𝐿𝑒𝑟𝑐ℎ𝑃ℎ𝑖}(x,1,3m/5)`$
where the LerchPhi function is defined as follows:

$$LerchPhi(z,a,v)=\underset{n=0}{\overset{\mathrm{}}{}}\frac{z^n}{(v+n)^a}$$
(43)
888This definition is valid for $`abs(z)<1`$. By analytic continuation, it is extended to the whole complex z-plane, for each value of a. If the coefficients of the series representation of a hypergeometric function are rational functions of the summation indices, then the hypergeometric function can be expressed as a linear sum of Lerch Phi functions. Reference: A. Erdelyi, 1953, Higher Transcendetal Functions, Volume 1, chapter 1, section 11.

If $`K=0`$, Eq. (41) reduces to Eq. (10) of V2000. As stressed by V2000, some factors give rise to an higher value of $`E/M`$ with respect the case of the late-formation value. The $`m/(m1)`$ value which accounts for the effect of early infall. The $`1/m`$ value in the square bracket of Eq. (41) which accounts for the cessation of cluster formation when $`t>>t_\mathrm{\Omega }`$. Finally in Eq. (41) a new term is present, which comes from the tidal interaction. In order to obtain an expression for the kinetic energy, starting from $`E/M`$, I use the virial theorem with the surface pressure term correction as in V2000 and utilizing the usual relation:

$$K=\frac{3\beta MkT}{2\mu m_\mathrm{p}}$$
(44)
(AC), where $`k`$ is the Boltzmann constant, $`\mu =0.59`$ is the mean molecular weight, $`m_\mathrm{p}`$ the proton mass and $`\beta =\frac{\sigma _\mathrm{v}^2}{kT/\mu m_\mathrm{p}}`$, being $`\sigma _\mathrm{v}`$ the mass-weighted mean velocity dispersion of dark matter particles. In this way, I finally get:

$$kT=\frac{2}{5}a\frac{\mu m_p}{2\beta }\frac{m}{m1}\left(\frac{2\pi G}{t_\mathrm{\Omega }}\right)^{\frac{2}{3}}M^{\frac{2}{3}}\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t}\right)^{\frac{2}{3}}+\frac{K(m,x)}{M_0^{8/3}}\right]$$
(45)
where $`a=\frac{\overline{\rho }}{2\rho (r_{\mathrm{vir}})\overline{\rho }}`$ is the ratio between kinetic and total energy (V2000) 999This term comes from the fact I used a modified version of the virial theorem in order to include a surface pressure term (see also V2000, AC). This correction is due to the fact that at the virial radius $`r_{\mathrm{vir}}`$ the density is non-zero and this requires a surface pressure term to be included in the virial theorem (Carlberg, Yee & Ellingson 1997) (the existence of this confining pressure is usually not accounted for in the top-hat collapse model). Using the relation $`\mathrm{\Delta }_{\mathrm{vir}}=\frac{8\pi ^2}{Ht^2}`$ (see V2000), and in the early-time limit: ($`t<<t_\mathrm{\Omega }`$), Eq. (45), reduces to:

$$kT=\frac{2}{5}\frac{m}{m1}a\frac{\mu m_p}{2\beta }GM^{\frac{2}{3}}\left(\frac{4\pi }{3}\rho _\mathrm{b}\mathrm{\Delta }_{\mathrm{vir}}\right)^{1/3}$$
(46)
which, in the case $`n2`$, $`a2`$ is identical to the late-formation formula, described in V2000 (see their Eq. (8)). Normalizing Eq. (45) similarly to V2000, I get:

$$kT8keV\left(\frac{M^{\frac{2}{3}}}{10^{15}h^1M_{}}\right)\frac{\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t}\right)^{\frac{2}{3}}+\frac{K(m,x)}{M^{8/3}}\right]}{\left[\frac{1}{m}+\left(\frac{t_\mathrm{\Omega }}{t_0}\right)^{\frac{2}{3}}+\frac{K_0(m,x)}{M_0^{8/3}}\right]}$$
(47)
where $`K_0(m,x)`$ indicates that $`K(m,x)`$ must be calculated assuming $`t=t_0`$