A QUANTITATIVE EXPLANATION OF THE OBSERVED POPULATION OF MILKY WAY SATELLITE GALAXIES

Sergey E. Koposov; Jaiyul Yoo; Hans-Walter Rix; David H. Weinberg; Andrea V. Macciò; Jordi Miralda Escudé

doi:10.1088/0004-637X/696/2/2179

1. INTRODUCTION

The inflationary cold dark matter (CDM) scenario predicts an initial fluctuation spectrum with power that continues down to small scales, and in consequence it predicts a mass function of dark matter (DM) halos that rises steeply toward low masses. A significant fraction of these halos survive as gravitationally self-bound units long after falling into more massive halos. As pointed out forcefully by Klypin et al. (1999) and Moore et al. (1999), the predicted number of subhalos within a Milky Way (MW) like galaxy halo greatly exceeded the then known numbers of MW or Local Group dwarf satellites, when subhalos and observed dwarfs were matched based on velocity dispersion or corresponding circular velocity (see also Kauffmann et al. 1993). This discrepancy between predicted and observed numbers has become known as the "missing satellite problem."

Proposed solutions fall into three general categories. The first modifies the properties of DM or the primordial fluctuations from inflation in a way that eliminates the low-mass DM subhalos themselves (e.g., Kamionkowski & Liddle 2000; Spergel & Steinhardt 2000; Bode et al. 2001; Zentner & Bullock 2003). The second appeals to astrophysical mechanisms that suppress star formation in low-mass halos so that they do not become observable dwarf satellites; photoheating by the metagalactic UV background is an attractive mechanism because it naturally introduces a cutoff at approximately the correct velocity scale (Bullock et al. 2000; Somerville 2002; Kravtsov et al. 2004). The third possibility, arguably a variant of the second, is that the numerous dwarf companions of the MW actually exist but have been missed by observational searches.

In this paper, we revisit the "missing satellite problem" with particular emphasis on the role of the new dwarf companions discovered in imaging data from the Sloan Digital Sky Survey (SDSS; York et al. 2000; Adelman-McCarthy et al. 2008). There are now about a dozen of these (Willman et al. 2005; Belokurov et al. 2006, 2007; Zucker et al. 2006; Irwin et al. 2007; Koposov et al. 2007; Walsh et al. 2007; a couple of systems still have ambiguous status), most of them at least an order of magnitude less luminous than the faintest of the previously known "classical" satellites.⁷ Spectroscopic follow-up (e.g., Martin et al. 2007; Simon & Geha 2007; Geha et al. 2009) for many of them indicates that they are indeed DM-dominated systems, even though most are fainter than typical globular clusters, as low as only ∼1000 L_☉ (e.g., Belokurov et al. 2007; Martin et al. 2008). Remarkably, almost all of the newly found faint satellite galaxies have stellar velocity dispersions in the range 3–10 km s⁻¹, though their luminosities vary widely. Similarly, the total masses within the inner 300 pc span less than an order of magnitude (Strigari et al. 2008).

Since the SDSS imaging in which these satellites have been discovered covers only ∼20% of the sky, a naive accounting would increase the estimated number of MW companions by 5 × 12 = 60, in addition to the 10 classical satellites. However, Koposov et al. (2008) use a well-defined identification algorithm to show that the SDSS dwarfs are also subject to strong radial selection effects. Most of the newly discovered objects could only have been found within distances of 50–100 kpc, much smaller than the inferred virial radius of the MW's DM halo (∼280 kpc for $\rho _{\rm vir}/\bar{\rho } = 340$ ; Xue et al. 2008). The faintest SDSS dwarfs are detectable over only 1/1000 of the halo virial volume (including the factor of 5 for sky coverage). Walsh et al. (2009) have recently reached similar conclusions based on an independent identification algorithm and independent Monte Carlo tests.

Such analyses are the basis for "volume corrections" for the faint MW satellite population. With proper volume corrections applied, the luminosity function of faint MW satellite galaxies turns out to be a rather shallow power law in the range −15 < M_V < −3 (Koposov et al. 2008). These results in turn imply that the number of satellites brighter than M_V = −3 is ∼80 or more, and the number above M_V = 0 could be a few hundred. Tollerud et al. (2008) reached a similar conclusion, adopting a radial satellite distribution based on the Via Lactea simulation of Diemand et al. (2007). Even this census counts only dwarfs that are above the effective surface brightness threshold for SDSS detection. With the Koposov et al. (2008) detection algorithm, this threshold is approximately 30 mag arcsec⁻² (V band), with a weak dependence on luminosity and distance. The dwarfs found in SDSS have surface brightnesses that range from 24 to 30 mag arcsec⁻².

Studies of the high-redshift Lyα forest indicate that the small-scale power expected in the standard ΛCDM scenario (inflationary CDM with a cosmological constant) is indeed present in the primordial fluctuation spectrum (Narayanan et al. 2000; Viel et al. 2005; Abazajian 2006; Seljak et al. 2006). Astrophysical suppression of star formation, and photoionization suppression in particular, has emerged as the most plausible and hence popular solution to the "missing satellite" conundrum. Within this category, there have been different proposals about what subhalos host the observed dwarf satellites. Bullock et al. (2000) suggested that the observed dwarfs are those whose subhalos assembled a substantial fraction of their mass before reionization, and thus before the onset of photoionization suppression. Stoehr et al. (2002) suggested that the measured stellar velocity dispersions are well below the virial velocity dispersions of the DM subhalos, and that the observed dwarfs occupy subhalos that are still above the velocity threshold where star formation suppression occurs. Kravtsov et al. (2004) used N-body simulations to show that roughly 10% of subhalos lose a large fraction (∼90%) of their mass during dynamical evolution without being completely disrupted; they suggested that the observed dwarfs occupy subhalos that were above the suppression threshold at the time they became satellites but have suffered extensive mass loss since then. These papers and others (e.g., Somerville 2002; Strigari et al. 2007; Orban et al. 2008) focus on explaining the classical (pre-SDSS) dwarf spheroidal population, with luminosities in the range −8 < M_V < −15 (excluding the Magellanic Clouds) and stellar velocity dispersions in the range 8 km s⁻¹ < σ_* < 25 km s⁻¹. The recently discovered SDSS dwarfs have much lower luminosities (−8 < M_V < −1.5), lower surface brightness, and somewhat lower velocity dispersion (σ_* ∼ 5 km s⁻¹), so they could have a distinct formation mechanism, or they could form a continuum with the classical dwarf spheroidals.

The new SDSS discoveries and their quantified detectability are the basis for the model-data comparison in this paper. We construct and test models of the MW dwarf satellite population that incorporate Monte Carlo realizations of merger trees for 10¹² M_☉ (main galaxy) halos, a detailed analytic model for the dynamical evolution and disruption of subhalos, and a variety of recipes for assigning stellar masses to these subhalos motivated by ideas in the existing literature. For most of our models, we assume that a subhalo can only accrete gas to form stars (1) before the epoch of reionization or (2) after reionization if its virial velocity exceeds a critical threshold before it enters the MW halo and becomes a satellite. The spirit of the exercise is similar to that of Bullock et al. (2000), but the dynamical modeling of subhalos is more sophisticated, and we are now in a position to directly include the (strong) constraints imposed by the SDSS dwarfs accounting for the radial selection function found by Koposov et al. (2008). In contrast to most previous studies, we treat the luminosity distribution as the primary test of models, rather than the stellar velocity dispersions or central masses (Strigari et al. 2007, 2008; Li et al. 2008; Macciò et al. 2009), or the inferred but unobservable subhalo circular velocities. This emphasis is motivated by the fact that the luminosity is the foremost quantity that matters for the observational selection. We consider stellar velocity dispersions and central masses as an additional test, but their interpretations are affected by the uncertainty in the DM profiles of the subhalos associated with observed dwarfs.

2. THE POPULATION OF DM SUBHALOS IN THE MILKY WAY

Our model for the MW satellites is based on the CDM scenario, with each satellite forming initially in a separate DM halo that at some point falls into the MW's DM halo. We refer to the bound DM satellites orbiting in the MW halo as subhalos. A subhalo may or may not correspond to a dwarf satellite galaxy, depending on whether it contains an observable number of stars. In this section, we describe our model for computing the dynamical evolution of subhalos.

We use the dynamical DM-only model of subhalos developed by Yoo et al. (2007) to compute the subhalo population and its orbital distribution. This model is described in detail in Yoo et al. (2007), where a much larger halo of 10¹⁵ M_☉ was considered as a model of a massive cluster of galaxies. Here, we consider instead a final halo of 10¹² M_☉ at the present time as a representation of the MW galaxy. Despite the change in the final halo mass, the model remains basically the same as described in Yoo et al. (2007), so here we make only a brief summary of its description.

The model uses the extended Press–Schechter formalism to generate a Monte Carlo merger tree of the parent halo at the present time (Press & Schechter 1974; Bond et al. 1991). We follow the dynamical evolution of all the subhalos with masses M_h>10⁶ M_☉ until they merge with the MW and lose their mass below M_h = 10⁵ M_☉. All halos start as isolated objects, and they grow in mass by accretion and mergers for as long as they remain isolated. At some redshift, z_sat, they merge into a larger halo (either the MW or another object that will become a MW subhalo). After this merger, the object has become a satellite or subhalo, and it stops growing in mass. It can subsequently lose mass by tidal stripping when it passes near the center of its parent halo or undergoes encounters with other subhalos. The subhalo is subject to dynamical friction, which tends to shrink its orbit, and to random encounters with other subhalos, which on average expand the orbit. The orbital eccentricity is also subject to random variations. The model allows for the presence of subhalos within other subhalos. When a subhalo is disrupted, any subhalos it contained are dispersed into the new, larger parent halo. This simple analytic model is able to reproduce the subhalo mass function, in reasonably good agreement with that found in numerical N-body simulations (Zentner et al. 2005; Shaw et al. 2006; Yoo et al. 2007). For the present purpose, this approach has the advantage (over N-body) of easily affording the required mass resolution and multiple halo realizations.

We adopt a flat ΛCDM cosmology with matter density Ω_m = 0.24, baryon density Ω_b = 0.04, power spectrum normalization σ₈ = 0.8, Hubble constant h = 0.7, and a primordial spectral index n_s = 0.95, consistent with recent measurements (Tegmark et al. 2006; Spergel et al. 2007). The matter power spectrum is computed by using the transfer function of Eisenstein & Hu (1999). We generate six Monte Carlo merger trees of a MW-sized halo. Each realization provides the subhalo mass function, their orbital elements, and density profiles at the present time. Our statistical results are the average of the six different realizations.

The dynamical model of Yoo et al. (2007) uses the Jaffe profile and its velocity dispersion to model subhalos and their dynamical interactions, for reasons of numerical simplicity and because large galaxies that are tidally limited satellites of a larger halo are reasonably well modeled by a Jaffe sphere for their baryon plus DM density profiles. However, the very low mass dwarf satellites tend to be dominated by DM even in their inner parts. We therefore make an adjustment to better connect our Monte Carlo simulation results to the observed MW dwarf galaxies: we use the subhalo masses and orbital elements, which are the quantities most robustly computed in the Yoo et al. (2007) model, but we calculate the density profiles and velocity dispersions of subhalos assuming that they have an NFW (Navarro, Frenk, & White) profile (Navarro et al. 1997). Using the standard spherical collapse model, the virial radius of an isolated halo is assigned as

$\begin{equation} R_{\rm vir}= \left[ \frac{3\, M_{\rm halo}}{4\, \pi \,\Delta _c \, \bar{\rho }_m(z)} \right]^\frac{1}{3}, \end{equation} \tag{ 1 }$

where Δ_c = (18 π² + 82 x − 39 x²)/(1 + x) (Bryan & Norman 1998), x = −(1 − Ω_m)/(Ω_m(1 + z)³ + 1 − Ω_m), and the mean cosmic density is $\bar{\rho }_m(z)=\Omega _m\rho _c(1+z)^3$ . The halo concentration c is computed using the relation from Bullock et al. (2001), scaled to σ₈ = 0.8 according to Macciò et al. (2007), with c = 0.8 × 9 × (M_h/10¹³ h⁻¹ M_☉)^−0.13/(1 + z).

For the model in this paper, we use in particular the subhalo masses at two different special epochs: M_rei ≡ M_tot(z = z_rei) when the universe reionizes and the photoionization background starts to suppress the star formation efficiency in low-mass halos, and M_sat ≡ M_tot(z_sat) at the epoch when a halo merges into a larger halo and we presume that subsequent star formation and gas accretion are halted in the subhalo. We shall also use below the halo circular velocity V_circ, which is the virial circular velocity, V_circ ≡ [(G M_tot)/R_vir]^1/2 (here M_tot refers to the total mass including DM and a universal fraction of baryons).

In Figure 1 we show the distribution of M_subhalo at redshifts z = 0 (left panel), z = z_sat (middle panel), and z = 8, 11, and 14 (right panel). As expected, the mass distribution is close to a power law, except near the resolution limit of our simulations. Also, in Figure 2 we show the accretion history of the MW subhalo population by plotting the halo masses of the present-day MW subhalos at the time of accretion versus the redshift at which they became satellites of larger halos. We see that most of the MW subhalos became satellites at z < 2. Most of the accreted satellites have small circular velocities V_circ < 20 km s⁻¹, so they lie in a range where gas accretion and star formation are likely to be suppressed after the epoch of reionization (Quinn et al. 1996; Thoul & Weinberg 1996; Bullock et al. 2000).

**Figure 2.** Epochs when (sub-)halos were accreted into larger halos, and masses at that time. This figure illustrates the results from one Monte Carlo realization of the semianalytic model, with each point showing the redshift z_sat at which a subhalo first became a satellite in a larger halo against its total mass M_sat at that epoch. The small panel on the right shows the distribution of z_sat. The solid and dashed lines show the locus of halos with V_circ(z_sat) = 40 and 20 km s⁻¹, respectively.
Download figure:
Standard image High-resolution image

3. POPULATING THE DM HALOS WITH STARS

3.1. Recipes to Assign Stellar Masses to Subhalos

To make direct observational predictions from these models, we populate each subhalo in a given Monte Carlo realization with stars according to a sequence of recipes, and then test how many of these satellites could have been found within the SDSS. Some of these recipes are mathematically simple illustrations, while others are motivated by the expected effects of ionization and cooling physics as discussed in the Introduction. For reference, the nomenclature of the recipes is summarized in Table 1. In all cases we calculate the stellar mass based on the subhalo mass (DM plus baryons in the universal fraction) at the accretion epoch z_sat, denoted as M_sat. We implicitly assume that satellites do not accrete new material to form additional stars and that tidal stripping of the DM does not affect the stellar content of the satellite if it survives to the present day. Simulations suggest that these assumptions are reasonable but not perfect approximations (Peñarrubia et al. 2008; Simha et al. 2008).

Table 1. List of Models Used

Model Name	Present-Epoch Stellar Mass
1A	M_* = f_* × M_sat
1B	M_* = f_* × `min`((M_sat/M₀)^α, 1) × M_sat
2
3A	$M_{} = \frac{f_{} \times (M_{\rm {sat}}-M_{\rm {rei}})}{(1+0.26\,(V_{\rm {crit}}/V_{\rm {circ}})^3)^3} + f_{*}\times M_{\rm {rei}}$
	Same as 3A for halos with V_circ(z_rei)>V_crit,r,
3B	For halos with V_circ(z_rei) < V_crit,r
	$M_{} = \frac{f_{} \times M_{\rm {sat}}}{(1+0.26\,(V_{\rm {crit}}/V_{\rm {circ}})^3)^3}$

Download table as: ASCII Typeset image

We begin with the simplest model (denoted by Model 1A) that the stellar mass is a constant fraction of the subhalo mass at the time of accretion into the main halo:

$\begin{equation} M_*=f_*\times M_{\rm {sat}}. \end{equation} \tag{ 2 }$

The arguments of Klypin et al. (1999) and Moore et al. (1999) suggest that this model will fail badly, and we show that it does indeed fail despite the new satellite discoveries and the radial selection biases that affect them. There is ample evidence that the efficiency of star formation declines rapidly toward low masses even well above the dwarf satellite regime (e.g., van den Bosch et al. 2007). In Model 1B, we allow the stellar fraction to vary as a power law of M_sat below a threshold M₀:

$\begin{equation} M_* = f_* \times {\tt min}\left(\left(\frac{M_{\rm {sat}}}{M_0}\right)^\alpha,1\right) \times M_{\rm {sat}}. \end{equation} \tag{ 3 }$

Our second approach to modeling stellar masses includes the effects of a pervasive energetic radiation field after the epoch of reionization, which heats gas and hence keeps it from accumulating at the centers of low-mass halos. Calculations by Quinn et al. (1996) and Thoul & Weinberg (1996) showed that gas accretion in halos with the circular velocities below V_circ ∼ 30–40 km s⁻¹ is strongly suppressed, while substantially larger halos are minimally affected (see also Weinberg et al. 1997; Gnedin 2000). In this spirit, we assume that halos below a critical circular velocity form no stars after reionization, and we thus assign stellar masses

$\begin{eqnarray} &&M_*= f_*\times M_{\rm {sat}} \mbox{if }V_{\rm {circ}}(z_{\rm {sat}})>V_{\rm {crit}} \end{eqnarray} \tag{ 4 }$

$\begin{eqnarray} &&\hspace{25pt} f_*\times M_{\rm {rei}} \mbox{if }V_{\rm {circ}}(z_{\rm {sat}})<V_{\rm {crit}}. \end{eqnarray} \tag{ 5 }$

This model (Model 2) has three adjustable parameters–f_*, V_crit, and z_rei—with expectations that V_crit ∼ 20–40 km s⁻¹ and z_rei ∼ 11 (e.g., Weinmann et al. 2007; Dunkley et al. 2009). The approach is similar to that of Bullock et al. (2000), except that we treat V_crit as free parameter, and the stellar mass formed before the epoch of reionization is assigned using M_* = f_* × M_rei, instead of simply dividing galaxies into "observable" or "unobservable" classes based on the fraction of the mass accreted by z_rei.

Our third class of models is similar to the second, but it replaces the sharp threshold of Equation (4) with the continuous transition found in numerical simulations by Gnedin (2000), Hoeft et al. (2006), and Okamoto et al. (2008). The numerical results in these papers can be described fairly well by a formula similar to that in Gnedin (2000), with the fraction of baryons that cool in low-mass halos suppressed by a factor of [1 + 0.26(V_crit/V_circ)³]⁻³; well after the reionization redshift, the critical velocity is found to be approximately independent of redshift. Gnedin (2000) found V_crit ∼ 40 km s⁻¹, but these results were artificially affected by numerical resolution (N. Gnedin 2008, private communication). Hoeft et al. (2006) and Okamoto et al. (2008) find V_crit ∼ 25–30 km s⁻¹. Including the pre-reionization contribution to M_*, this model (Model 3A) becomes

$\begin{equation} M_* = \frac{f_* \times (M_{\rm {sat}}-M_{\rm {rei}})}{(1+0.26\,(V_{\rm {crit}}/V_{\rm {circ}}(z_{\rm {sat}}))^3)^3}+f_*\times M_{\rm {rei}}. \end{equation} \tag{ 6 }$

The assumption that all halos can form stars before z_rei may not be justified because in halos with virial temperature T_vir ≲ 10⁴ K (V_circ ≲ 10 km s⁻¹) the gas does not get hot enough to cool by atomic processes, and simulations that include molecular cooling suggest that gas cooling and star formation are very inefficient in such halos (Haiman et al. 1997; Barkana & Loeb 1999; Machacek et al. 2001; Wise & Abel 2007; Bovill & Ricotti 2009; O'Shea & Norman 2008). We will therefore consider variant models (Model 3B) that eliminate stellar mass in pre-reionization halos below a critical threshold V_crit,r ∼ 10 km s⁻¹.⁸ In Model 3B, halos with V_circ(z_rei) < V_crit,r have stellar mass

$\begin{equation} M_* = \frac{f_* \times M_{\rm {sat}}}{(1+0.26\,(V_{\rm {crit}}/V_{\rm {circ}}(z_{\rm {sat}}))^3)^3}, \end{equation} \tag{ 7 }$

while halos with V_circ(z_rei)>V_crit,r have mass given by Equation (6).

To determine very roughly the plausible range of values for the stellar mass fraction f_*, we can refer to the results of Strigari et al. (2007), who derived M(<r_tidal)/L) = 30–800 M_☉/L_☉ for the classical dwarfs, and Simon & Geha (2007), who measured velocity dispersions for SDSS dwarfs and inferred total mass-to-light ratios of 140–1800 M_☉/L_☉. For a stellar mass-to-light ratio M_*/L_V = 1 M_☉/L_☉, we infer plausible values of f_* ∼ 10⁻⁴–10⁻², though these are very uncertain because all the dynamical mass-to-light ratio determinations suffer from the fact that the stars in luminous bodies of the dwarf spheroidals (dSphs) probe only the inner parts of the DM potential wells. Another line of argument comes from matching the mean space density of DM halos to that of observed field dwarfs: Tinker & Conroy (2009) find f_* ≈ 10^−3.6 at absolute magnitude M_r ≈ −10. In the rest of the paper, we will frequently refer to the stellar mass fraction normalized by the universal baryon fraction:

$\begin{equation} F_*\equiv \frac{f_*}{\Omega _b/\Omega _m} = 6.25 f_*. \end{equation} \tag{ 8 }$

Note that f_* and F_* refer to stellar fractions in halos where the efficiency is not suppressed, i.e., V_circ(z_sat)>V_crit. We will frequently refer to the quantity (M_*/M_sat) × (Ω_m/Ω_b) as the "star formation efficiency," by which we mean the efficiency with which the halo converted the baryons available to it at z_sat (for a universal baryon fraction) into stars observable at z = 0.

3.2. Detectability and Observable Properties for the Simulated Satellites

Color–magnitude diagrams for the faint dwarf spheroidal galaxies in the MW halo show that the stellar populations are predominantly "old" (older than several Gyr) and metal poor ([Fe/H] ≲ −1). To convert stellar masses to luminosities, we assume that all of our model dwarfs have a stellar mass-to-light ratio M_*/L_V ≈ 1 M_☉/L_☉ appropriate to an old, metal-poor population (Bruzual & Charlot 2003; Martin et al. 2008). The light of the lowest luminosity dwarfs can be dominated by a handful of bright stars and thus subject to stochastic variations. We ignore this complication; our "luminosities" are simply scaled stellar masses: L_V/L_☉ = M_*/M_☉. This seems appropriate, since the luminosities of the dwarfs galaxies are usually measured either by integrating over the luminosity function of old stellar population matched to the observed luminosity function of stars in dwarfs (Belokurov et al. 2006) or by averaging over possible stochastic variations of galaxy luminosity (Martin et al. 2008).

The detectability of a faint stellar MW satellite galaxy in an SDSS-like search depends on its luminosity and its distance from the Sun, as quantified by Koposov et al. (2008) (see also Walsh et al. 2009). On the basis of these results (Figure 12 of Koposov et al. 2008), we model the detectability of each simulated satellite as a binary decision using the criterion

$\begin{equation} \log _{10}(D_\odot /1\;{\rm kpc})<1.1\hbox{--}0.228 M_V. \end{equation} \tag{ 9 }$

Our simulations provide the current Galactocentric distance and orbital apocenter and pericenter for each subhalo, but not the orientation of the orbit. We therefore assign the heliocentric distance of the satellites

$\begin{equation} D_\odot =\sqrt{8.5^2+D_{\rm GC}^2-2\times 8.5\times D_{\rm GC}\,{\rm cos}(\phi)}, \end{equation} \tag{ 10 }$

where D_GC is the Galactocentric distance (in kpc) from the simulations and cos(ϕ) is a random variable uniformly distributed between −1 and 1 (ϕ is the angle between radial vectors from the GC to the Sun and to the subhalo). This method assumes that the satellite orbits are isotropically distributed across the sky (see Tollerud et al. 2008, for discussion of the validity of this approximation). As expected from Koposov et al. (2008), accounting for the detectability of satellites causes the "observable" population to differ strongly from the "simulated" one; only the brightest satellites are observable throughout the virial volume.

Not surprisingly, the Koposov et al. (2008) analysis also reveals a surface brightness threshold for dwarf detection, which is approximately 30 mag arcsec⁻² with little dependence on distance. We assume that any model dwarf that passes the luminosity threshold also passes the surface brightness threshold. Many recent SDSS satellite discoveries do lie near that survey's surface brightness limit; this assumption can therefore only be tested with the next generation of sky surveys. We discuss implications of this assumption in Section 5.

With a model that assigns stellar luminosities to each satellite halo, we can predict the expected stellar velocity dispersions for comparison with those measured for MW satellites by Walker et al. (2007), Simon & Geha (2007), and Martin et al. (2007). This can be done straightforwardly if we assume that the stars are test particles—an assumption supported by the observed (M/L)_dyn(<R_eff) ≫ (M/L)_*(<R_eff)—orbiting in an NFW potential with an isotropic velocity dispersion. Then, we can use the Jeans equation (Jeans 1919) to derive the velocity dispersion profile of stars:

$\begin{equation} \frac{d(\nu (r)\,\sigma ^2(r))}{\textit{d\,r}} + \nu (r) \, \frac{\textit{G\,M}(r)}{r^2}=0, \end{equation} \tag{ 11 }$

where ν is the density distribution of stars (see Strigari et al. 2007 for a more detailed treatment). Here we assume that the density of stars follows a Plummer profile ν ∝ [1 + (r/r_p)²]⁻² (Plummer 1911), which seems to fit observed density profiles reasonably well (Wilkinson et al. 2002; Belokurov et al. 2007). The mass profile M(r) used here is computed based on the virial radii and concentrations at the redshift z_sat of subhalo accretion. While the outer parts of the subhalos are tidally stripped, Peñarrubia et al. (2008) show that the stars and the inner part of the DM subhalo are stripped only at a very late stage, when the subhalo is close to complete disruption. They also show that the velocity dispersion in subhalos is a function of the total DM mass remaining bound inside the luminous body and therefore remains nearly constant until this late stage.

After numerically solving the Jeans equation, we compute the expected light-weighted velocity dispersion within the optical radius as

$\begin{equation} \sigma _*=\frac{\int \nu (r)\sigma (r) \,dx\,dy\,dz}{\int \nu (r) \,dx \,dy \,dz}, \end{equation} \tag{ 12 }$

where the integration is done over a cylinder within a radius $R = \sqrt{x^2+y^2}$ equal to the Plummer radius of the galaxy; the integral extends over ±∞ in z. The stellar velocity dispersion depends on the radial extent of the stellar tracers, which cannot be predicted within our simple modeling context (see also Benson et al. 2002). We therefore use the observed properties of the faint MW satellites to choose stellar radii, based on Martin et al. (2008). Specifically, we adopt Plummer radii r_p = 150 pc for M_V < −5, and for fainter dwarfs we adopt a linear relation between log r_p and M_V with r_p rising from 20 pc at M_V = 0 to 150 pc at M_V = −5.

The additional important component of the detectability is the tidal disruption of the satellite galaxies. Although our semianalytic model of DM subhalo evolution properly accounts for the tidal disruption of subhalos, it does not allow for the possibility that stars have been dispersed in a tidal stream while a small core of the subhalo survives. Here we simply classify a subhalo as unobservable if its current tidal radius is less than the expected Plummer radius of the stellar body, which would imply substantial tidal disruption of the stellar component. We also presume that a satellite is unobservable if its host subhalo has lost more than 99% of its original mass to tidal stripping.

4. RESULTS

4.1. Stellar Mass Function of the Full Satellite Populations

Figure 3 shows the predicted distribution of the stellar masses of satellites within R_virial = 280 kpc, assuming 4π sky coverage and complete satellite detectability. In the left panel, the solid curve shows Model 1A with a constant F_* = 10⁻³, making the stellar mass function a scaled version of the DM subhalo mass function. Introducing mass-dependent suppression, Model 1B with α = 1 (dashed) and α = 2 (dotted) lowers the low-mass end of the stellar mass function as expected. Since this model also adopts F_* = 10⁻³ = const. above M_sat = M₀ = 10¹⁰ M_☉, the high-mass end of the mass function is unchanged.

The middle panel of Figure 3 shows Model 2, with post-reionization suppression of star formation in halos below a sharp circular velocity threshold, either V_crit = 40 km s⁻¹ (solid) or V_crit = 20 km s⁻¹ (dashed), where we have adopted F_* = 10⁻³ and a reionization redshift z_rei = 11. The resulting stellar mass functions for the satellite galaxies are strongly bimodal, with the low-mass portion corresponding to dwarfs in which all stars formed before reionization and the high-mass portion corresponding to halos that exceeded the critical velocity threshold before becoming satellites, V_circ(z_sat)>V_crit. The low-mass portion is just a scaled version of the subhalo mass function at z = z_rei. Above M_* ≈ 10^6.5 M_☉ the host halos are all massive enough to have star formation after z_rei, and the mass function is the same as that of Model 1. If the velocity threshold is lowered to V_crit = 20 km s⁻¹, the high-mass peak in the distribution of satellite stellar masses extends to lower values before photoionization suppression cuts it off.

The bimodal appearance of the middle panel of Figure 3 is a direct consequence of the sharp V_circ threshold for photoionization suppression. The right-hand panel shows predictions for several variants of Models 3A and 3B, with the Gnedin (2000) formula (Equation (6)) used to describe photoionization suppression. With this smooth suppression, the "pre-reionization" and "post-reionization" portions of the mass function join to make a smooth overall mass function. The low-mass end of the mass function is now a mix of satellites that formed their stars before reionization and satellites with V_circ(z_sat) < V_crit whose post-reionization star formation was strongly suppressed but not completely eliminated. Lowering the assumed reionization redshift from z_rei = 11 to z_rei = 8 boosts the stellar mass function below M_* = 10⁴ M_☉. Conversely, if we eliminate pre-reionization SF in dwarfs with V_circ(z_rei) < V_crit,r = 10 km s⁻¹ (thick solid line, Model 3B), the number of satellites with M_* ⩽ 10³ M_☉ drops by a large factor, while at higher masses the stellar mass function is unaffected. The difference between the thin and thick solid lines is the contribution of satellites that formed stars primarily before reionization in halos with V_circ(z_rei) < 10 km s⁻¹, for z_rei = 11 and V_crit = 40 km s⁻¹.

4.2. Distribution of Observed Dwarf Satellite Luminosities, N(M_V)

Figure 4 illustrates the impact of selection effects on the observable satellite population. For one realization of Model 3B (with parameters that yield a good match to observations), filled circles show satellites that would be detectable in an all-sky, SDSS-like survey (Koposov et al. 2008), and open circles show undetectable satellites. The low end of the luminosity distribution, with M_V ≳ −5, is strongly affected by the radial selection bias.

For direct comparison with observations, we therefore select only those model satellites whose combination of luminosity and distance would make them detectable. At the bright end, M_V < −11, we assume that existing photographic surveys are complete to D_☉ = 280 kpc, and we thus compare the total number of dwarfs across the whole sky to the total population of satellites within the virial radius in the simulation. For M_V ⩾ −11, we randomly select 1/5 of the model galaxies to mimic the 20% sky coverage of SDSS DR5, and we count only those satellites that would be detectable according to the criteria of Koposov et al. (2008). We focus our data-model comparison on N_obs(M_V), the luminosity distribution of known MW satellites. We look at additional tests against stellar velocity dispersions, central masses, and the heliocentric radial distribution in Section 4.3.

The luminosities, distances, and velocity dispersions of the observed MW satellites that we use in all subsequent model-data comparisons were taken from various authors (Mateo 1998; Metz & Kroupa 2007; Martin et al. 2008) and are compiled in Table 2. The sample of SDSS satellites used here consists of those systems above the 50% completeness limits of Koposov et al. (2008). We do not include two systems, BooII and LeoV (Walsh et al. 2007; Belokurov et al. 2008), which do not formally satisfy the very conservative selection limits from Koposov et al. (2008). These limits were chosen to avoid the issue of significant "false positive" detections, at the expense of leaving out two objects that deeper follow-up found to be "real." For the analysis presented here, it is most important that the same selection criteria are applied to the mock satellite observations and the SDSS data. As our analysis subsequently shows, such a small difference in sample size is smaller than the model halo-to-halo variation of the number of galaxies. Therefore the inclusion of omission of these two objects does not affect our results significantly.

Table 2. Satellites Used for the Analysis and Parameters Adopted

Galaxy Name	M_V (mag)	σ_* (km s⁻¹)	D_☉ (kpc)
Bootes	−6.3	6.6	60
Canes Venatici II	−4.9	4.6	150
Carina	−9.4	6.8	100
Coma	−4.1	4.6	45
Canes Venatici I	−8.6	7.6	220
Draco	−8.75	10.0	80
Fornax	−13.2	10.5	138
Hercules	−6.6	5.1	130
Leo I	−11.5	8.8	250
Leo II	−9.6	6.7	205
Leo IV	−5.0	3.3	160
LMC	−18.6	⋅⋅⋅	49
Sagittarius	−12.1	11.4	24
Sculptor	−11.1	6.6	80
Sextans	−9.5	6.6	86
Segue 1	−1.5	4.3	23
SMC	−17.2	⋅⋅⋅	58
Ursa Minor	−9.0	9.3	66
Ursa Major I	−5.5	7.6	100
Ursa Major II	−4.2	6.7	30
Willman I	−2.7	4.3	40

Download table as: ASCII Typeset image

Anyway, as we will see later, the halo-to-halo variation of the number of galaxies in our models is noticeable, so we believe that the fact that we do not include two galaxies should not significantly affect our results.

The left panel of Figure 5 compares our simplest model (M_* ∝ M_sat, Model 1A) to the observed satellite counts, now including the satellite galaxy selection effects in the model. We randomly sample each of the six Monte Carlo halo simulations five times (choosing 1/5 of the faint satellites but always keeping the full set for M_V < −11), compute the mean model prediction as the mean of these 30 samplings, and compute the rms dispersion among these 30 in each absolute magnitude bin. Despite the selection bias against low-luminosity satellites, this model fails drastically for any choice of F_*, predicting a much steeper luminosity function than observed. For example, the model with F_* = 10⁻⁴ matches the observed counts near M_V = −9 but predicts far too many satellites fainter than M_V = −6. Selection effects and newly discovered satellites have not altered this basic discrepancy, first emphasized by Klypin et al. (1999) and Moore et al. (1999). The green band shows the 1σ dispersion in predicted counts, and it is clear that statistical fluctuations will not resolve the discrepancy either.

In the right panel we apply our purely empirical modification, M_*/M_sat ∝ M^α below a halo mass M_sat = M₀ = 10¹⁰ M_☉ (Model 1B). With F_* = 10⁻³ and α = 2, this model achieves reasonable agreement with the observed N_obs(M_V) over the full range 0 ⩾ M_V ⩾ −15. The agreement can be further improved by adjusting F_* and M₀, so it appears that this level of mass-dependent suppression is approximately what is needed to explain the observed shape of N_obs(M_V). Linear suppression (α = 1, green band) is not sufficient, predicting an excess of faint dwarfs when normalized to the bright dwarfs. All of our models fail to match the brightest bin (comprising the Small and Large Magellanic Cloud (SMC and LMC, respectively)); we defer discussion of this discrepancy to the end of this section.

Figure 6 shows the expected N_obs(M_V) distributions for Model 2, which has a sharp V_crit threshold for the suppression of SF after reionization in small halos. As in Figure 3, the predicted N_obs(M_V) is bimodal, with a bright peak corresponding to halos that exceeded V_crit before z_sat and a faint peak corresponding to stars formed before reionization. Raising the stellar fraction F_* with other parameters fixed (orange versus blue) shifts both peaks horizontally to higher M_V; the faint peak also increases in height because the brighter (though still faint) satellites can be seen over a larger fraction of the MW virial volume. Lowering V_crit with other parameters fixed (red versus blue) has no impact on the faint peak, but the bright peak extends to fainter magnitudes and grows in height because lower mass halos can now be populated with stars after reionization. Raising z_rei (green versus blue) with other parameters fixed has no impact on the bright peak, but it shifts the faint peak downward in amplitude and slightly downward in location because halos have accreted less mass by this higher redshift. While photoionization suppression reduces the discrepancy with the number of faint satellites seen in Model 1A, these sharp threshold models predict a gap between the faint and bright satellites that is clearly at odds with the data.

Figure 7 compares the Model 2 predictions with those of Model 3A, which uses the Gnedin (2000) formula to incorporate a smoothly increasing suppression of the stellar mass fraction in halos with V_circ(z_sat) ≲ V_crit. In both cases we use parameters F_* = 10⁻³, V_crit = 35 km s⁻¹, z_rei = 11. Model 3A is more physically realistic than Model 2, with a mass-dependent suppression that is calibrated on numerical simulations (and is approximately consistent with three independent numerical studies). Galaxies formed in halos with V_circ(z_sat) ≲ V_crit now fill the gap that was present in Model 2, producing a luminosity distribution that rises continuously from M_V = −14 down to M_V = −2, before radial selection effects finally cut it off. With these parameter choices, pre-reionization dwarfs dominate the counts (and exceed the observations) for M_V ⩽ −4, but suppressed post-reionization dwarfs dominate the counts at all brighter magnitudes.

Since Model 3 is both more physically realistic and more empirically successful than Models 1 and 2, we focus on it for the remainder of the paper, including Model 3B in which pre-reionization star formation is suppressed below a circular velocity threshold. Figure 8 systematically explores the impact of parameter variations in Models 3A and 3B. In the first three panels, the green band shows the Model 3A predictions for a fiducial set of parameter choices, F_* = 10⁻³, V_crit = 35 km s⁻¹, and z_rei = 11. Changing F_* (top left) shifts the predicted distribution horizontally to higher or lower luminosities, with some change in shape at the faint end because of the luminosity dependence of radial selection effects. Changing V_crit alters the predicted counts at intermediate luminosities, −4>M_V> − 11, while having little effect at the faint end (where pre-reionization dwarfs dominate) or at the bright end (where most galaxies exceed the highest threshold considered here). Changing z_rei alters the height of the pre-reionization peak at faint luminosities but has minimal impact for M_V < −7.

**Figure 8.** Predicted N_obs(*M_V*) for Models 3A and 3B with a variety of parameter choices, in the same format as Figure 5. In the first three panels, green bands show Model 3A predictions for a reference parameter set F_* = 10⁻³, V_crit = 35 km s⁻¹, z_rei = 11. Red and blue curves show the impact of changing the stellar mass fraction to F_* = 10⁻² or 10⁻⁴ (top left), the critical velocity threshold to V_crit = 45 or 25 km s⁻¹ (top right), or the reionization redshift to z_rei = 8 or 14 (lower left). The lower right panel compares the prediction of this reference model (now shown by the red curve) to predictions of Model 3B with a pre-reionization critical threshold V_crit,r = 6 km s⁻¹ (green band) or 10 km s⁻¹ (blue curve).
Download figure:
Standard image High-resolution image

With our fiducial parameter choices, Model 3A substantially overpredicts the number of satellites with M_V ≈ −3. Raising the reionization redshift to z_rei = 14 erases this discrepancy, but this value of z_rei seems implausible given the strong and rapidly evolving opacity of the intergalactic medium at z ≈ 6 seen in quasar spectra (Fan et al. 2006), and it is only marginally consistent with the Wilkinson Microwave Anisotropy Probe 5 (WMAP5) results. In the lower right panel, we return to z_rei = 11 but suppress pre-reionization star formation in halos with V_circ(z_rei) < 6 km s⁻¹ (green) or 10 km s⁻¹ (blue), motivated by the inefficient gas cooling expected below the threshold for atomic line excitation (Model 3B). The V_crit,r = 10 km s⁻¹ model yields acceptable agreement with the observed number counts over the full range 0 ⩾ M_V ⩾ −15. The V_crit,r = 6 km s⁻¹ model still yields an excess of faint satellites; results for V_crit,r = 8 km s⁻¹ (not shown) are nearly identical to those for 10 km s⁻¹, indicating that an 8 km s⁻¹ threshold is already sufficient to essentially eliminate the contribution of pre-reionization dwarfs. This pre-reionization suppression appears to be critical to explain the number of dwarfs observed by the SDSS.

Within Model 3B, there is strong degeneracy between the values of F_* and V_crit. Figure 9 shows that the parameter combinations (F_*, V_crit) = (3 × 10⁻³, 45 km s⁻¹), (10⁻³, 35 km s⁻¹), and (3 × 10⁻⁴, 25 km s⁻¹) all yield similar predictions and acceptable agreement with the observed number counts. The lower values of V_crit are favored by the numerical studies of Hoeft et al. (2006) and Okamoto et al. (2008). For the remainder of the paper we will adopt (F_*, V_crit, z_rei, V_crit,r) = (10⁻³, 35 km s⁻¹, 11, 10 km s⁻¹) as the fiducial parameter values for Model 3B.

For this fiducial model, Figure 10 illustrates in more detail the relative importance of stars formed before and after reionization. For systems with V_circ(z_rei)>V_crit,r, filled circles show the fraction of their stars that formed before reionization. For systems with V_circ(z_rei) < V_crit,r, open circles show the fraction of stars that would have formed before reionization, but because of the V_crit,r threshold these galaxies have no pre-reionization stars in this model. At every satellite luminosity, the average fraction of pre-reionization stars is small, or even zero, but albeit for different reasons at high and low luminosities. The host halos for the brighter, "classical" dwarf satellites were typically massive enough at z_rei to exceed V_crit, but that initial population of stars was subsequently swamped by the much larger post-reionization population. In contrast, the halos that now host the very faintest known satellites (M_V> − 4) did not exceed V_crit,r at z_rei and hence—in Model 3B—did not form any stars before z_rei. A small fraction of the satellites with M_V ≈ −5 have large populations of pre-reionization stars; these are subhalos that just exceeded V_crit,r at z_rei but have low enough values of V_circ(z_sat) that their post-reionization star formation was strongly suppressed. If the pre-reionization threshold at V_crit,r were smooth rather than sharp, then some additional fainter systems might have significant fractions of pre-reionization stars. However, the general conclusion that pre-reionization star formation should be a small fractional contribution at all satellite luminosities seems fairly robust, provided this star formation is suppressed in halos below the atomic cooling threshold, as seems to be required to match the observed luminosity distribution.

Figure 11 shows the complete stellar luminosity function of MW satellites inside 400 kpc, in the absence of any selection effects or incompleteness, again for the fiducial model (we choose 400 kpc for ease of comparison to Tollerud et al. 2008). In contrast to other figures, it shows the luminosity function for the whole sky (4π sr) and in terms of dN/dM_V (i.e., in bins of 1 mag). In the absence of selection effects, the luminosity function continues to rise toward faint magnitudes (as noted by Koposov et al. 2008), contrary to the almost flat luminosity distribution of observed dwarfs. The total number of satellites within 400 kpc brighter than M_V = 0 expected for the fiducial Model 3B is 230 ± 35. This value is somewhat lower than the 400 derived by Tollerud et al. (2008), but since both estimates extrapolate the number of known dwarfs by a factor of ∼10, we do not place much weight on this difference.

None of the models shown in Figures 5–9 reproduce the brightest observed bin—i.e., they all fail to produce satellites as bright as the SMC and the LMC. Our successful models have low stellar mass fractions, F_* ∼ 10⁻³, even well above the photoionization threshold V_crit. The most massive subhalos in our Monte Carlo realizations have typical mass M_sat ∼ 10¹¹ M_☉ (ranging from 10^10.5 M_☉ to 10^11.4 M_☉), with second-ranked halos that are 0.2–0.4 dex less massive. Reproducing the ∼10⁹ M_☉ stellar masses of the Magellanic Clouds then requires much higher stellar mass fractions F_* ∼ 0.05. To reproduce the full satellite population, the efficiency of gas accretion and star formation must continue to rise with halo mass above V_crit, or at least it must be higher for the SMC and LMC hosts. Since the number of bright SMC and LMC-like objects in our model is determined mainly by one parameter F_* (because these objects are not suppressed by the photoionization), that rise of star formation efficiency cannot be accommodated with our simple model without introducing additional parameters.

4.3. Velocity Dispersions, Central Masses, and Radial Distributions

As discussed in Section 3.2, predicting stellar velocity dispersions requires assumptions beyond those needed to compute N_obs(M_V). In particular, we assume that the satellites' host subhalos have NFW profiles with concentration given by the theoretically expected mean c(M) relation at z_sat, and that subsequent dynamical evolution (e.g., tidal stripping) does not alter the mass distribution of the inner parts of the subhalo probed by the stars. We also take the observed stellar radii (20–150 pc; see Section 3.2 for details) as input rather than predicting them from a physical model. With these assumptions, the right panel of Figure 12 shows the predicted distribution of stellar velocity dispersions for Model 3B with our fiducial parameter choices. The characteristic value and narrow spread of velocity dispersions for the newly discovered SDSS dwarfs arise quite naturally from these models, despite the large range of stellar luminosities and host subhalo masses. The predicted distribution is more sharply peaked than the observed one, probably because we did not include scatter in the halo concentration–mass relation and did not include observational uncertainties in the dispersion measurements. The mean value of σ_* differs by <20% between data and model, but we consider this small discrepancy is not worrisome, given the simplicity of our dynamical modeling.

**Figure 12.** Predictions of Model 3B with the fiducial parameters (F_*, V_crit, z_rei, V_crit,r) = (10⁻³, 35 km s⁻¹, 11, 10 km s⁻¹) compared to the observed distributions of absolute magnitude (left) and stellar velocity dispersions (right). The format of the left panel is the same as Figure 5. The right panel shows predicted and observed velocity dispersions only for the SDSS dwarfs—i.e., those with *M_V*> − 11—with data taken from Simon & Geha (2007).
Download figure:
Standard image High-resolution image

The total masses of dwarf satellites are difficult to determine observationally because of the small extent of the stellar distributions relative to the expected extent of the DM subhalo. However, Strigari et al. (2008) show that the total mass (principally DM) within a radius of 300 pc, M₃₀₀, can be inferred robustly from observations for nearly all of the known satellites. The top panel of Figure 13 compares the fiducial model predictions of M₃₀₀ to the Strigari et al. (2008) measurements. The model (red diamonds) naturally reproduces the key result of Strigari et al. (2008): over an enormous range of luminosities, the satellites have a narrow range of M₃₀₀, tightly concentrated around 10⁷ M_☉. The theoretical prediction is artificially tight because we have not included scatter in halo concentrations, which would produce roughly 0.15 dex (rms) of scatter in M₃₀₀ (see Macciò et al. 2009, Figure 1). The model predicts a weak trend of M₃₀₀ with luminosity, which is not evident in the data (but is similar to that predicted by Macciò et al. 2009).

**Figure 13.** Masses of the DM subhalos within the central 300 pc (top), their total present-day masses (middle) and their masses at the time of accretion into larger halos (bottom). We only show halos hosting observable satellites within the MW virial radius, as a function of satellite luminosity. Red diamonds show all the observable galaxies from six realizations of the fiducial Model 3B with (F_*, V_crit, z_rei, V_crit,r) = (10⁻³, 35 km s⁻¹, 11, 10 km s⁻¹). Blue filled circles show the predictions of Model 3A, which includes pre-reionization dwarfs (or, equivalently, has V_crit,r = 0). Error bars show the estimates of M₃₀₀ for observed MW satellites from (Strigari et al. 2008). The solid lines in the bottom panel show, from top to bottom, M_*/M_sat = 10⁻⁵, 10⁻⁴, and 10⁻³. Our models do not incorporate scatter in the concentration–mass relation; adding the theoretically expected scatter would add roughly 0.15 dex of rms scatter to the M₃₀₀ predictions.
Download figure:
Standard image High-resolution image

While the M₃₀₀ range of the satellites is low, the range of total subhalo masses (at z = 0) is more than three orders of magnitude, as shown in the middle panel of Figure 13. The trend of total mass with luminosity is much stronger than the trend for M₃₀₀, though there is a large scatter in mass at fixed luminosity because of tidal stripping. The near constancy of M₃₀₀ is a consequence of the density profiles of CDM halos: NFW halos with the theoretically predicted c(M) relation have only a weak dependence of M₃₀₀ on total mass over the range ∼10⁷–10¹⁰ M_☉ that hosts observed MW satellites (see Macciò et al. 2009 for further discussion). Thus, our models and the models of Macciò et al. (2009) are able to reproduce the narrow observed range of M₃₀₀ without much difficulty (see also Li et al. 2008, who examine M₆₀₀ rather than M₃₀₀). We note, however, that if we also allow satellites to form stars with efficiency F_* = 10⁻³ before reionization (Model 3A), then the M₃₀₀ range for the lowest luminosity dwarfs, with M_V> − 3, extends downward to M₃₀₀ ∼ 10^6.5 M_☉ (blue circles in Figure 13). Thus, careful dynamical measurements for the faintest dwarfs could in principle distinguish whether they arise mainly from pre-reionization star formation or from highly suppressed post-reionization star formation in more massive halos. It is noticeable that our model as well as the models of Macciò et al. (2009) and Li et al. (2008) predicts that M₃₀₀ or M₆₀₀ should slightly increase with galaxy luminosity contradicting the observations, where there is no correlation at all of M₃₀₀ versus luminosity (Strigari et al. 2008). The reason of this disagreement is yet to be understood. It can be caused either by some problems with the data (selection effects or systematics in M₃₀₀ measurements) or by some astrophysical effects. For example, Macciò et al. (2009) eliminate the correlation of M₃₀₀ versus luminosity by assuming that the inner profile of the halos with low concentration (i.e., massive halos) is modified during the process of tidal stripping (Kazantzidis et al. 2004).

Comparing the middle and upper panels shows that a small number of objects have M(z = 0) lower than M₃₀₀, which is possible because we calculate M₃₀₀ based on the subhalo profile at accretion. The tidal radii of these systems are <300 pc, but they are all faint satellites for which the stellar Plummer radii are small. While their true M₃₀₀ values should be M(z = 0), the values calculated in the upper panel are probably more directly comparable to the quantities estimated by Strigari et al. (2008), who extrapolate to 300 pc for the faintest systems assuming that they are not tidally truncated within this radius. To minimize the tidal effects one may also compute the masses within 100 pc instead of 300 pc. For our simulated galaxies we also derive M₁₀₀, which are in the range 1 × 10⁶–4 × 10⁶ M_☉ and are also consistent with the M₁₀₀ ≈ 1 × 10⁶–3 × 10⁶ M_☉ measurements from Strigari et al. (2008, supplementary information).

The bottom panel of Figure 13 shows the value of M_sat as a function of luminosity. The relation obviously reflects the underlying formula used to assign stellar masses to the DM halos (Equation (6)), and the scatter caused by the range of accretion redshifts (which affects the M_sat − V_circ mapping) is small. Even the faintest observable dwarfs have M_sat ∼ 10^8.5 M_☉, but they have star formation efficiencies of only ∼10⁻⁵. The difference between the middle and bottom panels illustrates the effect of tidal stripping. Nearly all the spread of M(z = 0) at fixed M_V comes from different degree of tidal stripping.

Figure 14 compares the distribution of heliocentric distances of the MW satellites found in the SDSS to the predicted distribution for M_V> − 11 satellites from our fiducial model. We show one distribution for each of the six Monte Carlo halo realizations. There are significant halo-to-halo variations in the predicted distributions, and the observed distribution follows the lower envelope of the predictions. The distance distribution is strongly influenced by the radial selection effects (the model predictions would be very different if we did not include them), but it also depends on the radial profile of subhalos and the dependence of this profile on M_sat and z_sat, so matching the observed distribution is a significant additional success of the model.

5. CONCLUSIONS

The satellite discoveries in the SDSS (Willman et al. 2005; Belokurov et al. 2006, 2007; Zucker et al. 2006; Irwin et al. 2007; Koposov et al. 2007; Walsh et al. 2007) have transformed our understanding of the MW's dwarf satellite population, extending the luminosity range by two orders of magnitude and the implied number of systems by a factor of 20. Careful quantification of the SDSS satellite detection efficiency (Koposov et al. 2008; Walsh et al. 2009) allows models that specify the relation between DM subhalos and their stellar content to be tested quantitatively against the observations. We have shown that CDM-based models incorporating previously advocated, physically plausible mechanisms for suppressing the stellar content of low-mass halos can reproduce the observed properties of the known satellite population, including their numbers, luminosity distribution, stellar velocity dispersions, central masses, and heliocentric radius distribution. However, parameters of these models are tightly constrained, and alternative assumptions lead to conflict with the data. In summarizing our results, it is useful to review both what works and what does not.

What works is a model in which the photoionizing background suppresses gas accretion onto halos with V_circ(z_sat) < V_crit ≈ 35 km s⁻¹ (Quinn et al. 1996; Thoul & Weinberg 1996; Bullock et al. 2000), with the smooth mass-dependent suppression suggested by numerical simulations (Equation (6); Gnedin 2000; Hoeft et al. 2006; Okamoto et al. 2008), and inefficient molecular cooling (and/or stellar feedback) drastically reduces the efficiency of star formation in pre-reionization halos below the hydrogen atomic line cooling threshold V_crit,r ≈ 10 km s⁻¹ (Haiman et al. 1997; Barkana & Loeb 1999; Machacek et al. 2001; Wise & Abel 2007; Bovill & Ricotti 2009; O'Shea & Norman 2008). There is some degeneracy between this model's two main parameters, V_crit and F_*, as shown in Figure 9, but with either parameter fixed the other is fairly well constrained (Figure 8). The other two parameters, z_rei and V_crit,r, just need to be in a range that keeps pre-reionization star formation too low to affect the observable luminosity function. For the values V_crit = 25–35 km s⁻¹ favored by numerical simulations, F_* must be ≲10⁻³, so even subhalos above the V_crit threshold have star formation efficiency far lower than the values F_* ≈ 0.1–0.4 found for bright galaxies (e.g., Pizagno et al. 2005; Mandelbaum et al. 2006; Dutton et al. 2007; Gnedin et al. 2007; Xue et al. 2008).

If we assign stellar extents based on observations, and make the reasonable dynamical assumptions discussed in Section 3.2, then our fiducial model naturally explains the characteristic value and narrow spread of stellar velocity dispersions found for SDSS dwarfs by Simon & Geha (2007). It also explains the characteristic value and narrow range of M₃₀₀ values found by Strigari et al. (2008). The M₃₀₀ values do not depend on the assumed stellar extent, and their narrow range arises from the theoretically predicted structure of CDM halos, which have a weak dependence of M₃₀₀ on total halo mass over the range M_halo ∼ 10⁸–10¹¹ M_☉. Thus any CDM-based model that prevents formation of observable dwarfs in halos below ∼10⁷ M_☉ should qualitatively reproduce the Strigari et al. (2007, 2008) results (e.g., Li et al. 2008; Macciò et al. 2009). Tempering this success, however, is the fact that the total z = 0 subhalo masses in our model span three orders of magnitude; some of this range is a consequence of tidal stripping, but the span of M_sat values is only slightly narrower. The model, in combination with the radial selection biases found by Koposov et al. (2008), also explains the observed heliocentric radius distribution of the SDSS dwarfs, which tests the predicted Galactocentric radius distribution of subhalos and its dependence on mass and accretion redshift.

Many alternative models fail badly in reproducing the observed luminosity distribution. Models with constant M_*/M_sat predict far too many faint satellites relative to bright satellites. The SDSS discoveries and luminosity-dependent selection biases do not in themselves resolve the "missing satellite" discrepancy highlighted by Klypin et al. (1999) and Moore et al. (1999); strong mass-dependent suppression of star formation efficiency is still required to reconcile CDM predictions with observations. A simple model in which M_*/M_sat = 10⁻³(Ω_b/Ω_m)(M_sat/10¹⁰ M_☉)² for M_sat < 10¹⁰ M_☉ is reasonably successful at matching the observations. This successful "empirical" model has a mass dependence of star formation efficiency roughly like that of the successful, physically motivated photoionization model (Equation (5); note that M_sat ∝ V³_circ at fixed z_sat).

Models with sharp suppression of star formation below the photoionization threshold V_crit fail at intermediate luminosities, M_V ∼ −8. Pre-reionization star formation can provide the population of faint dwarfs in such a model, but there is an unacceptable gap between the faint and bright populations (or, for parameter choices that fill the gap, there is an excess of dwarfs at other luminosities). It is striking, therefore, that the form of the mass-dependent photoionization suppression found in numerical simulations is just that required to match the shape of the observed luminosity distribution. However, the conversion of accreted baryons to stars must be very inefficient for our fiducial model to work, and it is not obvious why this conversion efficiency should be mass independent.

The most interesting of our "negative" conclusions is that star formation in halos before reionization must be extremely inefficient to avoid producing too many satellites in the range 0 ≳ M_V ≳ −6. Examination of Figure 8 suggests that the upper limit on the fraction of halo baryons converted to stars is a few ×10⁻⁴ for z_rei = 11, or 10⁻³ if reionization is pushed back to z_rei = 14. Madau et al. (2008) have reached exactly the same conclusion, with a similar numerical value for the efficiency limit, using the Via Lactea II simulation instead of a semianalytic method to predict the model subhalo population. Suppression of star formation in halos below the hydrogen atomic line cooling threshold is physically plausible, as the metallicity is low and molecular cooling should be inefficient. For agreement with N_obs(M_V), we require pre-reionization suppression in halos with V_circ(z_rei) < V_crit,r ≈ 10 km s⁻¹.

There are several caveats to these conclusions. First, as discussed in Section 4.2, reproducing the Magellanic Clouds requires that the most massive subhalos have M_*/M_sat ∼ 0.05–0.1, well above the F_* ∼ 10⁻³ of our fiducial model. Thus, the photoionization suppression described by Equation (5) must join onto a continuing increase of star formation efficiency with subhalo mass above V_crit, an increase that is presumably driven by other physical mechanisms. Indeed, there is nothing about our results that necessarily picks out photoionization as the suppression mechanism in low-mass subhalos, but it is a mechanism that comes in naturally (one might argue inevitably) at the desired scale (Bullock et al. 2000), and the numerically calibrated form yields a good match to the observed luminosity distribution.

In our fiducial model, even the faintest SDSS dwarfs form most of their stars after reionization, but they have V_circ(z_sat) far enough below V_crit that their star formation is highly suppressed according to Equation (6). The SDSS dwarfs are physically a continuum with the classical dwarfs, and their much lower luminosities are a consequence of the highly nonlinear relation between star formation efficiency and halo mass below V_crit. Halos with V_circ(z_rei)>V_crit,r form pre-reionization stars, but in nearly all cases they grow large enough by z_sat that the post-reionization population dominates by a large factor. A small number of systems with M_V ≈ −5 could have large fractions of pre-reionization stars, but at any luminosity such systems are rare. These conclusions are robust within our framework, but if we allowed for departures from our adopted prescriptions—in particular if photoionization suppression for V_circ ≪ V_crit were more aggressive than Equation (6) implies and pre-reionization suppression weaker than we have assumed—then it might be possible to construct models in which many dwarfs with M_V ≳ −6 are pre-reionization "fossils." The efficiency of converting halo baryons to stars in these systems must still be ∼10⁻⁴ or less to avoid producing too many faint satellites. Bovill & Ricotti (2009) and Salvadori & Ferrara (2009) have argued that halos cooling by H₂ before reionization naturally give rise to the physical and chemical properties of the SDSS dwarfs. However, even the low star formation efficiencies ∼0.5%–2% found by Salvadori & Ferrara (2009) appear far too high to be consistent with the observed number counts. On the other hand, Busha et al. (2009) propose a model in which post-reionization suppression of star formation is highly efficient (a sharp threshold) but the star formation efficiency in pre-reionization halos is strongly mass dependent, effectively spreading the low-luminosity peak evident in our Figure 6 up toward higher luminosities so that it fills out the entire faint end of the luminosity function.

A third caveat is that we do not explain the origin of the observed stellar extents; we just show that once the observed extents are adopted as inputs, then the observed stellar velocity dispersions emerge naturally. One possible explanation is that the baryons in low-mass halos condense until they reach a scale at which the velocity dispersion is a few km s⁻¹, and that this minimum dispersion provides the conditions necessary for star formation. We also have not attempted to explain the chemical abundance distributions or star formation histories of the satellites (see, e.g., Orban et al. 2008; Salvadori et al. 2008; Salvadori & Ferrara 2009).

A final caveat is that we have assumed that all dwarfs luminous enough to be found in the SDSS also lie above the surface brightness threshold for detection, which is about 30 mag arcsec⁻² (Koposov et al. 2008). Since some of the known satellites approach this threshold, it is possible that others fall below it. A large population of lower surface brightness dwarfs would change the number counts that our model reproduces. Note also that a large population of pre-reionization dwarfs would be observationally allowed if they lie below the surface brightness threshold; however, even in this scenario the pre-reionization dwarfs do not account for the presently known satellites. Deeper large area imaging surveys, such as PanSTARRS, the Dark Energy Survey, and Large Synoptic Survey Telescope (LSST), will show whether the MW satellite population includes a significant number of lower surface brightness systems.

Our model makes several predictions that can be tested by these upcoming surveys or by further follow-up studies of known dwarfs. Deeper surveys should reveal many more satellites, more than 200 with M_V < 0 and D_☉ < 400 kpc over the full sky, with the luminosity function shown in Figure 11. Deep imaging of Andromeda and other nearby galaxies can show whether they have similar satellite systems, though these searches will not reach the extremely low luminosities that can be probed in the MW. Most satellites in our model have stellar extents that are substantially smaller than the present-day tidal radius of their host halo. Tidal tails and tidal disruption should be rare, an implication that may be challenged by photometric evidence on the profiles and shapes of the ultra-faint galaxies, which have been interpreted as signs of tidal distortion or disruption (e.g., Martin et al. 2008). Measurements of the total subhalo masses of known dwarfs would provide a powerful test of the model predictions in Figure 13, but the small stellar extents may make such measurements impossible. Our models predict that satellites continue to form stars down to z_sat or below, and many observable systems should have z_sat = 1–2 (see Figure 2). These predictions may be testable with detailed stellar population modeling.

Our results greatly strengthen the argument (Bullock et al. 2000; Benson et al. 2002; Somerville 2002; Kravtsov et al. 2004) that photoionization naturally reconciles the CDM-predicted subhalo population with the observed dwarf spheroidal population, thus solving the "missing satellite problem" highlighted by Klypin et al. (1999) and Moore et al. (1999). The fiducial model presented here offers a detailed, quantitative resolution of this problem in light of new, greatly improved observational constraints, while relying on previously postulated and physically reasonable mechanisms to suppress star formation in low-mass halos. The MW satellites provide a fabulous laboratory for studying galaxy formation at the lowest mass scales, and much remains to be understood about gas cooling, star formation, feedback, and chemical enrichment in these systems. These issues provide challenging targets for numerical simulations and semianalytic models, whose predictions can be tested against detailed studies of the dynamics and stellar populations of the known dwarf satellites and of the many new satellites that will be revealed by the next generation of sky surveys.

S.K. was supported by the DFG through SFB 439 and by a EARA-EST Marie Curie Visiting fellowship. J.Y. is supported by the Harvard College Observatory under the Donald H. Menzel fund. D.W. acknowledges support from NSF grant AST-0707985 and the hospitality of the Institut d'Astrophysique de Paris during part of this work. S.K. acknowledges hospitality from the Kavli Institute for Theoretical Physics (KITP) Santa Barbara during the workshop "Building the Milky Way." We thank James Bullock for his helpful comments on the paper and the anonymous referee for prompt review and constructive comments.

This paper relies heavily on data from the SDSS. Funding for the SDSS and SDSS-II was provided by the Alfred P. Sloan Foundation, the Participating Institutions, the National Science Foundation, the U.S. Department of Energy, the National Aeronautics and Space Administration, the Japanese Monbukagakusho, the Max Planck Society, and the Higher Education Funding Council for England. The SDSS was managed by the Astrophysical Research Consortium for the Participating Institutions, which are listed at the SDSS Web site, http://www.sdss.org/.

A QUANTITATIVE EXPLANATION OF THE OBSERVED POPULATION OF MILKY WAY SATELLITE GALAXIES

Article metrics

Permissions

Author e-mails

Author affiliations

Dates

ABSTRACT

1. INTRODUCTION

2. THE POPULATION OF DM SUBHALOS IN THE MILKY WAY

3. POPULATING THE DM HALOS WITH STARS

3.1. Recipes to Assign Stellar Masses to Subhalos

3.2. Detectability and Observable Properties for the Simulated Satellites

4. RESULTS

4.1. Stellar Mass Function of the Full Satellite Populations

4.2. Distribution of Observed Dwarf Satellite Luminosities, N(M_V)

4.3. Velocity Dispersions, Central Masses, and Radial Distributions

5. CONCLUSIONS

Footnotes

A QUANTITATIVE EXPLANATION OF THE OBSERVED POPULATION OF MILKY WAY SATELLITE GALAXIES

Article metrics

Permissions

Share this article

Author e-mails

Author affiliations

Dates

ABSTRACT

1. INTRODUCTION

2. THE POPULATION OF DM SUBHALOS IN THE MILKY WAY

3. POPULATING THE DM HALOS WITH STARS

3.1. Recipes to Assign Stellar Masses to Subhalos

3.2. Detectability and Observable Properties for the Simulated Satellites

4. RESULTS

4.1. Stellar Mass Function of the Full Satellite Populations

4.2. Distribution of Observed Dwarf Satellite Luminosities, N(MV)

4.3. Velocity Dispersions, Central Masses, and Radial Distributions

5. CONCLUSIONS

Footnotes

4.2. Distribution of Observed Dwarf Satellite Luminosities, N(M_V)