ON THE ASSEMBLY OF THE MILKY WAY DWARF SATELLITES AND THEIR COMMON MASS SCALE

A total of 1.6 million tagged particles end up in the main host at redshift z = 0 (0.35% of all particles within r₂₀₀): most are tidally stripped from their subhalos after infall, and only 52% of them remain bound to the 1925 self-gravitating clumps that partially survive the accretion event. The total mass of the stellar halo today is ${\cal M}_*=1.6 \times 10^9\,{M_\odot }$ within 402 kpc (of which 2.3 × 10⁸ M_☉ is in bound substructures), and Table 1 lists the contribution to ${\cal M}_*$ from each snapshot: the peak contribution comes from redshift 2. The mass in the stripped component between 1 kpc and 40 kpc is 5 × 10⁸ M_☉, comparable to recent estimates for the MW halo mass within the same range (Bell et al. 2008).

In terms of its mass assembly history, the stellar component of the VLII halo is already in place by redshift z = 1.5, while the dark component continues to grow from mergers with smaller subunits down to the present epoch. The faster assembly of the stellar component follows from the very steep mass-star formation efficiency relation adopted in our model. A similar behavior was seen in the assembly history of the six Aquarius stellar halos by Cooper et al. (2010). About 50% of today's stellar halo can be traced back to a subhalo of M_* = 8.8 × 10⁸ M_☉ that fell in at z = 2. The three next most massive contributors, with M_*/ M_☉ = (4.2 × 10⁸, 1.2 × 10⁸, 4.5 × 10⁷), were accreted at redshifts z = 7.8, 3.2, and 1.6, respectively. Combined, these four progenitor systems account for about 88% of the total stellar halo.

Figure 2 shows an image of the projected VLII dark matter and stellar halo mass density at three different redshifts. The disruption of the stellar component of an individual satellite is depicted in the right panels, where particles initially tagged in one subhalo at infall (the system reaches its M_h at z = 2) are violently stripped by tides and form shell- and stream-like structures. The bottom left and middle panels make it clear how the adopted steep stellar-mass–halo-mass relation offers an empirical solution to the "missing satellite problem." While the dark halo exhibits a vast amount of surviving self-bound substructure, only the most massive subhalos form enough stellar mass at infall to still be visible in the stellar halo today.

Zoom In Zoom Out Reset image size

Surviving satellites in the z = 0 VLII snapshot were identified by considering all particles that fall within the tidal radius of a bound subclump (identified in the dark matter density field by the 6DFOF group finder; see Diemand et al. 2006). We then performed a peculiar velocity cut of V_p ⩽ 2V_max in order to remove particles belonging to the main halo that just happen to be flying through the given satellite at that time. The total luminosity of each satellite was computed by summing up the contributions of all stellar particles linked to it and its metallicity determined from their mass-weighted average. Satellites of absolute magnitude M_V located closer than the SDSS completeness limit r_max (Tollerud et al. 2008) from the main halo center,

$$\begin{equation} r_{\rm max} = \left(\frac{3}{4\pi f_{\rm DR5}}\right)^{1/3} \times 10^{(-0.6M_V-5.23)/3}\;{\rm Mpc}, \end{equation} \tag{ 3 }$$

were flagged as observable in current surveys (f_DR5 = 0.194 here is the sky covering fraction of the SDSS DR5). Only 65 of the 1925 surviving satellites present today in the VLII stellar halo are bright enough to be detected by an all-sky SDSS-like survey, out to a galactocentric distance of 280 kpc, corresponding to the MW virial radius (Xue et al. 2008). Note that, as the SDSS detection efficiency is nearly independent of central surface brightness for μ_{0, V} ≲ 30 mag arcsec⁻² (Koposov et al. 2008), no surface brightness threshold has been applied to our sample. We have checked that all but 13 of the observable VLII satellites have a central surface brightness of at least 30 mag arcsec⁻². The 13 low surface brightness systems are all faint (M_V ≳ −5) and have half-light radii r_eff close to the gravitational softening scale of the simulation, i.e., they have been artificially puffed up. We have estimated that, without such a numerical effect, these subhalos would actually lie above the SDSS surface brightness limit. Our prescriptions do not therefore produce a large population of low surface brightness "stealth" satellites awaiting discovery in the MW halo, as suggested by Bullock et al. (2010).

The stellar masses today of the observable VLII satellites range from 2000 M_☉ to 10⁸ M_☉, and their dark masses from 4 × 10⁶ M_☉ to 3 × 10⁹ M_☉. They have mass-to-light ratios that range from about 20 to 10⁵, and infall redshifts mostly between 1.5 and 8 (with a few falling in later, between redshifts 0.6 and 1.1). Their mass-to-light ratios have decreased by a factor ranging between a few and 20, as their extended dark matter halos have been tidally stripped more than their compact luminous components. The 10 most luminous satellites (L > 10⁶ L_☉) in the simulation are hosted by subhalos with peak circular velocities today in the range V_max = 10–40  km s⁻¹ that have shed between 80% and 99% of their dark mass after being accreted at redshifts 1.7 < z < 4.6. Their V_max values have decreased by a factor of 2.5 on the average (and up to 6 in individual cases) after accretion, so their present-day V_max is not a good indicator of the subhalo maximum mass.

The differential satellite luminosity function is shown in the left panel of Figure 3, plotted against the luminosity function of observed MW dwarfs (see also works by Koposov et al. 2009; Macciò et al. 2010). As in Koposov et al. (2009), we use different scales for the number of satellites brighter and fainter than M_V = −11, introduced to produce a continuous luminosity function despite the limited SDSS DR5 sky coverage of f_DR5 = 0.194. The shaded region to the left of M_V = −11 spans the 1σ extent in the number of faint, but detectable satellites that fall within the SDSS DR5 footprint in each of 2000 randomly oriented sky projections of the VLII stellar halo. In the case that the satellites are isotropically distributed, our exercise should produce a result equivalent to a scaling of the number counts of observed faint satellites by the sky coverage fraction of the SDSS survey. Satellites to the right of the dashed line, corresponding to M_V = −11, are always detectable. The shaded region to the right of M_V = −11 covers the 1σ Poisson errors in each bin there. The lack of Magellanic Cloud equivalents in VLII is responsible for the discrepant bright end of the luminosity function, as the simulation does not contain satellites more luminous than M_V ≲ −14. Studies suggest that satellites this bright are rare around MW-like galaxies (Boylan-Kolchin et al. 2010; Liu et al. 2011). In the right panel of the same figure, we plot the luminosity–metallicity relation of surviving satellites compared to data for MW dSphs (Kirby et al. 2008, and references therein) and their best-fit model to the data. The relation assumed at infall (Equation (2)) has been modified by tidal mass losses and by the aging of stellar populations.

Zoom In Zoom Out Reset image size

In Figure 4, the cumulative galactocentric distance distribution of VLII luminous subhalos is compared to that of MW satellites fainter than M_V ⩾ −12. Even though the final spatial distribution of dwarfs depends on their accretion and orbital history, the two distributions track each other reasonably well, though VLII shows an excess of close-in satellites compared to the observations. The figure also shows (right top panel) the present-day stellar-mass–subhalo-mass relation. Much like the metallicity–luminosity relation in Figure 3, it has been modified from its infall value by tidal mass losses. As accreted systems lose a disproportionate amount of their dark mass, the apparent (inferred today) star formation efficiency of dSphs may exceed the true one by as much as a factor of 100. We also explore (bottom left panel of Figure 4) how the average one-dimensional (1D) line-of-sight velocity dispersion σ_los of dSphs (an observable quantity) compares to the peak (three-dimensional, 3D) circular velocity V_max of the system's host halo (a quantity of theoretical interest). Early attempts at this conversion assumed a proportionality factor of $\sqrt{3}$ motivated by the isotropy of these systems (Moore et al. 1999). Our model suggests that the baryonic tracers are deeply embedded in their dark matter potentials today, resulting in higher values of V_max, proportional to 2.2σ_los.⁴ The difference is even more pronounced when comparing to the highest V_max values these systems ever achieved. Shown in the figure as downward pointing triangles, these values indicate the amount of tidal mass loss that the subhalos hosting dSphs have undergone between infall and today. In Figure 5, we re-iterate this point by plotting the time evolution of V_max as a function of satellites' present-day luminosities.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Finally, we show how VLII dSphs compare to the Tully–Fisher relation (data from Pizagno et al. 2007) for disk galaxies (Figure 4, bottom right panel). Just like the MW dSph population, our star formation prescription fails to produce dwarfs that lie on the extrapolated Tully–Fisher relation. The flattening of the distribution of spheroidal galaxies to a constant velocity value in comparison to field spiral galaxies is reminiscent of the "common mass scale" seen in the M₃₀₀ mass function (see Figure 9), in which all currently detectable dSphs cluster around the value M₃₀₀ ∼ 10⁷  M_☉. As we shall see below, in our model the "common mass scale" is the result of a luminosity bias, and future detection of fainter dSph satellites would break it and reveal lighter dark matter subhalos hosting stellar systems.

We now turn to the structural properties of the surviving VLII dwarfs. For all satellites observable in an all-sky SDSS-like survey, we calculate the circularly averaged surface brightness seen by a galactocentric observer and fit it with a two-dimensional (2D) projected King profile (King 1962),

$$\begin{equation} I_*(R)=k\left[\left(1+\frac{R^2}{r_k^2}\right)^{-1/2}-\left(1+\frac{r_t^2}{r_k^2}\right)^{-1/2}\right]^2, \end{equation} \tag{ 4 }$$

where I_*(R) is the surface brightness at radius R from the satellite center, k is the profile normalization, r_k is the core radius of the satellite, and r_t is its tidal (cutoff) radius. This formula provides a good fit to the surface brightness profiles of MW dSphs (Mateo 1998) and has the advantage of having an analytical de-projected 3D form in spherical geometry. Our use of the King profile as a fit to the projected surface brightness of VLII satellites instead of, say an exponential profile, is motivated by its standard use in published studies of dSph mass estimation. We find that the luminous component of the observable simulated dSph is well fit by a King profile across the full range of observable luminosities. Note that, while the tidal radii in the luminous component do not correspond to the physical truncation scale in the satellites' dark matter extent, the two scales are indeed comparable. On average, the ratio between the tidal radius in the luminous component and the tidal radius in the dark matter extent varies between 0.5 and 1 for satellites more luminous than 10⁶ L_☉, while it is between 0.3 and 0.8 for fainter systems.

A comparison between the surface brightness and line-of-sight velocity dispersion profiles of two dSph satellites (Sculptor and Carina) and individual VLII counterparts is shown in Figure 6. Both surface brightness and kinematics should be compared only at radii larger than twice the VLII softening length of 40 pc, as the artificial "heating" of particles inside this region from numerical resolution effects would give a poor match to the data. In the surface brightness plot, we show the diverse light profiles of all VLII dwarfs that have a luminosity within ±80% of the corresponding MW satellite. Equally diverse is the span of velocity dispersion profiles these dwarfs have, with values ranging between 4 km s⁻¹ and 15 km s⁻¹, pointing to the differences between the potential wells in which these systems reside. The insets in the figure show the evolution of the luminous systems' parameters (King profile core and tidal radii and half-light radii) from "tagging" (infall) until today. The behavior seen here is representative of the other SDSS detectable dwarfs in VLII—the core radii increase with time (though they remain close to the simulation softening length), while the half-light and tidal radii decrease as the satellites lose stars from their outskirts to the main halo. In the velocity dispersion plots, we try to reproduce the "noisy" observed profile by sparsely sampling the simulation data. Using the same number of stellar particles in each corresponding bin as the one in the observed profile (M. Walker 2011, private communication), we note that the scatter in the simulated profile increases significantly.

Zoom In Zoom Out Reset image size

In a systematic analysis, Strigari et al. (2010) have shown that at least one simultaneously structural and kinematic counterpart to each of the five classical dSph with good kinematic and photometric data can be found in the Aquarius simulations. Our ability to obtain simultaneous brightness and kinematic fits for dSphs like Sculptor and Carina can be seen as a confirmation of their results and a success of our model, and is a prerequisite for testing the "common mass scale" of MW satellites. Structural parameters of the VLII satellite population (half-light radii and velocity dispersions measured directly from the particle data, and central surface brightnesses from the King model fits as a function of absolute magnitude) are compared to the observed properties of MW dSphs in Figure 7. The agreement is reasonably good. Reproducing the properties of the luminous component of dSph galaxies is a necessary condition for generating realistic "mock" data to test the mass measurement technique.

Zoom In Zoom Out Reset image size

The Jeans equation is the first moment over velocity of the collisionless Boltzmann equation. It provides a link between the underlying potential of a pressure-supported system and the distribution function of a tracer population in equilibrium with it. In spherical geometry its radial component is

$$\begin{equation} r\frac{d\left(\nu _*\sigma _r^2\right)}{dr} = -\nu _*V_c^2-2\beta \nu _*\sigma _r^2, \end{equation} \tag{ 5 }$$

where ν_*(r) is the luminosity density profile of the tracer stars, σ_r is their radial velocity dispersion, and β = 1 − σ²_t/σ_r² is the velocity anisotropy, i.e., the difference between the radial, σ_r, and the tangential, σ²_t = (σ_θ² + σ²_ϕ)/2, components of the velocity dispersion. Here, V_c(r) is the circular velocity profile of the system, a measure of the underlying gravitational potential Φ,

$$\begin{equation} V_c(r)^2 = G\frac{M({<}r)}{r} = -r\frac{d\Phi }{dr}. \end{equation} \tag{ 6 }$$

Since 3D kinematic data for MW dSph stars are not available, one has to use a 1D projection of the spherical Jeans equation along the line of sight between the observer and the studied system. This introduces an additional geometric effect involving the velocity anisotropy (Binney & Tremaine 2008):

$$\begin{equation} \sigma _{\rm los}^2(R) = \frac{2}{I_*(R)}\int _R^{\infty }\left(1-\beta \frac{R^2}{r^2}\right)\frac{\nu _*\sigma _r^2 r}{\sqrt{r^2 - R^2}}dr, \end{equation} \tag{ 7 }$$

where I_*(R) is the projected surface brightness profile of the tracer population. Of all quantities involved in the projected Jeans equation, σ_los and I_*(R) are readily determined observationally, and ν_*(r) can be obtained through a de-projection of the surface brightness profile, assuming it follows some analytical form—either a spherical King (1962) or a Plummer (1911) profile. Both of these profiles adequately fit the light distribution of MW satellites and the choice of one or the other does not affect the results of our kinematic study (Strigari et al. 2008). The remaining parameters, namely, β(r) and V_c(r), are not constrained observationally and have to be modeled appropriately.

We adopt a generalized NFW profile (Navarro et al. 1996) for the mass distribution of subhalos:

$$\begin{equation} \rho (r) = \frac{\rho _s}{(r/r_s)^{\gamma }(1+r/r_s)^{\delta }}, \end{equation} \tag{ 8 }$$

where ρ_s gives the density normalization, r_s is the scale radius, and γ and δ define the inner and outer slopes of the density falloff. The fiducial NFW form has γ = 1 and δ = 2. To allow for a range of inner density cusps as well as outer tidally modified profiles, we allow these two parameters to vary in the intervals 0.7 ≲ γ ≲ 1.2 and 2 ≲ δ ≲ 3 (Strigari et al. 2007). We have verified that the upper limit on the outer slope δ is high enough to accommodate the profiles of even the most massive VLII subhalos, whose density profiles fall off faster than those of field halos in the tidally stripped outer regions (Diemand et al. 2008).

The velocity anisotropy parameter β is largely unconstrained by current observations. According to its definition, β can vary between −∞ and 1, with β = −∞ when σ_tot = σ_t and β = 1 when σ_tot = σ_r. Available data are often in the form of 1D line-of-sight velocity dispersion, insufficient to constrain σ_r and σ_t individually. These two quantities are well determined only at the very center and very edge of a system, since σ_los(0) = σ_r and σ_los(r_t) = σ_t. In between these two extremes, the mass of a dSph enclosed within a radius close to the system's half-light radius can be calculated such that the uncertainty from the velocity anisotropy term is minimized (Wolf et al. 2010). At any other intermediate radius, the β parameter will more strongly affect the mass measurement and must be estimated using a theoretical model. In principle, β can be a function of radius within a particular system (Strigari et al. 2007). Previous studies have shown that marginalizing over a set of parameters that define a general anisotropy profile does not constrain the mass distribution better than assuming β to be constant with radius (Walker et al. 2009). The latter case, however, greatly simplifies the solution to the spherical Jeans equation (Mamon & Łokas 2005):

$$\begin{equation} \nu _*\sigma _r^2(R) = GR^{-2\beta } \int _R^{\infty }r^{2\beta -2}\nu _*(r)M(r)dr. \end{equation} \tag{ 9 }$$

The standard procedure in mass modeling of pressure-supported systems is to generate a theoretical velocity dispersion profile that fits the observed one, and to marginalize the solution over all free parameters but the density profile normalization. While most model parameters are largely unconstrained by this fitting procedure, the density profile normalization (defined as the total amount of mass enclosed in a fixed aperture radius, typically 300 pc or 600 pc) is very well determined by the data (Strigari et al. 2008; Walker et al. 2009). It is this analysis that yields the so-called common mass scale of MW dSph satellites.

We apply the mass measurement technique above to the simulated data using the following procedure.

• 1.  
  Obtain σ_los(R) and I_*(R) for a satellite of interest.
• 2.  
  Choose the density profile and anisotropy parameters ρ_s, r_s, γ, δ, and β.
• 3.  
  Integrate the Jeans equation (5) to solve for σ_r(r).
• 4.  
  Project σ_r(r) along the line of sight (Equation (7)).
• 5.  
  Compare the solution to the measured σ_los(R).
• 6.  
  Iterate (2)–(5) to find the model that best matches the kinematic data.

Once the best-fit mass density profile is determined, we integrate it to a given aperture radius, either 300 pc or 600 pc from the satellite center to obtain the total enclosed mass, M₃₀₀ or M₆₀₀, and compare it directly to the "true" value in the simulation. To obtain the "mock" observations of VLII dSphs in step 1, we place "hypothetical observers" in the simulation volume and reduce computational time by "observing" each of the 65 detectable satellites only 11 times. One of the "observers" is always in the center of the VLII main host, while the other 10 are randomly distributed at the same distance from the satellite as the galactocentric "observer." In this way we attempt to capture the anisotropy of subhalos' physical shapes and velocity ellipsoids.

Since integrating the Jeans equation requires a boundary condition on σ_r(r_t), which is not known a priori from the data, the line-of-sight projected solution cannot be trusted across the full extent of the surface brightness fit. The observed line-of-sight velocity dispersion profile is therefore calculated in 10 linear distance bins out to 60% of each satellite's tidal radius r_t. The simulated data used to determine each subhalo's kinematic profile do not have an associated observational uncertainty. To make a fair comparison with the derived mass estimates, "observational" errors should be consistent with the ones in the current literature (Walker et al. 2009). This is done by randomly assigning error values σ_i ranging between 1 and 2  km s⁻¹ to the velocity dispersion data when evaluating the goodness of fit in Equation (10) below.

Sample line-of-sight velocity dispersion profiles for two VLII simulated satellites and their assigned "observational" errors are shown in the top left panel of Figure 8. The two satellites have masses of M₃₀₀ = 7.22 × 10⁶ M_☉ and M₃₀₀ = 2.65 × 10⁷ M_☉, respectively. For each satellite, we plot all 11 best fit σ_los(R) profiles. The two sets of data points (in black and gray) correspond to the kinematic profile fits highlighted in darker colors, which were chosen to span the range of M₃₀₀ for which the fit's likelihood falls to 10% of its peak value (see the right top panel of Figure 8). The average M₃₀₀ values inferred from these 11 fits are 6.65 × 10⁶ M_☉ and 2.44 × 10⁷ M_☉, respectively, within 8% of the true values.

Zoom In Zoom Out Reset image size

To get the above estimates, we fit a model to each "observed" kinematic profile using a Monte Carlo Markov Chain (MCMC). A theoretical velocity dispersion profile σ_{los, th}(R) is computed for each parameter set in step 2 and the likelihood of the fit is evaluated as

$$\begin{equation} \mathcal {L} = \prod _{i} \frac{1}{\sqrt{2\pi \sigma _i^2}}\hbox{exp}\left[-\frac{1}{2}\frac{(\sigma _{\rm los}(R_i)-\sigma _{\rm los,th}(R_i))^2}{\sigma _i^2}\right], \end{equation} \tag{ 10 }$$

where σ_los(R_i) is the line-of-sight velocity dispersion measured at radius R_i and σ_i is its associated error. Each MCMC chain is initialized to maximize the likelihood $\mathcal {L}$ on a coarse model grid that spans the full parameter space. It is then left to explore the parameter space following the Metropolis–Hastings algorithm. Two concurrent chains are run on each data set and the MCMC is stopped after 10,000 attempts of each chain. This number is large enough to provide satisfactory fits to the kinematic data (Walker et al. 2009). The 11 resulting likelihood functions for M₃₀₀ of the two satellites above are shown in the right top panel in Figure 8.

Once the best-fit model for each data set has been found, it provides the full mass density profile that can be integrated to any distance from the dSph center, with 300 pc or 600 pc being of particular interest. There are 11 likelihood distributions for the estimates of M₃₀₀ and M₆₀₀ of each VLII satellite, and these can be compared to the true values. The results are shown in the bottom panels of Figure 8, in which we plot the average of the 11 best mass estimates as a function of the true dark matter mass. The nature of the MCMC data allows us to make two estimates for the error in each case. The smaller error bars (in black) in Figure 8 reflect the span of the 11 best-fit values for each satellite. The larger errors (in gray) span the full range in values where the likelihood of the fit has fallen to 10% of its maximum value. For the particular two satellites in the top panels of Figure 8, we thus have the following mass estimates: M₃₀₀ = 6.65^{+0.41 + 4.40} _{− 0.64 − 2.95} × 10⁶M_☉ for the lower-mass satellite and M₃₀₀ = 2.44^{+0.48 + 1.48} _{− 0.30 − 1.07} × 10⁷M_☉ for the higher-mass one, where the first set of errors corresponds to the best-fit value range and the second set of errors corresponds to the 10% likelihood value range.

We infer that the spherical Jeans analysis is a robust method to estimate the central masses of dSph satellite galaxies. The mean best-fit M₃₀₀ values determined through the Jeans analysis are on average within 12% of the true VLII values. The normalized scatter among the 11 best fits about the mean value is ^+0.16 _{− 0.15} on average and is shown as solid lines in the bottom panels of Figure 8. The normalized full range scatter about the mean where the likelihood of the fit is greater than 10% of its peak value is ^+0.78 _{− 0.53} on average (dashed lines). No significant difference to these results arises when just the galactocentric observer's estimates are considered: the best-fit galactocentric estimate of M₃₀₀ across the full mass range is on average within 16% of the true VLII value.

When evaluating M₆₀₀, we exclude systems that have M₆₀₀ ≲ 1.5 × 10⁷ M_☉. It is apparent in the bottom right panel of Figure 8 that these exhibit a larger scatter than the systems with M₆₀₀ ⩾ 1.5 × 10⁷ M_☉. The smallest dSphs host stellar systems smaller than 0.6 kpc in radial extent. The estimate of M₆₀₀ for them is therefore an extrapolation to a region in parameter space where tracers are simply not present. Among the dSphs with kinematic data available all the way to r = 0.6 kpc, the mean best-fit M₆₀₀ value determined with the Jeans analysis is on average within 17% of the true VLII value. The normalized scatter among the 11 best fits about the mean value is ^+0.29 _{− 0.26} on average. The normalized full range scatter about the mean where the likelihood of the fit is greater than 10% of its peak value is ^+1.08 _{− 0.55} on average. The VLII galactocentric observer would report a M₆₀₀ value within 23% of the true VLII value.

The range in σ_los(R) and I_*(R) when measured from different viewpoints is due to the anisotropic velocity tensor of the dSph and any deviation from perfect sphericity of its shape. The mass estimate uncertainties among the best-fit values presented here should therefore be regarded as intrinsic, cosmologically induced ones that will eventually be minimized when 3D kinematic data are collected and better models of the nature of dSph satellites in high-resolution, fully cosmological hydrodynamical simulations are carried out.

In the left panel of Figure 9, we plot the mass enclosed in the inner 300 pc versus total luminosity for all satellites that would be visible according to the SDSS detection limit. Overplotted are data from Strigari et al. (2008), suggesting a "common mass scale." Observable VLII dSphs exhibit a behavior similar to the MW dSph population, though a linear fit suggests a shallow but significant correlation between M₃₀₀ and the total luminosity: M₃₀₀∝L^{0.088 ± 0.024}, where the error on the slope reflects the 1σ deviation from bootstrap re-sampling of the sample. This is steeper than the L^{0.03 ± 0.03} relation derived from observational data (Strigari et al. 2008). Lowering the derived M₃₀₀ values of some of the ultrafaint dSphs, as suggested by Adén et al. (2009) for the case of the Hercules dwarf, might be sufficient to reconcile the observations with the slope predicted from our model.

Zoom In Zoom Out Reset image size

In our model, the nearly constant M₃₀₀ is explained by the combination of the steep adopted stellar-mass–halo-mass relation (Equation (1)) and the SDSS detectability magnitude limit. The remaining ∼2000 VLII surviving subhalos to which we assigned stars (about 200 of which are more luminous than M_V = 0) are just too faint to be detected by an SDSS-like survey. Some of them have M₃₀₀ significantly lower than 10⁷ M_☉ (see the right panel of Figure 9).

Boylan-Kolchin et al. (2011) have recently pointed out that state-of-the-art collisionless simulations of Galactic structure, both VLII and Aquarius, predict a population of massive subhalos that appear to be too dense (i.e., have too high a V_max for a given r_vmax) to host any of the known MW dSph. While a failure of the CDM model is one interpretation of this result, another possibility is that star formation becomes stochastic below some mass, in which case the MW halo would harbor several massive but dark satellites. It is also quite plausible that baryonic feedback processes in the more massive subhalos redistribute some of their dark matter, thereby lowering their central density (see e.g., Governato et al. 2010). For reference we have marked the six VLII subhalos that are too dense according to the Boylan-Kolchin et al. (2011) analysis in the left panel of Figure 9. Of these six, five have a luminosity greater than M_V = −11 and would have been detected as luminous dSph satellites. In our model, it is thus not possible for these overdense subhalos to have evaded detection because of an effectively low star formation efficiency resulting, for example, from their accretion or orbital history. As expected, these six halos have the highest values of M₃₀₀ for their luminosity among all VLII satellites. A reduction in M₃₀₀ by about a factor of two would bring them into agreement with the classical MW dSph sample. Our study of the standard Jeans analysis demonstrates that the "observed" M₃₀₀ values cannot be changed by that amount, i.e., that errors in the mass estimates are not responsible for the discrepancy between the models and the observations. Future, high-resolution hydrodynamical simulations of MW analogs that include an accurate treatment of baryonic processes and their effects on the central densities of luminous satellites will be crucial to the settling of this outstanding issue.