USING GAPS IN N-BODY TIDAL STREAMS TO PROBE MISSING SATELLITES

The host galaxy is modeled as a dark matter halo, as well as a set of subhalos which orbit around the halo's potential. A Milky-Way-sized halo is modeled with a static spherical Navarro–Frenk–White (NFW) profile (Navarro et al. 1997) with v_max = 210 km s⁻¹ located at r_max = 30 kpc. Each individual subhalo is modeled by a spherical Hernquist profile (Hernquist 1990)

$$\begin{equation} \Phi _i(r) = \frac{GM_i}{h_i+r} \end{equation} \tag{ 1 }$$

for simplicity, compared to Einasto profiles which produce better fits in simulations but are more complicated to compute (Springel et al. 2008). We use the formula found in Carlberg (2009), which approximates the results from both Springel et al. (2008) and Neto et al. (2007), where h_i is independent of galactocentric position, and is related to M_i by

$$\begin{equation} h(M) = 6\,\mathrm{kpc} \times \left(\frac{M}{10^{10}\ M_\odot } \right)^{0.43}. \end{equation} \tag{ 2 }$$

We use the mass and spatial distributions of the subhalos from the results of the Aquarius simulations (Springel et al. 2008), where the mass function is independent of the spatial distribution. The mass function is a power law

$$\begin{equation} \frac{dN}{dM} = 3.26\times 10^{-5}\ M_\odot ^{-1} \left(\frac{M}{2.52\times 10^{7}\ M_\odot } \right)^{-1.9}, \end{equation} \tag{ 3 }$$

and the spatial distribution follows an Einasto profile

$$\begin{equation} n(r) \propto \exp \left\lbrace -2.95 \left[ \left(\frac{r}{199\,\mathrm{kpc}} \right)^{0.678} - 1 \right] \right\rbrace. \end{equation} \tag{ 4 }$$

The subhalos' velocities are initialized with a Gaussian distribution where the velocity dispersion is the solution to the isotropic Jean's equation (Binney & Tremaine 2008) using the halo's potential. The subhalos orbit around this potential as test masses.

The progenitor of the stream, which is an approximation to a globular star cluster, is initialized using 10⁶ particles of equal mass as a King model with parameters w = 4.91, total mass 4.29 × 10⁴ M_☉, and a core radius of 0.01 kpc. This results in a zero-density radius of 0.103 kpc. Each N-body particle in the system interacts with the dark matter halo's and subhalos' potentials. With the Galactic center at the origin, the satellite is initially put at (x, y, z) = (30, 0, 0) kpc and velocity (v_x, v_y, v_z) = (0, 140, 0) km s⁻¹. The resulting orbit is confined on the xy-plane with eccentricity 0.33, peri- and apogalacticon at r_p = 15 kpc and r_a = 30 kpc, respectively. The azimuthal and radial periods are about 0.70 Gyr and 0.47 Gyr, respectively.

We use Gadget-2 (Springel 2005) for our N-body simulations. Since the public distribution¹ does not have functionality for external potentials, we modify the code such that in every time step, an external acceleration term which accounts for the potentials of the halo and all the subhalos is added to the accelerations of all the particles after their N-body interactions are computed.

Each of our simulations lasts 10 Gyr, and we impose a maximum time step of 1 Myr. The particle softening is 5 pc. Each simulation produces 500 snapshots, one every 20 Myr, and each consists of the positions and velocities of the N-body particles and subhalos.

The star cluster is modeled as an N-body system which forms a stream as the cluster is disrupted by the tidal field of the massive host. When the cluster is isolated, the energies of the individual particles are conserved to a few percent over 10 Gyr. Using the softening as minimal impact distance, the relaxation timescale in the core is ≳ 110 Gyr, which is much greater than the orbital period at ≲ 1 Gyr. Figure 1 shows the mass enclosed inside 0.103 kpc of the star cluster's center it is orbiting in the absence of subhalos. Because the stream is repeatedly stretched and compressed longitudinally along the eccentric orbit, the mass enclosed in a fixed radius is not always decreasing in time. The bottom panel of Figure 1 shows that the mass loss is driven purely by bulge shocking (Binney & Tremaine 2008), as the periodic bursts have exactly the same period as the radial period of the orbit.

Zoom In Zoom Out Reset image size

Figure 2 shows the velocity dispersion and traverse FWHM of our simulated stream without any subhalos. The stream is chosen so that its properties are on the same orders of magnitude as Pal-5 (Dehnen et al. 2004; Odenkirchen et al. 2009) and GD-1 (Koposov et al. 2010; Willett et al. 2009). In the derivation in Carlberg (2012), the gap formation rate is expressed as a function of galactocentric distance of the orbit and width of the stream. For our simulated stream, we adapt average values of 22 kpc and 0.3 kpc, respectively, over the entire stream.

Zoom In Zoom Out Reset image size

In a time independent and spherical potential Φ(r), three interesting conserved quantities are the radial, azimuthal, and latitudinal actions where

$$\begin{eqnarray} &&J_r = \frac{1}{\pi }\int _{r_p}^{r_a} dr \sqrt{2E-2\Phi (r)-\frac{L^2}{r}} \end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray} &&J_\phi = L_z \end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray} &&J_\theta = L-|L_z| \end{eqnarray} \tag{ 7 }$$

respectively (Binney & Tremaine 2008), where r_a and r_p are the apo- and pericentric distances of the orbit, respectively. Since our stream's orbital plane is the xy-plane, J_θ ≈ 0 (though not exactly 0 because the stream has a finite thickness), so we only consider J_r and L_z. For simplicity, we ignore the progenitor and only consider the particles which have already escaped from the cluster, so when computing J_r we assume that the potential due to the progenitor's potential is negligible.

J_r and L_z are useful since their dispersions are the origins of the stream's average width. For example, in the epicyclic approximation (Binney & Tremaine 2008), where κ and a are the epicyclic frequency and amplitude, respectively, the radial motion can be written as

$$\begin{equation} r(t)=a \cos (\kappa t + \psi), \end{equation} \tag{ 8 }$$

where ψ is an arbitrary phase angle. Then it can be shown that

$$\begin{eqnarray} &&J_r = \frac{1}{2\pi } \oint p_r dr \propto \int \dot{r}\,dr \propto \kappa a^2 \end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray} &&L_z = R_g^2\Omega, \end{eqnarray} \tag{ 10 }$$

where Ω and R_g are the orbital frequency and radius of the guiding center, respectively. Clearly, dispersions in both a and R_g can affect the width of the stream. Therefore, conserved quantities J_r and L_z are especially valuable in understanding the width of the stream. Figure 3 shows the distributions in (J_r, L_z) when the stream is 8 Gyr old. The two lobes at higher and lower L_z are the trailing and leading branches of the stream, respectively. The absolute dispersions in J_r and L_z are on the same order of magnitude, in rough agreement with the formula ΔJ_r/ΔL ∼ (r_a − r_p)/πr_p (Eyre & Binney 2011). In Section 4.1, we will show how the spread in actions affects the morphology of stream gaps.

Zoom In Zoom Out Reset image size

Plots similar to Figure 3 can be found in Eyre & Binney (2011) where they used angle-action variables extensively to study the relation between the stream and the orbit. A similar plot can also be found in Yoon et al. (2011), but in scaled energy and angular momentum which were first used by Johnston (1998) to describe the dynamics of tidal streams.

Yoon et al. (2011) divided their subhalo mass spectrum into separate mass ranges in order to resolve the contributions from each mass range. Using the same approach, we divide the subhalos from 6.5 × 10⁴ to 10⁸ M_☉ into 13 mass ranges (Table 1). Each mass range contains incrementally more subhalos, starting from the higher mass end toward the lower mass end. In the higher mass end, the mass ranges are chosen such that the increase in subhalo masses are roughly the same. In the lower mass end, the mass ranges are chosen such that the increase in subhalo numbers are roughly the same.

In each set of subhalos, we reduce their numbers by eliminating those whose orbits are always inside the perigalacticon and outside the apogalacticon of the progenitor's orbit. The largest subhalo in our simulations has a length scale of ∼1 kpc (Equation (2)), so all the subhalos with perigalacticon (apogalacticon) larger (smaller) than that of the progenitor's orbit by ∼2 kpc will interact minimally with the stream. This allows us to safely eliminate the subhalos with perigalacticon larger than 32 kpc, and apogalacticon smaller than 13 kpc.

We run 14 simulations—1 "smooth stream" without any subhalos, and 13 "ΛCDM streams" containing the subhalos in the mass ranges in Table 1—with identical initial conditions and dark matter halo potential. This allows us to resolve the effects of the lower mass subhalos whose existence are in question.

Gaps are manifested as local minima in the linear density along the stream. To obtain the linear density along a stream in an eccentric orbit, we first fit the stream with two degree-6 polynomials—one each for the leading and trailing streams—in polar coordinates centered at the Galactic center. The points along each line are spaced at 0.002 radians apart. Between each pair of adjacent points a cylinder of radius 1 kpc is drawn which lies lengthwise along the pair of points. The linear density is then the number of particles inside this cylinder divided by the length of the cylinder. This spacing is chosen so that the gaps as wide as the stream are well resolved.

The method used to find gaps in stream densities is inspired by the technique first used by Carlberg et al. (2012) to find gaps in observations. They used matched filters of the estimated shape of a density gap at various length scales to look for positions in the stream which potentially contain gap signals. The filter consists of a local minimum which is the underdensity of stars, and two local maxima on both sides of the minimum due to conservation of mass (Figure 4). This method is similar to the wavelets approach, where the integral of the filter function is constructed to vanish inside a certain domain. The potential gap signals are then easily identified as local maxima in the convolution between the filter and the signal.

Zoom In Zoom Out Reset image size

To obtain the significance of each potential gap signal against noise, Carlberg et al. (2012) produced bootstrap samples from the sky background. With the simulations in this study, we can estimate noise levels using the smooth stream. Note that the "smooth stream" itself is not totally smooth. As we will show in Section 4.2, there are large density fluctuations near the progenitor due to the coherent in epicyclic motion of the particles, as first explained by Küpper et al. (2008). When the particles become unbound from the progenitor, they pile up near the base of their cycloid trajectories, creating epicyclic overdensities (hereafter EO) along the stream. Although this intrinsic process to mimic gaps can be confused with gaps caused by subhalos, EOs are only apparent within ≲ 5 kpc away from the progenitor in our streams (Figures 7 and 8). The details of the dynamics of EOs is beyond the scope of this study, but this effect can be understood in terms of orbital actions. EOs occur due to coherent epicyclic motions of the particles, which nevertheless have finite dispersions in orbital actions (Figure 3) and are not perfectly coherent. Therefore, although the escaping particles' orbits stay roughly coherent in the first few clumps, their orbits eventually drift out of phase as they travel along the stream. This explains why the density peaks of EOs further downstream are not as apparent as the peaks closer to the progenitor (Just et al. 2009; Küpper et al. 2010).

After masking 10 kpc of the smooth stream centered at the progenitor, the rest of the smooth stream is simply noise. Our method to find gaps in a given stream can be summarized as follows.

• 1.  
  Compute

  $$\begin{equation} C_s(x) = \frac{1}{s} \int _{x-1.5s}^{x+1.5s} [ \rho (x^\prime)-\bar{\rho }_s(x) ]\, f\left(\frac{2(x-x^\prime)}{s} \right) dx^\prime , \end{equation} \tag{ 11 }$$

  where f(t) = (t⁶ − 1)exp (− 1.2321t²) is a matched filter function (Carlberg et al. 2012; Carlberg & Grillmair 2013), and s is the filter scale. $\bar{\rho }_s(x)$ is the mean of ρ(x) inside [x − 1.5s, x + 1.5s], the domain in which the integral of f(2x/s) itself vanishes. Each potential gap signal would appear as a local maximum in C_s(x). This convolution is computed at 12 logarithmically spaced filter scales from 0.1 to 5 kpc (Figure 4), and then all the local maxima of each stream are sorted by C_s.
• 2.  
  Repeat the above step using the smooth stream, but with ±5 kpc from the progenitor masked along the stream. The set of local maxima in $\tilde{C}_s(x)$ from this convolution is the noise, which are also sorted by $\tilde{C}_s$.
• 3.  
  Each local maximum in the signal set are compared against the noise set. A signal element at any position along either branch of the stream that ranks higher than 99% in the noise set is identified as a gap.

Inevitably, this method may detect the same gap at 99% confidence at different scales but in very close proximity. To avoid over-counting, we employ the following scheme to eliminate overlapping gaps. First, we define an overlap as two gaps whose C_s local maxima are identified at C_s1(x₁) and C_s2(x₂) that are within s₁ away from each other along the stream, where s₂ < s₁. When this occurs, the gap with higher C_s eliminates the lower.

Our gap detection method requires no prior knowledge whether a given gap is an EO or a subhalo perturbation, both of which can be identified as a series of over- and under-densities. When we count the number of gaps in the end, EOs will be included. One key result of our study is that gaps due to EOs are distributed very differently in lengths compared to gaps due to subhalo perturbation.

According to Yoon et al. (2011), gaps in general are diagonal and not perpendicular to the stream due to a gradient in angular momentum (hence a gradient in orbital velocities) across the width of the stream, which can shear a gap longitudinally. Figure 3 allows us to estimate the shearing effect in our streams using the distributions in angular momenta. For each branch of the stream, the FWHM spread in angular momentum is about ΔL ∼ ΔL_z ∼ 10 kpc km s⁻¹. For a narrow stream at r = 22 kpc, the spread in velocity is Δv = ΔL/r ∼ 0.5 km s⁻¹. Therefore, a gap that spans the width of the stream will be sheared by less than 1 kpc per Gyr.

Figures 5 and 6 show the time evolution of the smooth and a ΛCDM stream, respectively, from 7 to 8 Gyr. The EOs near the progenitor appear to shear by different amounts at different times, but this is due to the radial oscillation in the orbit where the radial period is ∼0.5 Gyr. Upon closer inspection of Figure 6, we also note that not only do subhalo gaps have complicated morphologies, but their orientations flip back and forth in a radial period due to the spread in J_r. Nevertheless, comparing panels of the same radial phase at one radial period apart, the end points of each gap across the width of the stream do not shift by any appreciable amount. Rather, the morphologies of the subhalo gaps are already apparent as each gap first appears.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

If the linear density of a stream is calculated by integrating the entire thickness of the stream in traversing slices along the stream, then the contrast of the gap will be reduced. This is because the edges of the gaps are not perfectly straight across the width of the stream, so dividing the stream into slices will smear out the density contrast. To investigate how much the smearing will affect gap detection, we calculate the linear densities in two ways. (1) Integrating cylindrical slices of radius 1 kpc along the stream, hereafter the "whole width," where 1 kpc was chosen to cover the entire thickness of the whole stream. (2) Integrating only the cylindrical slices of radius 0.04 kpc centered along the best fit line of each branch of the stream, hereafter the "central width," where 0.04 kpc is chosen to mimic the GD-1's observed width of 0.08 kpc (Carlberg & Grillmair 2013). This central width then encloses about 30%–40% of the mass of the whole width, depending on its orbital phase where, for example, the stream is radially compressed during pericentric passage.

Figures 7–10 show the densities along the whole and central widths of the smooth stream and a ΛCDM stream from 3 to 10 Gyr. The streams younger than 3 Gyr are not shown as the stream is ≲ 10 kpc long at those ages, so the gaps are dominated by very prominent EOs. Moreover, the stream itself does not yet have a large enough cross-section to produce enough gaps for meaningful statistics. In each panel, the shaded columns represent the gaps that are found on the scale of the columns' widths. Although these gaps are identified as being 99% significant, the density contrasts of the gaps have not been quantified. For the rest of this paper, we assume that all gaps identified at 99% significance can be observed. Note that because the gap finding process is applied independently to each snapshot, the shaded columns do not necessarily represent the time evolutions of individual gaps. Instead, the shaded columns show the general distributions of gaps—both in space and in gap lengths.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Our gap finding method has a number of problems. In the smooth stream (Figures 7 and 8), our method by construction identifies 1% of the noise as gaps. This is why there can be spurious gaps detected well beyond 5 kpc away from the progenitor, even though EOs tend to form very close to the progenitor. Also, the overall profile of the stream density can sometimes be confused as a gap as well. One example is a gap at 4 kpc at 5 Gyr in Figure 9, where a smooth density gradient from 3 to 8 kpc is mistaken as the right half of a long gap. At the 95% confidence limit, both kinds of false positives are quite common and can often be identified by eye. For the results below, we show gaps that are 99% significant, which minimizes the occurrence false positives.

Following the idealized experiments in Carlberg (2012), Carlberg & Grillmair (2013) derived an analytical relation between gap formation rate R_∪ which is the cumulative number of gaps longer than length l per unit stream length per unit time as a function of gap length (hereafter the "gap spectrum") such that

$$\begin{equation} R_\cup ^{\mathrm{ideal}} = 0.060\, \hat{r}^{0.44}\, l^{-1.16}\, \mathrm{kpc}^{-1}\,\mathrm{Gyr}^{-1} , \end{equation} \tag{ 12 }$$

where $\hat{r}\equiv r/30\,\mathrm{kpc}$, and we adapt r = 22 kpc for the average galactocentric radius of the stream. In this section, we aim to study the validity of Equation (12) in our self-consistent stream. We set the "length" of each gap as the scale s of the matched filter which identified the gap (Section 3.2)

4.3.1. Smooth Stream without Subhalos

4.3.2. ΛCDM Stream with Independent Sets of Subhalos

As an ideal case, Equation (12) ignores the visibility of gaps when the same position of a stream suffers impacts by multiple subhalos at different times. For instance, after one major impact by a massive subhalo which results in a long and high contrast density gap at an early time, subsequent impacts by less massive subhalos in that same region at a later time may not be visible.

Gap overlapping can be minimized by the following experiment. We run 13 separate simulations with the same initial conditions as the star cluster, but the subhalo masses are selected differentially from Table 1. This allows each stream to interact with an independent set of subhalos of a very small range of masses. Overlapping can still occur within the same simulation for each set of subhalos (hence a small number of gaps can still be eliminated), but to a much lesser extent than using integrated mass ranges.

Figure 11 shows the measured gap spectrum from the gaps collected from all 13 simulations using independent sets of subhalos. In the top panel, which includes all gaps, the measured gap spectrum matches Equation (12) reasonably well. However, this is a coincidence as the gaps contain EOs which are not described by Equation (12). In an attempt to eliminate EOs, in the bottom panel of Figure 11, the gaps that are located within 5 kpc away from the progenitor are eliminated. When computing the gap formation rates in these cases, the number of gaps are divided by a stream length which is reduced by 10 kpc and a stream age which is reduced by 2 Gyr (i.e., the age of the stream when it is 10 kpc long). This allows us to facilitate a fair comparison of gap spectra against the cases which include all gaps in the entire stream. Comparing the two panels in Figure 11, we can see the masking of the 10 kpc around the progenitor reduces the abundance of shorter gaps. This is expected since that region of the stream contains mostly EOs which occur at scales ≲ 1 kpc. In general it is difficult to tell whether a given gap within 5 kpc is due to EOs or subhalos, so in the process some subhalo gaps near the progenitor may have been eliminated as well.

Zoom In Zoom Out Reset image size

Whether Equation (12) is a good description of the gap spectrum in a general stream likely requires more simulations with varying orbital parameters. Nevertheless, it is apparent that the gap spectrum does not depend strongly on the age or the integrating width for the linear density of the stream. However, it is still an ideal case since a stream realistically interacts with all the subhalos at the same time. As we show in the following section, gap overlapping can significantly alter the gap spectrum.

4.3.3. ΛCDM Stream with all Subhalos

We now consider the validity of Equation (12) for stream gaps in the presence of all subhalos in each cumulative mass range in Table 1. Figure 12 compares the measured gap spectra of both the whole and central streams for three mass ranges, with and without the gaps within 5 kpc away from the progenitor. Clearly, in all cases, the ideal gap spectrum over estimates the measured spectrum by nearly an order of magnitude due to gap overlapping.

Zoom In Zoom Out Reset image size

Similar to the simulations with independent sets of subhalos (Figure 11), the gap spectra produced by full sets of subhalos do not have strong dependence on stream age and integrating width. The only exception is the youngest stream shown at 5 Gyr which consistently has higher R_∪ than the older streams. However, when the gaps near the progenitor are eliminated, the numbers of gaps at 5 Gyr in all cases decrease significantly, where the gap spectra are dominated by a single gap at 2–3 kpc, and a number of extremely short gaps. This is likely because the stream is still young, and the effective length of the stream (after masking 10 kpc centered at the progenitor) is only ∼15 kpc. This eliminates a significant part of the stream, making its stream statistics unreliable.

The weak dependence of the gap spectrum on the integrating width for linear density is also worth noting. Figure 6 shows that gaps in general have much more complicated morphologies than straight edges across the width of the stream. The explanation for these morphologies requires detailed understanding of how subhalo perturbations manifest in a self-consistent stream, which is beyond the scope of this study. While Carlberg (2013) studied the dynamics of subhalo perturbations for an idealized stream, we defer the self-consistent case to a future study.

Perhaps the most surprising result is that the gap spectra do not show obvious dependence on subhalo masses. The spectra are difficult to distinguish between the mass range of subhalos which causes the gaps. This is in disagreement with Carlberg (2012), which derived a relation between the length of a gap and the mass of the subhalo that caused it such that

$$\begin{equation} l(M) = 8.3 \left(\frac{r}{30\,\mathrm{kpc}} \right)^{0.37} \left(\frac{M}{10^8\ M_\odot } \right)^{0.41}\,\mathrm{kpc}. \end{equation} \tag{ 13 }$$

From this formula, it is reasonable to expect the inclusion of lower-mass subhalos to show more gaps at the shorter end. In their ideal simulations, however, Carlberg (2012) did not account for the time evolution of each gap. An example can be seen in Figures 9 and 10. The gap located at about −7 kpc at 4 Gyr evolves into a much longer gap centered at about −9 kpc at 10 Gyr. Evidently Equation (13) requires revision for self-consistent streams before it can be used to understand the relation between gap spectra and subhalo masses.

We now consider the issues when interpreting gap spectra from observations. A gap spectrum for GD-1 has been observed by Carlberg & Grillmair (2013), but we emphasize that the gap spectra from our simulated streams in this study should not be directly compared to the one in Carlberg & Grillmair (2013) because our models for both the star cluster and the galaxy halo are chosen in favor of a simple interpretation, and may be missing some complications discussed in Section 4.6.

We first project each stream onto sky coordinates. For simplicity, we put the hypothetical observer at the center of the galaxy, and then project each particle onto the azimuthal and altitudinal plane in galactocentric coordinates. Since the our stream progenitor is orbiting along the xy-plane in a spherical potential, the smooth stream appears as a straight line along the azimuthal direction, and each ΛCDM stream appears only a few degrees off the azimuthal plane due to subhalo perturbations.

The density along the stream is simply the number of particles in bins of 01 in the azimuthal direction. The match filter approach to detect gaps remain the same as the analysis above, but the 12 filter scales (Figure 4) are now logarithmically spaced in angular units from 034 to 14°, and the noise levels are obtained from the regions at >10°, rather than >5 kpc, away from the progenitor. The choice of bin size and filter scales are on the same orders of magnitude as Carlberg & Grillmair (2013), but putting the hypothetical observer at the Galactic center may affect angular sizes by factors of ∼2.

To ensure that the behaviors of the simulated streams are typical, we simulate each ΛCDM stream 10 times with the same initial conditions for the star cluster, but different realizations of the same subhalo distributions. At the end we take the median numbers of gaps of the 10 streams to avoid outliers.

4.4.1. Orbital Phase

One surprising result from Section 4.3.3 is that the gap spectrum has little to no dependence on age and subhalo masses. To investigate what this means when interpreting observations, the top panel of Figure 13 shows the cumulative numbers of gaps longer than 034 (i.e., all gaps detected in the entire stream) as functions of time. At t > 5 Gyr, the numbers of gaps due to subhalos vary according to the orbital phase of the stream progenitor. The "bursts" in numbers of gaps in the ΛCDM streams occur when the streams are stretched as they pass through the pericenters of their orbits. At t < 5 Gyr, on the other hand, this correlation does not exist for two reasons. First, our detection method (Section 3.2) uses the parts of the stream that are >5 kpc (before sky projection) or >10° (after sky projection) away from the progenitor in order to estimate noise. At t < 5 Gyr, the length of the stream varies between a few to 20 kpc, which may not be long enough to estimate noise. Second, in only 5 Gyr the stream does not yet have enough time and to interact with subhalos.

Zoom In Zoom Out Reset image size

Both the total number of gaps and the dynamical age of the stream are difficult to measure, as the some parts of a stream may not be observable. We define a more useful quantity S_∪ which is the cumulative number of gaps longer than a given gap length per unit stream length. In other words, leaving the age of the stream as an unknown, S_∪ differs from R_∪ in Section 4.3 by a normalization by age, and that S_∪ is after sky projection. In the next section, we show that gaps due to EOs and subhalos have very different S_∪(l) distributions which are directly observable.

4.4.2. Signal-to-noise Ratios

Another important distinction between simulations and observations is the signal-to-noise ratio (S/N). At 5 and 8 Gyr our smooth stream is represented by about 60,000 and 80,000 particles in total, respectively (Figure 1). Koposov et al. (2010) estimated that the 60° visible segment of GD-1 consists of 3000 stars. At an average distance of ∼10 kpc, the visible segment is ∼10 kpc long. In our simulations, after 5 and 8 Gyr, the average stream lengths are about 20 and 40 kpc, respectively (Figure 1). This means that our simulated stream should be represented about eight times fewer particles in order to be comparable to observations. With the progenitor masked, we reduce the number of particles in the stream by randomly sampling the stream using two, four, and eight times fewer particles than the original stream. The particle reduction applies to both the stream of interest and the smooth stream which is the source for estimating noise. This allows us to investigate the importance of high S/N.

Note that in our simulation each particle is equivalent to about 0.043 M_☉, which is less massive than the typical stars that are detected in observations. Our simulations are not meant to be physical models of the real stream. In this section, we are only concerned about matching the numbers of particles in the simulations to the numbers of stars in the observation. As the stars escape from the progenitor, the stream's self-gravity becomes negligible (Johnston 1998), and the particles' masses are no longer important.

Figure 14 shows the density profiles of the ΛCDM stream with subhalo masses 1.5 × 10⁵ < M/M_☉ < 10⁸ at 8 Gyr projected onto the sky. The panels show the gaps detected in the same stream after three levels of particle reduction. Even after reduction by a factor of eight, the stream appears to have retained most of its gaps despite a lower S/N.

Zoom In Zoom Out Reset image size

In Figure 15, each line shows the median of ten gap spectra from the same stream but with 10 realizations of the same subhalo distribution. In each panel, the solid (dashed) lines represent the times when the progenitor is at the pericenter (apocenter) of its orbit. When the stream is compressed and stretched as it oscillates radially (see Figure 1), its length can differ by up to a factor of two. Careful inspection of Figures 13 and 15 shows that during pericentric passages, the numbers of gaps are at maximum, but S_∪ is at minimum because the stream length is also at maximum. For a ΛCDM stream at high S/N (upper left panel in Figure 15), the gap spectra are not sensitive to this oscillation, except with an excess of shorter gaps and fewer longer gaps, which are expected as the stream, including its longitudinal structure, is compressed during apocentric passage. At low S/N (lower left panel in Figure 15), however, the gap spectrum during apocentric passage is consistently higher than that during pericentric passage. This is because the length of the stream is insensitive to the S/N, but the number of gaps is not. Therefore, high S/N data for the stream is important when studying stream gaps. Otherwise, the spectrum may be over- or underestimated depending on the orbital phase.

Zoom In Zoom Out Reset image size

The right panels in Figure 15 show the gap spectra of the smooth stream. They are also somewhat sensitive to S/N, but the most obvious difference from the spectra of ΛCDM stream is the shape of the spectra. This is especially obvious during pericentric passage where the gap spectra rapidly drop off to zero for gaps longer than ∼1°. Even during apocentric passage, the gap spectra remain flat at gap lengths ≳ 1°. If the gaps originated from subhalo perturbations, then the gap spectrum should be steep and extend well beyond 1°.

4.4.3. The Case of GD-1

We consider the effects of subhalos as massive as 10⁸ M_☉, since the effects by more massive subhalos are not relevant to us. Gaps caused by subhalos at these masses produce long gaps with high-density contrasts. For example, an obvious gap located at −7 kpc at 4 Gyr shown in Figures 9 and 10 are caused by a 4.5 × 10⁷ M_☉ subhalo. In fact, the perturbation by M ≳ 10⁸ M_☉ (h ≳ 1 kpc) subhalos can be so catastrophic that the stream is warped and divided into segments. As a result, a stream which originated from one progenitor can be observed as a few separate streams. Observations of Pal-5's and GD-1's gaps, on the other hand, show small-scale density fluctuations in a long, narrow stream, so these two streams are not sensitive to subhalos above 10⁸ M_☉. By coincidence, this upper limit approximately coincides with the upper limit beyond which the models of warm dark matter can be no longer be distinguished from CDM. Therefore, 10⁸ M_☉ is a reasonable upper limit where our simulations can be useful.

In the low-mass end, we only consider the effects of subhalos down to ∼6 × 10⁴ M_☉. From a separate simulation of the same stream but with only the subhalos with masses 6.5 × 10⁴ < M/M_☉ < 7.5 × 10⁴, the density profile is indistinguishable from the smooth stream, and the gap statistics are identical. Furthermore, Section 4.3.3 shows that the gap spectra have very little dependence on mass. Changing the mass lower limit from 7.5 × 10⁴ M_☉ to 6.5 × 10⁴ M_☉ produced indistinguishable gap spectra. This means that subhalos less massive than ∼10⁵ M_☉, even though they are much more abundant than those of higher masses (Equation (3)), have minimal effects on our stream.

The Milky Way has about 160 known GCs (Harris 1996), and a few hundred DGs brighter than L ≳ 10³L_☉ after bias corrections (see Bullock (2010) for a review). It is interesting to ask whether these known satellite systems, rather than the truly "missing" satellites, can contribute to the observed stream gaps. Typical GCs have masses ∼10⁵ M_☉, which correspond to the low end of our mass spectrum of subhalos. In the same mass range, though, there are orders of magnitudes more subhalos (e.g., ∼10⁵ subhalos at 10⁵ < M/M_☉ < 10⁶) than GCs, so GCs are unlikely to contribute significantly to observed stream gaps. On the other hand, DGs are commonly found at ≳ 10⁷ M_☉ (Strigari et al. 2008) which is the high end of our mass spectrum of subhalos. At that mass range (∼2000 subhalos at M > 10⁷ M_☉), the number of known DGs are only one order of magnitude below the number of subhalos, so DGs may contribute to some observed gaps. However, a common limitation in understanding the contributions from both GCs and DGs is their orbits, especially when the kinematics of these satellites are not well constrained. As done in our simulations (Section 3.1), subhalos that do not approach the stream's orbit will interact minimally with the stream. Table 1 shows that in our realizations of subhalos, only ∼15% of them would approach to within 2 kpc of GD-1's orbit. This means that most known satellites may never interact with a GD-1-like stream, and that stream gaps, if they were indeed due to satellites and were not EOs, are more likely due to satellites that have never been observed.

In order to keep our results simple, the galaxy is modeled as stationary, spherical NFW potential, the subhalos as static, test masses, and the satellite as a collisionless King model. These models ignore a number of known dynamical complications.

Two-body Relaxation. The star cluster is modeled as a collisionless system with relaxation timescale of ∼110 Gyr. Globular clusters typically have relaxation timescales of ≲ 10 Gyr (Harris 1996, 2010 Edition), so mass loss should originate from dynamical evaporation, in addition to tidal disruption. As a result, the star cluster should be disrupted even faster than we measured in Figure 1. This may have an important effect on the formation of gaps, which depends on the details of the dynamics of a stream (Carlberg 2013). The relation between gaps and mass loss mechanism will be investigated in a future study.

Dynamical Friction (DF). Both the star cluster and subhalos should suffer from DF as they orbit around the dark matter halo. Comparing the magnitudes of the accelerations due to DF and due to the orbit, a_DF/a_orbit ∼ 10⁻⁸ln Λ for both the star cluster at 22 kpc and a 10⁶ M_☉ subhalo at 100 kpc, where ln Λ ≡ ln (b_max/b_min) ≈ 10 is the log of the ratio of the maximum and minimum impact distances (Binney & Tremaine 2008). Therefore, DF is negligible throughout our model.

Disk Shocking. Dehnen et al. (2004) found that the evolution of Pal-5 is driven by the tidal shocks when crossing the Galactic disk, which is not modeled in our simulations. The orbit of Pal-5 in Dehnen et al. (2004) has peri- and apogalacticon at 5.5 kpc and 19 kpc, respectively, whereas our smooth stream has peri- and apogalacticon at 15 kpc and 30 kpc, respectively. Being farther away from the Galactic center, if our simulations contained a disk, its effect should be less severe for our simulated stream than Pal-5. Moreover, Dehnen et al. (2004) concluded that disk shocking is not responsible for the observed structure in Pal-5, while Küpper et al. (2010) concluded that EOs persist even under the influence of disk shocks, so the absence of a disk should not significantly change our conclusion. The gap formation rate with and without subhalos in the presence of a disk is beyond the scope of this study.

Halo Shape and Collapse History. Siegal-Gaskins & Valluri (2008) found that the shape of the halo potential can have a larger effect than subhalos have on the over all structure of a stream. However, their simulations focused on streams which originated from DGs at ∼10⁹ M_☉, as well as subhalos at ≳ 10⁷ M_☉ which is the high end of our mass spectrum of subhalos. The small gaps from a stream originating from a GC in a non-spherical halo has yet to be studied. In fact, since the initial collapse of the entire halo, the potential cannot be stationary throughout a Hubble time, which is the timescale of our simulations. In the future, we aim to repeat a similar study using potentials which resulted directly from high-resolution simulations such as Madau et al. (2008) and Springel et al. (2008). The self-consistent halo and subhalo potentials from those simulations can eliminate the idealized models in Section 2.1.