Streams in the Aquarius stellar haloes Free

Abstract

1 INTRODUCTION

In the last decade, the characterization of the phase-space distribution of stars in the vicinity of the Sun has become a subject of great interest. It is now well established that by studying the phase-space distribution function of the Milky Way we may be able to unveil its formation history (Freeman & Bland-Hawthorn 2002; Helmi 2008).

According to the current paradigm of galaxy formation, the stellar haloes of large galaxies such as our own are mostly formed through the accretion and mergers of smaller objects (Searle & Zinn 1978), rather than by a rapid collapse of a large and pristine gas cloud (Eggen, Lynden-Bell & Sandage 1962). During the process of accretion, these systems are tidally disrupted and contribute gas and stars to the final object. Whereas the gas component rapidly forgets its dynamical origin due to its dissipative nature, stars should retain this information for much longer time-scales because their density in phase space is conserved. This memory is especially easily accessible in the outer regions of galaxies, or outer stellar haloes, where the dynamical mixing time-scales are very long. A large amount of observed stellar streams in the outer regions of the Milky Way supports this picture for the formation of its stellar halo. The Sagittarius stream (Ibata, Gilmore & Irwin 1994; Ibata et al. 2001) and the Orphan stream (Belokurov et al. 2007) are just two examples among many other overdensities recently discovered (see also Newberg et al. 2002; Ibata et al. 2003; Yanny et al. 2003; Belokurov et al. 2006; Grillmair 2006; Starkenburg et al. 2009).

This formation model also predicts that for Milky Way like galaxies, the portion of the accreted stellar halo that dominates in the solar neighbourhood is dynamically old (e.g. Cooper et al. 2010; Tissera, White & Scannapieco 2012). Therefore, fossil signatures of the most ancient accretion events that a galaxy has experienced should be buried in their inner regions or, in the case of the Milky Way, close to or in the solar neighbourhood. However, due to the very short dynamical time-scale characteristic of these regions, the debris from accretion events is expected to be spatially well mixed and therefore rather difficult to detect. Furthermore, in cold dark matter models, dark matter haloes are expected to be strongly triaxial (see, e.g. Allgood et al. 2006; Vera-Ciro et al. 2011). As a consequence, chaotic orbits may be quite important and result in much shorter mixing time-scales (see, e.g. Vogelsberger et al. 2008), which leads to a phase-space distribution that is effectively smooth. On the other hand, the dissipative condensation of baryons can induce a transformation of the structure of dark haloes to a more oblate and axisymmetric shape in the inner regions while preserving their triaxial shape in the outer parts (Debattista et al. 2008; Abadi et al. 2010; Valluri et al. 2010; Bryan et al. 2012).

The first attempt to quantify the amount of substructure in the form of stellar streams expected in the solar neighbourhood was carried out by Helmi & White (1999). They suggested that ∼300–500 kinematically-coherent substructures should be present in a local volume around the Sun. Although they considered a fixed Galactic potential, their results were later confirmed by Helmi, White & Springel (2003), who analysed a full cosmological simulation of the formation of a cluster dark matter halo scaled down to galaxy size. They also showed that the debris from accretion events appeared to mix on time-scales comparable to those expected for integrable potentials. However, their results were based on a single N-body simulation and were limited by its numerical resolution (despite the 66 million particles used).

Stellar streams in the inner regions of the Galactic halo cannot be detected simply by looking for overdensities in the distribution of stars on the sky. Thus, several different spaces have been proposed and explored to identify substructure. Examples are the velocity space (Helmi & White 1999; Helmi et al. 1999; Williams et al. 2011), metallicity, colour and chemical abundance spaces (e.g. Font et al. 2006; Schlaufman et al. 2011; Monachesi et al. 2013; Valluri et al. 2013) or spaces defined by integrals of motion and their associated orbital frequencies (e.g. Helmi & de Zeeuw 2000; Bullock & Johnston 2005; Knebe et al. 2005; McMillan & Binney 2008; Klement et al. 2009; Morrison et al. 2009; Gómez & Helmi 2010; Gómez et al. 2010). At first sight, these studies as well as that by Helmi et al. (2003) might be in tension with the results reported by Valluri et al. (2013). These authors found no substructure in the integral of motion space of halo stars formed in a fully cosmological hydrodynamical simulation. They attributed the lack of substructure to the strong chaotic mixing present in this simulation. Their results were based on samples of ≈10⁴ stellar particles distributed across the entire halo. This likely imposes a strong limitation on the ability to identify and quantify substructure. If, as discussed by Helmi et al. (2003), there are 300–500 streams crossing the solar neighbourhood, then a resolution of at least a few thousand stars in a local small volume (of a few kpc across) is required. Hence, a much larger number across the whole halo. Therefore, samples of a total of 10 000 halo stars would really be too small to find or characterize substructure. We confirm this intuition below, where we show that substructure is clearly visible when the numerical resolution (particle number) is large enough.

Thanks to current surveys such as RAVE (Steinmetz et al. 2006) and SEGUE (Yanny et al. 2009), as well as the forthcoming astrometric satellite Gaia (Perryman et al. 2001), model predictions are becoming testable. Although to date only a handful of stellar streams have been detected in the solar neighbourhood (e.g. Helmi et al. 1999; Klement et al. 2009; Smith et al. 2009), with the full 6D phase-space catalogues that are already becoming available (Bond et al. 2010; Breddels et al. 2010; Carollo et al. 2010; Zwitter et al. 2010; Burnett et al. 2011; Beers et al. 2012) it may soon be possible to start deciphering the formation history of the Milky Way.

In this work, we revise the predictions made by Helmi et al. (2003). Our goal is to characterize and quantify the amount of substructure inside solar neighbourhood like volumes obtained from the fully cosmological very high resolution simulation of galactic stellar haloes modelled by Cooper et al. (2010). In Section 2, we briefly describe these models. We characterize the phase-space distribution of stellar particles in ‘solar neighbourhood’ spheres in Section 3. In Section 4, we measure the number of stellar streams, and explore the relevance of orbital chaos for halo stars near the Sun, as well as the impact of particle resolution on the quantification of substructure. Finally, in Section 5 we summarize our results.

Throughout this work, we use the term substructure only to refer to substructure in the form of stellar streams.

2 THE SIMULATIONS

In this work, we analyse the suite of high-resolution N-body simulations from the Aquarius Project (Springel et al. 2008a,b), coupled with the galform semi-analytic model (Cole et al. 1994, 2000; Bower et al. 2006) as described by Cooper et al. (2010).

2.1 N-body simulations

The Aquarius Project has simulated the formation of six Milky Way like dark matter haloes in a Λ cold dark matter (ΛCDM) cosmology. The simulations were carried out using the parallel Tree-PM code gadget-3 (an upgraded version of gadget-2, Springel 2005). Each halo was first identified in a lower resolution version of the Millennium II Simulation (Boylan-Kolchin et al. 2009) which was carried out within a periodic box of side 125 h⁻¹ Mpc in a cosmology with parameters Ω_m = 0.25, Ω_Λ = 0.75, σ₈ = 0.9, n_s = 1 and Hubble constant  H₀ = 100 h km s⁻¹ Mpc⁻¹ = 73 km s⁻¹ Mpc⁻¹. The haloes were selected to have masses comparable to that of the Milky Way and to be relatively isolated at z = 0. By applying a multi-mass particle ‘zoom-in’ technique, each halo was resimulated at a series of progressively higher resolutions. The results presented in this work are based on the simulations with the second highest resolution (Aq-2) available, with a Plummer equivalent softening length of 65.8 pc. In all cases, 128 outputs, starting from redshift z ≈ 45,¹ were stored. In each snapshot, dark matter haloes were identified using a Friends-of-Friends (Davis et al. 1985) algorithm, while subhaloes were identified with SUBFIND (Springel et al. 2001). For more details about the simulations, we refer the reader to Springel et al. (2008a,b). The main properties of the resulting haloes are summarized in Table 1. Note that we exclude from our analysis the halo Aq-F. This is due to its assembly history; more than 95 per cent of its accreted stellar halo mass comes from a single galaxy that was accreted very late, at z ≈ 0.7.

2.2 Semi-analytic model

To study the growth and properties of stellar haloes, Cooper et al. (2010) coupled the semi-analytic model galform to the Aquarius dark matter N-body simulations. The basic idea behind semi-analytic techniques is to model the evolution of the baryonic components of galaxies through a set of observationally and/or theoretically motivated analytic prescriptions. Critical physical processes that govern galaxy formation, such as gas cooling, star formation feedback, supernovae, winds of massive stars and active galactic nuclei, etc., are taken into account. The parameters that control these processes were fixed on the basis of an implementation of galform on the Millennium Simulation (Springel et al. 2005), as this set of choices allows one to successfully reproduce the properties of galaxies on large scales as well as those of the Milky Way (Bower et al. 2006; Cooper et al. 2010).

2.3 Building up stellar haloes

Within the context of CDM, stellar haloes may be formed via hierarchical aggregation of smaller objects, each one of them imprinting their own chemical and dynamical signatures on the final halo's properties (Freeman & Bland-Hawthorn 2002). As previously explained, the galform model provides a detailed description of the composite stellar population of every galaxy in the simulation at any given time. However, it does not follow the dynamics of these systems fully. To study the dynamical properties of the resulting galactic stellar haloes, Cooper et al. (2010) assumed that the most strongly bound dark matter particles in progenitor satellites could be used to trace the phase-space evolution of their stars. In every snapshot of the simulation, a fixed fraction of the most bound dark matter particles were selected to trace any newly formed stellar population in each galaxy in the simulation. This implies that each tagged particle has a different final stellar mass associated with it. The fraction of selected bound particles is set such that properties of the satellite population at z = 0, like their luminosities, surface brightness, half-light radii and velocity dispersions, are consistent with those observed for the Milky Way and M31 satellites. For a detailed description of this procedure, we refer the reader to section 3 of Cooper et al. (2010). The properties of the resulting stellar haloes obtained in each simulation, such as total mass, half-light radius and the root-mean-square (rms) scatter in the logarithm of the stellar mass assigned to the individual particles, are listed in Table 1. Note that, although the method introduced in Cooper et al. (2010) leads to a successful match of various observables regarding the structure and characteristics of the Galactic stellar halo and its satellite population, the dynamical evolution of the baryonic components of galaxies is much simplified. This likely has an effect on, e.g. the efficiency of satellite mass-loss due to tidal stripping, or the satellite's internal structural changes due to adiabatic contraction, and possibly even on its final radial distribution (Libeskind et al. 2010; Romano-Díaz et al. 2010; Schewtschenko & Macciò 2011; Geen, Slyz & Devriendt 2013).

3 CHARACTERIZATION OF SUBSTRUCTURE IN SOLAR VOLUMES

In this section, we characterize the phase-space distribution of stellar particles inside a ‘solar neighbourhood’ sphere located at 8 kpc from the galactic centre of each Aquarius stellar halo. Following Gómez et al. (2010), we chose for the spheres a radius of 2.5 kpc as this is approximately the distance within which the astrometric satellite Gaia (Perryman et al. 2001) will be able to provide extremely accurate 6D phase-space measurements for an unprecedentedly large number of stars. The final configuration of the host dark matter haloes is, in all cases, strongly triaxial. Therefore, to allow a direct comparison between the ‘solar neighbourhood’ spheres from different haloes, we have rotated each halo to its principal axis and placed the sphere along the direction of the major axis. In all cases, the ratios and directions of the principal axis were computed using dark matter particles located within 6–12 kpc.

Table 2 shows that 90 per cent of the total stellar mass enclosed in these spheres comes, in all cases, from three to five significant contributors, |$N_{\rm ms}^{\it sn}$|⁠. This is in agreement with De Lucia & Helmi (2008) and Cooper et al. (2010), who find that stellar haloes are predominantly built from fewer than five satellites with masses comparable to the brightest classical dwarf spheroidals of the Milky Way. Fig. 1 shows the stellar mass fractions contributed by the five most significant contributors to the total stellar mass enclosed in each sphere.

Figure 1.

Cumulative stellar mass fractions obtained from our solar neighbourhood like spheres, centred at 8 kpc from the galactic centre. The step-like cumulative fractions show the contribution from the five most massive contributors to each sphere in rank order. The most massive contributor is ranked as number 1. The red dashed line shows the 90 per cent level.

Open in new tabDownload slide

In Fig. 2, we present two different projections of velocity space. The different colours indicate particles coming from different satellites that have contributed with at least five particles (while those from satellites contributing fewer particles are shown in black). It is interesting to observe how these distributions in velocity space vary from halo to halo. Essentially, less massive (dark) haloes have smaller velocity dispersions and thus the distribution of particles in velocity space is more compact. The values of the velocity dispersions of these distributions are listed in Table 2. A comparison with the estimated values of the velocity ellipsoid of the local stellar halo (see Chiba & Beers 2000), (σ_R, σ_ϕ, σ_Z) = (141 ± 11, 106 ± 9, 94 ± 8) km s⁻¹, shows that the ellipsoids of the haloes Aq-A-2, -C-2 and -D-2 have amplitudes comparable to those observed for the Milky Way. Note, however, that the dynamics of the stellar particles in these simulations are only determined by the underlying dark matter halo potential. If the mass associated with the disc and the bulge were to be included, the velocity dispersions would be significantly increased, since we may relate the velocity dispersion σ² to the circular velocity |$V_{\rm c}^2$| and |$V_{\rm c}^2 \propto M(<\!r) = M_{\rm \rm disc} + M_{\rm halo}$|⁠, and M_disc ∼ M_halo near the Sun (Binney & Tremaine 2008). In Table 2, we also show the local stellar halo average densities, ρ₀. We find that, except for the halo Aq-B-2, the values obtained are in reasonable agreement with the estimates for the solar neighbourhood, ρ₀ = 1.5 × 10⁴ M_⊙ kpc⁻³ (see, e.g. Fuchs & Jahreiß 1998). However, one should bear in mind that, if a disc component were to be included, this should lead to a contraction of the halo, and hence to a larger density, especially on the galactic disc plane (e.g. Debattista et al. 2008; Abadi et al. 2010). The highest local density is found for the halo Aq-B-2. This halo has a low stellar mass compared to the haloes Aq-D-2 and E-2, but is much more centrally concentrated than Aq-D-2. In comparison to Aq-E-2, which is also very centrally concentrated, Aq-B-2 has a strongly prolate shape, which explains its much higher value of ρ₀ on the major axis, where our ‘solar neighbourhood’ sphere is placed.

Figure 2.

Distribution in the V_R versus V_ϕ (top panels) and V_R versus V_Z (bottom panels) space of the stellar particles located inside a sphere of 2.5 kpc radius at 8 kpc from the galactic centre. The different colours represent different satellites. The corresponding host halo is indicated at the lower left-hand corner of each panel.

Open in new tabDownload slide

From Fig. 2 we can also appreciate the very large amount of substructure, coming from many different objects, that is present in these volumes. This substructure, in the form of stellar streams, can be seen as groups of unicoloured star particles. Recall that these stellar haloes are built up in a fully cosmological scenario. Therefore, effects such as violent variation of the host potential due to merger events, and chaotic orbital behaviour induced by the strongly triaxial dark matter haloes (Vera-Ciro et al. 2011), are naturally accounted for in these simulations. Nonetheless, these physical processes have not been efficient enough to completely erase the memory of the origin of the stellar halo particles. Note that, for this analysis, we have treated all tagged particles equally. As described in Section 2.3, the stellar-to-total mass ratio varies from particle to particle. Thus, the relative masses of the identified streams as judged by their particle number could be quite different from their relative stellar mass. We will explore this further in the following section.

Figure 3.

Distribution in E versus L_z space for the same set of particles as shown in Fig. 2. The corresponding host halo is indicated at the top left-hand corner of each panel.

Open in new tabDownload slide

Figure 4.

Distribution in V_az versus V_ΔE (top panels) and ν versus V_ΔE (bottom panels) space for the set of particles shown in previous figures. The corresponding host halo is indicated at the lower right-hand corner of each panel.

Open in new tabDownload slide

4 QUANTIFICATION OF SUBSTRUCTURE IN SOLAR VOLUMES

In the previous section, we found that the phase-space distribution of stellar particles inside a ‘solar neighbourhood’ sphere is very rich in substructure in all of our haloes. We will now quantify the number of stellar streams crossing the ‘solar neighbourhoods’ of our Aquarius haloes and characterize their evolution in time. Our goal is to compare the results of this analysis with previous studies (Helmi & White 1999; Helmi et al. 2003) that have analytically estimated this number and concluded that the merger history of the Milky Way could be recovered from the phase-space distribution of solar neighbourhood stars. Furthermore, we will address what fraction of the stellar particles that appear to be smoothly distributed in phase space are actually in streams that could not be resolved due to the high but limited particle resolution of our N-body simulations. We will also establish the importance of chaotic mixing for the quantification of substructure.

4.1 Resolved substructure

To quantify the number of streams inside our spheres, we first need to specify how to identify them. Our definition of a stream must ensure that particles in the same stream share the same orbital phase and progenitor. Following Helmi et al. (2003), a stream is identified when (i) two or more particles from the same parent satellite are found within one of our 2.5 kpc solar neighbourhood spheres at redshift z = 0 and (ii) they share at least one particle from the same parent satellite that has never been separated by more than a distance r_stream from either of them.² The value of r_stream is set to account for the stellar streams that, due to numerical resolution, are poorly resolved in the ‘solar neighbourhood’. In this work, we adopt r_stream = 0.8 × r_apo for each particle, where r_apo is its apocentre at z = 0, estimated from the outputs of the N-body simulations. Varying the value of r_stream between 0.6 × r_apo and 0.9 × r_apo did not significantly affect our results. This is required because a parent satellite can contribute with multiple stellar streams that may cross each other in configuration space (see, e.g. Helmi & White 1999). Note that using r_apo as a scale allows us to naturally adapt the algorithm to different orbital configurations.

In practice, we proceed as follows:

• For each particle inside our solar neighbourhood sphere of 2.5 kpc radius at redshift z = 0, we identify its parent satellite.

• We measure the time t_form (prior to the time of accretion) when the spatial extent of all particles associated with this parent takes its smallest value (as measured by the rms dispersion of distance from the centre of mass of the satellite).

• At t = t_form we identify all particles located in a sphere of 4 kpc radius, centred on the selected particle. This radius has to be large enough to initially contain all potential neighbouring members of this given particle. We have performed various tests and found the results to be robust to the sphere's extent.

• We follow these particles forward in time until z = 0, and discard those that, at any time, are more distant than r_stream from the selected particle.

• The selected particle is considered to be in a stream if, at the final time, another particle from the same parent satellite can be found within our solar neighbourhood sphere and, moreover, they have in common at least one of the original neighbouring particles lying within their respective r_stream.

Fig. 5 shows, as an example, the time evolution of one of the streams crossing the ‘solar neighbourhood’ of the halo Aq-A-2, at z = 0. The corresponding progenitor is the second most massive contributor to the stellar halo and, prior to accretion, it had a stellar mass of ≈1.25 × 10⁸ M_⊙. Its redshift of infall is z ≈ 2.3 and by z = 0 it has been fully disrupted. The stream is indicated with an arrow in Fig. 6. The black dots correspond to the particles that have always been neighbours (i.e. within r_stream) of a reference particle in this stream, which is indicated with a red dot. This figure clearly shows the full spatial extent of the stream, which probes regions far beyond the ‘solar neighbourhood’ volume.

Figure 5.

Time evolution of a stellar stream that crosses the ‘solar neighbourhood’ sphere of the Aquarius halo A-2 at redshift z = 0. The green dots show all the stellar particles from the corresponding parent satellite. With black dots we show the particles that have always been within r_stream from the reference particle (red dot), which is indicated with an arrow in Fig. 6.

Open in new tabDownload slide

Figure 6.

Distribution in V_R versus V_ϕ space of the stellar particles located inside a sphere of 2.5 kpc radius centred at 8 kpc from the galactic centre. The corresponding host halo is indicated at the lower left-hand corner of the top panels. The top panels show the distribution of particles that were linked to stellar streams (resolved component), whereas the bottom panels the remaining particles (unresolved component). As in Fig. 2, different colours indicate different satellites. The arrow in the top left-hand panel of this figure indicates the stream that is shown in Fig. 5.

Open in new tabDownload slide

The numbers of streams with at least two particles found inside the ‘solar neighbourhood’ spheres located at 8 kpc from the galactic centre are listed in Table 3. These numbers are well in the range of the ∼300–500 stellar streams around the Sun predicted by the models of Helmi & White (1999) (see also Gómez et al. 2010). The halo Aq-A-2 has the smallest number of resolved streams, as well as the lowest fraction of particles associated with them (see Table 3), and the largest associated particle fraction is found in the halo Aq-B-2. Table 3 also shows the total number of accreted particles found within each sphere as well as the total number of progenitor satellites that have contributed them. Note that the fraction of stellar mass in resolved streams is always larger than the fraction of star particles in resolved streams. Thus, the star particles with high stellar-to-total mass ratio, which come from the denser and inner parts of the original satellites, are more likely to be found in resolved streams.

From this table, we see that disrupted satellites contribute many more star particles at 8 kpc in Aq-B-2 than in the other haloes, and as a result the streams are much better resolved. The next best resolved ‘solar neighbourhood’ is found in Aq-D-2. For the other haloes, a significant fraction of the particles have not been associated with any substructure and therefore the total number of streams (which would include those unresolved) could be significantly larger. We will explore this further in Section 4.2.

In Fig. 6 we compare the velocity distributions of stellar particles in streams (top panels) to those not associated with any structure (bottom panels) according to our algorithm. Different columns correspond to distributions of stellar particles inside ‘solar neighbourhoods’ of different haloes. As in Fig. 2, the different colours represent different satellites. In addition to a very large number of inconspicuous substructures, we can see that the most prominent streams have been recovered in all cases. Note that the particles in streams tend to populate the wings of the distribution. In general, particles in the core of the velocity distribution will have been accreted earlier and tend to have shorter periods (and hence to mix faster). Therefore, streams will be more difficult to identify because of the limited resolution. In addition, and as explained below, it is possible that due to the triaxiality of the dark matter potential some of these orbits exhibit chaotic mixing. The associated streams would then be too convoluted to be observable. Note, however, that thanks to the much larger number of star particles and the smaller velocity dispersion (see Section 3), a significant fraction of the particles in the core of the distribution for Aq-B-2 are in resolved streams. Fig. 7 shows the distribution of stellar particles in streams in integrals of motion space. Substructure in this space is clearly visible. We emphasize that this is the case even in the presence of violent variation of the host potential due to merger events and chaotic orbital behaviour induced by the strongly triaxial dark matter haloes.

Figure 7.

Distribution in E versus L_z space of the stellar particles that were linked to stellar streams (resolved component) located inside a sphere of 2.5 kpc radius centred at 8 kpc from the galactic centre. The corresponding host halo is indicated at the lower left-hand corner of the top panels. The colour coding is the same as in the top panel of Fig. 6.

Open in new tabDownload slide

In the top panel of Fig. 8, we explore how the number of streams found inside spheres of 2.5 kpc radius changes as a function of galactocentric distance along the direction of the major axis. It is interesting to see that in all cases the number of streams decreases as a power law with radius. Note, however, that the haloes Aq-B-2 and E-2 have steeper profiles than the other haloes. As described by Cooper et al. (2010), the majority of the stellar mass in these two cases comes from one or two objects that deposit most of their mass in the inner regions of the haloes. Interestingly, these two haloes have the largest fraction of streams coming from the five most significant stellar mass contributors (see Table 3). As shown in the bottom panel of Fig. 8, the same behaviour is observed in the local density profiles of these haloes, ρ₀(r), measured along the major axis. Note that, due to the triaxial nature of the resulting haloes, varying azimuthally the location of our spheres results in local stellar densities that are, in general, an order of magnitude smaller than the observed value in the solar neighbourhood. For example, at a distance of 8 kpc along the intermediate axis, the halo Aq-A-2 has ρ₀ = 4.9 × 10³ M_⊙ kpc³ and a total number of streams N_streams = 37. An isodensity contour defined by the value of ρ₀ at 8 kpc along the major axis places the solar neighbourhood like sphere at a distance of ≈4 kpc on the intermediate axis. At this location N_streams raises to 68, indicating that azimuthal variations in the number of streams mainly reflect changes in local stellar density.

Figure 8.

Top panel: number of streams as a function of galactocentric radius found within spheres of 2.5 kpc radii. Bottom panel: spherically averaged stellar density profiles. The different colours indicate different haloes.

Open in new tabDownload slide

4.2 Smooth/unresolved component

In addition to the stellar particles associated with streams with at least two particles, our ‘solar neighbourhood’ spheres contain an important number of stellar particles that are not linked to any substructure, i.e. they appear to be smoothly distributed in phase space. The reason for this could be either of dynamical origin or, as posited before, simply due to the numerical resolution of the simulation. It is well known that phase-space regions of triaxial dark matter haloes may exist that are occupied by chaotic orbits (e.g. Vogelsberger et al. 2008). On such orbits, stellar particles that are initially nearby diverge in space exponentially with time. As a result, substructure that is present in localized volumes of phase space left by accretion events is rapidly erased. This process is known as chaotic mixing. On the other hand, initially nearby particles in phase space on regular orbits diverge in space as a power law in time. Therefore, the time-scale in which these particles fully phase-mix may be much larger.

Based on the previous discussion, we will now analyse the evolution in time of the local (stream) density around the particles found at the present day in each ‘solar neighbourhood’ sphere. Our goal is to assess the influence of chaotic mixing on the underlying phase-space distributions and to estimate the fraction of particles that, due to the high but nonetheless limited numerical resolution of our simulations, were not associated with any stream.

As in Section 4.1, we place a sphere of 4 kpc radius around each particle at t = t_form and tag all the surrounding neighbours. We track these particles forward in time until redshift z = 0 and those neighbours that depart from our particles by more than r_stream at any time are discarded. We compute the local spatial density of the selected particle, ρ(t), by counting the number of neighbours within a distance r_dens = 0.5 × r_apo. In Fig. 9, we show the time evolution of the spatial density of particles from two different satellites that contribute to the ‘solar neighbourhood’ of the halo Aq-A-2. From this figure, we can appreciate that the rate at which the spatial density of a stream decreases with time varies strongly (likely as a consequence of the different orbits). We can also observe that in many cases, when a particle becomes unbound from its parent satellite, the density of the newly formed stream exhibits initially a very rapid and quasi-exponential decrease. As shown by Gómez (2010) (see also Helmi & Gomez 2007), this initial transient is also present in regular orbits and does not necessarily imply a manifestation of chaotic behaviour. It is therefore necessary to follow the evolution of the stream's density for long time-scales (much longer than a crossing time) to disentangle regular from chaotic behaviour.

Figure 9.

Time evolution of the density of streams originating in two different satellites that contribute, at redshift z = 0, to the solar neighbourhood sphere obtained from the halo Aq-A-2. Each panel corresponds to a different satellite. The dashed lines show the power-law fits applied to the particles identified as stream members, as explained in Section 4.2, whereas the dot–dashed lines correspond to the remaining particles.

Open in new tabDownload slide

Fig. 9 shows the results of applying such a fit to the stream densities as a function of time. Since we are interested in the behaviour of the density on long time-scales, we do not include the initial quasi-exponential transient in the fits. We proceed as follows:

• We define the formation time of a stream, t_stream, as the time when its density decreases from one snapshot to the next one by 50 per cent.

• Starting from the snapshot associated with t_stream, we iteratively fit 10 times a power-law function to the stream's density by increasing in each iteration the starting time snapshot by snapshot. The fits are always normalized to the value of the density at the starting point. Furthermore, each data point is weighted according to the number of particles used to measure the corresponding density.

• We estimate the goodness of each fit by computing the rms error (RMSE), defined as

  \begin{equation} \displaystyle {\rm RMSE} = \sqrt{\frac{\sum _{i=1}^{m}\left(\hat{\rho }(t_{i}) - \rho (t_{i})\right)^{2}}{m-2}}, \end{equation}

  (5)

  where ρ(t_i) is the density estimated from the particle count and m is the number of data points used in the fit.

• Of all these experiments, we only keep the value of n obtained from the fit where the RMSE is smallest.

The results of applying this procedure to all the particles inside 2.5-kpc-radius spheres of different haloes are shown in Fig. 10. We explore spheres located at three different galactocentric radii, namely 8, 15 and 20 kpc. The black solid histograms show the distribution of n-values for all the particles in each ‘solar neighbourhood’ sphere. At R = 8 kpc (top panels), with the exception of the haloes Aq-B-2 and E-2, we find narrow distributions peaking at values slightly larger than n = 3. Note, however, that as we move outward from the galactic centre, the peak of all distributions tends to shift towards values much closer to n = 3. A clear example is the halo Aq-B-2, where the peak of the distribution shifts from n ≈ 7 at R = 8 kpc to n ≲ 3 at R = 20 kpc (bottom panel).

Figure 10.

The black solid histograms show the distribution of n-values found for particles in ‘solar neighbourhood’ spheres of 2.5 kpc radius. Each column corresponds to a different satellite. Different rows correspond to volumes centred at different galactocentric radii. The black (respectively grey) dashed histograms correspond to those particles that are (not) associated with resolved stellar streams. The vertical red dashed line is located at n = 3. For regular orbits, n ≤ 3 is expected. Each panel shows the median values |$\tilde{n}$|⁠, obtained from the distributions for particles in resolved and in unresolved streams.

Open in new tabDownload slide

Fig. 10 also shows in black dashed histograms the distribution n-values for those particles associated with ‘resolved’ stellar streams, while the grey-dashed histograms correspond to the remaining particles. From this figure, we can observe that although the two distributions overlap significantly, there is a clear offset between them. In each panel, we show the median values of n, |$\tilde{n}$|⁠. Stellar particles in streams (which are more likely to be on regular orbits) tend to have the smallest values of n and, in general, their distribution peaks at a value much closer to n = 3 than their unresolved counterparts. Let us recall that this is the value of n expected for the dependence of density on time for regular orbits in a 3D time-independent gravitational potential (see, e.g. Helmi & White 1999; Vogelsberger et al. 2008).

It is interesting to compare the n-distributions obtained from the haloes Aq-A-2 and B-2 at R = 8 kpc. As previously discussed, while the halo B-2 has the largest number of streams and the greatest fraction of accreted particles associated with them, the opposite is found for the halo Aq-A-2. The value of |$\tilde{n} \approx 7$| obtained from the resolved component in the halo B-2 suggests that these particles are undergoing chaotic mixing. Nonetheless, thanks to the high particle resolution, substructure can be efficiently identified. Note the large difference between the medians obtained from the resolved and unresolved components. For the halo Aq-A-2, the distribution of n-values associated with the two components shows a very similar value of |$\tilde{n} \approx 4.3$|⁠, considerably closer to what is expected for regular orbits. This suggests that a significant fraction of the particles that have not been associated with streams are on (nearly) regular orbits, and that because of the limited numerical resolution, they are not found in (massive) streams.

Interestingly, we find that, on average, particles in resolved streams in Aq-A-2 were released from their parent satellite³ 2.5 Gyr later and have more extended orbits than the particles in the unresolved component. The average apocentric distances for these two subsets are 32 and 12 kpc, respectively. Thus, the shorter orbital periods and longer times since release partially explain why the latter are not found to be part of (massive) streams. A lower limit to the fraction of accreted particles in resolved stellar streams, corrected by resolution effects, can be obtained by quantifying the fraction of particles in the unresolved component with n-values smaller than the median of the resolved component, |$\tilde{n}_{\rm res}$|⁠. We find in the halo Aq-A-2 a total of 411 particles that satisfy this condition. In addition to the 283 particles associated with the resolved component, we estimate that at least 50 per cent of all accreted particles in this halo should be part of a resolved stellar stream. Furthermore, the corresponding stellar mass fraction should at least be as large as 65.7 per cent.

In general, particles in resolved streams were released later (than those in the unresolved component) for all haloes. We find for the haloes Aq-B-2, C-2, D-2 and E-2 a time difference of 1.9, 1, 2.2, and 2 Gyr, respectively. However, differences in the mean apocentric distances are found only in the halo Aq-A-2. As for the halo Aq-A-2, we can obtain for the remaining haloes an estimate for the lower limit of particles in resolved stellar streams corrected by resolution effects. We find a total of (380, 348, 268, 350) particles with n-values smaller than |$\tilde{n}_{\rm res}$| in haloes (B-2, C-2, D-2, E-2), respectively. Thus, we estimate that the corresponding fractions of particles and stellar mass in resolved streams for these haloes should be, at least, as large as (55, 75, 62, 84 per cent) and (96, 80, 91, 84 per cent), respectively.

5 SUMMARY

In this work, we have characterized the phase-space distribution as well as the time evolution of debris for stars found inside ‘solar neighbourhood’ volumes. Our analysis is based on simulations of the formation of galactic stellar haloes, as described by Cooper et al. (2010). The haloes were obtained by combining the very high resolution fully cosmological ΛCDM simulations of the Aquarius Project with a semi-analytic model of galaxy formation.

We find that, for the haloes Aq-A-2, -C-2 and -D-2, the measured velocity ellipsoid at 8 kpc on the major axis of the halo is in good agreement with the estimate for the local Galactic stellar halo (Chiba & Beers 2000). Similarly, for all haloes but Aq-B-2, the local stellar halo average density on the major axis is in reasonable agreement with the observed value (Fuchs & Jahreiß 1998). However, we note that if the contribution of the disc and bulge were to be included, the resulting velocity dispersions as well as the local stellar halo average density would be significantly increased. This could suggest that these dark haloes are too massive to host the Milky Way (in agreement with Battaglia et al. 2005; Sales et al. 2007; Smith et al. 2007; Wang et al. 2012; Vera-Ciro et al. 2013). On the other hand, if located in a prolate low-mass halo such as Aq-B-2, the measured local stellar density would imply that the Sun should be on the intermediate or minor axis, i.e. the Galactic disc's angular momentum would be aligned with the major axis of the dark halo (Helmi 2004). On the minor axis, the haloes B and E have velocity ellipsoids of smaller amplitude and hence could perhaps be made consistent with those observed near the Sun, if the effects of the contraction induced by the disc were taken into account. However, for the haloes A, C and D, the dispersions on the minor axis are larger than those on the major axis at the equivalent of the solar circle.

In agreement with previous work, we find that 90 per cent of the stellar mass enclosed in our ‘solar neighbourhoods’ comes, in all cases, from three to five significant contributors. These ‘solar neighbourhood’ volumes contain a large amount of substructure in the form of stellar streams. In the five analysed stellar haloes, substructure can be easily identified as kinematically cold structures, as well as in well-defined lumps of stellar particles in spaces of pseudo-conserved quantities.

By applying a simple algorithm to follow the time evolution of local densities in the neighbourhood of reference particles, we have quantified the number of (massive) streams that can be resolved in solar neighbourhood like volumes. In haloes where solar neighbourhood like volumes are best resolved, we find that the number of resolved streams is in very good agreement with the predictions given by Helmi & White (1999). As we move outward from the galactic centre, the number of streams decreases as a power law, in a similar fashion to the density profiles.

To explore whether the halo-to-halo scatter in the total number of resolved streams is due to poor particle resolution or to chaotic mixing, we have estimated the rate at which local (stream) densities decrease with time. Interestingly, we find that in our best resolved solar neighbourhood like volume (Aq-B-2) a large amount of substructure can be identified, even in the presence of chaotic mixing. In this halo, up to 82 per cent of the star particles (and 92 per cent of the stellar mass) can be associated with a resolved stream. On the other hand, our poorest resolved solar neighbourhood like volume (Aq-A-2) shows a smaller amount of substructure but significantly longer mixing time-scales (a rate that is close to what is expected from regular orbits). In this case, only 20 per cent of the accreted star particles (and 31 per cent of the accreted stellar mass) can be linked to streams. A comparison of the orbital properties of resolved and unresolved stream populations shows that, in general, particles in resolved streams are released later and have longer orbital periods. Moreover, diffusion in both populations occurs at very similar rates. Thus, this analysis suggests that our strongest limitation to quantifying substructure is mass resolution rather than diffusion due to chaotic mixing. We have estimated a lower limit to the fraction of accreted particles in resolved stellar streams by correcting for resolution effects. We find that this fraction should be, at least, as large as 50 per cent by number and 65 per cent by stellar mass in all haloes.

The results presented in this work suggest that the phase-space structure in the vicinity of the Sun should be quite rich in substructure. Recall that our analysis is based on fully cosmological simulations and therefore chaos and violent changes in the gravitational potential that are so characteristic of a cold dark matter model are naturally included.

With the imminent launch of the astrometric satellite Gaia (Perryman et al. 2001), we will be able for the first time to robustly quantify the amount of substructure in the form of stellar streams present within ≈2 kpc from the Sun. This mission from the European Space Agency will provide accurate measurements of positions, proper motions, parallaxes and radial velocities of an unprecedentedly large number of stars.⁴ In addition, recent and ongoing surveys such as SEGUE (Yanny et al. 2009), RAVE (Steinmetz et al. 2006), LAMOST (Cui et al. 2012), HerMES (Barden et al. 2010) and Gaia-ESO (Gilmore et al. 2012) will allow us to complement these Gaia observations. A direct comparison of the degree of substructure found in the solar neighbourhood with that found in numerical models such as those presented here will allow us to establish directly, for example, how much of the stellar halo has been accreted and how much has been formed in situ.

The authors would like to thank the anonymous referee for his/her very useful comments and suggestions that helped improve the clarity of this paper. FAG would like to thank Brian W. O'Shea for providing useful comments and suggestions and Monica Valluri for valuable discussions. FAG was supported through the NSF Office of Cyberinfrastructure by grant PHY-0941373 and by the Michigan State University Institute for Cyber-Enabled Research (iCER). This work was initiated thanks to financial support from the Netherlands Organisation for Scientific Research (NWO) through a VIDI grant. AH acknowledges the European Research Council for ERC-StG grant GALACTICA-240271. APC acknowledges support from the Natural Science Foundation of China, grant No. 11250110509.

REFERENCES