ON THE PERSISTENCE OF TWO SMALL-SCALE PROBLEMS IN ΛCDM

The ΛCDM model of cosmology has long been known to be plagued by numerous "small-scale problems." These include the much larger predicted number of sub-halos compared to the number of observed satellite galaxies (missing satellites problem, Klypin et al. 1999; Moore et al. 1999), the discrepancy between a predicted density peak of dark matter at the centers of halos and observations indicating a constant-density core (core-cusp problem, Dubinski & Carlberg 1991; Walker & Peñarrubia 2011), the apparent over-abundance of backsplash galaxies in the Local Group (Pawlowski & McGaugh 2014a), and many other tensions between ΛCDM predictions and observations (summarized for example by Kroupa 2012; Famaey & McGaugh 2013; Weinberg et al. 2013; Walker & Loeb 2014; McGaugh 2015). Two problems that have recently emerged and are gaining increasing attention are the Too-big-to-fail (TBTF) and the Satellite-planes problems.

The TBTF problem is concerned with the internal dynamics of dwarf galaxies. When comparing the central masses of MW dwarf Spheroidal (dSph) satellites deduced from their kinematics with those of dark matter sub-halos in simulations, Boylan-Kolchin et al. (2011, 2012) found that simulations of MW equivalents each contain ≈10 (actually 5–40) denser sub-halos than those compatible with the most-luminous observed dSphs. This can be related to Figure 2 of Kroupa et al. (2010) where it was shown that the observationally deduced DM halo masses of the MW satellites show a significant overabundance of ${M}_{0.3\;\mathrm{kpc}}\simeq {10}^{7}{M}_{\odot }$ halos and a lack of both less and more massive values compared to the theoretically predicted distribution for luminous sub-halos. This would indicate that the most-massive sub-halos do not host the most-luminous dSphs but remain virtually dark, while the most-luminous dSphs live in sub-halos of only intermediate mass. This raises the question of what prevented the massive sub-halos from forming galaxies. The TBTF problem has been identified not only among the MW satellite galaxies, but also for the M31 satellite galaxies (Tollerud et al. 2014), within the Local Group (Garrison-Kimmel et al. 2014a; Kirby et al. 2014) and also appears to be present for field galaxies (Papastergis et al. 2015).

The TBTF problem is more difficult to resolve than the missing satellites problem since it requires either that luminous galaxies do not form in the sub-halos with the largest central dark matter density, or that one introduces a process that reduces the central dark matter density of these most massive sub-halos to values consistent with the observed velocities in dwarf galaxies. One suggested solution to this problem has been a "light" MW. If the MW halo is less massive than generally assumed (at least ≲8 × 10¹¹M_⊙ instead of (1–2) × 10¹²M_⊙) this translates into a lower number of massive sub-halos, thus alleviating the TBTF problem (Boylan-Kolchin et al. 2012; Wang et al. 2012; Vera-Ciro et al. 2013). However, in general, such a low MW mass is disfavored, for example, by the Local Group timing argument (Li & White 2008), the analysis of the positions, line-of-sight velocities and proper motions of the MW satellite population (Watkins et al. 2010), the Galactic escape speed as a function of radius (Piffl et al. 2014), and the modeling of stellar streams in the MW halo (Küpper et al. 2015). Furthermore, that the TBTF problem is also present for the M31 satellites (Tollerud et al. 2014), for more distant dwarfs in the Local Group (Garrison-Kimmel et al. 2014a; Kirby et al. 2014) and possibly even beyond (Papastergis et al. 2015) is a strong argument against purely local or environment-dependent solutions.

Not only the internal properties of dwarf galaxies are problematic for ΛCDM. The overall spatial distribution of satellite galaxies also appears to be at odds with cosmological predictions. Kunkel & Demers (1976) and Lynden-Bell (1976) already noted that the then-known six satellites of the MW are anisotropically distributed and align in a common plane with the Magellanic Stream. But it was Kroupa et al. (2005) who first argued that this flattened distribution of the 11 brightest ("classical") MW satellites is problematic for the cosmological model. Not only are these 11 classical satellites situated in this flattened structure, but also the later discovered faint and ultra-faint satellites define the same structure (Metz et al. 2009; Kroupa et al. 2010). Furthermore, Pawlowski et al. (2012b) have shown that this "Vast Polar Structure" (VPOS) also consists of globular clusters categorized as young halo objects (which are likely associated with satellite galaxies), and that about half of the stellar and gaseous streams in the MW halo align with the satellite plane. The recently discovered MW satellites in the southern hemisphere have also been found to align with the VPOS (Pawlowski et al. 2015). Finally, the potentially most important piece of information is provided by the proper motions measured for the 11 brightest satellites. These reveal that nine of the satellites are consistent with orbiting within the VPOS, eight of these even co-orbit in the same direction (Metz et al. 2009; Pawlowski & Kroupa 2013).

A similar plane of satellite galaxies consisting of about half of the M31 satellites has been discovered by Ibata et al. (2013). It is seen edge-on from the MW. Ibata et al. (2013) found that 13 out of 15 satellites in this Great Plane of Andromeda (GPoA) show coherent line-of-sight velocities, indicating that they might co-orbit M31, analogous to the satellites in the VPOS. In addition to the satellite galaxies, the isolated dwarf galaxies in the Local Group are also confined to two very narrow and highly symmetric planes (Pawlowski et al. 2013). Similarly, the dwarf galaxy members of the nearby NGC 3109 association are confined to a very narrow, linear distribution (Bellazzini et al. 2013; Pawlowski & McGaugh 2014a). Coherent alignments of satellite galaxies and streams around more distant hosts have also been discovered in recent years (e.g., summarized in Pawlowski & Kroupa 2014). Most recently, Tully et al. (2015) have reported two narrow dwarf galaxy planes around Centaurus A, and Müller et al. (2015) reported the discovery of 16 dwarf galaxy candidates distributed anisotropically around M83. Ibata et al. (2014a, 2015) have performed a systematic study of diametrically opposite satellite galaxies around hosts within the Sloan Digital Sky Survey (SDSS; York et al. 2000) and found that they preferentially display anti-correlated line-of-sight velocities, which is consistent with approximately 60% of all satellite galaxies being situated in thin, co-orbiting planes. This interpretation has been challenged by Phillips et al. (2015), who performed a similar analysis of SDSS satellite pairs and find a signal consistent with only 30% of all satellites to be in co-orbiting planes. Another related analysis was performed by Cautun et al. (2015), who suggest that the observed signal is sensitive to selection effects. Given this situation, one will have to wait for further observational data, increasing the number of opposite satellite pairs with known line-of-sight velocities, to put tighter constraints on the fraction of co-orbiting satellites in the universe. However, the existence of satellite-plane configurations in the Local Group is not subject to debate.

The observed narrowness of the satellite galaxy planes of the Local Group and their coherent kinematic signature have been found to be extremely unlikely among satellite galaxies drawn from a variety of cosmological simulations. Metz et al. (2007) found the spatial distribution of the then known MW satellites to be inconsistent with an isotropic or prolate dark matter substructure distribution at a 99.5% level. Pawlowski et al. (2012a) determined that the correlation and polar orientation of satellite orbits, seen among the eight satellites for which proper motion measurements were available at that time, were very rare (0.4%–2%) among dark matter sub-halos drawn from the Via Lactea (Diemand et al. 2007, 2008) and Aquarius simulations (Springel et al. 2008). A later re-analysis (Pawlowski & Kroupa 2013) used improved proper motion information and included three additional satellites for which proper motions were determined in the meantime. They found that the observed kinematic coherence of the 11 brightest MW satellites alone is present in only 0.6% of satellite systems found in the Via Lactea and Aquarius simulations, and is also very unlikely when compared to the ΛCDM predictions of Libeskind et al. (2009) and Deason et al. (2011). Ibata et al. (2014b) determined that finding a GPoA-like structure with a similar kinematic signature in the Millennium-II simulation has a chance of only 0.04%. This result was confirmed by Pawlowski et al. (2014), who also demonstrated that a satellite structure as flattened and coherently orbiting as the VPOS is found in only 0.06% of all cases in this simulation. This is a conservative upper limit since the analysis used only the spatial and kinematic alignment of the 11 brightest satellites. Most recently, Gillet et al. (2015) were unable to find a GPoA equivalent in their analysis of a simulation of a Local Group equivalent from the CLUES project. Their results also show that simulated satellites, which can be assigned to a planar distribution for a given snapshot of a simulation, display large velocities perpendicular to that plane. While not testable for the M31 case, this is inconsistent with the observed alignment of the orbital directions of the MW satellites in the VPOS. A similar analysis by Pawlowski & McGaugh (2014b) used the "Exploring the Local Volume in Simulations" (ELVIS) suite of simulations by Garrison-Kimmel et al. (2014b) and found that only one out of 4800 realizations (0.02%) was able to reproduce the observed flattening and orbital correlation of the 11 brightest MW satellites. They also conclude that host galaxies in Local Group environments are as unlikely to harbour satellite galaxy planes as isolated hosts. In contrast to these analyses, claims of consistency of the observed satellite galaxy planes have frequently been shown to ignore aspects such as the kinematic coherence of the observed structure, to be based on unrealistic assumptions, or to not properly account for survey selection effects (e.g., Metz et al. 2009; Pawlowski et al. 2012a; Ibata et al. 2014b; Pawlowski et al. 2014).

The small-scale problems are often not seen as immediate failures of the ΛCDM model, but rather as possibly caused by flaws in constructing ΛCDM predictions from non-dissipative "dark-matter-only" simulations (but, see, e.g., Kroupa 2012, 2015). On scales of galaxies, the effects of baryons such as cosmic reionization, star formation, and feedback from stars, supernovea or AGNs become important (Navarro et al. 1996; Bullock et al. 2000; Mashchenko et al. 2008; Bovill & Ricotti 2009; Governato et al. 2010; Sawala et al. 2010; Pontzen & Governato 2012; Del Popolo et al. 2014). Consequently, there are now numerous attempts to include baryonic effects in cosmological simulations by modeling the hydrodynamics of the gas, as well as star formation and feedback processes on a sub-grid level (e.g., Libeskind et al. 2010; Vogelsberger et al. 2014). Unfortunately, these models often suffer from unphysical modeling of the baryonic processes (see, for example, the discussions in Kroupa 2015; Schaye et al. 2015).

While the effect of baryons can be substantial for ΛCDM predictions and cause a re-evaluation of some of the problems (such as the missing satellites problem, Koposov et al. 2009; or the core-cusp problem, Pontzen & Governato 2012; but see Peñarrubia et al. 2012 why solutions to these two problems might be mutually exclusive), others are more independent of baryonic effects. This is particularly true for the satellite plane problem: the distribution of satellite dwarf galaxies and their kinematic coherence on scales of hundreds of kiloparsecs—where the gravitational effect of the dark matter is dominant—can be safely assumed to not be affected by the internal physics of individual dwarf galaxies.

In this respect, it can be seen as surprising that Sawala et al. (2014, S14 hereafter; the same arguments apply to Sawala et al. 2015, in which they present the same results without addressing the concerns made in our Sections 2 and 3) have claimed that the set of hydrodynamic cosmological simulations they analyze resolves not only the missing satellites and TBTF problems, but also the satellite plane problem. Their work is based on the "Evolution and Assembly of GaLaxies and their Environments" (EAGLE) simulations by the Virgo Consortium (Schaye et al. 2015). They analyze a set of hydrodynamical cosmological zoom-simulations leading to the formation of Local-Group-like galaxy pairs, which model a number of baryonic physics effects (see Schaye et al. 2015 and the supplementary information of S14 for more details).

In Section 2, we investigate whether the measure of flattening of a satellite system from the satellites' projected positions (Starkenburg et al. 2013, S14) provides a reliable test of similarity to the observed situation (Section 2.1), whether the distribution of orbital poles of simulated satellite systems shown by S14 is similarly correlated as those of the observed MW satellites (Section 2), and whether the range of flattening of satellite systems found in hydrodynamic cosmological simulations are indeed different to purely collision-less "dark-matter-only" simulations of Local Group equivalents (Section 2.3). We will then discuss the persistence of the TBTF problem in Section 3. Our results will be summarized and discussed in Section 4.

The satellite plane problem for the eleven classical MW satellite galaxies⁵ can be summarized as follows: the (three-dimensional) positions of the MW satellite galaxies are distributed such that they all lie close to one common plane, which has a root-mean-square (rms) height of 19.6 kpc and an rms minor-to-major axis ratio of c/a = 0.182. Most of these satellites have orbital planes that are closely aligned with the plane defined by their position and they also share the same orbital direction (they co-orbit). This is indicated by the close clustering of 8 out of 11 orbital poles (directions of angular momentum) close to one normal vector describing the orientation of the best-fitting plane. In addition, one of the remaining 3 of the 11 orbital poles is directed along the opposite normal vector (i.e., retrograde relative to the 8 others).

Therefore, the defining characteristics of the satellite plane problem are that

• 1.  
  the satellites are distributed in a highly flattened, planar structure in three-dimensional space,
• 2.  
  the majority of the satellites co-orbit in the same sense, and
• 3.  
  these satellites orbit within the plane, indicating that the plane is not just a transient alignment.

2.1. Positions

The surprising degree of flattening of the MW satellite system has commonly been described in either absolute terms, as the rms (minor-axis) height of, e.g., the 11 classical MW satellites (${r}_{\mathrm{per}}=19.6\;{\rm{kpc}},$ e.g., Kroupa et al. 2005; Zentner et al. 2005; Metz et al. 2007, 2009; Kroupa et al. 2010; Wang et al. 2012; Pawlowski et al. 2013, 2014), or as a relative flattening of the distribution's minor-to-major axis ratio ($(c/a{)}_{\mathrm{std}}^{\mathrm{MW}}=0.18,$ e.g., Metz et al. 2007; Libeskind et al. 2005, 2009 ⁶, Deason et al. 2011; Pawlowski et al. 2013, 2014; Wang et al. 2013; Pawlowski & McGaugh 2014b). Commonly, the major and minor axis directions are determined from the eigenvectors of the tensor of inertia (ToI) defined by the non-weighted satellite positions.

Starkenburg et al. (2013) and S14 chose to define flattening in a different and idiosyncratic way. They characterize the spatial anisotropy of a satellite system using the ratio of eigenvalues of the reduced ToI (Bailin & Steinmetz 2005). For the 11 classical MW satellites, this results in⁷ $(c/a{)}_{\mathrm{red}}^{\mathrm{MW}}=0.36.$ Their approach is equivalent to measuring the axis ratios of the satellite distribution after the satellites have been projected onto a unit sphere around the host galaxy's center. Such a measure does not characterize a planar flattening of the satellite system, but only the preference of satellites to lie close to one common great circle, or to cluster about two opposed directions on the sky.

No justification was provided by these authors for their non-standard way of characterizing the "flattening" of a satellite system. S14 refer to Bailin & Steinmetz (2005). That study used the reduced ToI to measure and compare the principal axes of the dark matter particle mass distribution among six concentric shells within dark matter halos. The radial distance normalization was introduced so that substructures in the outer part of a shell would not dominate the ToI. These arguments do not apply to the satellite plane problem, which is not concerned with the principal axes of the mass distribution but merely with the distribution of satellite positions (each satellite carries the same weight regardless of its mass); substructure is not an issue in the satellite distribution (which consists only of ∼10¹ objects), and no comparison between different radial shells is intended (only one axis ratio of the whole satellite system is desired).

One potential benefit of the reduced ToI is that it might be less sensitive to outliers at large distances. For the observed MW satellites, the most distant satellite Leo I is possibly not associated with the VPOS, because its proper motion results in a most-likely orbital pole not aligned with the satellite plane normal (see Section 2.2.1). Furthermore, Leo I is rapidly receding from the MW, such that—depending on the potential of the MW—it might not even be a bound satellite of the MW (Boylan-Kolchin et al. 2013).⁸ To test the sensitivity of the standard and reduced ToI on this most-distant satellite, we have repeated the ToI fits using only the 10 classical MW satellites excluding Leo I. As expected, the axis ratio of the reduced ToI fit is not strongly affected ($(c/a{)}_{\mathrm{std}}^{\mathrm{noLeoI}}=0.37,$ compared to 0.36 including Leo I), but neither is that of the standard ToI ($(c/a{)}_{\mathrm{std}}^{\mathrm{noLeoI}}=0.20,$ compared to 0.18 including Leo I). The orientation of the satellite plane derived using the reduced ToI changes by less than one degree, and that derived using the standard ToI by 64. For the latter, excluding Leo I results in a slightly better alignment with the average orbital poles in Figure 2.

However, the standard ToI fit is less affected by the satellites closest to the MW. For a non-infinitesimally thin, planar structure these nearby satellites are expected to be the most offset from the overall distribution in their angular positions. Furthermore, they reside in a region in which processes that have the potential to change their orbit, such as precession and satellite–satellite encounters, are more important. Thus, while maybe less sensitive to distant outliers, the reduced ToI is more sensitive to outliers at close distances. In the case of the MW, this can be demonstrated using Sagittarius, which has a Galactocentric distance of only 18 kpc but is orbiting perpendicular to the VPOS. Applying the two ToI methods to the 10 classical satellites, except Sagittarius, results in almost the same axis ratios for the standard ToI ($(c/a{)}_{\mathrm{std}}^{\mathrm{noSag}}=0.17,$ compared to 0.18, including Sagittarius), but in largely different parameters for the reduced ToI fit ($(c/a{)}_{\mathrm{red}}^{\mathrm{noSag}}=0.25,$ compared to 0.36 including Sagittarius). Similarly, the best-fit plane orientation is essentially unchanged for the standard ToI (the directions differ by only 20), but, for the reduced ToI, the orientation is changed by 95 if Sagittarius is excluded from the fit. Thus, compared to the standard ToI, the reduced ToI trades a somewhat decreased sensitivity to distant outliers for a higher sensitivity to nearby outliers.

In the following, we will show that the use of the reduced ToI ignores available information and is a less robust and less meaningful test of the flattening of a satellite system.

If the position vector of satellite i is ${{\boldsymbol{r}}}_{i}=({r}_{i}^{x},{r}_{i}^{y},{r}_{i}^{z}),$ the standard ToI is

$$\begin{eqnarray}{I}_{\mathrm{std}}({\boldsymbol{r}}) & = & \displaystyle \sum _{i=1}^{11}\;\left({{\boldsymbol{r}}}_{i}^{2}\cdot {\mathbb{1}}\right)-\left({{\boldsymbol{r}}}_{i}{{\boldsymbol{r}}}_{i}^{{\rm{T}}}\right)\\ & = & \displaystyle \sum _{i=1}^{11}\;\left(\begin{array}{rcl}{\left({r}_{i}^{y}\right)}^{2}+{\left({r}_{i}^{z}\right)}^{2} & -{r}_{i}^{x}{r}_{i}^{y} & -{r}_{i}^{x}{r}_{i}^{z}\\ -{r}_{i}^{x}{r}_{i}^{y} & {\left({r}_{i}^{x}\right)}^{2}+{\left({r}_{i}^{z}\right)}^{2} & -{r}_{i}^{y}{r}_{i}^{z}\\ -{r}_{i}^{x}{r}_{i}^{z} & -{r}_{i}^{y}{r}_{i}^{z} & {\left({r}_{i}^{x}\right)}^{2}+{\left({r}_{i}^{y}\right)}^{2}\end{array}\right)\end{eqnarray} \tag{ 1 }$$

where ${{\boldsymbol{r}}}_{i}^{{\rm{T}}}$ is the transposed position vector ${{\boldsymbol{r}}}_{i}.$ The elements of the reduced ToI are defined as (Bailin & Steinmetz 2005):

$$\begin{eqnarray}&&{I}_{\mathrm{red}}^{\alpha ,\beta }({\boldsymbol{r}})=\displaystyle \sum _{i=1}^{11}\;\displaystyle \frac{{r}_{i}^{\alpha }{r}_{i}^{\beta }}{{r}_{i}^{2}},\end{eqnarray} \tag{ 2 }$$

so the total reduced ToI is

$$\begin{eqnarray}{I}_{\mathrm{red}}({\boldsymbol{r}})=\displaystyle \sum _{i=1}^{11}\;\left(\begin{array}{rcl}{\left({r}_{i}^{x}\right)}^{2} & {r}_{i}^{x}{r}_{i}^{y} & {r}_{i}^{x}{r}_{i}^{z}\\ {r}_{i}^{x}{r}_{i}^{y} & {\left({r}_{i}^{y}\right)}^{2} & {r}_{i}^{y}{r}_{i}^{z}\\ {r}_{i}^{x}{r}_{i}^{z} & {r}_{i}^{y}{r}_{i}^{z} & {\left({r}_{i}^{z}\right)}^{2}\end{array}\right)/{{\boldsymbol{r}}}_{i}^{2}.\end{eqnarray} \tag{ 3 }$$

Therefore, if each position vector ${{\boldsymbol{r}}}_{i}$ is normalized to ${{\boldsymbol{r}}}_{i}^{\prime }=\displaystyle \frac{{{\boldsymbol{r}}}_{i}}{\left|{{\boldsymbol{r}}}_{i}\right|},$ then ${I}_{\mathrm{red}}({{\boldsymbol{r}}}^{\prime })={I}_{\mathrm{red}}({\boldsymbol{r}})$ and

$$\begin{eqnarray}{I}_{\mathrm{std}}({{\boldsymbol{r}}}^{\prime })+{I}_{\mathrm{red}}({{\boldsymbol{r}}}^{\prime }) & = & \displaystyle \sum _{i=1}^{11}\;\left({\left({r}_{i}^{x}\right)}^{2}+{\left({r}_{i}^{y}\right)}^{2}+{\left({r}_{i}^{z}\right)}^{2}\right){\mathbb{1}}\\ & = & \displaystyle \sum _{i=1}^{11}\;{\mathbb{1}}.\end{eqnarray} \tag{ 4 }$$

We can see that the standard ToI ${I}_{\mathrm{std}}$ constructed using the normalized satellite positions is equivalent to the negative of the reduced ToI, ${I}_{\mathrm{red}},$ except for an additive factor proportional to the unit matrix ${\mathbb{1}}.$ This has the important consequence that eigenvectors of ${I}_{\mathrm{std}}$ constructed from the normalized satellite positions are identical to those of ${I}_{\mathrm{red}}$ (every vector is an eigenvector to the unit matrix). Either ToI will therefore give the same axes of most and least flattening of a satellite system with normalized positions. Due to the negative sign in this relation, the eigenvector corresponding to the shortest axis of the satellite distribution belongs to the smallest eigenvalue of ${I}_{\mathrm{red}},$ but to the largest eigenvalue of ${I}_{\mathrm{std}}.$ In the following, we will use both the reduced ToI and the more commonly employed standard ToI to compute the respective axis ratios ${(c/a)}_{\mathrm{red}}$ and ${(c/a)}_{\mathrm{std}}$ of satellite systems.

In spite of this connection, however, the reduced ToI ${I}_{\mathrm{red}}$ contains less information than the standard ToI, because only two angular coordinates of the satellites are used, not the full positions. It can thus be expected that using ${(c/a)}_{\mathrm{red}}$ as a metric for flattening of satellite systems will be less robust than a test based on ${(c/a)}_{\mathrm{std}}.$ We now demonstrate this. In what follows, we focus on the 11 brightest, classical MW satellites, whose positions we take from the compilation of McConnachie (2012).

2.1.1. Random Satellite Positions from an Isotropic Distribution

To assess whether testing the flattening of a satellite distribution using the reduced ToI gives comparable results to using the standard ToI, we construct 10,000 randomized satellite systems drawn from an isotropic distribution. One random realization is constructed by assigning a random angular position on a sphere (centered on the Galactic center) to each of the 11 classical MW satellite galaxies, while keeping each satellite's Galactocentric distance fixed. This ensures that the satellite distribution follows the observed radial distribution of the MW satellites exactly. We then measure ${(c/a)}_{\mathrm{std}}$ and ${(c/a)}_{\mathrm{red}}$ for each realization and determine whether they fall below the observed value for the MW satellite system. The resulting cumulative distributions are plotted in the upper left panel of Figure 1.

Zoom In Zoom Out Reset image size

We find that ${(c/a)}_{\mathrm{std}}\leqslant (c/a{)}_{\mathrm{std}}^{\mathrm{MW}}$ in 1.05% of all cases. The spatial flattening of the MW system measured via the minor-to-major axis ratio of the full three-dimensional positions of the 11 brightest satellites is thus very unlikely under the hypothesis that the satellite positions are isotropically distributed. For the reduced ToI method, using normalized satellite distances, we find ${(c/a)}_{\mathrm{red}}\leqslant (c/a{)}_{\mathrm{red}}^{\mathrm{MW}}$ in 1.95% of all cases. Testing anisotropy using this method is roughly twice as likely to agree with an isotropic distribution. In this sense, the reduced ToI test is less discriminating.

Of the 195 realizations fulfilling the ${(c/a)}_{\mathrm{red}}$ criterion, only 38 simultaneously fulfill the ${(c/a)}_{\mathrm{std}}$ criterion. The majority, 81%, of the realizations that pass the ${(c/a)}_{\mathrm{red}}$ test are less flattened than the MW satellite system; that is, they have a larger axis ratio than observed. The upper right panel of Figure 1 demonstrates this point. The figure compares the axis ratios ${(c/a)}_{\mathrm{std}}$ and ${(c/a)}_{\mathrm{red}}$ for the set of 10,000 randomized satellite distributions. While a weak correlation between the two measures is present, there is considerable scatter, and most realizations, which have as low a ${(c/a)}_{\mathrm{red}}$ as the MW satellite system (left of the red dashed line), are not simultaneously as flattened in their full 3D distribution measured via ${(c/a)}_{\mathrm{std}}$ (below the black solid line).

We conclude that using the reduced ToI test as a statistic for the degree of spatial flattening produces more than 80% false positives when applied to the MW satellite system. The reduced ToI test is an insufficient criterion for the comparison with the observed three-dimensional satellite distribution.

2.1.2. Fixed Satellite Positions, but Reshuffled Distances

The spatial flattening of a satellite system implies a certain degree of correspondence between a satellite's distance, and its angular coordinate. The reduced ToI axis ratio ${(c/a)}_{\mathrm{red}}$ is completely insensitive to this correspondence since it considers only the normalized satellite positions. To demonstrate this, we generated 10,000 "reshuffled" satellite distributions, keeping the angular positions of the 11 classical MW satellites fixed, but selecting their distances randomly, without replacement, from the observed set of 11 distances. This procedure again conserves the observed radial distribution since each distance is selected exactly once per realization.

The lower left panel of Figure 1 shows the results. Not surprisingly, most of distance-reshuffled realizations (more than 92.5%⁹ ) have larger ${(c/a)}_{\mathrm{std}}$ than the observed MW satellite population and are therefore less strongly flattened. The MW satellite distribution is thus indeed plane-like in three-dimensional space and not merely anisotropic in the angular positions around the MW (in which case one would expect the axis ratios of shuffled realizations to more evenly distribute below and above the observed $(c/a{)}_{\mathrm{std}}^{\mathrm{MW}}$). The reduced ToI test, however, is by definition unable to discriminate between these different cases and results in the same ${(c/a)}_{\mathrm{red}}$, independent of the reshuffling.

The lower right panel of Figure 1 provides a more intuitive illustration of this. It shows the positions of the 11 classical MW satellites (blue squares) such that the best-fit plane (determined using the reduced ToI) is seen edge-on. In addition, the red plus signs mark the satellite positions in one of the reshuffled distributions: they clearly do not look like a flattened planar distribution. This impression is confirmed by their (standard ToI) axis ratio of ${(c/a)}_{\mathrm{std}}=0.46.$ However, by construction, the reduced ToI axis ratio is the same as the observed one. Similarly, the green crosses mark the satellite positions if they would all reside at the exact same Galactocentric distance of 200 kpc (i.e., their distances are normalized to 200 kpc). Again, their ${(c/a)}_{\mathrm{std}}$ axis ratio is much larger than that of the observed satellite positions, but the reduced ToI test would classify this system to be just as flattened as the observed MW satellite system.

As expected, using the reduced ToI is not a robust test of the spatial flattening of satellite galaxy systems. Ignoring the radial distances of the satellites when testing their spatial anisotropy biases toward finding agreement with the observed situation even though the full three-dimensional distribution can be markedly different.

2.2. Kinematics

Here, we suggest a simple, first-order test, to check whether an orbital pole distribution, for instance, such as that of the EAGLE model satellites shown in Figure 4 of S14, indicates a similarly unexpected clustering of orbital poles close to the normal direction of the best-fitting plane to the satellite positions as is present for the observed data.

2.2.1. The Observed Orbital Poles

For the MW satellites, we determine their orbital poles from their observed positions and proper motions as discussed in detail in Pawlowski & Kroupa (2013). If more than one proper motion measurement is available for a particular satellite, we use the uncertainty-weighted proper motion to determine the satellite's most-likely orbital pole direction. For Draco, we update the Hubble Space Telescope proper motion measurement from the preliminary value used in Pawlowski & Kroupa (2013) to the one published in Pryor et al. (2015). Compared to the orbital pole published in Pawlowski & Kroupa (2013) this improves the alignment of Draco's orbital pole with the VPOS plane normal (angle to the normal of the best-fitting plane to the 11 classical satellites ${\theta }_{\mathrm{VPOS}}^{\mathrm{class}}={14\buildrel{\circ}\over{.} 1}_{-0.1}^{+2.6}$) and results in an even tighter concentration of orbital poles. The eight most-concentrated poles now have a spherical standard distance of ${{\rm{\Delta }}}_{\mathrm{sph}}\;=$ 272 instead of 293. This is in line with the trend that better proper motion data results in closer alignments of the derived orbital poles with the satellite plane axis and a reduction of the scatter around it.

Figure 2 plots these observed orbital poles. The upper panel shows the orbital poles and the normal directions to the VPOS, the plane fitted to the young halo globular clusters (Pawlowski et al. 2012b) and the average stream normal (Pawlowski & Kroupa 2014) in Galactic coordinates. The lower panel shows the poles in a coordinate system in which the "north" and "south" poles of the all-sky-plot are pointing along the minor axis of the distribution of the 11 classical MW satellites, as determined using the reduced ToI method. This causes the projected satellite positions, plotted as blue squares, to lie close to the "equator" of this coordinate system.

Zoom In Zoom Out Reset image size

2.2.2. Binned Orbital Poles

To test whether a simulated satellite system has an orbital pole distribution that is as strongly clustered as the observed orbital poles of the MW satellites, we start by drawing 100,000 realizations of 11 satellite positions from an isotropic distribution with normalized distances. For each satellite, a velocity direction is selected at random from an isotropic distribution. Then the orbital pole of each of the satellites is constructed.

Using the reduced ToI method, we determine the minor axis of the satellite distribution for each realization of 11 satellites. We then transform the satellite positions and orbital pole directions of a given realization into a "satellite plane coordinate system": we align the minor axis with the poles of the coordinate system. Because the minor axis is an axial direction, we are left with two possibilities of how to determine the "northern" (upper) coordinate system pole. We decide to define it as that of the plane normal directions, which has the larger number of orbital poles aligned to within 60^◦ (one-quarter of the sky). If both directions have the same number of aligned poles, one is chosen at random.

Eight of the observed MW satellites have orbital poles pointing to within 60^◦ of the satellite plane normal. Out of the 100,000 random realizations, only 0.79% have at least this observed number of 8 orbital poles within 60° (bin 1) of the northern plane pole. Even finding only at least 7 poles in this area is rare and occurs in 4.1% of the realizations. This demonstrates that the observed orbital pole distribution for the MW is rare and significant at the 99% level, but note that we are using a simplified test that does not utilize the full information available for the MW satellite orbital pole distribution such that this underestimates the significance (for a better test, see Pawlowski & Kroupa 2013).

We define four equal-area bins in the satellite coordinate system. Bin 1 is defined as <60° from the "northern" pole, bin 2 as $\leqslant 30^\circ$ "north" of the equator, bin 3 as $\leqslant 30^\circ$ "south" of the equator and bin 4 covers the area $\lt 60^\circ$ from the "southern" pole. They contain 8, 1, 1, and 1 of the observed MW satellite orbital poles, respectively. This allows one to easily compare the observed orbital pole distribution with those of simulations. As an example, we use the simulation results shown in Figure 4 of S14. From that plot, we read 4, 2, 2, and 3 orbital poles in these bins.

A perfectly isotropic orbital pole distribution would yield an average of $\frac{11}{4}=2.75$ poles per bin (black dotted line in Figure 3). However, the coordinate system is defined by the angular positions of the satellites, which are by definition 90° away from the associated orbital pole (the direction of angular momentum around the origin is always perpendicular to the position vector). Thus, in a coordinate system such as ours, which is defined to have an equator 90° away from the reduced ToI minor axis, a completely isotropic distribution of positions and velocities (which we sampled with the 100,000 realizations discussed before) results in slightly more orbital poles close to the coordinate system's "north" and "south" poles, and thus bins 1 and 4, than would be the case for a completely randomly oriented coordinate system. The resulting distribution is plotted as the black solid line in Figure 3, the dashed lines indicate the standard deviation from this average (about ±1.2 poles per bin).

Zoom In Zoom Out Reset image size

The green dots in the figure illustrate the distribution of the most-likely positions of the observed orbital poles (Pawlowski & Kroupa 2013). The strong peak of eight poles within the "northernmost" (right) bin 1 differs by more than three standard deviations from the expected number of about four in the case of an isotropic distribution. In contrast to the observed orbital poles, the simulated orbital poles (red crosses, read from Figure 4 of S14), follow the expectation derived from the isotropic distributions of satellite velocities if the effect of defining the coordinate system by the satellite positions on the sky is taken into account. For this example, our test thus concludes that there is no indication that the simulation of S14 reproduces a similarly strongly rotationally supported disk of satellites like the one around the MW.

Zoom In Zoom Out Reset image size

To estimate how much the proper motion uncertainties of the observed MW satellites affect this comparison, we have generated 10,000 realizations of orbital pole distributions by randomly varying the proper motions within their uncertainties. For each realization, the counts in bins 1–4 were recorded. Figure 3 shows the resulting average counts (smaller open dots) and the standard deviations around them (error bars). Even though the orbital poles of Sextans and Carina have large uncertainties that can move them out of bin 1, on average, there are more than seven orbital poles in bin 1 and that number does not drop below 6. Thus, even though some of the observed orbital poles have considerable uncertainties, the tension with S14 remains significant when accounting for them.

An extreme kinematic coherence is the most important defining characteristic of the VPOS (and also the GPoA and other satellite galaxy planes). A simulation, which does not contain a similar number of closely concentrated orbital poles as the observed MW satellite system, can therefore not be claimed to resolve the disk of satellites problem. In particular, it is clear that the simulation analyzed and plotted by S14 does not have as strong of a correlation of orbital poles with the normal direction to the satellite plane as the observed MW satellite system.

2.3. Comparison to Collision-less Simulations

The TBTF problem can be formulated as follows. The measured densities at the (de-projected) half-light radii ${r}_{1/2}$ of local galaxies are too low compared to the most massive subhalos predicted in ΛCDM simulations (Boylan-Kolchin et al. 2012). The enclosed dynamical mass (and thus density) within r_1/2, or its corresponding circular velocity ${v}_{\mathrm{circ}},$ can be directly and accurately estimated from the observed line-of-sight velocity dispersion (Wolf et al. 2010).

In contrast, some theoretical papers, such as S14, chose to compare the maximum circular velocity ${v}_{\mathrm{max}}$ of galaxies formed in their simulation to that of the MW satellites. This velocity cannot be directly measured for observed dwarf galaxies. Inferring ${v}_{\mathrm{max}}$ from observational data requires assumptions about how the dark matter density profile extrapolates from ${r}_{1/2}\approx 500\;\mathrm{pc},$ where the kinematics are measured, to the radius ${r}_{\mathrm{max}}\approx 1-5\;\mathrm{kpc},$ where the circular velocity curve of the satellite's expected dark matter halo has its maximum. Assigning ${v}_{\mathrm{max}}$ values to observed satellites is therefore very uncertain and depends on the assumptions made. Consequently, comparing simulations to observations based on ${v}_{\mathrm{max}}$ is considerably less straightforward than the comparison of ${v}_{\mathrm{circ}}$ at ${r}_{1/2},$ which, for best accuracy, can use the simulation particle data directly (see, e.g., Boylan-Kolchin et al. 2012; Tollerud et al. 2014). Due to these possible systematics, counting satellites above a given ${v}_{\mathrm{max}}$ obscures the underlying issue of the TBTF problem and does not directly address it. We therefore refer to this formulation as the extrapolated TBTF problem.

A meaningful comparison with the observed MW satellites requires a careful assessment of the best available observational information, in particular, when the method of comparison can be severely affected by systematics. S14 only compare their simulated galaxies with the ${v}_{\mathrm{max}}$ values of Peñarrubia et al. (2008). To estimate ${v}_{\mathrm{max}}$ for the observed MW satellites, their study assumes that the dark matter halos follow an NFW profile (Navarro et al. 1996) and makes use of the relation between ${v}_{\mathrm{max}}$ and ${r}_{\mathrm{max}}$ as given by field halos in ΛCDM. Thus, already, the comparison of these velocities to the dark-matter-only simulations by S14 is problematic (the analysis is concerned with satellite sub-halos, not with isolated field halos), but even more so for the hydrodynamical simulations, which result in less concentrated sub-halos (having a lower central dark matter density) than dark-matter-only simulations (Zolotov et al. 2012). Furthermore, more recent estimates result in significantly lower values of ${v}_{\mathrm{max}}$ for the same MW satellites (Kuhlen 2010; Boylan-Kolchin et al. 2012).

To demonstrate the differences in the estimated ${v}_{\mathrm{max}},$ we plot the cumulative distribution of the number of sub-halos above a given ${v}_{\mathrm{max}}$ in Figure 5. The curves for the simulated galaxies (and sub-halos for the dark-matter-only runs) have been extracted from Figure 3 of S14 by carefully tracing them. In addition to the ${v}_{\mathrm{max}}$ values that Peñarrubia et al. (2008) report for the MW satellites (magenta crosses), we also plot the two more recent estimates (yellow plus signs and black dots). In particular, the ${v}_{\mathrm{max}}$ values by Boylan-Kolchin et al. (2012) are significantly lower than those of Peñarrubia et al. (2008). According to Boylan-Kolchin et al. (2012), this difference is due to their much more detailed analysis. They compare the observed properties of the MW satellites with a large number of sub-halos from the Aquarius simulations (Springel et al. 2008), use more recent dynamical constraints based on the mass within the half-light radius (which has smaller uncertainties than previous mass estimates, Wolf et al. 2010), and even correct for the effect of the gravitational softening length on the inner halo densities. Comparing the resulting ${v}_{\mathrm{max}}$ values of Boylan-Kolchin et al. (2012) instead of those of Peñarrubia et al. (2008) with the simulations of S14 shows that the hydrodynamical simulations merely alleviate but do not resolve the extrapolated TBTF problem. Below ${v}_{\mathrm{max}}\lesssim 25\;\mathrm{km}\;{{\rm{s}}}^{-1},$ the simulations contain about twice as many satellites with a given velocity than observed.

Zoom In Zoom Out Reset image size

Figure 5 illustrates that the value of ${v}_{\mathrm{max}}$ that is assigned to a satellite depends very sensitively on theoretical priors. It is not a direct measurement. None of these issues would arise if, to test whether the TBTF problem is resolved, the observed ${v}_{\mathrm{circ}}({r}_{1/2})$ would have been directly compared with the rotation curves of the simulated galaxies (as, for example, done by Boylan-Kolchin et al. 2012; Tollerud et al. 2014). However, even at the highest resolution level L1, S14 still have a softening length of $\epsilon =94\;\mathrm{pc}$, and for their extrapolated TBTF problem comparison they use resolution level L2, which has $\epsilon =216\;\mathrm{pc}$. According to Boylan-Kolchin et al. (2012), softening has the effect to reduce the density on scales of about $3\epsilon ,$ such that cumulative quantities like ${v}_{\mathrm{circ}}$ are underestimated. The de-projected half-light radii of several of the classical MW satellites are of the same order, five out of the nine dSphs analyzed in Boylan-Kolchin et al. (2012) have ${r}_{1/2}$ between 200 and 400 pc. It therefore appears unlikely that the simulations of S14 reliably resolve the radii relevant for the TBTF problem. They might even lack the resolution to conclusively investigate the problem at all.

Finally, Figure 5 (or Figure 3 of S14) also shows that there appears to be no problem at the high-mass end of the cumulative ${v}_{\mathrm{max}}$ distribution for the dark-matter-only simulations (gray lines). If using the Peñarrubia et al. (2008) velocities, the dark-matter-only distribution matches the satellite data above $30\;\mathrm{km}\;{{\rm{s}}}^{-1}.$ The dark-matter-only simulations begin to over-predict the ${v}_{\mathrm{max}}$ function only at lower masses; therefore, it looks as if the extrapolated TBTF problem was solved without the need for any hydrodynamics. One likely reason for the low number of simulated halos above ${v}_{\mathrm{max}}=30\;\mathrm{km}\;{{\rm{s}}}^{-1}$ is that the host halos are at the low-mass end of what is reasonable for the MW and M31. Already Boylan-Kolchin et al. (2012) have discussed that a significantly less massive MW halo could, in principle, avoid the TBTF problem, but their estimate of the required halo mass of about 0.5 × 10¹² ${M}_{\odot }$ appears unrealistically small (see the discussion in their Section 5.1, but also Wang et al. 2012; Vera-Ciro et al. 2013). According to the supplementary information of S14, the second-most-massive galaxies in their Local Group equivalent simulations (which can be identified as the MW equivalents) have a median mass of only 0.9 × 10¹² and a lower limit of 0.5 × 10¹² ${M}_{\odot }$. This means that half of the MW equivalents have halo masses between 0.5 and 0.9 × 10¹² ${M}_{\odot }$.

Similarly, the stellar masses of the MW and M31 analogs in the S14 hydrodynamical simulations lie in the range of $(1.5-5.5)\times {10}^{10}\;{M}_{\odot }.$ These are very low compared to the observed values. The lower limit of this range, $1.5\times {10}^{10}\;{M}_{\odot },$ is extremely low compared to the stellar mass of the MW ($(4.9-5.5)\times {10}^{10}\;{M}_{\odot };\;$Flynn et al. 2006) and even the upper end of the range, $5.5\times {10}^{10}\;{M}_{\odot },$ is very low compared to the stellar mass of M31 ($(10-15)\times {10}^{10}\;{M}_{\odot };\;$Tamm et al. 2012). This again suggests that the associated halos are of too low mass, resulting in fewer sub-halos for a given ${v}_{\mathrm{max}}$ and thus alleviating the TBTF problem before any hydrodynamical effects come into play.

The addition of baryonic physics to ΛCDM models of hierarchical structure formation can ameliorate certain discrepancies between the observed and simulated structures of galaxies. In this paper, we investigate recent claims that baryonic physics can do more, providing a solution even to the "satellite planes" problem, the nearly planar distribution of satellite galaxies around the Milky Way and other galaxies. Such claims are surprising given that the scale of the satellite systems (hundreds of kiloparsecs) is much greater than the scale over which baryonic physics can be expected to act. By analyzing dissipationless simulations in a consistent way to a previous analysis of dissipational simulations, we find no significant difference between the predicted distribution of satellite plane flattenings. Claims to the contrary were shown to result in part from a non-standard measure of plane thickness that ignores the radial positions of the satellites.

In addition to being thin, the Milky Way satellite system also exhibits coherent rotation. We introduced a simple test for kinematic coherence that compares the alignment of orbital poles with the satellite plane normal. We apply the new test to the results of a hydrodynamical simulation and show that the latter fail to reproduce the concentration of orbital poles observed in the Milky Way satellite system. We find that the distribution of orbital poles in the hydrodynamical simulation is consistent with that produced by random, isotropic velocities, as in the dark-matter-only simulations. Here again, we find no evidence that the addition of baryonic physics is useful in reconciling simulations with observations.

Yet another issue of ΛCDM is the TBTF problem. TBTF refers to the high densities of the inner parts of simulated sub-halos compared with the central densities of observed satellites; the latter computed using the circular velocity, v_1/2, measured at the half-light radius r_1/2. Recent claims that the TBTF problem is resolved in hydrodynamical simulations are based on a comparison of the maximum circular velocities, ${v}_{\mathrm{max}},$ of simulated and observed galaxies. Since v_max cannot be measured for the Milky Way satellites, this method relies on a model-dependent extrapolation of the mass distribution. When we account for the fact that the simulated satellites are situated in sub-halos, not field halos, we find that the TBTF problem is still present, even in its extrapolated form. Once again, confrontation of simulations with observational data shows no preference for hydrodynamical over dark-matter-only simulations.

We thank Shea Garrison-Kimmel and the ELVIS collaboration for making their simulations publicly available. The contribution of M.S.P. to this publication was made possible through the support of a grant from the John Templeton Foundation. D.M. was supported by the National Science Foundation under grant no. AST 1211602 and by the National Aeronautics and Space Administration under grant no. NNX13AG92G.