Subhalo Spreading of Thin Tidal Star Streams

Thin tidal streams are powerful probes of the shape and the internal substructure of the dark matter which dominates the mass of a galactic halo (Lynden-Bell & Lynden-Bell 1995; Morrison et al. 2000; Ibata et al. 2002; Johnston et al. 2002). The tidal field of the galaxy pulls stars away from the outskirts of a globular cluster, releasing the stars to orbit freely in the galaxy. Energy exchange between stars within a cluster (Chandrasekhar 1942; von Hoerner 1957; Hénon 1961) boosts inner cluster stars to orbit out to the cluster outskirts, where tidal heating at pericenter passages (Aguilar & White 1985; Binney & Tremaine 2008; Meiron et al. 2021) continues to drive mass loss from the cluster. The locations and velocities of the stars along the stream provides a snapshot of the potential along the stream's path from which the mean mass profile of the galactic halo and its large-scale asymmetries can be derived (Binney 2008).

The subhalo mass function in a galactic halo is a power law for cold dark matter (Klypin et al. 1999; Moore et al. 1999), N(>M) ∝ M^−0.9 (Springel et al. 2008). On the other hand, warm dark matter has a particle free-streaming cutoff scale (Bond & Szalay 1983), which depletes the numbers of subhalos below masses ≲ 10⁹ M_☉ (Angulo et al. 2013; Benson et al. 2013; Lovell et al. 2014), with details depending on the dark matter particle properties. The sizes and masses of the subhalos set the scale size and velocity of the perturbation to a stream. The volume density of the substructures in the galactic halo determines the rate of encounters with the stream (Ibata et al. 2002; Johnston et al. 2002; Carlberg 2012). Measurements of stream velocity dispersions then offer the prospect of a quantitative measurement of the subhalo numbers.

Finding tidal streams immersed in the Galactic background of unrelated stars requires well-calibrated wide-field photometric data (Irwin et al. 1994; Gunn et al. 1998; York et al. 2000) and color–magnitude diagram filtering techniques (Grillmair et al. 1995; Odenkirchen et al. 2001; Rockosi et al. 2002), leading to the discovery of a whole range of tidal streams from globular clusters and dwarf galaxies (Ibata et al. 1994; Belokurov et al. 2006; Grillmair & Dionatos 2006). The number of streams has now reached about 100 and continues to grow (Mateu 2022) with the Gaia Mission (Gaia Collaboration et al. 2016) and spectroscopic facilities providing proper motions and radial velocities (Cooper et al. 2022), along with colors and spectra that can be used to distinguish stream stars from the field. Individual streams have been examined for density variations, i.e., gaps, along a stream, which provide evidence that low-mass 10^6–8 M_☉ dark subhalos are present (Yoon et al. 2011; Carlberg 2012) particularly in the exceptionally long and thin GD-1 stream (Grillmair & Dionatos 2006; Koposov et al. 2010; Carlberg & Grillmair 2013; Bonaca et al. 2019). Measurement of the properties of a population of streams complements the detailed study of individual streams.

This paper analyzes simulations of sets of streams from star clusters in realistic halos to measure the evolution of the kinematics of streams and characterize their current-time properties. The elements of stream spreading have been studied in simulations of ever increasing complexity, from a subhalo population in a spherical static halo (Carlberg 2009), a subhalo population in an aspherical static potential (Ngan et al. 2015, 2016; Qian et al. 2022; Ferrone et al. 2023), and here the subhalos found in dynamically assembled halos which develop different subhalo numbers. The simulation setup is described in Section 2. Section 3 discusses the kinematic properties of the streams. The basic dynamics of the streams are reviewed in Section 4. The spreading of the streams with time is measured in Section 5. Section 6 examines the role of subhalos in the time evolution of individual streams. Predictions for the velocity distribution, both over the visible length and with distance along a set of streams, are given in Section 7. Section 8 summarizes the results and some implications for observational studies.

The simulations are started with the dark matter distribution used in Carlberg & Keating (2022), which was derived from the initial conditions of FIRE model m12i (Hopkins et al. 2014). The 70.514 million dark matter particles, each having a mass of 3.518 × 10⁴ M_☉, are assigned a softening of 40 pc. The early-time particle distribution is evolved to 0.8 Gyr at which time the dark matter halos in the simulation were identified using the Amiga Halo Finder code (Gill et al. 2004; Knollmann & Knebe 2009) as host dwarf galaxies for a population of high-redshift, low-metal-abundance globular clusters, which are the progenitors of most thin streams (Martin et al. 2022). Globular-cluster-like star particle clusters are constructed using a W = 7 King model (King 1966) with masses randomly drawn from the range (4–20) × 10⁴ M_☉ with an N(M) ∝ M⁻¹ distribution. The clusters are placed in randomly oriented, rotating, exponential disks in the inner 0.5 kpc of the dark matter halos, which have a typical radius of the peak circular velocity of a few kiloparsecs. The total mass of star clusters in a halo is limited to 10⁻⁴ times the mass of the halo (Hudson et al. 2014) for halos more massive than 5 × 10⁸ M_☉. The procedure generates 4.718 million star particles in 271 star clusters placed in 133 dark matter halos. The outer radii of the clusters are started with the tidal radius at the radial positions of the clusters. The clusters quickly evolve into equilibrium as they orbit in their initial local subhalo. The star particles have masses of 4 M_☉ with an assigned softening of 2 pc. The resulting star clusters have a typical half-mass radius, r_h , of 4 pc, Figure 1. The star particles are accurately followed from the half-mass radius through the tidal zone out into the stream. The softening allows reasonable time steps in the simulation and reduces the relaxation time into a regime that can be simulated accurately. The code includes a cluster star relaxation process using a Monte Carlo velocity perturbation scheme (Carlberg 2018). The internal heating algorithm is calibrated to precise Nbody6 (Aarseth 1999) runs. The mass-loss rate of a cluster depends on its mass, which increases with decreasing cluster mass (Figure 1) and the tidal field along its orbit, which on average is decreasing with time as accretion builds up the main halo.

Zoom In Zoom Out Reset image size

The largest dark matter halo at 4.8 Gyr is augmented with a central bulge and a growing disk potential to create a Milky Way–like main halo gravitational potential. A galactic disk potential depletes the subhalos in the inner halo (Garrison-Kimmel et al. 2017). The disk also introduces an axis of symmetry which reduces the triaxiality of the inner halo. The dark mass of the main dark matter halo assembles very quickly, with 70% of the mass of the main halo in place within 5 Gyr. The disk of the Milky Way grows at a much slower rate than the halo mass assembly, which we approximate as linear in time, both in mass and in cylindrical radius (Frankel et al. 2019). The final mass of the main halo in the simulation is 5.7 × 10¹¹ M_☉ within 100 kpc, slightly below the (6.1–6.9) × 10¹¹ M_☉ of recent measurements (Deason et al. 2021; Shen et al. 2022). We start the disk growing at time 4.8 Gyr, shortly after a major merger and the beginning of the phase of relatively slow accretion growth of the main halo. A particle of mass 0.05 × 10¹⁰ M_☉ with a softening of 0.5 kpc is placed at the center of the halo. This particle is a lighter version of the nuclear bulge (Launhardt et al. 2002), which is not important to the dynamics at a few kiloparsecs but avoids small time steps near the center. A Miyamoto & Nagai (1975) disk of mass 3.4 × 10¹⁰ f(t) M_☉, with parameters a = 2.0f(t) kpc and b = 0.28 kpc, is centered on the bulge particle. The growth of both mass and disk scale length with time is 0 before 5 Gyr, after which it is f(t) = (t − 5)/9.1 to a final time of 14.1. The final disk is 3.4 × 10¹⁰ M_☉, slightly more than the 3.2 × 10¹⁰ M_☉ (combined thin and thick disks) of Zhou et al. (2023), but half the MW2014 model disk (Bovy 2015). The baseline and double disk models give 8.1 kpc circular velocities of 195 and 232 km s⁻¹, respectively. Three simulations (Table 1) are run, one with no disk and one with the Bovy (2015) disk (double our baseline), to generate a range of subhalo numbers in the inner halo and disk influence on streams.

Figure 2 shows the mass functions of the subhalos in the main halo at the end of the simulation for several radial bins. The mass distributions are very near a power law of N(>M) ∝ M^−0.9 (Springel et al. 2008). The solid lines are for the baseline model with a disk. The dashed line is for a model started from the initial dark matter distribution, but with no added disk. The dotted lines are for a model with a disk twice the mass of the baseline model with the disk scale length boosted from 2 to 3 kpc to maintain a similar surface density in the two runs. Figure 3 shows the numbers of inner halo subhalos in the 10^7–8 M_☉ mass interval and between 10 and 60 kpc. The simulations with a disk deplete the numbers of inner subhalos relative to a disk-less simulation (Garrison-Kimmel et al. 2017). Accretion brings in new subhalos at a declining rate of addition in a noisy process.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The graininess of the n-body simulation causes the dark matter particle population to increase the velocity dispersion of the stream star particles slightly. The process is a random walk so the velocities increase approximately in proportion to the square root of time (Chandrasekhar 1942). Peñarrubia (2018, 2019) has shown that the tidal field of the nearby particles dominates the increase in velocity dispersion of small coherent structures such as star clusters and streams. The effects are coherent over the interparticle spacing of the dark matter particles, which is in the range of 0.1–1 kpc for the radial range of 3–100 kpc in our simulations. Streams respond to particle passages with coherent velocity changes over these length scales, which over several orbital times develop into density variations on similar scales in the streams (Banik & Bovy 2021). We undertake a special purpose simulation to measure the particle-noise-driven velocity rise for our simulation. The particles from the final time main halo in the baseline simulation are placed at their radial locations but randomly rotated on a sphere to symmetrize and smooth the mass distribution. A ring of 30,000 star particles on circular orbits is started at 20 kpc, in one simulation with a ring at the local circular velocity and no random velocities, and another with Gaussian perturbations in radial velocity. The results are shown in Figure 4. The simulations show that the angular momentum is conserved to about 1.1% over 10 Gyr, and the radial and tangential velocities rise from zero to about 3 km s⁻¹ over 10 Gyr. Increasing the softening of the dark matter particles to as much as 1 kpc provides only a modest reduction in the particle heating, as expected from the logarithmic behavior of two-body scattering. The particle-noise-induced velocity changes are below the perturbations from subhalos for these simulations but will add noise to streams as cold as 1 km s⁻¹.

Zoom In Zoom Out Reset image size

The Gadget-4 code (Springel et al. 2021) evolves the simulation. The code is augmented to include the slow heating of star particles in clusters with the algorithm of Carlberg (2018). The gravitational field of the growing galactic disk is added to the gravitational force. The small sizes of the star clusters require small time steps, typically 5000–10,000 years, so the simulation requires more than a million time steps to complete. Time is reported in computational gigayear time units, where 1 unit is 1.022 true gigayears.

The x–y projections of the stars released from the clusters in the last 4 Gyr are shown in Figure 5. The streams in the inner 60 kpc have orbital shapes ranging from nearly circular to nearly radial, with a tendency for more radial streams at larger radii. The set of streams at (−70, 22, −11) on the left middle of the plot are in a subhalo of 1.45 × 10¹⁰ M_☉ which passed within 15 kpc and is now at 74 kpc, continuing to move away at 150 km s⁻¹. We focus much of our analysis on the region closer than 60 kpc, partly to avoid this structure and its clusters, which orbit within its subhalo and are not yet centered on the main halo. Known nondwarf galaxy streams have median distances within 50 kpc (Mateu 2022).

Zoom In Zoom Out Reset image size

The streams with progenitor clusters currently between 10 and 150 kpc are plotted using local stream coordinates in Figure 6. The inner limit at 10 kpc is to exclude streams that encounter the galactic disk, which completely disperses about half the streams. Each stream of Figure 6 has been rotated into its instantaneous great-circle frame, which uses the velocity of the cluster (a reference particle at the center of each cluster in the simulation) to define the equator of the orbital plane perpendicular to the cluster's angular momentum, and the current cluster position defines the zero of the stream longitude. Cluster IDs are shown on the left and their current distances on the right of Figure 6. The great-circle approximation transforms the somewhat jumbled distribution of streams of Figure 5 into orderly, nearly linear distributions. The inner streams have greater angular lengths and are more widely dispersed than those at large distances from the center; however, there is a range of stream lengths and widths at all radii. In the simulation setup, clusters were put into host halos in order of halo mass, so there are sets of clusters that have remained close in radius, such as 72–78 in the middle panel. The streams generally have a thin core, surrounded by a wide component (Carlberg 2018; Malhan et al. 2019).

Zoom In Zoom Out Reset image size

The average sky density of star particles (each 4 M_☉ here) perpendicular to the stream axis is shown in Figure 7 for the 10–30 and 30–60 kpc distance ranges. Near the progenitor cluster the streams have an average half-width at half-maximum of 05. The thin cores of the streams have a non-Gaussian tail around them extending beyond 10° in width. The asymmetries of Figure 7 are likely due to the small numbers of streams on average, with the streams at 30–60 kpc showing a tilt in the opposite direction. On average, streams become significantly wider with distance along the stream, which lowers the peak surface density of stars in the stream. The average surface density declines by approximately a factor of 5 from 1° to 32° along the stream length. For an individual stream the density along the stream varies with orbital phase, with streams being stretched in longitude but narrowed in latitude at their orbital pericenter. The density along each stream varies with orbital phase, reaching a minimum at pericenter, but being at a smaller orbital radius is often observationally easier to find. On the average, at an angular distance down the streams of 10° the mean stream surface density is about 10% of the initial stream and the distribution becomes much broader; nevertheless, most streams do have a thin core over some fraction of their length.

Zoom In Zoom Out Reset image size

This section uses the basic dynamics of streams to estimate the expected effects of the measured subhalos on the streams. Tides pull stars away from the cluster into orbits that are offset outward and inward from the cluster by approximately the Jacobi radius (Binney & Tremaine 2008, their Equation (8.14)),

$$\begin{eqnarray}&&{r}_{{\rm{J}}}=r{\left(\displaystyle \frac{{m}_{c}}{2M(r)}\right)}^{1/3}=r{\left(\displaystyle \frac{{{Gm}}_{c}}{2{v}_{c}^{2}r}\right)}^{1/3},\end{eqnarray} \tag{ 1 }$$

where m_c is the mass of the cluster, r is its galactic radius, M(r) is the total galactic mass inside r, and v_c is the circular velocity at r. The time-varying gravitational field of the simulations pulls star particles away from the clusters with a range of angular momenta. For this calculation, we make the approximation that stars are added to the leading and trailing streams with the same velocity as the cluster but offset radially inward and outward by the tidal radius. The stars pulled inward have less angular momentum than the cluster and will orbit at a higher angular rate, developing the leading tidal stream, and vice versa for the trailing stream. The stars that escape at zero relative angular velocity move away from the progenitor cluster at velocity Δv_s ≃ Ωr_J, where Ω = v_c /r,

$$\begin{eqnarray}&&{\rm{\Delta }}v\simeq {\left(\displaystyle \frac{{{Gm}}_{c}{v}_{c}}{2r}\right)}^{1/3}.\end{eqnarray} \tag{ 2 }$$

Stars leaving the cluster have a spread of velocities, σ_r , comparable to Δv, because the tidally induced velocities are proportional to the tidal acceleration, which varies linearly with distance from a cluster's center. Consequently, the leading edge of the stream will move down the stream at ∼2Δv. For typical values of m_c = 10⁵ M_☉, v_c = 200 km s⁻¹, and an orbital pericenter of r = 10 kpc, 2Δv = 3.2 km s⁻¹, which over 10 Gyr develops into a tidal tail of 16 kpc, for a total length of leading and trailing of 32 kpc. The stream width is ≃σ_r /κ ≃ r_J, where κ is the radial oscillation frequency, equal to the epicyclic frequency for small eccentricity orbits. Most stars are stripped from the outskirts of a star cluster at its orbital pericenter. Streams are thinnest and most visible as they pass near their orbital pericenter. The angular length of a stream should then be about 32/10 ≃ 3.2 radians, or 180°, that is, streams should typically wrap halfway around the sky. The four longest streams attributed to star clusters are 100°–140° long (GD-1, M68, Murrumbidgee, and NGC 3201) with Orphan at 230° long possibly having a dwarf galaxy progenitor (Mateu 2022). The decline of surface density and increasing width of streams along their length seen in observed streams (Bernard et al. 2016; Bonaca et al. 2020) is also seen in the simulations. That is, the simulated streams generally have members around their entire orbit; however, the stars most distant from the progenitor will be observationally difficult to find as their surface density declines relative to unrelated field stars.

Subhalos crossing the stream at velocity v change the velocities of stream stars. An approximate analysis of the resulting stream heating helps to understand the process. The time for lower-mass subhalos to cross a stream is so short that the impact approximation is satisfied, δ v ≃ 2GM_h /(bv), which we approximate as 2GM_h /(r_m v), where r_m is the radius of the peak of the circular velocity curve for a subhalo. The velocity change for the more massive subhalos, ≳10⁷ M_☉, which have virial radii over 4 kpc, will have increased dependence on the orbital details. Fitting to the subhalo peak circular velocity radii in the simulation we find a subhalo changes the velocities of the nearest stream stars by,

$$\begin{eqnarray}&&\delta v=1.3{\left(\displaystyle \frac{{M}_{h}}{{10}^{8}{M}_{\odot }}\right)}^{0.49}\left(\displaystyle \frac{300}{v}\right)\,\mathrm{km}\,{{\rm{s}}}^{-1},\end{eqnarray} \tag{ 3 }$$

with a 10% error on the velocity scale and 6% on the power coefficient. The power-law fit to the radius where the halo circular velocity peaks is ${r}_{m}={(1.4\pm 0.2)({M}_{h}/{10}^{8}{M}_{\odot })}^{0.48}$ kpc.

The rate at which the velocity dispersion rise is then (see Carlberg 2013 for details),

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}\langle \delta {v}^{2}\rangle }{{{dM}}_{h}\,{dt}}=\int \int \displaystyle \frac{4{G}^{2}{M}_{h}^{2}}{{b}^{2}{v}^{2}}\displaystyle \frac{{dn}({M}_{h},r)}{{{dM}}_{h}}{vf}(v)\pi \,{db}\,{dv},\end{eqnarray} \tag{ 4 }$$

where f(v) is the velocity distribution function, ∫f(v)dv = 1, in the frame of stream, and n(M_h , r) is the radial density distribution of subhalos of mass M_h . The virial radius is 3–6 times r_m , so we integrate to an impact parameter b = 3r_m , which is close to the virial radius,

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}\langle \delta {v}^{2}\rangle }{{{dM}}_{h}\,{dt}}\simeq 20\pi \displaystyle \frac{{G}^{2}{M}_{h}^{2}}{{r}_{m}v}n({M}_{h},r).\end{eqnarray} \tag{ 5 }$$

Combining the estimate of δ v of Equation (3) with the measured halo mass function of approximately n(M) ∝ M^−1.9 leads to a heating rate proportional to M^0.14.

The heavier subhalos have a slight margin in heating stream stars; however, their effect diminishes as their encounters with a stream become less likely to occur. When there are fewer than two subhalo stream crossings in some time interval a single subhalo stream crossing leads to coherent stream changes, not heating. The 10–60 kpc radial distance range dominates stream heating, with a total heating up to and including the two-crossing mass of 10^7.6 M_☉ (for 4 Gyr) of $\sqrt{\langle \delta {v}^{2}\rangle }=0.3{(t/\mathrm{Gyr})}^{1/2}$ km s⁻¹. Figure 8 averages the heating and stream crossing rates over the last 4 Gyr at intervals of 0.2 Gyr. The plotted error bars give the rms spread in the heating rate. About half of the heating comes from subhalos between 10⁷ and 5 × 10⁷ M_☉. We note that the simulated lower-mass subhalos may be artificially large and too easily disrupted (Errani & Navarro 2021), so the heating from lower-mass halos is likely somewhat underestimated.

Zoom In Zoom Out Reset image size

The calculations summarized in Figure 8 give a total heating of ∼1 km s⁻¹ over 10 Gyr to the two-crossing subhalo mass of M ∼ 10^7.6 M_☉, beyond which subhalo stream crossings are rare and the random quadrature sum is not valid. Hence, subhalo heating alone cannot by itself explain all of stream spreading. However, in a triaxial or more complex potential, which is generally the case, the large-scale potential asymmetries along the streams cause yet more orbital spreading (Carlberg 2015; Ngan et al. 2016).

A 10^7.5 M_☉ subhalo has r_virial ∼ 5 kpc, so the dominant subhalos will lead to coherent velocity changes over 10 kpc sections of a stream. The velocity perturbations will lead to orbital changes which will add to these velocities and mix the stars along the length a stream. Exactly how the velocity perturbations play out to measurable changes in the streams will be examined in the next sections.

The time at which a star leaves the cluster is measured in a simulation as the first time at which it is beyond some minimum distance from the cluster center. We use 0.2 kpc as the minimum distance, which is typically twice the largest tidal radius for clusters in the main halo. The age of the star in the tidal stream is the difference between the current simulation time and the time at which the star crossed the minimum distance radius. Figure 9 plots the sky positions of the stream stars in four age ranges from clusters that are in the radial range of 10–60 kpc with a minimum angular momentum of $\left|L\right|\gt 1000$ kpc km s⁻¹. The selected radial range reduces the numbers of streams that pass through the disk and those orbiting in infalling subhalos. The stream age range of the stars is indicated in the legend for each panel, with the youngest stars in the upper left and the oldest in the lower right.

Zoom In Zoom Out Reset image size

The length on the sky of the stream stars increases with their age in streams, as expected for stars that are being pulled into the leading and trailing tidal tails. The older stars are in longer stream segments, but with a considerable overlap with the younger stars, a consequence of the spread in angular momentum of the stream stars. After about 4 Gyr the stars dramatically spread out from the streams, with streams largely blurred out on the sky after 8 Gyr. The orbital perturbations are larger than the approximate calculation of Equation (5) and larger than seen in a static symmetric halo or aspherical potentials with subhalos (Ngan et al. 2016), indicating that it is the combination of subhalo heating and orbit dispersion in the aspherical and time-variable galactic halo potential.

The angular momentum of a star stream particle is an approximately conserved quantity, useful for analyzing stream evolution. The angular momentum of an individual particle oscillates around its orbit in the aspherical potential. Figure 10 shows the stream particle angular momenta, which is the magnitude of the angular momenta vectors, as a function of the sky longitude for the same stars and age bins as Figure 9. The blurred structures at large angular momenta are streams orbiting within subhalos (dwarf galaxies) that are themselves orbiting within the main halo. The near vertical streams at low angular momenta are being tidally disrupted in the inner galaxy.

Zoom In Zoom Out Reset image size

Streams at low angular momenta are usually closer to the halo's center, with a higher subhalo density, so blur much quicker than those at large angular momenta. Consequently, the thin parts of inner streams are generally younger than those at large radii. Stars of a range of ages in a stream overlap in azimuth, as is most visible for streams at larger angular momenta. As stars leave a cluster to join a stream with a range of velocities, they have a range of angular momenta, meaning that at any given angular distance from the progenitor, there is a range in star ages, with the mean age of star particles increasing with distance from the cluster.

Figure 11 shows the angular momenta of stream particles versus their stream ages for those with clusters in the 10–60 kpc radial range. The tidal force broadens the angular momentum distribution of a stream at every pericenter passage, which leads to a spurt of new particles added to a stream, spread in angular momentum. The low-angular-momentum clusters have shorter radial orbital periods than higher-angular-momentum clusters, which sets the interval between spurts of particle loss. The range of angular momenta relative to the cluster means that stream particles that newly join a stream close together in time will spread to a range of angular distances along the stream.

Zoom In Zoom Out Reset image size

The panels of Figure 11 show the angular momentum history for particles with the amount of time they have been in streams. The top panel shows the angular momentum for stars leaving the clusters over the course of the entire simulation and the bottom panel measures the angular momentum at time 9.8 Gyr with the 4.3 Gyr younger ages at that time. The figure demonstrates that stars leaving the clusters at earlier times have a similar angular momentum spreading in time as those at later times. The earliest 3 Gyr particles show a large spread resulting from the clusters orbiting within their early-time subhalos.

Figure 12 shows how the rms spread of the angular momentum of the stream particles evolves depending on the angular momentum of the stream. The streams at low angular momentum spread quickly with the streams at higher angular momentum, ∼10⁴ kpc km s⁻¹ and above, requiring the entire age of the simulation to be substantially dispersed. Streams at lower angular momentum orbit at smaller radii where the higher subhalo density leads to higher encounter rates.

Zoom In Zoom Out Reset image size

The particle-by-particle history of the spreading of the angular momenta of stream stars is shown in Figure 13 for two sets of particles emerging over short intervals of time, 0.2 and 0.4 Gyr, from two clusters (numbers 208 and 228 in Figure 6). The orbits traverse the radial range of approximately 17–33 kpc (cluster 208) and 12–60 kpc (cluster 228). Subhalos pass near and through streams at speeds of 100–200 km s⁻¹, inducing nearly instantaneous velocity changes in the particles of the stream. The fractional changes in angular momentum are typically only 2%–3% in a single encounter, with particles ahead of the point of closest encounter losing angular momentum and those behind gaining. After several orbits the small angular momentum differences develop into significant positional changes, shown in Figure 14 for cluster 268. The orbits spread a factor of two in radius and azimuth in about 2.5 Gyr, suggesting that the perturbed orbits may have become chaotic and exponentially diverging. Both subhalo perturbations and orbital spreading in the aspherical, time-varying halo potential play a role in spreading thin streams (Ngan et al. 2016). The small coherent oscillations of the angular momenta are the result of the motion of the bulge particle at the center of the halo, which is used to define the center of momentum of the system.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The mass of the subhalo that induces the angular momentum spike at time 11.78 Gyr in the lower panel of Figure 13 is readily estimated from the velocity change. The stream segment is at 50 kpc, so the peak velocity change is about 3 km s⁻¹. Equation (3) requires a subhalo of mass 5 × 10⁸(v/300)^2.04 M_☉. Figure 15 shows the distances of subhalos near the stream segment with time. The subhalo closest to the stream segment is at 3.2 kpc, moving at 326 km s⁻¹, for which Equation (3) gives 3.8 × 10⁸ M_☉, whereas the directly measured mass of the subhalo found in the simulation is 3.6 × 10⁸ M_☉ with a peak circular velocity of 14.9 km s⁻¹ at 1.7 kpc in the subhalo. The rate of subhalo encounters rises at pericenter, a consequence of the increased angular length of the stream (Figure 14). The small wiggles in the angular momenta are the result of the motion of the bulge particle that defines the center of mass.

Zoom In Zoom Out Reset image size

Figure 16 shows the rise in velocity dispersion with time along streams in the 10–30 kpc range and the 30–60 kpc range. Stream particles are binned in ages of 0.2 Gyr with the velocity dispersion of individual streams calculated from the mean subtracted variance. The line is the unweighted average of the stream's velocity dispersions. There are an average of about 100 particles per time bin per stream, with a large spread in the numbers. Streams at smaller orbital radii have a quicker rise in velocity dispersion with time because they are interacting with a higher density of subhalos, as expected from the calculations of stream heating, and shown in Figure 8. The age of a star in the stream is not an observational quantity, although there is a predicted variation in the stellar mass function with stream age (Webb & Bovy 2022). The distance of a star along the stream is a proxy for stream age, but as shown above, at any given distance there are stars with a wide range of stream ages, which blur the distribution.

Zoom In Zoom Out Reset image size

Measurement of stream spreading is a potentially powerful test of whether the amount of subhalos predicted in a cold dark matter cosmology are present in the Milky Way. Viable observational tests of stream spreading require the use of practical observational measurements with no knowledge of how long a star has been in the stream. A straightforward approach to a practical measurement is to measure the velocity dispersion with angular distance along the stream from the progenitor. Another approach is to sum the velocity distribution along the entire length of the visible stream, which does not require any knowledge of the progenitor cluster's location. We will concentrate on the range of 10–60 kpc where currently known streams have been found (Mateu 2022) and observational velocity measurements with an accuracy of a 1 km s⁻¹ are practical. Spectroscopic observations are important to measure the chemical abundance patterns characteristic of stars from the same globular cluster because streams spread to large angles from the orbital equator of the stream. For potentially practical measurements we use stream stars within 1° or 5° of the stream equator for our illustrative measurements.

7.1. Velocities with Stream Longitude

The distribution of velocity dispersions of stream stars is shown as a function of angular distance for the leading (positive angles) and trailing tidal arms for the eight streams in the baseline simulation within 10–30 kpc and with L > 2000 kpc km s⁻¹ in Figure 17. Lower-angular-momentum streams come close to the disk which provide a large-scale tide that further disperses the streams. The velocities are measured in spherical coordinates in each stream's great-circle frame, where the cluster angular momentum defines the orbital plane and its position the zero of stream longitude and latitude. The plots show stars within 5° of the orbital equator. Using stars only within 1° of the stream equator suppresses the rise in velocity dispersion beyond angles 10°–20° along the stream. Including streams stars out to a width 20° gives an ≃50% increase in the velocity dispersion in all components. The mean and standard deviations of the velocities are made in bins of 4° along the streams. Measurements in half the width and half the length give very similar results, but with doubled errors. In these measurements all orbital phases of streams are taken and there is no tracing of the high-density centerline of the stream with length.

Zoom In Zoom Out Reset image size

Figure 17 shows that the radial velocity dispersion is generally larger than the tangential or perpendicular velocity dispersions. The tidal forces that pull stars away from the clusters are essentially radial, giving the stars a larger initial radial velocity dispersion than in the tangential or perpendicular directions. The primary heating action of subhalos is to pull stars on either side of a crossing point toward the subhalo, which causes the tangential and vertical motions to have a more straightforward dependence on subhalo stream heating. The velocity dispersions in the leading and trailing streams are sufficiently symmetric that they are folded together to increase the statistical significance of simple σ ∝ ϕ^α fit with angular distance ϕ.

The particle velocities along the streams are shown in the lower row of panels of Figure 17. The mean radial and tangential velocities along the stream, shown as the open circles, primarily reflect the orbital motion along the various streams. The mean vertical velocity along the stream (the latitude component of the velocity in the stream's great-circle frame) is an indicator that streams can be tilted with respect to the orbital equator and that there can be relatively large-scale velocity perturbations from the more massive orbiting subhalos (Erkal et al. 2019; Shipp et al. 2019). The velocity dispersions along the streams shown in Figure 17 have a lot of scatter, as should be expected for the essentially random process of subhalo encounters with streams. For the 29 streams identified in the 10–60 kpc region the average radial velocity dispersion rises from 10 to 50 km s⁻¹ in 20°–30°, the tangential velocity dispersion rises from 7 to 40 km s⁻¹ in 40°, and the vertical velocity dispersion from 7 to 25 km s⁻¹ in 40°.

Table 2 reports the number of subhalos of mass 10^7–8 M_☉ in the 10–60 kpc range in the second column and the number of streams in the same radial range in the third column. The double mass disk has a few more subhalos than the baseline model at the final moment measured, although the time average over the last 1.5 Gyr is lower than the baseline.

Table 2. Subhalo and Stream Counts

Model	Nhalos	Nstreams
Disk	10–60 kpc
Baseline	39 ± 4	29
No disk	54 ± 4	66
2× disk	31 ± 3	36

Download table as:  ASCIITypeset image

Figure 18 plots the coefficients of the power law to the σ–ϕ velocity dispersion longitude relation for the three components of the velocity. The fits are done as least squares to the linear logarithmic relation, $\mathrm{log}{\sigma }_{\mathrm{fit}}+{\alpha }_{\mathrm{fit}}\phi$. The plotted results show a weak dependence on the maximum longitude fitted, which ranges from 30° to 60°. Each simulation has a distinct set of σ–ϕ relations. The simulation without a disk has the most subhalos and has the largest velocity dispersion. The baseline and the double disk simulation have similar ϕ dependences with the double disk model, which has more subhalos at late times (Table 2), having larger velocity dispersions. The measurements of the velocity dispersions of a population of streams offer some prospect for estimating the numbers of subhalos within a galactic halo. The results are sensitive to which streams are selected and how the measurements are made, so will need comparison to a set of simulations. The most challenging aspect for these measurements is that the stream progenitor needs to be identified, which is not possible for many streams (Bonaca et al. 2021).

Zoom In Zoom Out Reset image size

7.2. Longitude Summed Velocity Distributions

The distribution of velocities along a set of streams does not require any knowledge of the location of the progenitor; instead it requires summing the velocities of stars of all longitudes along the visible length of the stream. Figure 19 shows a set of streams longer than 50° in the 10–60 kpc range. These simulated streams are approximate analogs for the well-studied GD-1 stream (Grillmair & Dionatos 2006). The stream particles are placed into a grid of cells (1° × 05 on the sky) and smoothed with a 2D Gaussian, with σ_ϕ = 10 and σ_θ = 075. Streams are required to have a standard deviation of the longitudes of their particles of more than 20°, which removes streams with many stars in a fuzzy distribution around the progenitor. The cells above a threshold of 2 M_☉/□° are fitted to a fourth-order polynomial. The velocities in each component are fitted with the same procedure. The panels in each row show the particle latitude offsets from the mean stream, radius (no streamline correction), and the offsets in radial velocity, tangential velocity, and perpendicular velocities with stream longitude on the sky. In the leftmost panel the logarithm of the filtered peak density along the centerline is inserted into the figure with a scaling to the axis values of $5{\mathrm{log}}_{10}({D}_{p})-15$.

Zoom In Zoom Out Reset image size

Figure 20 shows the star particle distributions in streamline fit-adjusted latitude, radial velocity, tangential velocity, and velocity in the latitude direction, from left to right, respectively, for the three streams longer than 50° plotted in Figure 19. Streams having an rms spread in longitude greater than 10° are analyzed. Stream ends are defined as where the surface density drops below 2 M_☉/□° or the stream path changes more than 2° or the velocity changes more than 30 km s⁻¹. The entire population of 12 streams longer than 10° has similar values. The velocity distributions all have comparably narrow cores of ∼2–5 km s⁻¹, with broad non-Gaussian wings that vary substantially from stream to stream. Streams longer than 50° in the baseline simulation have a mean radial velocity dispersion of 2.1 km s⁻¹ within ±1°. GD-1 has a FWHM width of 05 (Koposov et al. 2010), which is a 1σ width of 021, and a measured radial velocity dispersion of 2.1 ± 0.3 km s⁻¹ (Gialluca et al. 2021). The simulation without a disk has long streams that are more than twice the width and velocity dispersion of the baseline simulation. The simulations convincingly show that the greater subhalo numbers in a simulation without a galactic disk leads to significantly greater subhalo heating of the streams. The baseline simulation streams have widths and velocities comparable to GD-1, but a statistical test requires consideration of observational selection effects and should be done with more simulations.

Zoom In Zoom Out Reset image size

The stream fitting procedure used for the long streams of Figure 19 is applied to all streams longer than 10°. The 3σ and 6σ velocity dispersions are measured for each stream. The average 3σ and 6σ clipped velocity dispersions are displayed in Figure 21 with the error bars showing the standard deviations of the mean values. The number of streams on average is given in the legend of the plot. The differences between the streams are significant on average over the stream population. An important caveat is that these streams all had progenitors with masses above 4 × 10⁴ M_☉, with a typical cluster losing about 60% of its mass. Therefore even lower-mass clusters, which lose more mass earlier, are likely to have even higher velocity dispersions.

Zoom In Zoom Out Reset image size

Velocity dispersions measured within 1° of the stream are most comparable to current velocity measurements, which find a radial velocity dispersion of 2.1 ± 0.3 km s⁻¹ (Gialluca et al. 2021). The average radial velocity dispersions of streams longer than 50° in the simulations are 2.1, 3.0, and 4.2 km s⁻¹ for the models with our baseline disk, the double mass disk, and the no-disk model, respectively. All three of the baseline disk simulation have velocity dispersions below 2.55 km s⁻¹, whereas for the no-disk model four of the 12 long streams are below this limit. A single stream is indicative of the subhalo population, but confident conclusions require study of the entire population of streams.

Dark matter subhalos orbiting within a galactic halo perturb the velocities of stars along the streams, which increases their velocity dispersion. The result is a dramatic increase in the spread of stream stars on the sky with time since they left their progenitor cluster. Subhalos in approximately the 10^7–8 M_☉ range dominate stream heating, based on their numbers and the impact approximation scattering calculation. There is considerable stream-to-stream variation because there are typically only 30–50 subhalos of the dynamically important 10^7–8 M_☉ mass interval in the 10–60 kpc region where streams are currently found. The divergence of orbits in the aspherical potential adds to the stream spreading.

Our simulations with disks contain streams with FWHMs and velocity dispersions comparable to those observed, although a statistical comparison to the stream population allowing for observational selection effects and field star rejection procedures has yet to be done. Most notably there are thin streams longer than 50°. For those streams the average 3σ clipped radial velocity dispersion within ±1° is 2.1 km s⁻¹ in the baseline simulation, rising to 4.2 km s⁻¹ in the no-disk simulation. The 3σ and 6σ clipped velocity dispersions have ratios ranging around 2. The significantly non-Gaussian velocity and width distributions are signature indicators of the subhalo heating of streams over their lifetimes.

The angular distance of stars along the stream from the progenitor cluster is a proxy for the time a star has been in the stream, although there is a spread of ages at any given distance. We measure the velocities within ±5° because restricting the measurements to stars within ±1° of the stream equator misses many higher-velocity stars. We fit the streams in the 10–60 kpc region above a minimum angular momentum of 2000 kpc km s⁻¹ to avoid streams that come close to the disk. Figure 18 presents fits for the components of the velocity dispersion with distance along the streams. Average velocity dispersions rise from ≃3 to 5 km s⁻¹ near the clusters to ≃30–40 km s⁻¹ at 20°–30° away. The velocity dispersions depend on the numbers of subhalos in the same volume as the streams, with ∼50% higher velocity dispersion in a model with no disk, which leaves about 60% more subhalos in the inner 60 kpc (Table 2).

The distribution of stars in position and velocity within a strip close to a stream above a minimum density threshold are measurements that do not require locating the progenitor star cluster. The stream widths and velocity dispersions plotted in Figure 21 show that the no-disk simulation with the most subhalos has the highest velocity dispersion streams. The statistical significance of the differences are greater for measurements extending 5° from the stream. The lower panel of Figure 21 shows that the long streams in the baseline simulation have a radial velocity dispersion of 2.0 km s⁻¹ whereas the disk-less simulation streams have an average radial velocity dispersion of 4.2 km s⁻¹, twice that of measurements of GD-1 (Gialluca et al. 2021). However, four of the 12 no-disk long streams do have velocity dispersions comparable to GD-1. That is, a few streams with ∼2 km s⁻¹ radial velocity dispersions are present in all three simulations. These results underscore that unbiased velocity measurements around a set of Milky Way streams will be required to reach confident conclusions about the subhalo population.

Comments from an anonymous referee substantially improved this paper, as did comments from Ana Bonaca, Kathryn Johnston, and David Weinberg.

This research was supported by NSERC of Canada. Computations were performed on the Niagara Supercomputer at the SciNet HPC Consortium. SciNet is funded by the Canada Foundation for Innovation, the Government of Ontario, Ontario Research Fund—Research Excellence, and the University of Toronto.

Software: Gadget-4 (Springel et al. 2021), and Amiga Halo Finder (Gill et al. 2004; Knollmann & Knebe 2009).