format: text

text:
# Populating Dark Matter Haloes with Galaxies: Comparing the 2dFGRS with Mock Galaxy Redshift Surveys

## 1 Introduction

The distribution of galaxies contains important information about the large scale structure of the matter distribution. On large, linear scales the galaxy power spectrum is believed to be proportional to the matter power spectrum, therewith providing useful information regarding the initial conditions of structure formation, i.e., regarding the power spectrum of primordial density fluctuations. On smaller, non-linear scales the distribution and motion of galaxies is governed by the local gravitational potential, which is cosmology dependent. One of the main goals of large galaxy redshift surveys is therefore to map the distribution of galaxies as accurately as possible, over as large a volume as possible. The Sloan Digital Sky Survey (SDSS; York et al. 2000) and the 2 degree Field Galaxy Redshift Survey (2dFGRS; Colless et al. 2001) are two of the prime examples. These surveys, which are currently being completed, will greatly enhance and improve our knowledge of large-scale structure and will become the standard data sets against which to test our cosmological and galaxy formation models for the decade to come.

However, two effects complicate a straightforward interpretation of the data. First of all, the distribution of galaxies is likely to be biased with respect to the underlying mass density distribution. This bias is an imprint of various complicated physical processes related to galaxy formation such as gas cooling, star formation, merging, tidal stripping and heating, and a variety of feedback processes. In fact, it is expected that the bias depends on scale, redshift, galaxy type, galaxy luminosity, etc. (Kauffmann, Nusser & Steinmetz 1997; Jing, Mo & Börner 1998; Somerville et al. 2001; van den Bosch, Yang & Mo 2003). Therefore, in order to translate the observed clustering of galaxies into a measure for the clustering of (dark) matter, one needs to either understand galaxy formation in detail, or use an alternative method to describe the relationship between galaxies and dark matter (haloes). One of the main goals of this paper is to advocate one such method and to show its potential strength for advancing our understanding of large scale structure.

Secondly, because of the peculiar velocities of galaxies, the clustering of galaxies observed in redshift-space is distorted with respect to the real-space clustering (e.g., Davis & Peebles 1983; Kaiser 1987; Regos & Geller 1991; van de Weygaert & van Kampen 1993; Hamilton 1992). On small scales, the virialized motion of galaxies within dark matter haloes smears out structure along the line-of-sight (i.e., the so-called “finger-of-God” effect). On large scales, coherent flows induced by the gravitational action of large scale structure enhance structure along the line-of-sight. Both effects cause an anisotropy in the two-dimensional, two-point correlation function $`\xi (r_p,\pi )`$, with $`r_p`$ and $`\pi `$ the pair separations perpendicular and parallel to the line-of-sight, respectively. The large-scale flows compress the contours of $`\xi (r_p,\pi )`$ in the $`\pi `$ direction by an amount that depends on $`\beta \mathrm{\Omega }_m^{0.6}/b`$. The small-scale peculiar motions implies that $`\xi (r_p,\pi )`$ is convolved in the $`\pi `$-direction by the distribution of pairwise velocities, $`f(v_{12})`$. Thus, the detailed structure of $`\xi (r_p,\pi )`$ contains information regarding the Universal matter density $`\mathrm{\Omega }_m`$, the (linear) bias of galaxies $`b`$, and the pairwise velocity distribution $`f(v_{12})`$.

From the above discussion it is obvious that understanding galaxy bias is an integral part of understanding large scale structure. One way to address galaxy bias without a detailed theory for how galaxies form is to model halo occupation statistics. One simply specifies halo occupation numbers, $`N(M)`$, which describe how many galaxies on average occupy a halo of mass $`M`$. Many recent investigations have used such halo occupation models to study various aspects of galaxy clustering (Jing, Mo & Börner 1998; Peacock & Smith 2000; Seljak 2000; Scoccimarro et al. 2001; White 2001; Jing, Börner & Suto 2002; Bullock, Wechsler & Somerville 2002; Berlind & Weinberg 2002; Scranton 2002; Kang et al. 2002; Marinoni & Hudson 2002; Zheng et al. 2002; Kochanek et al. 2003). In two recent papers, Yang, Mo & van den Bosch (2003; hereafter Paper I) and van den Bosch, Yang & Mo (2003; hereafter Paper II) have taken this halo occupation approach one step further by considering the occupation as a function of galaxy luminosity and type. They introduced the conditional luminosity function (hereafter CLF) $`\mathrm{\Phi }(L|M)\mathrm{d}L`$, which gives the number of galaxies with luminosities in the range $`L\pm \mathrm{d}L/2`$ that reside in haloes of mass $`M`$. The advantage of this CLF over the halo occupation function $`N(M)`$ is that it allows one to address the clustering properties of galaxies as function of luminosity. In addition, the CLF yields a direct link between the halo mass function and the galaxy luminosity function, and allows a straightforward computation of the average luminosity of galaxies residing in a halo of given mass. Therefore, $`\mathrm{\Phi }(L|M)`$ is not only constrained by the clustering properties of galaxies, as is the case with $`N(M)`$, but also by the observed LFs and the halo mass-to-light ratios.

In Papers I and II we used the observed LFs and the luminosity- and type-dependence of the galaxy two-point correlation function to constrain the CLF in the standard $`\mathrm{\Lambda }`$CDM cosmology. In this paper, we use this CLF to populate dark matter haloes in high-resolution $`N`$-body simulations. The ‘virtual Universes’ thus obtained are used to construct mock galaxy redshift surveys with volumes and apparent magnitude limits similar to those in the 2dFGRS. This is the first time that realistic mock surveys are constructed that (i) associate galaxies with dark matter haloes, (ii) are independent of a model for how galaxies form, and (iii) automatically have the correct galaxy abundances and correlation lengths as function of galaxy luminosity and type. In the past, mock galaxy redshift surveys were constructed either by associating galaxies with dark matter particles (rather than haloes) using a completely ad hoc bias scheme (Cole et al. 1998), or by linking semi-analytical models for galaxy formation (with all their associated uncertainties) to the merger histories of dark matter haloes derived from numerical simulations (Kauffmann et al. 1999; Mathis et al. 2002).

We use our mock galaxy redshift survey to investigate a number of statistical measures of the large scale distribution of galaxies. In particular, we focus on the two-point correlation function in redshift space, its distortions on small and large scales, and the galaxy pairwise peculiar velocities. Where possible we compare our predictions with the 2dFGRS and we discuss the sensitivity of these clustering statistics to several details regarding the halo occupation statistics. We show that the halo occupation obtained analytically can reliably be implemented in $`N`$-body simulations. We find that the standard $`\mathrm{\Lambda }`$CDM model, together with the halo occupation we have obtained, can reproduce many of the observational results. However, we find significant discrepancy between the model predictions and observations on small scales. We show that to get consistent results on small scales, either the mass-to-light ratios for clusters of galaxies are significantly higher than normally assumed, or the linear power spectrum has an amplitude that is significantly lower than its ‘concordance’ value.

This paper is organized as follows. In Section 2 we review the CLF formalism developed in papers I and II. Section 3 introduces the $`N`$-body simulations and describes our method of populating dark matter haloes in these simulations with galaxies of different type and luminosity. Section 4 investigates several clustering statistics in real-space and focuses on the accuracy with which mock galaxy distributions can be constructed using our CLF formalism. In Section 5 we use these mock galaxy distributions to construct mock galaxy redshift surveys that are comparable in size with the 2dFGRS. We extract the redshift-space two-point correlation function from this mock redshift survey, investigate its anisotropies induced by the galaxy peculiar motions, and compare our results to those obtained from the 2dFGRS by Hawkins et al. (2003). In section 6 we discuss possible ways to alleviate the discrepancy between model and observations on small scales, and we summarize our results in Section 7.

## 2 The Conditional Luminosity Function

In Paper I we developed a formalism, based on the conditional luminosity function $`\mathrm{\Phi }(L|M)`$, to link the distribution of galaxies to that of dark matter haloes. We introduced a parameterized form for $`\mathrm{\Phi }(L|M)`$ which we constrained using the LF and the correlation lengths as function of luminosity. In Paper II we extended this formalism by constructing separate CLFs for the early- and late-type galaxies. In this paper we use these results to populate dark matter haloes, obtained from large numerical simulations, with both early- and late-type galaxies of different luminosities. For completeness, we briefly summarize here the main ingredients of the CLF formalism, and refer the reader to papers I and II for more details.

The conditional luminosity function is parameterized by a Schechter function:

$$\mathrm{\Phi }(L|M)\mathrm{d}L=\frac{\stackrel{~}{\mathrm{\Phi }}^{}}{\stackrel{~}{L}^{}}\left(\frac{L}{\stackrel{~}{L}^{}}\right)^{\stackrel{~}{\alpha }}\mathrm{exp}(L/\stackrel{~}{L}^{})\mathrm{d}L,$$
(1)
where $`\stackrel{~}{L}^{}=\stackrel{~}{L}^{}(M)`$, $`\stackrel{~}{\alpha }=\stackrel{~}{\alpha }(M)`$ and $`\stackrel{~}{\mathrm{\Phi }}^{}=\stackrel{~}{\mathrm{\Phi }}^{}(M)`$ are all functions of halo mass $`M`$111Halo masses are defined as the masses within the radius $`r_{180}`$ inside of which the average overdensity is $`180`$.. Following Papers I and II, we write the average total mass-to-light ratio of a halo of mass $`M`$ as

$$\frac{M}{L}(M)=\frac{1}{2}\left(\frac{M}{L}\right)_0\left[\left(\frac{M}{M_1}\right)^{\gamma _1}+\left(\frac{M}{M_1}\right)^{\gamma _2}\right],$$
(2)
which has four free parameters: a characteristic mass $`M_1`$, for which the mass-to-light ratio is equal to $`(M/L)_0`$, and two slopes, $`\gamma _1`$ and $`\gamma _2`$, that specify the behavior of $`M/L`$ at the low and high mass ends, respectively. A similar parameterization is adopted for the characteristic luminosity $`\stackrel{~}{L}^{}(M)`$:

$$\frac{M}{\stackrel{~}{L}^{}(M)}=\frac{1}{2}\left(\frac{M}{L}\right)_0f(\stackrel{~}{\alpha })\left[\left(\frac{M}{M_1}\right)^{\gamma _1}+\left(\frac{M}{M_2}\right)^{\gamma _3}\right],$$
(3)
with

$$f(\stackrel{~}{\alpha })=\frac{\mathrm{\Gamma }(\stackrel{~}{\alpha }+2)}{\mathrm{\Gamma }(\stackrel{~}{\alpha }+1,1)}.$$
(4)
Here $`\mathrm{\Gamma }(x)`$ is the Gamma function and $`\mathrm{\Gamma }(a,x)`$ the incomplete Gamma function. This parameterization has two additional free parameters: a characteristic mass $`M_2`$ and a power-law slope $`\gamma _3`$. For $`\stackrel{~}{\alpha }(M)`$ we adopt a simple linear function of $`\mathrm{log}(M)`$,

$$\stackrel{~}{\alpha }(M)=\alpha _{15}+\eta \mathrm{log}(M_{15}),$$
(5)
with $`M_{15}`$ the halo mass in units of $`10^{15}h^1\mathrm{M}_{}`$, $`\alpha _{15}=\stackrel{~}{\alpha }(M_{15}=1)`$, and $`\eta `$ describes the change of the faint-end slope $`\stackrel{~}{\alpha }`$ with halo mass. Note that once $`\stackrel{~}{\alpha }(M)`$ and $`\stackrel{~}{L}^{}(M)`$ are given, the normalization of the conditional LF, $`\stackrel{~}{\mathrm{\Phi }}^{}(M)`$, is obtained through equations (1) and (2), using the fact that the total (average) luminosity in a halo of mass $`M`$ is

$$L(M)=_0^{\mathrm{}}\mathrm{\Phi }(L|M)LdL=\stackrel{~}{\mathrm{\Phi }}^{}\stackrel{~}{L}^{}\mathrm{\Gamma }(\stackrel{~}{\alpha }+2).$$
(6)
Finally, we introduce the mass scale $`M_{\mathrm{min}}`$ below which we set the CLF to zero; i.e., we assume that no stars form inside haloes with $`M<M_{\mathrm{min}}`$. Motivated by reionization considerations (see Paper I for details) we adopt $`M_{\mathrm{min}}=10^9h^1\mathrm{M}_{}`$ throughout.

In order to split the galaxy population in early and late types, we follow Paper II and introduce the function $`f_{\mathrm{late}}(L,M)`$, which specifies the fraction of galaxies with luminosity $`L`$ in haloes of mass $`M`$ that are late-type. The CLFs of late- and early-type galaxies are then given by

$$\mathrm{\Phi }_{\mathrm{late}}(L|M)\mathrm{d}L=f_{\mathrm{late}}(L,M)\mathrm{\Phi }(L|M)\mathrm{d}L$$
(7)
and

$$\mathrm{\Phi }_{\mathrm{early}}(L|M)\mathrm{d}L=\left[1f_{\mathrm{late}}(L,M)\right]\mathrm{\Phi }(L|M)\mathrm{d}L.$$
(8)
As with the CLF for the entire population of galaxies, $`\mathrm{\Phi }_{\mathrm{late}}(L|M)`$ and $`\mathrm{\Phi }_{\mathrm{early}}(L|M)`$ are constrained by 2dFGRS measurements of the LFs and the correlation lengths as function of luminosity. We assume that $`f_{\mathrm{late}}(L,M)`$ has a quasi-separable form

$$f_{\mathrm{late}}(L,M)=g(L)h(M)q(L,M).$$
(9)
Here

$$q(L,M)=\{\begin{array}{ccc}1\hfill & \text{if }g(L)h(M)1\hfill & \\ \frac{1}{g(L)h(M)}\hfill & \text{if }g(L)h(M)>1\hfill & \end{array}$$
(10)
is to ensure that $`f_{\mathrm{late}}(L,M)1`$. We adopt

$$g(L)=\frac{\widehat{\mathrm{\Phi }}_{\mathrm{late}}(L)}{\widehat{\mathrm{\Phi }}(L)}\frac{_0^{\mathrm{}}\mathrm{\Phi }(L|M)n(M)dM}{_0^{\mathrm{}}\mathrm{\Phi }(L|M)h(M)n(M)dM}$$
(11)
where $`n(M)`$ is the halo mass function (Sheth & Tormen 1999; Sheth, Mo & Tormen 2001), $`\widehat{\mathrm{\Phi }}_{\mathrm{late}}(L)`$ and $`\widehat{\mathrm{\Phi }}(L)`$ correspond to the observed LFs of the late-type and entire galaxy samples, respectively, and

$$h(M)=\mathrm{max}(0,\mathrm{min}[1,\left(\frac{\mathrm{log}(M/M_a)}{\mathrm{log}(M_b/M_a)}\right)])$$
(12)
with $`M_a`$ and $`M_b`$ two additional free parameters, defined as the masses at which $`h(M)`$ takes on the values $`0`$ and $`1`$, respectively. As shown in Paper II, this parameterization allows the population of galaxies to be split in early- and late-types such that their respective LFs and clustering properties are well fitted.

In Papers I and II we presented a number of different CLFs for different cosmologies and different assumptions regarding the free parameters. In what follows we focus on the flat $`\mathrm{\Lambda }`$CDM cosmology with $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ and $`h=H_0/(100\mathrm{km}\mathrm{s}^1\mathrm{Mpc}^1)=0.7`$ and with initial density fluctuations described by a scale-invariant power spectrum with normalization $`\sigma _8=0.9`$. These cosmological parameters are in good agreement with a wide range of observations, including the recent WMAP results (Spergel et al. 2003), and in what follows we refer to it as the “concordance” cosmology. Finally, we adopt the CLF with the following parameters: $`M_1=10^{10.94}h^1\mathrm{M}_{}`$, $`M_2=10^{12.04}h^1\mathrm{M}_{}`$, $`M_a=10^{17.26}h^1\mathrm{M}_{}`$, $`M_b=10^{10.86}h^1\mathrm{M}_{}`$, $`(M/L)_0=124h(\mathrm{M}/\mathrm{L})_{}`$, $`\gamma _1=2.02`$, $`\gamma _2=0.30`$, $`\gamma _3=0.72`$, $`\eta =0.22`$ and $`\alpha _{15}=1.10`$. This model (referred to as model D in Paper II) yields excellent fits to the observed LFs and the observed correlation lengths as function of both luminosity and type222Note that the parameters listed here are slightly different from those given in the orignal version of Paper II, as they are based on a corrected version of the galaxy luminosity function. As shown in Paper I, a change in the overall amplitude of the luminosity function in the fitting has some effect on the best-fit values of the correlation lengths. This is due to the combination of the following two effects. First, our model assumes a fixed mass-to-light ratio for massive haloes and so a change in the amplitude of the luminosity function leads to a change in the relative number of galaxies in small/large haloes. Second, although the correlation length as a function of luminosity was used as input in our fitting of the conditional luminosity function, there is some freedom for the ‘best-fit’ values of the correlation lengths to change in the fitting, because the errorbars on the observed correlation lengths are quite large.

## 3 Populating Haloes with Galaxies

### 3.1 Numerical Simulations

The main goal of this paper is to use the CLF described in the previous Section to construct mock galaxy redshift surveys, and to study a number of statistical properties of these distributions that can be compared with observations from existing or forthcoming redshift surveys. The distribution of dark matter haloes is obtained from a set of large $`N`$-body simulations (dark matter only). The set consists of a total of six simulations with $`N=512^3`$ particles each, that have been carried out on the VPP5000 Fujitsu supercomputer of the National Astronomical Observatory of Japan with the vectorized-parallel P3M code (Jing & Suto 2002). Each simulation evolves the distribution of the dark matter from an initial redshift of $`z=72`$ down to $`z=0`$ in a $`\mathrm{\Lambda }`$CDM ‘concordance’ cosmology. All simulations consider boxes with periodic boundary conditions; in two cases $`L_{\mathrm{box}}=100h^1\mathrm{Mpc}`$ while the other four simulations all have $`L_{\mathrm{box}}=300h^1\mathrm{Mpc}`$. Different simulations with the same box size are completely independent realizations and are used to estimate errors due to cosmic variance. The particle masses are $`6.2\times 10^8h^1\mathrm{M}_{}`$ and $`1.7\times 10^{10}h^1\mathrm{M}_{}`$ for the small and large box simulations, respectively. One of the simulations with $`L_{\mathrm{box}}=100h^1\mathrm{Mpc}`$ has previously been used by Jing & Suto (2002) to derive a triaxial model for density profiles of CDM haloes, and we refer the reader to that paper for complementary information about the simulations. In what follows we refer to simulations with $`L_{\mathrm{box}}=100h^1\mathrm{Mpc}`$ and $`L_{\mathrm{box}}=300h^1\mathrm{Mpc}`$ as $`L_{100}`$ and $`L_{300}`$ simulations, respectively.

Dark matter haloes are identified using the standard friends-of-friends (FOF) algorithm (Davis et al. 1985) with a linking length of $`0.2`$ times the mean inter-particle separation. Haloes obtained with this linking length have a mean overdensity of $`180`$ (Porciani, Dekel & Hoffman 2002), in good agreement with the definition of halo masses used in our CLF analysis. For each individual simulation we construct a catalogue of haloes with $`10`$ particles or more, for which we store the mass (number of particles), the position of the most bound particle, and the halo’s mean velocity and velocity dispersion. Note that the FOF algorithm can sometimes select poor systems (those with small number of particles) that are spurious and have abnormally large velocity dispersions. We therefore have made a check to make sure that the particles assigned to a system according to the FOF algorithm are gravitationally bound. Our test showed that this correction is important only for low-mass haloes, and that it has almost no effect on our results. The left panel of Fig. 1 plots the $`z=0`$ halo mass functions for one of the $`L_{100}`$ simulations and for one of the $`L_{300}`$ simulations (histograms), with all spurious haloes erased. For comparison, we also plot (solid line) the analytical halo mass function given in Sheth & Tormen (1999) and Sheth, Mo & Tormen (2001)333This same mass function is used in the CLF analysis described in Section 2.. The agreement is remarkably good, both between the two simulations and between the simulation results and the theoretical prediction.

Note that our choice for box sizes of $`100h^1\mathrm{Mpc}`$ and $`300h^1\mathrm{Mpc}`$ is a compromise between high mass resolution and a sufficiently large volume to study the large-scale structure. The impact of mass resolution is apparent from considering the conditional probability function

$$P(M|L)\mathrm{d}M=\frac{\mathrm{\Phi }(L|M)}{\mathrm{\Phi }(L)}n(M)\mathrm{d}M,$$
(13)
(see Paper I), which gives the probability that a galaxy of luminosity $`L`$ resides in a halo with mass in the range $`M\pm \mathrm{d}M/2`$. The right panel of Fig. 1 plots this probability distribution obtained from the CLF given in Section 2 for four different luminosities: $`L=L^{}/100`$, $`L=L^{}/10`$, $`L=L^{}`$, and $`L=10L^{}`$. Whereas $`10L^{}`$ galaxies are typically found in haloes with $`10^{13}h^1M10^{15}h^1\mathrm{M}_{}`$, galaxies with $`L=L^{}/10010^8h^2\mathrm{L}_{}`$ typically reside in haloes of $`M5\times 10^{10}h^1\mathrm{M}_{}`$. Comparing these probability distributions with the halo mass functions in the left panel, we see that the $`L_{300}`$ simulations can only yield a complete galaxy distribution down to $`L0.4L^{}`$. The $`L_{100}`$ simulation, however, resolves dark matter haloes down to masses of $`10^{10}h^1\mathrm{M}_{}`$, which is sufficient to model the galaxy population down to $`L0.01L^{}`$. On the other hand, luminous galaxies may be under-represented in this small box simulation, because it contains fewer massive haloes than expected. Combining these two sets of simulations, however, will enable us to study the clustering properties of galaxies covering a sufficiently large volume and a sufficiently large range of luminosities.

### 3.2 Halo Occupation Numbers

When populating haloes with galaxies based on the CLF one first needs to choose a minimum luminosity. Based on the mass resolution of the simulations we adopt $`L_{\mathrm{min}}=10^7h^2\mathrm{L}_{}`$ throughout. The mean occupation number of galaxies with $`LL_{\mathrm{min}}`$ for a halo with mass $`M`$ then follows from the CLF according to:

$$N(M)=_{L_{\mathrm{min}}}^{\mathrm{}}\mathrm{\Phi }(L|M)dL.$$
(14)
In order to Monte-Carlo sample occupation numbers for individual haloes one requires the full probability distribution $`P(N|M)`$ (with $`N`$ an integer) of which $`N(M)`$ gives the mean, i.e.,

$$N(M)=\underset{N=0}{\overset{\mathrm{}}{}}NP(N|M)$$
(15)
As a simple model we adopt

$$P(N|M)=\{\begin{array}{ccc}N^{}+1N(M)\hfill & \text{if }N=N^{}\hfill & \\ N(M)N^{}\hfill & \text{if }N=N^{}+1\hfill & \\ 0\hfill & \text{otherwise}\hfill & \end{array}$$
(16)
Here $`N^{}`$ is the largest integer smaller than $`N(M)`$. Thus, the actual number of galaxies in a halo of mass $`M`$ is either $`N^{}`$ or $`N^{}+1`$. This particular model for the distribution of halo occupation numbers is supported by semi-analytical models and hydrodynamical simulations of structure formation (Benson et al. 2000; Berlind et al. 2003) which indicate that the halo occupation probability distribution is narrower than a Poisson distribution with the same mean. In addition, distribution (16) is successful in yielding power-law correlation functions, much more so than for example a Poisson distribution (Benson et al. 2000; Berlind & Weinberg 2002).

### 3.3 Assigning galaxies their luminosity and type

Since the CLF only gives the average number of galaxies with luminosities in the range $`L\pm \mathrm{d}L/2`$ in a halo of mass $`M`$, there are many different ways in which one can assign luminosities to the $`N_i`$ galaxies of halo $`i`$, and yet be consistent with the CLF. The simplest approach would be to simply draw $`N_i`$ luminosities (with $`L>L_{\mathrm{min}}`$) randomly from $`\mathrm{\Phi }(L|M)`$. Alternatively, one could use a more deterministic approach, and, for instance, always demand that the $`j^{\mathrm{th}}`$ brightest galaxy has a luminosity in the range $`[L_j,L_{j1}]`$. Here $`L_j`$ is defined such that a halo has on average $`j`$ galaxies with $`L>L_j`$, i.e.,

$$_{L_j}^{\mathrm{}}\mathrm{\Phi }(L|M)𝑑L=j.$$
(17)
We adopt an intermediate approach in most of our discussion, giving special treatment only to the one brightest galaxy per halo. The luminosity of this so-called “central” galaxy, $`L_c`$, is drawn from $`\mathrm{\Phi }(L|M)`$ with the restriction $`L>L_1`$ and thus has an expectation value of

$$L_c(M)=_{L_1}^{\mathrm{}}\mathrm{\Phi }(L|M)LdL=\stackrel{~}{\mathrm{\Phi }}^{}\stackrel{~}{L}^{}\mathrm{\Gamma }(\stackrel{~}{\alpha }+2,L_1/\stackrel{~}{L}^{}),$$
(18)
The remaining $`N_i1`$ galaxies are referred to as “satellite” galaxies and are assigned luminosities in the range $`L_{\mathrm{min}}<L<L_1`$, again drawn from the distribution function $`\mathrm{\Phi }(L|M)`$. In Section 4.2, we test the effect of luminosity sampling by comparing the results obtained from all the three approaches.

Finally, the galaxies are assigned morphological types as follows. For each galaxy with luminosity $`L`$ in a halo of mass $`M`$ we draw a random number $``$ in the range $`[0,1]`$. If $`<f_{\mathrm{late}}(L,M)`$ then the galaxy is a late-type, otherwise an early-type.

### 3.4 Assigning galaxies their phase-space coordinates

Once the population of galaxies has been assigned luminosities and types, they need to be assigned a position within their halo as well as a peculiar velocity. The central galaxy is assumed to be located at the “center” of the corresponding dark halo, which we associate with the position of the most bound particle, and its peculiar velocity is set equal to the mean halo velocity (cf. Yoshikawa, Jing & Börner 2003). For the satellite galaxies we follow two different approaches. In the first, we assign the $`N_i1`$ satellites the positions and peculiar velocities of $`N_i1`$ randomly selected dark matter particles that are part of the FOF halo under consideration. This thus corresponds to a scenario in which satellite galaxies are completely unbiased with respect to the density and velocity distribution of dark matter particles in FOF haloes. We refer to satellite galaxies populated this way as “FOF satellites”.

We also consider a more analytical model for the satellite distribution. This allows us first of all to assess whether a simple analytical description can be found to describe the population of satellite galaxies, and secondly, provides us with a simple framework to investigate the sensitivity of various clustering statistics to details regarding the density and velocity bias of satellite galaxies. We assume that the number density distribution of satellite galaxies follows a NFW density distribution (Navarro, Frenk & White 1997):

$$\rho (r)=\frac{\overline{\delta }\overline{\rho }}{(r/r_\mathrm{s})(1+r/r_\mathrm{s})^2},$$
(19)
where $`r_s`$ is a characteristic radius, $`\overline{\rho }`$ is the average density of the Universe, and $`\overline{\delta }`$ is a dimensionless amplitude which can be expressed in terms of the halo concentration parameter $`c=r_{180}/r_s`$ as

$$\overline{\delta }=\frac{180}{3}\frac{c^3}{\mathrm{ln}(1+c)c/(1+c)}.$$
(20)
Here $`r_{180}`$ is the radius inside of which the halo has an average overdensity of $`180`$. Numerical simulations show that halo concentration depends on halo mass, and we use the relation given by Bullock et al. (2001), converted to the $`c`$ appropriate for our definition of halo mass. The radial number density distribution of satellite galaxies is assumed to follow equation (19) with a concentration $`c_g=c`$, and the angular position is assumed to be random over the $`4\pi `$ solid angle. Peculiar velocities are assumed to be the sum of the peculiar (mean) velocity of the host halo plus a random velocity which is assumed to be distributed isotropically and to follow a Gaussian, one-dimensional velocity distribution:

$$f(v_j)=\frac{1}{\sqrt{2\pi }\sigma _{\mathrm{gal}}}\mathrm{exp}\left(\frac{v_j^2}{2\sigma _{\mathrm{gal}}^2}\right).$$
(21)
Here $`v_j`$ is the velocity relative to that of the central galaxy along axis $`j`$, and $`\sigma _{\mathrm{gal}}`$ is the one-dimensional velocity dispersion of the galaxies, which we set equal to that of the dark matter particles, $`\sigma _{\mathrm{dm}}`$, in the halo under consideration. We refer to satellite galaxies populated this way as “NFW satellites”.

## 4 Results in Real Space

Fig. 2 and 3 show slices of mock galaxy distributions (hereafter MGDs) constructed from $`L_{100}`$ and $`L_{300}`$ simulations, respectively. Satellite galaxies are assigned positions and velocities using the NFW scheme outlined above. Results are shown for all galaxies (upper right panels), and separately for early types (lower right panels) and late types (lower left panels). For comparison, we also show the distribution of dark matter particles in the upper left panels. Note how the large scale structure in the dark matter distribution is delineated by the distribution of galaxies, and that early-type galaxies are more strongly clustered than late-type galaxies.

In this section we discuss the general, real-space properties of these MGDs. In Section 5 below we construct mock galaxy redshift surveys to investigate the impact of redshift distortions. The main goal of this section, however, is to investigate with what accuracy the combination of numerical simulations and our CLF analysis can be used to construct self-consistent mock galaxy distributions. In particular, we want to examine to what accuracy these MGDs can recover the input used to constrain the CLFs. Note that this is not a trivial question. The CLF modeling is based on the halo model, which only yields an approximate description of the dark matter distribution in the non-linear regime (see discussions in Cooray & Sheth 2002 and Huffenberger & Seljak 2003). In addition, as described in Section 3, the CLF alone does not yield sufficient information to construct MGDs, and we had to make additional assumptions regarding the distribution of galaxies within individual haloes. A further goal of this section, therefore, is to investigate how these assumptions impact on the clustering statistics.

### 4.1 The luminosity function

The CLFs used to construct the MGDs shown in Fig. 2 and 3 are constrained by the 2dFGRS luminosity functions for early- and late-type galaxies obtained by Madgwick et al. (2002). Therefore, as long as the halo mass function is well sampled by the simulations, the LFs of our MGDs should match those of Madgwick et al. (2002). Fig. 4 shows a comparison between the 2dFGRS LFs (symbols with errorbars) and the ones recovered from the MGDs (solid lines). To emphasize the level of agreement between the recovered LFs and the input LFs, Fig. 5 plots the ratio between the two. Over a large range of luminosities, the recovered LFs match the observational input extremely well. In the $`L_{300}`$ simulation, however, the LFs are under-estimated for $`L3\times 10^9h^2\mathrm{L}_{}`$ ($`M_{b_J}5\mathrm{log}h18.4`$). This owes to the absence of haloes with $`M2\times 10^{11}h^1M_{}`$ (see Fig. 1). Note how this discrepancy sets in at higher $`L`$ for the late-type galaxies than for the early-types, because the latter are preferentially located in more massive haloes. For the early-types the $`L_{300}`$ mock is virtually complete down to $`M_{b_J}5\mathrm{log}h17`$ (see Fig. 10 of Paper II), reflecting the fact that only a very small fraction of the early-type galaxies brighter than this magnitude reside in haloes below the mass resolution limit. In the $`L_{100}`$ simulations, on the other hand, the LFs accurately match the data down to the faintest luminosities, but here the MGD underestimates the LFs at the bright end ($`M_{b_J}5\mathrm{log}h22`$). This owes to the limited boxsize, which causes the number of massive haloes (the main hosts of the brightest galaxies) to be underestimated (cf. Fig. 1). Note that even the LFs of the $`L_{300}`$ simulations underestimate the observed number of bright galaxies. This, reflects a small inaccuracy of our CLF to accurately match the observed bright end of the LFs (see paper II).

### 4.2 The real-space correlation function

In addition to the LFs of early- and late-type galaxies, the CLFs used here to construct our MGDs are also constrained by the luminosity and type dependence of the correlation lengths as measured from the 2dFGRS by Norberg et al. (2002a). Here we check to what degree this “input” is recovered from the MGDs.

The left panel of Fig. 6 plots the real-space two-point correlation functions (2PCFs) for dark matter particles in the $`L_{100}`$ (dashed line) and $`L_{300}`$ (dotted line) simulations. The solid line corresponds to the evolved, non-linear dark matter correlation function of Smith et al. (2003) and is shown for comparison444In fitting the CLF we have used this function to compute the correlation length of the dark matter (see Paper II).. As one can see, on large scales ($`r6h^1\mathrm{Mpc}`$) the correlation amplitude obtained from the $`L_{100}`$ simulations is systematically lower than both that obtained from the $`L_{300}`$ simulations and that obtained from the fitting formula of Smith et al. , suggesting that the box-size effect is non-negligible in the $`L_{100}`$ simulations. Note also that the large scale correlation amplitude given by the $`L_{300}`$ simulations is slightly higher than Smith et al. ’s model. It is unclear if this discrepancy is due to the inaccuracy of the fitting formula, or due to cosmic variance in the present simulations. As we will see below, this discrepancy limits the accuracy of model predictions.

The right-hand panel of Fig. 6 plots the 2PCFs for the galaxies in the $`L_{100}`$ (dashed line) and $`L_{300}`$ (dotted line) MGDs. Note how the galaxies reveal the same trend on large scales as the dark matter particles, with larger correlations in the $`L_{300}`$ than in the $`L_{100}`$ MGD.

Fig. 7 shows the correlation lengths $`r_0`$ as function of luminosity for all (upper panel), early-type (middle panel) and late-type (lower panel) galaxies. These have been obtained by fitting $`\xi (r)`$ with a power law relation of the form $`\xi (r)=(r/r_0)^\gamma `$ over the same range of scales as used by Norberg et al. (2002a). Solid squares and open stars correspond to correlation lengths obtained from the $`L_{300}`$ and $`L_{100}`$ MGDs, respectively. Note that the errobars on the predicted correlation lengths are based on the scatter among independent simulations boxes. They are significantly smaller than the errorbars on the observational data, because the model predictions are based on real-space correlation functions, while the observational results are based on projected correlation functions in redshift space. The agreement with the data (open circles) is reasonable, even though several systematic trends are apparent. In particular, the correlation lengths obtained from the $`L_{300}`$ simulation are slightly higher than the observations while the opposite applies to the $`L_{100}`$ simulation. These discrepancies are due to two effects. First of all, as shown in Fig. 6 the dark matter on large scales is more strongly clustered in the $`L_{300}`$ simulations than in the $`L_{100}`$ simulations. That this can account for most of the differences between the scale-lengths obtained from the $`L_{300}`$ and $`L_{100}`$ simulations, is illustrated by the dotted and solid horizontal lines, which indicate the correlation lengths of the dark matter particles in the $`L_{300}`$ and $`L_{100}`$ simulations, respectively. Secondly, the measured correlation lengths correspond to a non-zero, median redshift which is larger for the more luminous galaxies. In determining the best-fit parameters for the CLF this redshift effect is taken into account (see Papers I and II). However, in the construction of our MGDs, we only use the dark matter distribution at $`z=0`$. As discussed in Paper I, this can over-estimate the correlation length by about $`10\%`$. Given these sources of systematic errors, one should be careful not to over-interpret any discrepancy between the correlation lengths in the mock survey and those obtained from real redshift distributions.

In order to investigate the sensitivity of the 2PCF in the MGDs to the way we assign luminosities and phase-space coordinates to the galaxies within the dark matter haloes, we construct MGDs using one of the $`L_{300}`$ simulations with different models for the luminosity assignment and spatial distribution of satellite galaxies within haloes. We have confirmed that using one of the $`L_{100}`$ simulations instead yields the same results. We first test the impact of the luminosity assignment. Here, instead of the fiducial model for the luminosity assignment (the intermediate approach discussed in Section 3.3), we use both the deterministic and random assignments (see Section 3.3 for definitions) to construct the MGDs. In Fig. 8 we shown the ratios between the correlation functions obtained from these MGDs and those obtained from the fiducial MGD. For bright galaxies, the deterministic model gives the lowest amplitudes on small scales ($`r1h^1\mathrm{Mpc}`$), while the random model gives the highest amplitudes. This is expected. The mean number of bright galaxies in a typical halo is not much larger than 1 and so not many close pairs of bright galaxies are expected in the deterministic model. More such pairs are expected in the random model because more than one galaxies in a typical halo can be assigned a large luminosity due to random fluctuations. The dashed lines in Fig. 8 correspond to a MGD with FOF satellites (see Section 3.4). The agreement of the 2PCFs between this MGD with ‘FOF satellites’ and our fiducial MGD indicate that the spherical NFW model is a good approximation of the average density distribution of dark matter haloes. We have also tested the impact of changing the concentration of galaxies, $`c_g`$; increasing (decreasing) $`c_g`$ with respect to the dark matter halo concentration, $`c`$, increases (decreases) the 2PCFs on small scales ($`r1h^1\mathrm{Mpc}`$). However, even when changing the ratio $`c_g/c`$ by a factor of two, the amplitude of this change is smaller than the differences resulting from changing the luminosity assignment.

All in all, changes in the way we assign luminosities and phase-space coordinates to the galaxies only have a mild impact on the 2PCFs, and only at small scales $`1h^1\mathrm{Mpc}`$. This is in good agreement with Berlind & Weinberg (2002) who have shown that these effects are much smaller than changes in the second moment of the halo occupation distributions. For example, assuming a Poissonian $`P(N|M)`$, rather than equation (16) has a much larger impact on the 2PCFs than any of the changes investigated above. As we show in Section 5 below, with the $`P(N|M)`$ of equation (16) we obtain correlation functions that are in better agreement with observations, providing empirical support for this particular occupation number distribution.

It is interesting to note that although small changes in the way we assign luminosities and phase-space coordinate do not have a big impact on the statistical measurements we are considering here, such changes can lead to quite different results for other statistical measures. As shown in van den Bosch et al. (2004), various statistics of satellite galaxies around bright galaxies can be used to distinguish models that make similar predictions about the clustering on large scales.

### 4.3 Pairwise velocities

The peculiar velocities of galaxies are determined by the action of the gravitational field, and so are directly related to the matter distribution in the Universe. Observationally, the properties of galaxy peculiar velocities are inferred from distortions in the correlation function. We defer this discussion to Section 5. Here we derive statistical quantities directly from the simulated peculiar velocities of galaxies.

We define the pairwise peculiar velocity of a galaxy pair as

$$v_{12}(r)[𝐯(𝐱+𝐫)𝐯(𝐱)]\widehat{𝐫},$$
(22)
with $`𝐯(𝐱)`$ the peculiar velocity of a galaxy at $`𝐱`$. The mean pairwise peculiar velocity and the pairwise peculiar velocity dispersion (PVD) are

$$v_{12}(r)\text{and}\sigma _{12}(r)[v_{12}(r)v_{12}(r)]^2^{1/2},$$
(23)
where $`\mathrm{}`$ denotes an average over all pairs of separation $`r`$.

In order to gain insight, we compute $`v_{12}(r)`$ and $`\sigma _{12}(r)`$ from the $`L_{300}`$ simulations for both dark matter particles and for galaxies with $`M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h18.4`$ (which corresponds to the completeness limit of these simulations, see Fig. 4).

Results are shown in Fig. 9. The upper left panel compares the mean pairwise peculiar velocities of the dark matter particles (solid circles) with those of two realizations of the galaxies: one with ‘NFW satellites’ (open circles) and the other with ‘FOF satellites’ (stars). At sufficiently small separations, one probes the virialized regions of dark matter haloes, and one thus finds $`v_{12}=0`$. At larger separations, one starts to probe the infall regions around the virialized haloes, yielding negative values for $`v_{12}(r)`$. Finally, at sufficiently large separations $`v_{12}(r)0`$ due to the large scale homogeneity and isotropy of the Universe.

Both the dark matter particles and the galaxies from our MGDs indeed reveal such a behavior, with $`v_{12}(r)`$ peaking at $`3h^1\mathrm{Mpc}`$. However, there is a markedly strong difference between the $`v_{12}(r)`$ of galaxies in the MGD with NFW satellites and that of the dark matter. In this particular MGD, the galaxies experience significantly smaller infall velocities than the dark matter particles. However, this difference between dark matter and galaxies is almost absent in the MGD with FOF satellites. This is due to the fact that in the NFW model, we populate satellites with isotropic velocity dispersions within a sphere of radius $`r_{180}`$. We are thus assuming that the entire region out to $`r_{180}`$ is virialized in that there is no net infall. However, simple collapse models predict that for our concordance cosmology only the region out to $`r_{340}`$ (i.e., the radius inside of which the average overdensity is 340) is virialized (Bryan & Norman 1998). The difference between the MGDs with NFW satellites and FOF satellites indicates that the regions between $`r_{340}`$ and $`r_{180}`$ are still infalling, resulting in non-zero $`v_{12}`$.

In the lower-left panel, we compare the PVDs for galaxies and dark matter particles. Here the MGDs with FOF satellites and NFW satellites are fairly similar, and significantly lower than for the dark matter. This can be understood as follows. At small separations, the PVD is a pair weighted measure for the potential well in which dark matter particles (galaxies) reside. For the galaxies in our MGDs the halo occupation number per unit mass, $`N/M`$, decreases with the mass of dark matter haloes (see Paper II). Therefore, the massive haloes (with larger velocity dispersions) contribute relatively less to the PVDs of galaxies. Although the difference between the $`\sigma _{12}(r)`$ of the MGDs with FOF and NFW satellites shows that the PVDs have some dependence on the details regarding the infall regions around virialized haloes, these effects are typically small.

The upper-right and lower-right panels of Fig. 9 show how $`v_{12}(r)`$ and $`\sigma _{12}(r)`$ depend on galaxy type. Results are shown for the MGD based on NFW satellites. The mean velocities for early-type galaxies are larger than those for late-type galaxies on large scales, but smaller on small scales. In addition, the PVD of early-type galaxies is higher than that of late-type galaxies on all scales. All these differences are easily understood as a reflection of the fact that early-type galaxies are preferentially located in the larger, more massive haloes which have larger velocity dispersions and larger infall velocities.

Fig. 10 shows the pairwise velocity distributions for four different separations $`r`$, within a logarithmic interval of $`\mathrm{\Delta }\mathrm{log}r=0.125`$. On small scales, the distribution is well fit by an exponential for both dark matter particles and galaxies. This validates the assumption made in earlier analyses about this distribution (Davis & Peebles 1983; Mo, Jing & Börner 1993; Fisher et al. 1994; Marzke et al. 1995). It is also consistent with earlier results obtained from theoretical models and numerical simulations based on dark matter particles (Diaferio & Geller 1996; Sheth 1996; Mo, Jing & Börner 1997; Seto & Yokoyama 1998; Efstathiou et al. 1988; Magira, Jing & Suto 2000). For larger separations $`f(v_{12})`$ is skewed towards negative values of $`v_{12}`$, because galaxies tend to approach each other due to gravitational infall. Clearly, a single exponential function is no longer a good approximation to the pairwise peculiar velocity distribution at large separations. Although for $`v_{12}<0`$ (infall) the exponential remains remarkably accurate, for $`v_{12}>0`$ the pairwise velocity distribution reveals a more Gaussian behavior. This may have important implications for the derivation of PVDs (especially at large separations), which typically is based on the assumption of a purely exponential $`f(v_{12})`$. We shall return to this issue in more detail in Section 5.2.

## 5 Results in Redshift Space

The statistical quantities of galaxy clustering discussed in the previous section are based on real distances between galaxies in our MGDs. However, because of the peculiar velocities of galaxies, such quantities cannot be obtained directly from a galaxy redshift survey. On small scales the virialized motion of galaxies within dark matter haloes cause a reduction of the correlation power, while on larger scales the correlations are boosted due to the infall motion of galaxies towards overdensity regions (Kaiser 1987; Hamilton 1992). As discussed in the introduction, these distortions contain useful information about the Universal density parameter, the bias of galaxies on large (linear) scales, and the pairwise velocities of galaxies.

In this section, we use the MGDs presented above to construct large mock galaxy redshift surveys (hereafter MGRSs). The main goals are to compare various clustering statistics from these mock surveys with observational data from the 2dFGRS, and to investigate how the details about the CLF and the distribution of galaxies within haloes impact on these statistics. For the model-data comparison we use the large scale structure analysis of Hawkins et al. (2003; hereafter H03), which is based on a subsample of the 2dFGRS consisting of all galaxies located in the North Galactic Pole (NGP) and South Galactic Pole (SGP) survey strips with redshift $`0.01z0.20`$ and apparent magnitude $`b_J<19.3`$. This sample consists of $`166,000`$ galaxies covering an area on the sky of $`1090\mathrm{deg}^2`$.

In order to carry out a proper comparison between model and observation, we aim to construct MGRSs that have the same selections as the 2dFGRS. First of all, the survey depth of $`z_{\mathrm{max}}=0.2`$ implies that we need to cover a volume with a depth of $`600h^1\mathrm{Mpc}`$, i.e., twice that of our big $`L_{300}`$ simulations. In principle, we could stack $`4\times 4\times 4`$ identical $`L_{300}`$ boxes (which have periodic boundary conditions), so that a depth of $`600h^1\mathrm{Mpc}`$ can be achieved in all directions for an observer located at the center of the stack. However, there is one problem with this set-up; as we have shown in Figs. 1 and 4 the $`L_{300}`$ MGD is only complete down to $`M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h18.4`$. Taking account of the apparent magnitude limit of the survey, this implies that our MGRSs would be incomplete out to a distance of $`350h^1\mathrm{Mpc}`$. We can overcome this problem by using the higher resolution $`L_{100}`$ simulation, which is complete down to $`M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h14`$. We therefore replace the central $`2\times 2\times 2`$ $`L_{300}`$ boxes with a stack of $`6\times 6\times 6`$ $`L_{100}`$ boxes. The final lay-out of our virtual universe is illustrated in Figure 11. Unless stated otherwise, satellite galaxies are assigned to dark matter haloes based on our standard NFW method described in Section 3.4.

Observational selection effects, which are modelled according to the final public data release of the 2dFGRS (see also Norberg et al. 2002b), are taken into account using the following steps:

1. We place a virtual observer at the center of the stack of boxes (the solid dot in Figure 11), define a $`(\alpha ,\delta )`$-coordinate frame, and remove all galaxies that are not located in the areas equivalent to the NGP and SGP regions of the 2dFGRS.
2. Next, for each galaxy we compute the redshift as ‘seen’ by the virtual observer. We take the observational velocity uncertainties into account by adding a random velocity drawn from a Gaussian distribution with dispersion $`85\mathrm{km}\mathrm{s}^1`$ (Colless et al. 2001).
3. We compute the apparent magnitude of each galaxy according to its luminosity and distance. Since galaxies in the 2dFGRS were pruned by apparent magnitude before a k-correction was applied, we proceed as follows: We first apply a negative k-correction, then select galaxies according to the position-dependent magnitude limit (obtained using the apparant magnitude limit masks provided by the 2dFGRS team), and finally k-correct the magnitudes back to their rest-frame $`b_J`$-band. Throughout we use the type-dependent k-corrections given in Madgwick et al. (2002).
4. To mimic the position- and magnitude-dependent completeness of the 2dFGRS, we randomly sample each galaxy using the completeness masks provided by the 2dFGRS team. The incompleteness of the 2dFGRS parent sample is taking into account by randomly discarding $`9\%`$ of all mock galaxies (Norberg et al. 2002b).
5. Finally, we mimic the actual selection criteria of the 2dFGRS sample used in H03 by restricting the sample to galaxies within the redshift range $`0.01z0.20`$ and with completeness $`0.7`$.

Each MGRS thus constructed contains, on average, $`169000`$ galaxies, with a dispersion of $`5000`$ due to cosmic variance. The number of galaxies in our mock catalogues are consistent with the observations at the $`1\sigma `$ level. Note that the correlation functions presented by H03 have been corrected for the observational bias due to fiber collisions, and we therefore do not mimic these effects in our MGRSs.

Since we have 2 $`L_{100}`$ simulations and 4 $`L_{300}`$ simulations, we construct $`2\times 4=8`$ mock catalogues with different combinations of small- and big-box simulations. In what follows, we refer to this set of mock catalogues as MSBs (for Mock Small/Big). As an example, Fig. 12 shows the distribution of a sub-set of galaxies in one of these mock catalogues. Although each of our MSB catalogs covers an extremely large volume, and should thus not be very sensitive to cosmic variance, it is constructed using simulations with box sizes of $`100`$ and $`300h^1\mathrm{Mpc}`$ only. If, for instance, the $`L_{100}`$ simulation contains a big cluster, the $`6\times 6\times 6`$ reproduction of this box in our MGRSs might introduce some unrealistic features. Furthermore, as shown in Section 4.2 the $`L_{100}`$ box underestimates the amount of clustering on large scales. Therefore, this set of MGRSs, which replicate this box 27 times, might underestimate the clustering on large scales as well. In order to test the sensitivity of our results to these potential problems, and to have a better handle on the impact of cosmic variance in our mock surveys, we construct four alternative MGRSs. Each consists of a $`4\times 4\times 4`$ stack of one of the four $`L_{300}`$ simulations (i.e., we replace the $`6\times 6\times 6`$ stack of $`L_{100}`$ boxes by a $`2\times 2\times 2`$ stack of $`L_{300}`$ boxes). In what follows we refer to this set of mock catalogues as MBs (for Mock Big). These MGRSs, although incomplete for $`M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h>18.4`$, should not suffer from the lack of clustering power on large scales. The MSB set, on the other hand, does not suffer from incompleteness, but instead lacks some large scale power. As we will see below, both the MSB and MB mocks give similar results on large scales, suggesting that the box-size effect does not have a significant influence on our results.

### 5.1 Two-Point Correlation Functions

From our MGRSs we compute $`\xi (r_p,\pi )`$ using the estimator (Hamilton 1993)

$$\xi (r_p,\pi )=\frac{RRDD}{DR^2}1$$
(24)
with $`DD`$, $`RR`$, and $`DR`$ the number of galaxy-galaxy, random-random, and galaxy-random pairs with separation $`(r_p,\pi )`$. Here $`r_p`$ and $`\pi `$ are the pair separations perpendicular and parallel to the line-of-sight, respectively. Explicitly, for a pair $`(𝐬_\mathrm{𝟏},𝐬_\mathrm{𝟐})`$, with $`𝐬_𝐢=cz_i\widehat{𝐫}_𝐢/H_0`$, we define

$$\pi =\frac{𝐬𝐥}{|𝐥|},r_p=\sqrt{𝐬𝐬\pi ^2}$$
(25)
Here $`𝐥=\frac{1}{2}(𝐬_\mathrm{𝟏}+𝐬_\mathrm{𝟐})`$ is the line of sight intersecting the pair, and $`𝐬=𝐬_\mathrm{𝟏}𝐬_\mathrm{𝟐}`$. Random samples are constructed using two different methods. The first uses the mean galaxy number density at redshift $`z`$ calculated from the 2dFGRS LF. The second randomizes the coordinates of all mock galaxies within the simulation box. Both methods yield indistinguishable estimates of $`\xi (r_p,\pi )`$ and in what follows we only use the former. Following H03 each galaxy in a pair with redshift separation $`s`$ is weighted by the factor

$$w_i=\frac{1}{1+4\pi n(z_i)J_3(s)}$$
(26)
with $`n(z)`$ the number density distribution as function of redshift and $`J_3(s)=_0^s\xi (s^{})s^2ds^{}`$. Hence each galaxy-galaxy, random-random, and galaxy-random pair is given a weight $`w_iw_j`$. We substitute $`\xi (s^{})`$ with a power law using the same parameters as in H03. This redshift dependent weighting scheme is designed to minimize the variance on the estimated correlation function (Davis & Huchra 1982; Hamilton 1993).

Since the redshift-space distortions only affect $`\pi `$, the projection of $`\xi (r_p,\pi )`$ along the $`\pi `$ axis can get rid of these distortions and give a function that is more closely related to the real-space correlation function. In fact, this projected 2PCF is related to the real-space 2PCF through a simple Abel transform

$$w_p(r_p)=_{\mathrm{}}^{\mathrm{}}\xi (r_p,\pi )d\pi =2_{r_p}^{\mathrm{}}\xi (r)\frac{r\mathrm{d}r}{\sqrt{r^2r_p^2}}$$
(27)
(Davis & Peebles 1983). Therefore, if the real-space 2PCF is a power-law, $`\xi (r)=(r_0/r)^\gamma `$, the projected 2PCF $`w(r_p)`$ can be written as

$$w_p(r_p)=\sqrt{\pi }\frac{\mathrm{\Gamma }(\gamma /21/2)}{\mathrm{\Gamma }(\gamma /2)}\left(\frac{r_0}{r_p}\right)^\gamma r_p.$$
(28)
We start our investigation of the redshift-space clustering properties by computing $`w_p(r_p)`$ for a number of luminosity bins and for early- and late-type galaxies separately. To compare these projected correlation functions with the 2dFGRS results from Norberg et al. (2002a), we estimate $`w_p(r_p)`$ using volume-limited samples with the same redshift and magnitude selection criteria as those adopted by Norberg et al. (2002a). For the MSB mocks (which use a stack of $`6\times 6\times 6`$ $`L_{100}`$ boxes), however, these $`w_p(r_p)`$ reveal a systematic ‘break’ at $`r_p10h^1\mathrm{Mpc}`$. As we have shown in Section 4.2, this owes to the fact that, because of the small box-size of the $`L_{100}`$ simulation, the 2PCF is too small on large scales (see Figs 6 and 7). We can circumvent this problem by using MGRSs from the the MB set, in which the stack of $`6\times 6\times 6`$ $`L_{100}`$ boxes is replaced by a stack of $`2\times 2\times 2`$ $`L_{300}`$ boxes. However, these MGRSs are only complete down to $`M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h18.4`$ and can therefore only be used for galaxies brighter than this.

The upper panels of Fig. 13 plot $`w_p(r_p)`$ for different magnitude bins and for early- and late-type galaxies separately. Except for the faintest magnitude bin, these projected correlation functions are obtained from MGRS in the MB set. Results for the magnitude bin with $`17.5>M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h>18.5`$ (solid lines) are obtained from the MSB set. As discussed in Paper II, the projection significantly washes out the features in the real-space 2PCFs at $`2h^1\mathrm{Mpc}`$, and the projected 2PCFs better resemble a power-law. The exception is the $`w_p(r_p)`$ for the faintest subsample of galaxies, where the ‘break’ mentioned above is clearly visible. To highlight the luminosity and type dependence of $`w_p(r_p)`$, the lower panels of Fig. 13 plot the ratios of $`w_p(r_p)`$ to that of a reference sample defined as all (early-type plus late-type) galaxies with $`19.5>M_{b_J}5\mathrm{l}\mathrm{o}\mathrm{g}h>20.5`$. For a given luminosity, the correlation amplitude is higher, and the slope is steeper, for early-type galaxies than for late-type galaxies. Significant changes in the slope (and thus deviations from a perfect power-law) occur at separations $`r_p2h^1\mathrm{Mpc}`$, which is at least qualitatively in agreement with recent results from the SDSS (Zehavi et al. 2003).

In order to facilitate a more direct comparison with the 2dFGRS data, we fit a single power-law relation of the form (28) to these $`w_p(r_p)`$ over the range $`2h^1\mathrm{Mpc}<r_p<15h^1\mathrm{Mpc}`$. This range is also adopted by Norberg et al. (2002a) when fitting the projected 2PCFs obtained from the 2dFGRS. Fig. 14 plots the real-space correlation lengths $`r_0`$ and the slopes $`\gamma `$ thus obtained as function of luminosity and galaxy type. The agreement between our MGRSs and the 2dFGRS is acceptable. The slight but systematic overestimate of $`r_0`$ is due to the effects discussed in Section 4.2.

We now turn to a comparison of the projected correlation function for the entire, flux limited surveys. The upper-left panel of Fig. 15 compares the $`w_p(r_p)`$ obtained from our 8 MSB and 4 MB MGRSs with that of the 2dFGRS obtained by H03. The projected correlation functions from our MSBs and MBs agree well with each other (i.e., the 1-$`\sigma `$ errorbars overlap), and, at $`r_p3h^1\mathrm{Mpc}`$, with the 2dFGRS results. Note that at $`r_p10h^1\mathrm{Mpc}`$ the $`w_p(r_p)`$ obtained from the MB mocks is slightly larger than that obtained from the MSB mocks, again due to the effects discussed in Section 4.2.

At large scales, $`w_p(r_p)`$ is predominantly sensitive to the halo occupation numbers $`N(M)`$ and virtually independent of the second moment of $`P(N|M)`$ or of details regarding the spatial distribution of satellite galaxies. The good agreement at large scales among different MGRSs and with the observations, therefore strongly supports our CLF and it shows that any ‘cosmic variance’ among the different MGRSs has only a relatively small impact on $`w_p(r_p)`$. On small scales, however, the MGRSs reveal more correlation power (by about a factor 2) than observed. On such scales, $`w_p(r_p)`$ is sensitive to our assumptions about the second moment of $`P(N|M)`$ and, to a lesser degree, the spatial distribution of satellite galaxies. We shall return to this small-scale mismatch and its implications in Section 6 below.

Rather than projecting $`\xi (r_p,\pi )`$, one may also average $`\xi (r_p,\pi )`$ along constant $`s=\sqrt{r_p^2+\pi ^2}`$, yielding the redshift-space 2PCFs $`\xi (s)`$. The upper-right panel of Fig. 15 plots $`\xi (s)`$ obtained from our MGRSs, compared to the 2dFGRS results from H03. We find a similar behavior as with the projected correlation function; the 8 MSBs and 4 MBs agree quite well with each other and with the observations at $`s6h^1\mathrm{Mpc}`$. At smaller redshift-space separations, however, the MGRSs slightly overpredict the correlation power. Note that the MB samples predict higher $`\xi (s)`$ on small scales than the MSB samples. This difference comes from the fact that the MB samples are incomplete for galaxies fainter than $`M_{b_J}5\mathrm{log}h=18.4`$. To test this we construct a mock survey from the MSB sample, but only accepting galaxies brighter than this. This yields a $`\xi (s)`$ in excellent agreement with that of the MB samples over all scales. Thus, although the use of only large-box simulations can results in systematic errors on small scale, the use of small-box simulations in the MSB samples does not cause any signicifant, systematic errors on large scale.

### 5.2 Redshift Space Distortions

We now turn to a comparison of the detailed shape of $`\xi (r_p,\pi )`$. In particular, we focus on the distortions with respect to the real-space correlation function $`\xi (r)`$ induced by the peculiar velocities of galaxies.

The two-dimensional correlation function $`\xi (r_p,\pi )`$ is often modeled as a convolution of the real-space 2PCF $`\xi (r)`$ and the conditional distribution function $`f(v_{12}|r)`$:

$`1+\xi (r_p,\pi )=`$ (29)
$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}\left[1+\xi (\sqrt{r_p^2+(\pi v_{12}/H_0)^2})\right]f(v_{12}|r)dv_{12}`$
(Peebles 1980). Here $`v_{12}`$ corresponds to the pairwise peculiar velocity along the line-of-sight and $`r`$ corresponds to the real-space separation. It is standard practice to assume an exponential form for $`f(v_{12}|r)`$ and to ignore its dependence on separation $`r`$ (cf., Davis & Peebles 1983; Mo, Jing & Börner 1993, 1997; Fisher et al. 1994; Marzke et al. 1995; Guzzo et al. 1997; Jing, Börner & Suto 2002; Zehavi et al. 2002). However, as we have shown in Section 4, the exponential form is only adequate at small separations, and the PVD varies quite strongly with separation. Furthermore, equation (29) is only valid for an isotropic velocity field in the limit where the probability of a real-space pair separation $`r`$ is independent of the probability of an associated relative velocity $`v_{12}`$. Although perhaps a reasonable approximation on small, highly non-linear, scales, it is certainly not valid in linear theory where the velocity and density fields are tightly coupled. In an attempt to partially correct for this, one often assumes that $`f(v_{12})`$ is the probability distribution for the relative velocity about the mean. Using the self-similar infall model, this mean pairwise peculiar velocity, $`v_{12}`$, is modeled as

$$v_{12}(r)=H_0F\left(\frac{y}{1+(r/r_0)^2}\right)$$
(30)
(Davis & Peebles 1977) with $`y=|\pi v_{12}/H_0|`$ the separation in real-space along the line-of-sight. $`F=0`$ corresponds to a Universe without any flow other than the Hubble expansion, while $`F=1`$ corresponds to stable clustering. Given the fairly ad hoc nature of this model, and the strong sensitivity to the uncertain value of $`F`$ (Davis & Peebles 1983), great care is required when interpreting any results based on this model.

A more robust model is based on linear theory and directly modeling the infall velocities around density perturbations. Following Kaiser (1987) and Hamilton (1992) one can write the observed correlation function on linear scales as

$$\xi _{\mathrm{lin}}(r_p,\pi )=\xi _0(s)𝒫_0(\mu )+\xi _2(s)𝒫_2(\mu )+\xi _4(s)𝒫_4(\mu ).$$
(31)
Here $`𝒫_l(\mu )`$ is the $`l^{th}`$ Legendre polynomial, and $`\mu `$ is the cosine of the angle between the line-of-sight and the redshift-space separation $`𝐬`$. According to linear perturbation theory the angular moments can be written as

$$\xi _0(s)=\left(1+\frac{2\beta }{3}+\frac{\beta ^2}{5}\right)\xi (r),$$
(32)
$$\xi _2(s)=\left(\frac{4\beta }{3}+\frac{4\beta ^2}{7}\right)\left[\xi (r)\overline{\xi }(r)\right],$$
(33)
$$\xi _4(s)=\frac{8\beta ^2}{35}\left[\xi (r)+\frac{5}{2}\overline{\xi }(r)\frac{7}{2}\widehat{\xi }(r)\right],$$
(34)
with

$$\overline{\xi }(r)=\frac{3}{r^3}_0^r\xi (r^{})r^2dr^{},$$
(35)
and

$$\widehat{\xi }(r)=\frac{5}{r^5}_0^r\xi (r^{})r^4dr^{}$$
(36)
Given a value for $`\beta `$ and the real-space correlation function, which can be obtained from $`\xi (r_p,\pi )`$ via the projected correlation function $`w_p(r_p)`$, equation (31) yields a model for $`\xi (r_p,\pi )`$ on linear scales that takes proper account of the coupling between the density and velocity fields. To model the non-linear virialized motions of galaxies within dark matter haloes, one convolves this $`\xi _{\mathrm{lin}}(r_p,\pi )`$ with the distribution function of pairwise peculiar velocities $`f(v_{12}|r)`$.

$`1+\xi (r_p,\pi )=`$ (37)
$`{\displaystyle _{\mathrm{}}^{\mathrm{}}}\left[1+\xi _{\mathrm{lin}}(r_p,\pi v_{12}/H_0)\right]f(v_{12}|r)dv_{12}`$
Thus, by modeling $`\xi (r_p,\pi )`$ one can hope to get both an estimate of $`\beta `$ as well as information regarding the pairwise peculiar velocity distribution. We follow H03, and assume that the real-space 2PCF is a pure power-law, $`\xi (r)=(r/r_0)^\gamma `$, and that $`f(v_{12}|r)`$ is an exponential that is independent of the real-space separation $`r`$:

$$f(v_{12}|r)=f(v_{12})=\frac{1}{\sqrt{2}\sigma _{12}}\mathrm{exp}\left(\frac{\sqrt{2}|v_{12}|}{\sigma _{12}}\right)$$
(38)
Using a simple $`\chi ^2`$ minimization technique, we fit these models, described by the four parameters $`\beta `$, $`\sigma _{12}`$, $`r_0`$, and $`\gamma `$, to the $`\xi (r_p,\pi )`$ in each of our 8 MSB and 4 MB MGRSs. The $`\chi ^2`$ is defined as

$$\chi ^2=\left(\frac{\mathrm{log}[1+\xi ]_{\mathrm{model}}\mathrm{log}[1+\xi ]_{\mathrm{data}}}{\mathrm{log}[1+\xi +\mathrm{\Delta }\xi ]_{\mathrm{data}}\mathrm{log}[1+\xi \mathrm{\Delta }\xi ]_{\mathrm{data}}}\right)^2,$$
(39)
where the summation is over the $`\xi (r_p,\pi )`$ data grid with the restriction $`8h^1\mathrm{Mpc}<s<20h^1\mathrm{Mpc}`$ (see H03) and $`\mathrm{\Delta }\xi (r_p,\pi )`$ is the rms of $`\xi (r_p,\pi )`$ determined from each of our 8 MSB MGRSs (or of our 4 MB MGRSs). The averages (over the 8 or the 4 MGRSs) of the best-fit values for $`\beta `$, $`\sigma _{12}`$, $`r_0`$, and $`\gamma `$, along with the variances among different samples, are listed in the first two lines of Table 1. These should be compared with the values listed in the last line, which correspond to the best-fit values obtained from the 2dFGRS by H03 using exactly the same method. As one can see, the best-fit values for $`\beta `$ and correlation lengths of the MSB and the MB sample agree with each other and with the 2dFGRS value at better than the $`1\sigma `$ level. On the other hand, the discrepancies regarding $`\gamma `$ and $`\sigma _{12}`$ are significant, both of which are significantly higher in our MGRSs than in the 2dFGRS.

In order to investigate these discrepancies in more detail we compute two statistics of the redshift space distortions which we compare to the 2dFGRS. As above, we take great care in using exactly the same method and assumptions as H03. Therefore, even if some aspects of the model are questionable, this allows a meaningful comparison of our results with those obtained by H03. First of all, we compute the modified quadrupole-to-monopole ratio

$$q(s)\frac{\xi _2(s)}{\frac{3}{s^3}_0^s\xi _0(s^{})s^2ds^{}\xi _0(s)}.$$
(40)
where $`\xi _l(s)`$ is given by

$$\xi _l(s)=\frac{2l+1}{2}_1^1\xi (r_p,\pi )𝒫_l(\mu )d\mu .$$
(41)
The lower-left panel of Fig. 15 plots $`q(s)`$ for MSBs and MBs together with the 2dFGRS results (open circles with errorbars). Although the MGRSs reveal the same overall behavior as the 2dFGRS data, and are mutually consistent, they systematically overpredict $`q(s)`$. On small scales, where random peculiar velocities cause a rapid increase of $`q(s)`$, this indicates that the virialized motions in our mock surveys are larger than observed (see also below). On large scales, where $`q(s)`$ asymptotes to the linear theory value of

$$q(s)=\frac{\frac{4}{3}\beta \frac{4}{7}\beta ^2}{1+\frac{2}{3}\beta +\frac{1}{5}\beta ^2},$$
(42)
this might indicate that the value of $`\beta `$ inherent to our MGRSs is too small compared to the real Universe. On the other hand, Cole, Fisher & Weinberg (1994) have shown that non-linear, small scale power can affect $`q(s)`$ out to fairly large separations. Therefore the systematic overestimate of $`q`$ at large $`s`$ may simply be a reflection of the random peculiar velocities being too large, rather than an inconsistency regarding the value of $`\beta `$.

The second statistic that we use to compare the redshift space distortions in our MGRSs with those of the 2dFGRS are the PVDs, $`\sigma _{12}(r_p)`$, as a function of projected radius, $`r_p`$. Following H03, we keep $`r_0`$, $`\gamma `$ and $`\beta `$ fixed at the ‘global’ values listed in Table 1 and determine $`\sigma _{12}(r_p)`$ by minimizing $`\chi ^2`$ in a number of independent $`r_p`$ bins555Note that the PVDs thus obtained are a kind of average of the true PVD along the line-of-sight. Therefore, these PVDs should not be compared directly to the true PVD shown in Fig. 9.. The results are shown in the lower-right panel of Fig. 15. Whereas the 2dFGRS reveals a $`\sigma _{12}(r_p)`$ that is almost constant with radius at about $`500`$$`600\mathrm{km}\mathrm{s}^1`$, our MGRSs reveal a strong increase from $`\sigma _{12}600\mathrm{km}\mathrm{s}^1`$ at $`r_p=0.1h^1\mathrm{Mpc}`$ to $`\sigma _{12}900\mathrm{km}\mathrm{s}^1`$ at $`r_p=1.0h^1\mathrm{Mpc}`$, followed by a decrease to $`\sigma _{12}500\mathrm{km}\mathrm{s}^1`$ at $`r_p=10h^1\mathrm{Mpc}`$. Thus, at around $`1h^1\mathrm{Mpc}`$, our MGRSs dramatically overestimate the PVD. Although there is a non-negligible amount of scatter among the different mock surveys, reflecting the extreme sensitivity of the PVDs to the few richest systems in the survey, the variance among the 8 (4) MGRSs is small compared to the discrepancy.

As shown by Peacock et al. (2001) the best-fit values of $`\sigma _{12}`$ and $`\beta `$ are highly degenerate. We have tested the impact of this degeneracy on our $`\sigma _{12}(r_p)`$ by repeating the same exercise using a value for $`\beta `$ that is $`0.1`$ larger (smaller) than the values listed in Table 1. This leads to an increase (decrease) of $`\sigma _{12}(r_p)`$ of the order of 5 percent (20 percent) at projected radii of $`1h^1\mathrm{Mpc}`$ ($`10h^1\mathrm{Mpc}`$). Given that our MGRSs overpredict the PVD at $`r_p=1h^1\mathrm{Mpc}`$ by about 70 percent, it is clear that this discrepancy is not a reflection of the $`\beta `$$`\sigma _{12}`$ degeneracy. Thus, the standard $`\mathrm{\Lambda }`$CDM model seems to have a severe problem in matching the observed PVDs on intermediate scales.

## 6 Towards a Self-Consistent Model for Large Scale Structure

Our MGRSs, based on a flat $`\mathrm{\Lambda }`$CDM concordance cosmology with $`\mathrm{\Omega }_m=0.3`$ and $`\sigma _8=0.9`$, and on a CLF that is required to yield cluster mass-to-light ratios of $`(M/L)_{\mathrm{cl}}=500h(\mathrm{M}/\mathrm{L})_{}`$, reveals clustering statistics that are overall in reasonable agreement with the data from the 2dFGRS. Nevertheless, two discrepancies have come to light: the MGRSs predict too much power on small scales and PVDs that are too high. We now investigate possible ways to alleviate these discrepancies.

### 6.1 Halo occupation models

The discrepancies between our MGRS and the 2dFGRS results might indicate a problem with our halo occupation models. Although the CLF is fairly well constrained by the observed luminosity function and the observed luminosity dependence of the correlation lengths (see Papers I and II), we have made additional assumptions regarding the second moments of the halo occupation number distributions and regarding the distribution of galaxies within individual dark matter haloes.

As we have shown in Section 4.2 the real space correlation function depends only very weakly on our method of distributing satellite galaxies within dark matter haloes (cf. Fig. 8). We have verified, using a number of tests, that modifications of the spatial distribution of satellite galaxies within dark matter haloes have no significant influence on $`w_p(r_p)`$ or on $`\sigma _{12}(r_p)`$. Therefore, none of the discrepancies mentioned above can be attributed to errors in our satellite model.

Our results are more susceptible to changes in the second moment of the halo occupation number distributions. On small scales, $`\xi (r)`$ scales with the average number of galaxy pairs in individual haloes $`\frac{1}{2}N(N1)`$. Therefore, one can decrease the power on small scales, to bring our $`w_p(r_p)`$ in better agreement with observations, by lowering the second moments of our halo occupation distributions. However, our distribution (16) is already the narrowest distribution possible, and modifying the second moment of the halo occupation distributions can therefore only aggravate the discrepancies on small scales.

### 6.2 Velocity Bias

A seemingly obvious explanation for the too high PVDs is that the peculiar velocities of galaxies are biased with respect to the dark matter. We define the velocity bias (sometimes called ‘dynamical’ bias) as $`b_{\mathrm{vel}}=\sigma _{\mathrm{gal}}/\sigma _{\mathrm{DM}}`$, with $`\sigma _{\mathrm{gal}}`$ and $`\sigma _{\mathrm{DM}}`$ the peculiar velocity dispersions of (satellite) galaxies and dark matter particles in a given halo, respectively. Note that in our fiducial MGRSs we adopt $`b_{\mathrm{vel}}=1`$ (i.e. no velocity bias). Fig. 16 shows $`w_p(r_p)`$, $`\xi (s)`$, $`q(s)`$ and $`\sigma _{12}(r_p)`$ for MGRSs (with the MSB configuration of simulation boxes) in which $`b_{\mathrm{vel}}=0.6`$; i.e., the velocity dispersion of satellite galaxies is only 60 percent of that of the dark matter particles in the same halo (dashed lines). With such a pronounced velocity bias, both $`q(s)`$ and the PVDs, as well as $`\beta `$ and the global value of $`\sigma _{12}`$ listed in Table 1, are all consistent with the 2dFGRS results.

In two recent papers, Berlind et al. (2003) and Yoshikawa et al. (2003) measured the velocity bias of ‘galaxies’ in a smoothed particle hydrodynamics (SPH) simulation and found that $`b_{\mathrm{vel}}`$ decreases from $`b_{\mathrm{vel}}0.9`$$`1.0`$ for haloes with $`M3\times 10^{14}h^1\mathrm{M}_{}`$ to $`b_{\mathrm{vel}}0.6`$$`0.8`$ for haloes with $`M3\times 10^{12}h^1\mathrm{M}_{}`$. Thus, although simulations predict that low-mass haloes might have values for the velocity bias as low as $`b_{\mathrm{vel}}0.6`$, required here to bring our PVDs in agreement with observations, the PVD is dominated by galaxies in massive haloes for which these same simulations apparently predict close to $`b_{\mathrm{vel}}=1`$ (i.e., no velocity bias). Furthermore, the introduction of velocity bias cannot solve the excess power on small scales. After all, the real-space correlation function is independent of $`b_{\mathrm{vel}}`$ such that the discrepancies regarding $`w_p(r_p)`$ on small scales remain (see upper-left panel of Fig. 16). To make matters worse, reducing the peculiar velocities of galaxies inside dark matter haloes, increases the small scale power in redshift space. This means, that the redshift space correlation function $`\xi (s)`$ becomes actually more discrepant with the 2dFGRS data (see upper-right panel of Fig. 16). Therefore, although some amount of velocity bias might be expected, we do not consider it a viable solution for the problems mentioned above.

### 6.3 Cluster mass-to-light ratios

Because the PVD is a pair weighted statistic it is extremely sensitive to the few richest systems in the sample (i.e., Mo, Jing & Börner 1993, 1997; Zurek et al. 1994; Marzke et al. 1995; Somerville, Primack & Nolthenius 1997). The fact that the PVDs in our MGRSs are too large compared with observations therefore might indicate that either there are too many clusters of galaxies in our mock surveys (see Section 6.4 below), or that these clusters contain too many galaxies.

Our CLF was constructed under the constraint that the average mass-to-light ratio of haloes with $`M10^{14}h^1\mathrm{M}_{}`$ is equal to $`(M/L)_{\mathrm{cl}}=500h(\mathrm{M}/\mathrm{L})_{}`$ (in the photometric $`b_J`$ band). This value is motivated by the average mass-to-light ratio of clusters obtained by Fukugita, Hogan & Peebles (1998). To reduce the number of galaxies per cluster we now set $`(M/L)_{cl}=1000h(\mathrm{M}/\mathrm{L})_{}`$ and repeat the entire exercise: we first use the method described in Section 2 to compute the parameters of the new conditional luminosity function. This CLF is used to construct new MGRSs (using the same configuration of simulation boxes as in MSB), from which we determine the same statistics as before.

The results are listed in Table 1 and shown as dotted lines in Fig. 16. Clearly, increasing the mass-to-light ratio of clusters lowers $`q(s)`$ and $`\sigma _{12}(r_p)`$, bringing them in better agreement with the 2dFGRS results. Although the PVDs are still somewhat too high, especially at around $`0.4h^1\mathrm{Mpc}`$, the extent of this discrepancy is similar to its 1-$`\sigma `$ variance of the 8 MGRSs, indicating that this remaining difference is consistent with ‘cosmic variance’. As can be seen from Table 1, both $`\beta `$ and $`\gamma `$ are now in much better agreement with the 2dFGRS. In addition, the reduction of the number of galaxies in clusters significantly reduces $`w_p(r_p)`$ at small projected separations (see Fig. 16), bringing it in good agreement with the observations 666Note that the correlation amplitude predicted with $`(M/L)_{cl}=1000h(\mathrm{M}/\mathrm{L})_{}`$ is slightly lower than the observed amplitude, because in this model more galaxies are assigned to small haloes in order to match the observed luminosity function. Since the errorbars on the observed correlation lengths are relatively large, the model tends to compromise the accuracy of the fit to the correlation lengths.. A similar reduction of small scale power is also evident in $`\xi (s)`$. Thus, these particular MGRSs have clustering characteristics that are overall in good agreement with the 2dFGRS results. The question is therefore whether or not such a high mass-to-light ratio for clusters of galaxies is compatible with observations.

The cluster mass-to-light ratios quoted by Fukugita et al. (1998) are $`(450\pm 100)h(\mathrm{M}/\mathrm{L})_{}`$ in the $`B`$-band based on X-ray and velocity-dispersion data. Taking these numbers at face value, a cluster mass-to-light ratio of $`(M/L)_{\mathrm{cl}}=1000h(\mathrm{M}/\mathrm{L})_{}`$ is ruled out at the $`5\sigma `$ level. Using a variety of methods to estimate cluster masses, Bahcall et al. (2000) obtained $`(M/L)_B=(330\pm 77)h(\mathrm{M}/\mathrm{L})_{}`$, which is consistent with the results of Carlberg et al. (1996), $`(M/L)_B=(363\pm 65)h(\mathrm{M}/\mathrm{L})_{}`$ based on galaxy kinematics in clusters. Taking the average of these two measurements yields $`M/L_B_{\mathrm{cl}}=(350\pm 70)h(\mathrm{M}/\mathrm{L})_{}`$, which rules out the cluster mass-to-light ratio required to match the clustering power on small scales at more than $`7\sigma `$. Thus, unless the cluster mass-to-light ratios obtained from current observations are seriously in error, increasing the average cluster mass-to-light ratio to $`(M/L)_{\mathrm{cl}}1000h(\mathrm{M}/\mathrm{L})_{}`$ does not seem a viable solution for the problems at hand.

### 6.4 Power spectrum normalization

Rather than lowering the average number of galaxies per cluster, we may also hope to lower the PVDs by reducing the actual number of clusters. As we have shown in Section 5.2, our results, and thus the number density of rich clusters, is robust against cosmic variance. Therefore, a lower number density of clusters implies a different cosmological model. It is well known that the abundance of (rich) clusters is extremely sensitive to the power-spectrum normalization parameter $`\sigma _8`$. The too high PVDs could thus be indicative of a too high value for $`\sigma _8`$.

We therefore wish to compute the PVDs in a $`\mathrm{\Lambda }`$CDM cosmology with identical cosmological parameters as before, except that $`\sigma _8=0.75`$ rather than $`0.9`$. Note that the choice of $`\sigma _8=0.75`$ is somewhat arbitrary, but it does represent a compromise between the constraints on the value of $`\sigma _8`$ from various observations and the low value required by our results on the PVDs (see below). In principle, constructing new MGRSs for a different cosmology requires new $`N`$-body simulations of the dark matter distribution. This, however, is computationally too expensive, which is why we use an approximate method instead. First, we compute the new best-fit parameters of the CLF for this cosmology, again demanding that $`(M/L)_{\mathrm{cl}}=500h(\mathrm{M}/\mathrm{L})_{}`$. Next we populate the dark matter haloes in our $`\sigma _8=0.9`$ simulation boxes with galaxies according to this new CLF. Finally, we construct a new sample of 8 MSB MGRSs, in which we weigh each galaxy in a halo of mass $`M`$ by

$$w=\frac{n(M|\sigma _8=0.75)}{n(M|\sigma _8=0.9)}$$
(43)
with $`n(M)`$ the number density of dark matter haloes of mass $`M`$. This, to first order, mimics the effect of lowering $`\sigma _8`$ on the halo mass function, and so should be a reasonable approximation on small scales where the clustering properties are determined by the galaxy distribution in individual haloes777Note that for haloes with a given mass, the concentration parameters are smaller in the $`\sigma _8=0.75`$ model than they are in the standard $`\mathrm{\Lambda }`$CDM model. This change of concentration is taken into account in our analyses, eventhough its effect is almost negligible.. The results for these MGRSs are shown as dot-dashed lines in Fig. 16. As one can see, the agreement with observational results in this model is much better than in the standard $`\mathrm{\Lambda }`$CDM model. These results suggest that a $`\mathrm{\Lambda }`$CDM model with $`\sigma _80.75`$ may match all the observational results obtained from the 2dFGRS. Unfortunately, in the absence of proper $`N`$-body simulations for this model, it is impossible to make a more detailed comparison with observation.

The question is, of course, whether such a low $`\sigma _8`$ is compatible with other independent observations. Currently, the value of $`\sigma _8`$ is constrained mainly by three types of observations: weak lensing surveys, cluster mass functions, and anisotropy in the cosmic microwave background. Recent cluster abundance analyses give values of $`\sigma _8`$ (assuming $`\mathrm{\Omega }_m=0.3`$) in a wide range, from $`0.6`$ to $`1`$ (e.g. Borgani et al. 2001; Seljak 2002; Viana, Nichol & Liddle 2002; Pen 1998; Fan & Bahcall 1998; Reiprich & Böhringer 2002). Results from weak lensing surveys are equally uncertain, with $`\sigma _8`$ spanning the range $`0.7`$ to $`1`$ (e.g., Jarvis et al. 2003; Hoekstra, Yee & Gladders 2002; Refregier et al. 2002; Bacon et al. 2003). Thus, our preferred value, $`\sigma _8=0.75`$, is consistent with these observations. At the moment, the most stringent constraint on the value of $`\sigma _8`$ is from WMAP (Spergel et al. 2003): $`\sigma _8=0.84\pm 0.04`$ ($`1\sigma `$ error). Even taking this result at face value, one cannot rule out $`\sigma _8=0.75`$ with any high confidence. Thus, there is no strong observational evidence to argue against a $`\mathrm{\Lambda }`$CDM model with $`\sigma _8=0.75`$. Furthermore, as discussed in van den Bosch, Mo & Yang (2003), a value of $`\sigma _8`$ as low as $`0.75`$ can also help to alleviate several problems in current models of galaxies formation, such as the ones in connection to the Tully-Fisher relation and to the rotation curve shapes of low-surface brightness galaxies. Our results presented here give additional support for a relatively low value of $`\sigma _8`$.

## 7 Conclusions

In this paper, we have used realistic halo occupation distributions, obtained using the conditional luminosity function technique introduced by Yang, Mo & van den Bosch (2003), to populate dark matter haloes in high-resolution simulations of the $`\mathrm{\Lambda }`$CDM ‘concordance’ cosmology. The simulations follow the evolution of $`512^3`$ dark matter particles in periodic boxes of $`100h^1\mathrm{Mpc}`$ and $`300h^1\mathrm{Mpc}`$ on a side. Subsequently, the dark matter haloes identified in these simulations are populated with galaxies of different luminosity and different morphological type.

We have shown that the luminosity functions and the correlation lengths as function of luminosity, both for the early- and the late-type galaxies, are in good agreement with observations. Since these same observations were used to constrain the conditional luminosity functions, which in turn were used to populate the dark matter haloes, this agreement shows that the halo occupation statistics obtained analytically can be implemented reliably in $`N`$-body simulations to construct realistic, self-consistent, mock galaxy distributions. We have demonstrated that the details of the spatial distribution of galaxies within individual dark matter haloes have only a very mild effect on the two-point correlation function, and only at real-space separations $`r0.3h^1\mathrm{Mpc}`$.

The mean pairwise peculiar velocities, $`v_{12}`$, however, depend rather strongly on whether satellite galaxies (any galaxy in a dark matter halo other than the most luminous, central galaxy) are associated with random dark matter particles of the friends-of-friends (FOF) group, or whether they are assigned peculiar velocities assuming a spherical, isotropic velocity distribution around the central galaxy. In the former case, $`v_{12}`$, which indicates the amount of infall around overdensity regions, is similar to that of the dark matter. In the latter case, $`v_{12}`$ is significantly suppressed with respect to the dark matter. This difference indicates that the outer parts of the FOF-groups are not yet virialized.

The pairwise velocity dispersions (PVDs) of the galaxies are found to be significantly smaller than those of dark matter particles. Since the PVD is a pair weighted measure for the potential well in which dark matter particles (galaxies) reside, this can be understood as long as the average number of galaxies per unit halo mass, $`N/M`$, decreases with $`M`$ (Jing et al. 1998). Indeed, the halo occupation numbers inferred from our conditional luminosity function indicate that $`N/MM^a`$ with $`a0.2`$.

Stacking a number of $`100h^1\mathrm{Mpc}`$ and $`300h^1\mathrm{Mpc}`$ simulation boxes allows us to construct mock galaxy redshift surveys (MGRSs) that are comparable to the 2dFGRS in terms of sky coverage, depth, and magnitude limit. For each of these MGRSs we estimate the two-point correlation functions $`\xi (r_p,\pi )`$. These are used to derive a number of statistics about the large scale distribution of galaxies, which we compare directly with the 2dFGRS results. In particular, we calculate the projected 2PCFs $`w_p(r_p)`$ as function of luminosity and type. The best-fit power-law slope and correlation lengths of these projected correlation functions are found in good agreement with the 2dFGRS results obtained by Norberg et al. (2002a). In addition, we also compute $`w_p(r_p)`$ and the redshift space correlation function $`\xi (s)`$ for the entire MGRSs. These are compared to the 2dFGRS results obtained by Hawkins et al. (2003).

Although the agreement with the 2dFGRS data is excellent on scales larger than $`3h^1\mathrm{Mpc}`$, on smaller scales $`w_p(r_p)`$ is about a factor two larger than observed. To investigate this in more detail, we analyzed the redshift-space distortions present in $`\xi (r_p,\pi )`$ by computing the quadrupole-to-monopole ratios $`q(s)`$ and the pairwise velocity dispersions $`\sigma _{12}(r_p)`$. A comparison with the results of Hawkins et al. (2003) shows that the standard $`\mathrm{\Lambda }`$CDM model over-predicts the clustering power on small scales by a factor of about two, and the PVDs by about $`350\mathrm{km}\mathrm{s}^1`$. After examining a variety of possibilities, we find that the only viable solution to these problems is to reduce the power spectrum amplitude, $`\sigma _8`$, from $`0.9`$ to $`0.75`$.

No doubt, in the coming years, new results from the 2dFGRS and the SDSS will significantly improve the data on the large scale distribution of galaxies. The analysis presented here, based on the conditional luminosity function, will hopefully prove a useful tool to further constrain both galaxy formation and cosmology. In this respect, our results regarding constraints on $`\sigma _8`$ are an important illustration of the potential power of this approach.

## Acknowledgement

We thank Ed Hawkins and the 2dFGRS team for providing us with their data in electronic format, the anonymous referee for insightful comments that helped to improve the paper, and Gerhard Börner, Guinevere Kauffmann, Simon White and Saleem Zaroubi for useful discussions. XHY is supported by NSFC (No.10243005) and by USTCQN. YPJ is supported in part by NKBRSF (G19990754) and by NSFC (No.10125314). Numerical simulations presented in this paper were carried out at the Astronomical Data Analysis Center (ADAC) of the National Astronomical Observatory, Japan.


text:
# Accurate Internal Proper Motions of Globular Clusters

## 1. Proper Motions in the Field of M4

Some of us have recently developed a method (Anderson and King, 2000, 2003), based on multi-epoch HST images, which allows us to derive very accurate proper motions. The key to our success is the removal of systematic errors of many types from our astrometry. Presently, we are able to measure the position of a well exposed star in a single image with a precision of 0.02 pixel, i.e. 2 mas in the WFs, 1 mas in the PC. Multiple observations allow to reduce this uncertainty of a factor $`\sqrt{N}`$, where N is the number of frames.

### 1.1. The Internal Proper Motion of M4 and its Distance

The internal proper motions of M4 is $``$0.50 mas/yr. This means that with our 8 frames in each of two epochs, and with a time baseline of $``$5 yrs, we can reach a precision of $``$0.2 mas/yr. We can then subtract, in quadrature, the distribution of the errors to the observed distribution, and have a precise measure of the internal dipersion of proper motions. Figure 1 shows some preliminary results. We plan to furtherly improve the proper motions of M4 stars by the reducing the third epoch WF@ACS images which overlap with our fields.

We have recently undertaken a program aimed at obtaining radial velocities for stars located in the same radial interval and in the same magnitude range of the proper motion data in M4 (and other clusters). A comparison of the proper motion dispersion with the radial velocity dispersion will lead to the direct measurement of the GC distance. The major source of uncertainty being the sampling error (in first approximation $`1/\sqrt{2N_{tot}}`$, where $`N_{tot}`$ is the total number of measured stars), which gives an error 1.3 $`\%`$ for the 3,000 stars we plan to observe in each cluster with the multi-fiber high resolution facility FLAMES@VLT. The statistical limitation will be set by the number of radial velocities, as we have accurate proper motions for a large number of stars.

### 1.2. Measure of the Galactic Constant $`\mathrm{\Theta }_0/R_0`$

An obvious by-product of our proper motion measurements is a decontamination of the cluster sample from the field stars (King et al. 1998, Bedin et al. 2001). There is a number of studies which can be undertaken with the motions of the fields objects that we get for free. We give a brief description of an example of such applications.

In Bedin et al. (2001), two $`HST`$ deep observations of the globular cluster M4, separated by a time base line of $``$ 5 yrs, allowed us to obtain a pure sample of the cluster main sequence stars. The identification of an extra-Galactic point source enables us to use the proper motions of field stars (which were junk in Bedin et al. 2001) to measure a fundamental parameter of the Galaxy. M4 is a globular cluster projected on the edge of the Galactic bulge ($`\mathrm{}9^{}`$, $`b16^{}`$). We expect only a small number of foreground disk stars in our fields, but in the background we look through the edge of the bulge at a height of $`2`$ kpc. Although at such heights the density of the bulge is rather low, the volume we are probing is sizable, so that we see a large number of bulge members.

Their absolute proper motion (pm) is just the reflection of the Sun’s angular velocity with respect to that point; from that pm we can derive the value of the angular velocity of the local standard of rest (LSR) with respect to the Galactic center, which is the fundamental Galactic-rotation constant $`AB=\mathrm{\Theta }_0/R_0`$ (cf. Kerr & Lynden-Bell 1986).

To derive this value we need to: (1) find the mean distance of the bulge stars whose motion we are observing, (2) correct the observed pm for the velocity of the Sun with respect to the LSR, and (3) relate the corrected pm to the angular velocity of the LSR with respect to the Galactic center.

For the distance of the bulge stars that we are observing, we assume the following working hypotheses: (1) Disk and halo stars are a negligible component of the field stars in our M4 images, i.e., the field stars are mainly bulge members. (2) The bulge stars on our line of sight are part of a spherical spatial distribution around the Galactic center. (3) Our observations go deep enough that we do not loose stars on the far side of the bulge. From these assumptions, it follows that we can express the distance of the centroid of the bulge stars along our line of sight (we will refer to it as the bulge) as a geometrical constant$`\times `$the distance of the Sun from the Galactic center. This distance is $`R=R_0\mathrm{cos}\mathrm{}\mathrm{cos}b.`$ If we take $`R_0=8.0`$ kpc, then $`R=7.6`$ kpc.
The absolute motion of the bulge –corrected for the Sun’s peculiar motion– allows us to get a direct estimate of the Oort-constant difference $`(AB)`$, which is related to the transverse velocity of the LSR ($`\mathrm{\Theta }_0`$) and its Galactocentric distance ($`R_0`$), according to equations presented in Bedin et al. (2003a; in Bedin et al. 2003b we will give more details and a more concise formulation leading to identical results), we get: $`(AB)\pm \sigma _{(AB)}`$ $`=`$ $`\mathrm{\Theta }_0/R_0\pm \sigma _{\mathrm{\Theta }_0/R_0}`$ $`=`$ $`27.6\pm 1.7`$ $`\mathrm{km}/\mathrm{s}\mathrm{kpc}`$. The quoted error is internal and corresponds to an uncertainty of 7$`\%`$.

### 1.3. Galactic Model

Finally, we mention that the field star data can be used to probe in depth the Galactic structure and stellar content. We are presently collaborating with the research group at the University of Pisa on a project aimed at the definition of a Galactic model. An example on how the field star data from proper motions as those in King et al. (1998) and Bedin et al. (2001) can be used in this kind of investigation is in Castellani et al. (2001).

## References

Anderson, J. & King, I. R 2000, PASP, 112, 1360

Anderson, J. & King, I. R 2003, PASP, 115, 113

Bedin, L. R., Anderson, J., King, I. R. & Piotto, G. 2001 ApJ, 560, L75

Bedin, L. R., Piotto, G., Anderson, J. & King, I. R. 2003a, Mem. Soc. Astron. Italiana, 10-12 April 2002, 74, in press

Bedin, L. R., Piotto, G., Anderson, J. & King, I. R. 2003b, AJ, submitted

Castellani, V., Degli’Innocenti, S., Petroni, S., Piotto, G. 2001, MNRAS, 324, 167

Kerr, F. J. & Lynden-Bell, D. 1986, MNRAS, 221, 1023

King, I. R., Anderson, J., Cool, A. & Piotto, G. 1998 ApJ, 492L, 37