SIMULATING TIDAL STREAMS IN A HIGH-RESOLUTION DARK MATTER HALO

Originally developed to compute collisionless N-body dynamics for galaxies (Hernquist & Ostriker 1992), the SCF method solves the gravitational Poisson equation by basis decomposition. This method is optimized for dark matter halos since its lowest-order basis function is already an idealized profile, and the higher-order basis functions can be used to describe the radial and angular deviations from this idealized profile. The basis functions are bi-orthonormal and complete, so they can model any matter distribution as the decomposition order increases.

The approach to model the gravitational potentials in the snapshots of existing simulations of dark matter halos has been studied extensively by Lowing et al. (2011, hereafter L11). Even with moderate decomposition orders, L11 was able to recover much of the detailed dynamics inside the halos of the Aquarius simulations (Springel et al. 2008). We follow their approach to compute the forces in the halo of the VL-2 simulation at redshift zero. In this section we briefly summarize the method.

The density $\rho ({\boldsymbol{r}} )$ and gravitational potential ${\Phi }({\boldsymbol{r}} )$ are related by the Poisson equation ${{\nabla }^{2}}{\Phi }=4\pi G\rho$. Given a simulation snapshot that contains ρ, SCF solves the Poisson equation by decomposing Φ and ρ such that

$$\begin{eqnarray}&&{\Phi }\left( {\boldsymbol{r}} \right)=\mathop{\sum }\limits_{nlm}{{A}_{nlm}}{{{\Phi }}_{nlm}}\left( {\boldsymbol{r}} \right)\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\rho \left( {\boldsymbol{r}} \right)=\mathop{\sum }\limits_{nlm}{{A}_{nlm}}{{\rho }_{nlm}}\left( {\boldsymbol{r}} \right),\end{eqnarray} \tag{ 2 }$$

where ${\boldsymbol{r}} \equiv (r,\phi ,\theta )$ is expressed in spherical coordinates for convenience. We follow the derivation of Hernquist & Ostriker (1992), which used the Hernquist profile (Hernquist 1990), where ${{{\Phi }}_{000}}\equiv -1/(1+r)$ and ${{\rho }_{000}}\equiv 1/[2\pi r{{\left( 1+r \right)}^{3}}]$ are the zeroth-order basis functions. The general basis functions can be written as

$$\begin{eqnarray}&&{{{\Phi }}_{nlm}}\propto -\frac{{{r}^{l}}}{{{\left( 1+r \right)}^{2l+1}}}C_{n}^{\left( 2l+3/2 \right)}\left( \xi \right){{Y}_{\operatorname{lm}}}\left( \theta ,\phi \right)\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{{\rho }_{nlm}}\propto \frac{{{r}^{l}}}{r{{\left( 1+r \right)}^{2l+3}}}C_{n}^{\left( 2l+3/2 \right)}\left( \xi \right){{Y}_{\operatorname{lm}}}\left( \theta ,\phi \right),\end{eqnarray} \tag{ 4 }$$

where Y_lm are the spherical harmonics, $C_{n}^{\left( \alpha \right)}(\xi )$ are the Gegenbauer polynomials, and $\xi \equiv (r-1)/(r+1)$. Taking advantage of the bi-orthonormality of the basis functions, each A_nlm can be computed based on a given $\rho ({\boldsymbol{r}} )$ by

$$\begin{eqnarray}&&{{A}_{nlm}}\propto \int \rho ({\boldsymbol{r}} ){\Phi }_{nlm}^{*}({\boldsymbol{r}} )d{\boldsymbol{r}} .\end{eqnarray} \tag{ 5 }$$

The accelerations $-\nabla {\Phi }$ can be obtained by differentiating Equation (1) analytically.

In practice, a snapshot from an N-body simulation contains $\rho ({\boldsymbol{r}} )$, which is represented by discrete particles, and the computations are much more easily done in real space. We refer the reader to Hernquist & Ostriker (1992) for a straightforward recipe to compute accelerations given a list of particle positions and masses. A noteworthy feature of the SCF method is that all accelerations are proportional to the terms

$$\begin{eqnarray}&&\mathop{\sum }\limits_{k=1}^{N}{{m}_{k}}{{\tilde{{\Phi }}}_{nl}}\left( {{r}_{k}} \right){{P}_{\operatorname{lm}}}\left( {\rm cos} {{\theta }_{k}} \right)\left[ \begin{array}{ccccccccccccccc} {\rm cos} (m{{\phi }_{k}}) \\ {\rm sin} (m{{\phi }_{k}}) \\ \end{array} \right],\end{eqnarray} \tag{ 6 }$$

where m_k and $({{r}_{k}},{{\theta }_{k}},{{\phi }_{k}})$ are the mass and position of the kth particle, respectively, ${{\tilde{{\Phi }}}_{nl}}$ is the radial part of Equation (3), and P_lm are the Legendre polynomials. These linear summations over all N particles can be easily parallelized and precomputed given $(n,l,m)$. The decomposition is truncated at orders n_max and l_max (and m goes from 0 to l in real space), so both terms in Equation (6) are evaluated

$$\begin{eqnarray}&&\frac{{{n}_{{\rm max} }}({{l}_{{\rm max} }}+1)({{l}_{{\rm max} }}+2)}{2}\end{eqnarray} \tag{ 7 }$$

times. After Equation (6) is precomputed, the accelerations can be obtained analytically.

A detailed strategy to parallelize the SCF method has been proposed by Hernquist et al. (1995) to compute N-body dynamics using the SCF method itself. Note that our problem is even simpler because we are only extracting the accelerations from an existing N-body snapshot, and the accelerations do not need to be propagated back to the bodies. The extracted accelerations serve as "external forces" for each particle of our stream simulations.

SCF offers a convenient way to account for redshift dependence. Since the model of one static potential is entirely represented by the set of coefficients A_nlm in Equation (1), these coefficients can simply be interpolated in redshift. L11 showed that the time variation of the first few coefficients at the lowest orders are relatively small, but higher-order variations may be more difficult to capture. However, this is unlikely to be a problem for us because we account for subhalos separately from the smooth halo. We defer the construction of a redshift-dependent potential to future studies.

We use the Amiga Halo Finder ⁵ (AHF) to identify bound structures in the high-resolution region of VL-2. The most massive halo, hereafter the "main halo," is the halo that hosts the subhalos and the tidal streams of interest in this study. Many particles in the main halo are also bound into additional levels of substructures. In the first level there are 11,523 "subhalos" that consist of at least 150 particles. The minimum halo mass threshold is chosen such that the average peak circular velocity for those halos (${{V}_{{\rm max} }}\sim 2.5$ km s⁻¹) approximately corresponds to the peak circular velocity below which the subhalo abundance is suppressed by numerical limitations (Diemand et al. 2008; Diemand & Moore 2011). This threshold choice allows for a robust sampling of the subhalo mass function.

In addition to each subhalo's particle membership, AHF also automatically computes their intrinsic properties such as positions at the minimum potential, bulk velocities, scale radii, virial radii, etc. Table 1 summarizes the facts about our main halo, as well as the most and least massive subhalos. Note that the total number of particles and total mass of the main halo include all the particles of the subhalos, as returned by AHF. As we discuss in the following section, we remove the subhalo particles from the main halo, so the mass of the main halo in our simulations is slightly lower than the value stated in Table 1.

Table 1.  Main Halo and Subhalos in the z = 0 Snapshot of VL-2

Halo	N	M	rs	rvir
		(${{M}_{\odot }}$)	(kpc)	(kpc)
Main halo	$4.46\;\times \;{{10}^{8}}$	$1.89\;\times \;{{10}^{12}}$	60.4	400
Most massive subhalo	$1.15\;\times \;{{10}^{6}}$	$4.69\;\times \;{{10}^{9}}$	6.43	54.2
Least massive subhalo	150	$6.15\;\times \;{{10}^{5}}$	0.411	2.75

Note. The total number of particles (N), total masses (M), scale radii (r_s), and virial radii (r_vir) of three nominal halos extracted by the halo finder, which found a total of 11,523 subhalos inside the virial radius of the main halo.

Download table as:  ASCIITypeset image

In order to obtain a smooth halo from VL-2 that is a ΛCDM simulation, we remove all the subhalo particles from the main halo, leaving "voids" that were previously occupied by subhalos. When we construct the smooth halo's potential, we use a low-order decomposition that only captures the overall shape of the main halo, but not the lumpiness due to the removed subhalos.

For convenience, for the smooth halo we only consider the cases where ${{n}_{{\rm max} }}={{l}_{{\rm max} }}$, which we will refer to as the decomposition "order" collectively. Figures 1 and 2 show the accuracies of constructing the force field of the smooth halo using orders 10 and 20 by comparing against the force field that is obtained by directly summing up the force contributions from all the particles in the smooth halo. The force fields using both SCF and direct summation are computed in 94,747 randomly distributed positions inside a spherical shell of $r=[15,40]$ kpc in thickness where our stream will be orbiting.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

There are 660 and 4620 terms in the decompositions using orders 10 and 20, respectively (Equation (7)). Figure 3 shows that the improvements in accuracy by increasing the order from 10 to 20 are marginal. For the purpose of this study, it is not necessary to model VL-2 with the highest possible accuracy. Our goal is to study the effects that a dark matter halo has on tidal streams in general, but not to model any particular objects that have been observed.

Zoom In Zoom Out Reset image size

One alternative way to construct a smooth halo is to perform a decomposition with similar orders on all the main halo particles without removing the subhalo particles. We ran a similar analysis to Figures 1 and 2 using all particles, and we found that compared to having removed the subhalo particles, using all particles gave similar error levels except at positions near the subhalos, where the errors would typically be a few percent higher. This was expected because decompositions at orders 10 or 20 cannot model the forces inside the subhalos that were contributing to the directly summed forces. Removing the subhalo particles allowed us to confirm that our force errors were within 1% inside the smooth halo alone.

For the lumpy halo, we now add the forces by the subhalos back into the smooth halo. We decompose each subhalo using ${{n}_{{\rm max} }}=4,{{l}_{{\rm max} }}=0$. In addition to being simpler to compute, we prefer not to decompose subhalos at higher orders because some low-mass subhalos can have only a few hundred particles, so a low n_max and spherically symmetric decomposition can avoid unphysical clumpiness in each subhalo. Also, subhalos' effects on tidal streams are expected to be of short duration, so it is not necessary to model individual subhalos in detail.

Adding subhalos separately to the smooth halo provides two advantages. First of all, we can control the mass range of subhalos that are present in the galaxy. This can be used to test theories such as warm dark matter that suppress the formation of low-mass subhalos. Also, each individual subhalo can orbit freely around the smooth halo, so their encounters with tidal streams can be modeled realistically.

Adding subhalos separately introduces extra mass to the smooth halo when comparing lumpy and smooth halos. Our full range of subhalos have masses $1.067\;\times \;{{10}^{11}}{{M}_{\odot }}$ in total, which is about 6% of the mass of the smooth halo. In our simulations, we only use a subset of subhalos that (a) have pericentric distances of less than 40 kpc and (b) are less massive than ${{10}^{8}}{{M}_{\odot }}$. The reason for requirement (a) is that subhalos that do not approach the stream will have minimal effects on our streams with apocentric distances at 30 kpc in their orbits. The reason for requirement (b) is that encounters with massive subhalos do not leave interesting signatures. As shown in Yoon et al. (2011), massive subhalos are rarer than low-mass subhalos, so massive subhalos mostly influence a stream by distant encounters that do not cause gaps in the stream. On the other hand, when a massive subhalo makes a close encounter with a stream, the effect is likely catastrophic to the stream as opposed to leaving small gaps that have been observed (Carlberg et al. 2012; Carlberg & Grillmair 2013). A separate study for the effects of massive subhalos on streams in a realistic potential is currently under way. For this study, after imposing requirements (a) and (b), the lumpy halo has 3808 subhalos in our simulations, which have a total mass of $1.607\;\times \;{{10}^{10}}{{M}_{\odot }}$, or about 1% of the mass of the smooth halo. Therefore, the streams simulated with and without subhalos will only have slightly different orbits, and the difference can be safely ignored for our comparisons.

The streams' progenitor is a star cluster following a King profile with core size 0.01 kpc and w = 4.91, where w can be thought of as the ratio between the depth of the potential and central velocity dispersion of the cluster. Its total mass is $4.29\;\times \;{{10}^{4}}{{M}_{\odot }}$, central velocity dispersion is 2.4 km s⁻¹, and tidal radius is 0.103 kpc (identical to the one used in NC14). We use $N=1,000,000$ particles for the the main results of this study. Furthermore, as argued in NC14, the number of particles in the stream may be higher than the number of stars in the observed streams. Whether the effects we present in the following sections can be found in existing data requires more careful investigation, which is beyond the scope of this study.

We simulate the star cluster with two orbits shown in Figure 4—both start at $(x,y,z)=(30,0,0)$ kpc, but one orbit with $({{v}_{x}},{{v}_{y}},{{v}_{z}})=(0,120,0)$ km s⁻¹ ("Orbit 1"), and the other with $(0,41.0,113)$ km s⁻¹ ("Orbit 2"), which is the former with an initial velocity inclined at 70° from the xy-plane. The VL-2 main halo has not been aligned with any axis, so the directions of the initial velocities are arbitrarily chosen to explore the main halo as much as possible. Although neither orbit is confined in any orbital planes, both are tube orbits with apo- and pericentric distances at roughly ∼30 and ∼17 kpc, which are chosen to be similar to the inferred values of the GD-1 stream (Grillmair & Dionatos 2006; Willett et al. 2009).

Zoom In Zoom Out Reset image size

The SCF decomposition is performed using the procedure exactly described in Hernquist & Ostriker (1992). Since we consider decomposition orders that are much lower than the number of particles in VL-2, it is not necessary to apply a smoothening kernel to the VL-2 particles when computing the decomposition coefficients because a low-order decomposition, by construction, cannot capture the granularity of the VL-2 particles.

Our N-body simulations of streams are computed using Gadget-2 (Springel 2005), which is available to the public.⁶ We modify Gadget-2 so it uses the precomputed SCF decomposition coefficients to construct the external accelerations, which are added to the particles in our stream simulations after the stream particles' N-body forces have been computed. We impose a maximum time step of 1 Myr and softening of 5 pc for the particles so that they are essentially collisionless.

4.1.1. Orbital Convergence

The basis functions in the SCF method are complete, which means that the higher the decomposition order, the more accurate the constructed model will be. It is crucial to understand the decomposition order necessary for our purposes. Neither VL-2 nor our stream simulations are meant to be physical models of any observations, so it is not necessary for us to model the VL-2 main halo to the highest possible accuracy with a high decomposition order. In fact, since we have removed subhalos from the main halo, modeling these "voids" in the smooth halo may result in spurious forces. On the other hand, if the decomposition order is too low, then it compromises the value in using VL-2, which features shapes and profiles missing in idealized profiles. We require a decomposition order that is high enough to capture the shape and profile of the VL-2 smooth halo, but not high enough to capture the granularity of its substructures (the contribution by substructures will be added in after the smooth halo has been modeled).

The radial coordinates of the two stream orbits as functions of time for various decomposition orders are shown in Figure 5. Both orbits have converged to within a few percent everywhere for 10 Gyr. This is expected from Figures 1 and 2, which show that the forces almost everywhere are accurate to $\sim 1\%$ even at only order 10. Therefore, a decomposition order 10 is sufficient for our purposes. Using order 10 has the advantage that subhalo "voids" will not be captured. There are more than 10,000 subhalos found within the virial radius of the main halo, and the granularities of these subhalos cannot be modeled by polynomials of only the tenth degree in either the radial or angular components of SCF basis functions.

Zoom In Zoom Out Reset image size

L11 showed that decompositions lower than order 10 give significant errors in the region less than a few kiloparsecs away from the halo's center, so L11 adapted order 20 to suppress those errors. That region is not relevant for our streams because they orbit at galactocentric distances well beyond a few kiloparsecs. Nevertheless, our convergence results agree with L11, which showed that the acceleration errors compared to direct summation are at $\sim 1\%$ beyond 17 kpc, and that increasing from order 9 to 19 provides minimal gain at those regions.

4.1.2. Stream Density Convergence

As shown in NC14, signatures of subhalos can be detected in the linear density, or number of particles per unit length, along the stream. Even though the previous sections show that the accelerations and orbits have converged to within a few percent, we now investigate how well the stream density converges.

Figure 6 shows the linear density along four streams using 100,000 particles each and traced along a sky projection as seen by an observer at the galactic center (so the unit along the length of stream is angular) at 6.7 Gyr. The four streams started with identical progenitors in Orbit 1, but in the smooth halos that are constructed using orders 10, 15, 20, and 25. Although the overall profiles of the streams' densities on large scales are similar among the four orders, the detailed fluctuations are quite different. To investigate the significance of these differences, we simulated a few additional streams with different random realizations of the progenitor, but with the same number of particles, physical parameters, and orbit, in the same smooth halo potential using order 10. We find that the small-scale differences between these additional streams are indistinguishable from the differences between the four panels in Figure 6. Therefore, the differences shown in the figure are simply due to stochastic noise of the particles, and our decomposition has sufficiently converged at order 10.

Zoom In Zoom Out Reset image size

The convergence of our decomposition can also be seen in Figure 7, which shows the mass loss for each stream in Figure 6. Every large spike, or major mass loss event, in each stream is associated with a pericentric passage shown in the top panel of Figure 5. However, the minor mass-loss events do not occur at the same rate among the four orders. Our additional streams with different progenitor realizations as explained above show, again, that the small-scale differences for each decomposition order shown in Figure 7 are indistinguishable from stochastic noise.

Zoom In Zoom Out Reset image size

We adapt order 10 for the smooth halo for the rest of our results. As we show in the next sections, gaps caused by subhalos are longer than the spurious gaps caused by changing the decomposition order as shown in Figure 6. This is similar to the result of NC14, where the gap size distribution can be used to distinguish "intrinsic" gaps from subhalo gaps. In the future, we aim to repeat similar analyses to NC14 using various SCF decomposition orders for the realistic halo.

4.2.1. Stream 1 in the Smooth Halo

We first examine the intrinsic features in the stream in Orbit 1 (hereafter "Stream 1") in the smooth halo without any subhalos. Figure 8 shows the surface density of the stream projected on the sky from 8 to 8.5 Gyr as seen by an observer situated at the galactic center, and Figure 9 shows the linear density along the stream at the same times. This time frame roughly corresponds to one radial oscillation (though it is not a "radial period"; see Figure 5), so it allows us to investigate the streams' features as the streams stretch and compress during that oscillation.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The resulting stream is long and narrow, with a length-to-width ratio of about 50 (which is slightly lower than the observed values in GD-1 (Grillmair & Dionatos 2006)), so it can simply be represented by its linear density. The most prominent feature in the linear density is the spikes located within 10° away from either side of the progenitor. Upon inspection, we find that new spikes develop immediately after pericentric passages, and each spike originates from the progenitor and migrates toward the ends of the stream. One possible cause of the density spikes is EOs, which are caused by the "piling up" of particles in their epicyclic orbits as they escape from the progenitor. In the original derivations in Küpper et al. (2008, 2010), the spacings between the EOs are calculated assuming an axisymmetric potential. Although the VL-2 halo has no spatial symmetry, we can approximate it by spherically averaging our potential model (i.e., taking ${{l}_{{\rm max} }}=0$ in Equation (1)). Using a representative value of R for the angular frequency Ω and epicyclic frequency κ, we can estimate the spacing

$$\begin{eqnarray}&&|{{y}_{C}}|=\frac{4\pi {\Omega }}{\kappa }\left( 1-\frac{4{{{\Omega }}^{2}}}{{{\kappa }^{2}}} \right){{x}_{L}},\end{eqnarray} \tag{ 8 }$$

where ${{x}_{L}}={\rm GM}/(4{{{\Omega }}^{2}}-{{\kappa }^{2}})$, between the first EOs on either branch of the stream (Küpper et al. 2008). For R = 24 kpc (roughly the galactocentric radius of the progenitor at 8.3 Gyr), we obtain ${\Omega }\simeq 8$ and $\kappa \simeq 12\;{\rm Gy}{{{\rm r}}^{-1}}$. This gives $|{{y}_{C}}|\simeq 0.77$kpc, or 18 when projected onto the sky, which is almost a factor of 2 smaller than the measured spacing at $\sim 3{}^\circ$ between the first two density spikes at 8.3 Gyr. Note that as the stream compresses and stretches as it oscillates radially, the spacing between the spikes can easily differ by factors of a few (e.g., compare 8.2 and 8.4 Gyr in Figure 9).

The density profiles along the stream shown in Figure 9 indicate that the density spikes near the progenitor are very prominent inside the VL-2 potential and are similar to those in NC14 in a spherical potential—they are regularly spaced and appear very close to the progenitor. As we shall see in the next sections, gaps caused by subhalos are randomly spaced and can occur anywhere along the stream.

4.2.2. Stream 1 in the Lumpy Halo

Figures 10 and 11 show the sky projection and linear density of Stream 1 in a lumpy halo with 3807 subhalos. Compared to the same stream in the smooth halo, the density clearly shows gaps as local minima that are more pronounced and occur everywhere along the stream. The density spikes near the progenitor are still present with similar spacing to the stream in the smooth halo. This is not surprising, since the extra subhalos contribute only $\sim 1\%$ of the total mass of the smooth halo, so any effects due to the stream's orbit would not be different from the smooth halo case.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Recall from Section 2.4 that our simulations in the lumpy halo contains 3808 subhalos. An interesting but rare event that occurred in the simulation with 3808 subhalos was that the progenitor had a close encounter at less than one half-mass radius from a $1.8\;\times \;{{10}^{7}}{{M}_{\odot }}$ subhalo at an early time, which caused a sudden burst in mass loss for the progenitor. This resulted in a stream that was dominated by two large but smooth spikes, one in each branch of the stream. When this subhalo was eliminated, we were able to resolve the density fluctuations caused by close encounters between the tidal tails and subhalos as shown in Figures 10 and 11, rather than between the progenitor and subhalos. In the future we aim to study the probability that a progenitor or its tidal stream becomes catastrophically perturbed by subhalos in the VL-2 halo.

4.2.3. Stream 2 in the Smooth Halo

Figure 12 shows the sky-projected density of the stream in Orbit 2 (hereafter "Stream 2") from 6.0 to 6.7 Gyr in the smooth halo. Similar to our investigation for Stream 1, this time frame also covers roughly one radial oscillation even though a well-defined radial period does not exist in the VL-2 potential (Figure 5). Compared to Stream 1, Stream 2 becomes "fluffy" at the far ends of its tidal tails, with widths sometimes comparable to the length of the narrow part of the stream. This is very different from Stream 1 shown in Figure 8, which has a high length-to-width ratio throughout.

Zoom In Zoom Out Reset image size

An important implication is that only the narrow segment of the stream may be dense enough to be observable, and the rest of the stream may be too diffuse. Simulations in Pearson et al. (2015) in triaxial potentials also showed similarly diffuse features dubbed "stream-fanning" in the tips of tidal tails, which are not found in observations. The physical origin of stream-fanning has yet to be understood and is beyond the scope of our study, so it is not clear whether our Stream 2 is exhibiting stream-fanning.

Because of the diffuse tails, it is no longer appropriate to analyze the entire stream simply by the linear density. Nevertheless, we can trace the linear density along the narrow part of the stream out to about 10°–20° away from the progenitor. Figure 13 shows the densities along the stream at each time shown in Figure 12. The density profiles show prominent density spikes next to the progenitor, and the spacings between those spikes are in rough agreement with the derivation in Section 4.2.1. Other than this, the density profiles fall off very quickly and smoothly toward the diffuse ends of the stream.

Zoom In Zoom Out Reset image size

4.2.4. Stream 2 in the Lumpy Halo

Figure 14 shows the sky-projected density of Stream 2 from 6.0 to 6.7 Gyr in a lumpy halo with 3808 subhalos. The stream resembles the case in the smooth halo (Figure 12) with very fluffy ends, but in a lumpy halo they are even more prominent. In the fluffy regions the width can exceed the length of the narrow regions, which span only 10°–20°.

Zoom In Zoom Out Reset image size

The linear density along the stream can be plotted, but similar to the case in the smooth halo, it is only appropriate to the narrow regions near the progenitor. As seen in Figure 15, the density profiles in those regions are also smooth with the exception of the first or second density spikes. Therefore, for this stream, there are no distinguishable signatures in the linear density that traces the existence of subhalos.

Zoom In Zoom Out Reset image size

Perhaps a more interesting aspect of this stream is the tips of its tidal tails. The wide and diffuse tails shown in Figure 14 are reminiscent of "shell"-like features found in satellite systems that plunge near the galactic center on very eccentric orbits. With apo- and pericentric distances at ${{r}_{a}}\simeq 30$ and ${{r}_{p}}\simeq 17$ kpc, respectively, the stream's eccentricity is $e=({{r}_{a}}-{{r}_{p}})/({{r}_{a}}+{{r}_{p}})\simeq 0.3$. This implies that satellites do not need to be on eccentric orbits to have significant parts of their tidal tails disrupted such that their surface brightnesses may be even lower than that of the narrow stream itself, which is already difficult to observe. In NC14, the same progenitor orbiting a spherical and idealized halo with similar eccentricity to Streams 1 and 2 here could only produce a narrow stream. This indicates that using a realistic halo is important when studying the effects that subhalos have on GD-1-like streams. In a future study, we aim to simulate the tidal disruptions of globular-cluster-type satellites in many orbits in both the smooth and lumpy halos in order to study the survival rate of tidal tails.