Resonant stirring of dwarf galaxies

# The resonant nature of tidal stirring of disky dwarf galaxies orbiting the Milky Way

## Abstract

Using -body simulations we study the tidal evolution of initially disky dwarf galaxies orbiting a Milky Way-like host, a process known to lead to the formation of dwarf spheroidal galaxies. We focus on the effect of the orientation of the dwarf galaxy disk’s angular momentum with respect to the orbital one and find very strong dependence of the evolution on this parameter. We consider four different orientations: the exactly prograde, the exactly retrograde and two intermediate ones. Tidal evolution is strongest for the exactly prograde and weakest for the exactly retrograde orbit. In the prograde case the stellar component forms a strong bar and remains prolate until the end of the simulation, while its rotation is very quickly replaced by random motions of the stars. In the retrograde case the dwarf remains oblate, does not form a bar and loses rotation very slowly. In the two cases of intermediate orientation of the disk, the evolution is between the two extremes, suggesting a monotonic dependence on the inclination. We interpret the results in terms of the resonance between the angular velocity of the stars in the dwarf and its orbital motion by comparing the measurements from simulations to semi-analytic predictions. We conclude that resonant effects are the most important mechanism underlying the tidal evolution of disky dwarf galaxies.

galaxies: dwarf — galaxies: fundamental parameters — galaxies: kinematics and dynamics — galaxies: structure — Local Group

## 1. Introduction

The formation of dwarf spheroidal (dSph) galaxies in the Local Group remains an open question but one of the most promising scenarios for their origin is via the tidal interaction of their disky progenitors with more massive hosts like the Milky Way. The scenario, proposed by Mayer et al. (2001) explains the morphology-density relation observed among the dwarfs of the Local Group and accounts for the non-sphericity of the dSph objects.

The efficiency of the mechanism and its observational predictions have been investigated in detail by Klimentowski et al. (2009), Kazantzidis et al. (2011) and Łokas et al. (2011, 2012). These studies explored the dependence of the process on a large number of orbital and structural parameters of the dwarf. The general picture that emerged from these studies is that a disky dwarf progenitor, once accreted by a massive host, undergoes strong tidal stirring and mass stripping if the orbit is tight enough.

Typically, at the first pericenter passage, the disk transforms into a tidally induced bar (for a detailed description of the properties of such a bar see Łokas et al. 2014a). The bar becomes thicker and shorter in time leading in the end to the formation of a spheroidal stellar component. The morphological transformation is accompanied by strong changes in the kinematics as quantified by the amount of ordered to random motion. The latter starts to dominate at some point and at the end the galaxy is pressure supported.

Among the parameters expected to have a strong impact on the evolution is the inclination between the angular momentum of the dwarf galaxy disk and its orbital angular momentum. However, in the studies mentioned above only a narrow range of inclinations was studied in detail, namely those with values , and . Because of this range, and the way the properties of the dwarf were measured, no clear evidence for the dependence on this parameter was found. On the other hand, the difference between the prograde and retrograde galaxy encounters has been recognized for a long time and known to lead to very different outcomes (e.g. Holmberg 1941; Hénon 1970; Kozlov et al. 1972; Toomre & Toomre 1972; Keenan & Innanen 1975). The issue has been recently addressed again by D’Onghia et al. (2009, 2010) using an improved version of the impulse approximation applied to rotating systems. Although this approximation is not directly applicable to our simulations, we attempt a comparison between these results and the numerical ones.

In this paper we aim at clarifying the issue of the dependence of the results of tidal encounters between dwarfs and their hosts on the inclination of the dwarf’s disk. For this purpose we performed four simulations of tidal evolution of a dwarf galaxy orbiting a Milky Way-like host with disk inclinations , , and . The angles of and correspond to exactly prograde and exactly retrograde orientations of the dwarf’s disk. We also measured the properties of the dwarf galaxy in a different way that enables clear comparisons between different runs. Preliminary results of this study, using lower resolution simulations, were discussed in Łokas & Semczuk (2014).

The paper is organized as follows. In section 2 we present the simulations used in this study. In section 3 we discuss the properties of the dwarf galaxies as they evolve in time, focusing on their kinematics, morphology and density profiles. Section 4 compares the results of the simulations to the predictions of semi-analytic models. The conclusions follow in section 5.

## 2. The simulations

The simulations used in this study were similar to those described in detail in Łokas et al. (2014a). Here we therefore provide only a short summary. The initial conditions for the simulations consisted of -body realizations of two galaxies: the Milky Way-like host and the dwarf galaxy, generated via procedures described in Widrow & Dubinski (2005) and Widrow et al. (2008). Both galaxies contained exponential disks embedded in NFW (Navarro et al. 1997) dark matter haloes, each made of particles ( total). We note that the results presented here differ only slightly (are less noisy) from those of lower resolution simulations in Łokas & Semczuk (2014) where a smaller number of particles was used (by a factor of five). We are therefore confident that our present resolution is sufficient to grasp all the essential features of the evolution.

The dwarf galaxy model had a dark halo of mass M and concentration . Its disk had a mass M, an exponential scale-length kpc and thickness . The model is stable against formation of the bar in isolation for the time scales of interest here. The host galaxy was similar to the model MWb of Widrow & Dubinski (2005). It had a dark matter halo of mass M and concentration . The disk of the host had a mass M, the scale-length kpc and thickness kpc. The disk was also stable against bar formation to avoid strong variations of the host potential in time. The disk of the Milky Way was coplanar with the orbit of the dwarf. Although this may seem contrary to observational constraints where most of satellite orbits are found to be polar (e.g. Pawlowski & Kroupa 2014), the choice was motivated by the necessity to avoid any additional variability which may be due to the passages through the plane of the Milky Way disk. However, we have performed an additional simulation to verify that for the orbits used here the evolution of the dwarf depends very weakly on the orientation of the orbit with respect to the Milky Way disk.

The dwarf galaxy was initially placed at an apocenter of a typical, eccentric orbit around the Milky Way with apo- to pericenter distance ratio of kpc. The initial position was at the coordinates kpc of the simulation box and the velocity vector of the dwarf was toward the negative direction. We performed four simulations with different dwarf disk orientations with respect to the orbit: two coplanar with the orbit (prograde and retrograde) and two perpendicular to the orbit with angular momenta in the same and opposite direction to the dwarf’s orbital velocity. The different initial conditions, in particular the components of the unit angular momentum vector, are listed in Table 1. We will refer to the simulations by names indicating the initial inclination of the disk, I0, I90, I180 and I270, where the inclination is measured as the rotation angle around the axis of the simulation box.

The evolution of the system in each simulation was followed for 10 Gyr using the GADGET-2 -body code (Springel et al. 2001; Springel 2005) with outputs saved every 0.05 Gyr. The adopted softening scales were kpc and kpc for the disk and halo of the dwarf while kpc and kpc for the disk and halo of the host, respectively.

## 3. Evolution of the dwarfs

In this section we look at the inner properties of the dwarf galaxies as they are transformed by the tidal forces from the Milky Way. All measurements discussed below were made for stars (and dark matter particles) within the radius of 0.5 kpc from the center of the dwarf.

### 3.1. Mass content

We begin the analysis by measuring the mass inside this radius. Figure 1 compares the mass of stars (upper panel) and dark matter (lower panel) in the four simulations. As expected, the dark mass content decreases systematically, most significantly at pericenters that occur at , 3.3, 5.5, 7.6, 9.7 Gyr from the start of the simulation and there is no dependence on the inclination because the halos were spherical and isotropic in all cases.

On the other hand, for the stellar component (upper panel) there is a significant difference between dwarfs with varying initial inclination. However, the dependence may seem surprising because we expect the prograde disk to be more stripped, while the opposite is seen in the Figure: the stellar mass for the I0 simulation is even increased after the first pericenter passage, while it is decreased in the remaining cases. While the tidal stripping is indeed stronger in the outer parts, in the inner region we probe by this measurement, the stellar content increases as a result of the significant change in the structure of the stellar component due to the formation of a tidally induced bar.

### 3.2. Kinematics

The main criteria usually applied in order to verify if a dwarf galaxy transformed from a disky object to a dSph are based on kinematics and shape of the stellar component. A dSph galaxy is supposed to be characterized by the dominance of random motions of the stars over the amount of rotation and its shape should be sufficiently close to spherical. In this and the following subsection we look in detail at the evolution of these quantities.

In order to measure these properties, for each simulation output we determine the directions of the principal axes of the stellar component from stars within the radius of 0.5 kpc using the inertia tensor and rotate these stars so that the new coordinate system is aligned with the principal axes (the axis is along the longest, the axis along the intermediate and the axis along the shortest axis of the stellar component). We then introduce a standard spherical coordinate system such that measures the angle around the axis and the angle from the axis towards the plane.

The kinematic properties were estimated using these coordinates. In the top panel of Figure 2 we plot the rotation around the shortest axis as a function of time. We note that there is no other significant streaming motion along the other spherical coordinates or around the two other principal axes (see Łokas et al. 2014b for a brief discussion of this issue). Clearly, the rotation is decreasing most strongly for the exactly prograde inclination of the disk (I0), a significantly smaller decrease is seen for the perpendicular orientations (I90 and I270) and the effect is the weakest for the retrograde case (I180). Note that the decrease of rotation is not steadily monotonic even in the exactly prograde case (I0). This is due to the tidal torques acting on the bar at pericenter passages that can speed up or slow down the bar depending on its particular orientation at this moment (see Łokas et al. 2014a for details).

The behavior of the velocity dispersion is exactly the opposite. In the middle panel of Figure 2 we show the 1D velocity dispersion calculated as . The increase of at the first pericenter passage is strongest for the I0 case, intermediate for I90 and I270 and for I180 remains constant in time or even slightly decreases due to mass loss. The ratio shown in the lower panel of Figure 2 decreases for all simulations, but reaches a value significantly below unity only for the prograde case. For the three remaining cases a substantial amount of rotation is retained, although the hierarchy of lower for more prograde cases is preserved.

In Figure 3 we show the evolution of the different velocity dispersions (upper panel), (second panel) and (third panel) as a function of time. In all cases the increase of a given dispersion at the first pericenter passage and its level at later times is highest for the I0 case, intermediate for I90 and I270 and non-existent for I180. In addition, this increase is the most abrupt for the prograde dwarf, while for the intermediate inclinations I90 and I270 the increase occurs much more slowly in time and takes about half the orbital period between the first and second pericenter passage. Interestingly, significant increase is seen in all dispersions in spite of the fact that, due to the formation of the bar, one could expect the radial to increase much more significantly.

We end the analysis of the kinematics by plotting in the lower panel of Figure 3 the value of the anisotropy parameter where is the tangential second moment including rotation. A systematic dependence of this parameter on the inclination of the disk is also present: the stellar orbits are most radial for the I0 case with almost constantly at the level of 0.5, corresponding to mildly radial orbits of the stars, characteristic of the bar. For the intermediate cases I90 and I270 the values stay between 0 and 0.4, while for the retrograde case I180 remains negative due to the dominant presence of rotation. Interestingly, at the end of the evolution all three non-prograde cases I90, I180 and I270 have almost isotropic orbits () although they reach this special value via different evolutionary paths.

### 3.3. Shapes

Figure 4 illustrates the evolution of the shape of the stellar component of the dwarfs in time. In the first and second panels from top we plot the axis ratios (intermediate to longest) and (shortest to longest). The thickening of the dwarf, as quantified by the increasing value of is similar in all cases, although for the prograde case decreases after the first pericenter due to the formation of the bar. Much more significant differences are seen in the evolution of . This value decreases most strongly for the prograde I0 case signifying the prolate shape characteristic of the bar. The intermediate inclinations I90 and I270 lead to less prolate shapes, while I180 remains disky for the whole evolution, as indicated by remaining constantly close to unity.

The shape can be also quantified in terms of the triaxiality parameter which is shown in the third panel of Figure 4. The values of the parameter at all times confirm the prolate shape due to the bar in the case of I0. For the intermediate cases I90 and I270 we get indicating a triaxial shape. For the retrograde case I180 the parameter remains low, , at all times as is characteristic of a disk.

The presence of a bar is usually detected by measuring the bar mode of the Fourier decomposition of the stars projected along the shortest axis of the stellar distribution (see a more detailed discussion in Łokas et al. 2014a). Usually, is considered as high enough to be interpreted as a bar. As we can see in the lower panel of Figure 4 this is always the case for simulation I0 after the first pericenter passage and also for some significant periods of time for the intermediate cases I90 and I270, which means that the bar also forms there, but it is much weaker. The exactly retrograde disk does not form a bar as its at all times. Slight temporary increases of this value are due to stretching of the dwarf at the pericenters.

These measurements are confirmed by the maps of the surface density distribution of the stars in the four dwarfs plotted in Figure 5. The distributions are shown in projection along the intermediate axis of the stellar component so in their most non-spherical appearance. The snapshots were selected for the time Gyr after the start of the simulations, corresponding to the last apocenter passage. In all cases, except for the retrograde one I180, the remnant of the bar formed after the first pericenter passage is still visible in the inner parts. In addition, the distribution of the stars in the prograde case I0 is much less diffuse in the outer parts of the maps. In spite of the fact that this dwarf evolves most strongly, its neighborhood is not uniformly filled with stripped debris because the lost stars form well-defined, narrow tidal streams.

### 3.4. Density profiles

The dependence of the properties of the dwarfs on the initial inclination of the disk also manifests itself in the evolved density profiles. In Figure 6 we show examples of the stellar density profiles (different colors) and the dark matter density profile (black line) for different simulations considered here, measured at the fourth apocenter. The transition from the bound component to the tidal tails is visible as the break in the slope of the density profiles. This transition is however only well defined for the exactly prograde case I0 (and the dark matter profile). In this case the transition from the steeper to the shallower profile (where the slope is around ) occurs at around 3 kpc. For other simulations no such clear break radius is seen.

As discussed in Łokas et al. (2013) using similar (but only mildly prograde) simulation setups, the break radii can be interpreted as the tidal radii. In this case, the dependence on the orbit of the star within the satellite is expected (Keenan & Innanen 1975; Read et al. 2006) and we will attempt a detailed comparison in a follow-up paper. Since the stellar density profiles of our simulated dwarfs do not show clear signatures of the break radius (except for the exactly prograde case) here we propose a comparison in terms of density. The stellar profiles shown in Figure 6 demonstrate clear hierarchy: at the outer radii (larger than 1 kpc) the stellar density profile of I180 (green line) is above all the other profiles, the one of I0 (red lines) is the lowest, and the ones of I90 and I270 fall exactly on top of each other and between the other two. This means that we can quantify the amount of tidal stripping in these different cases by measuring the density of the stars, rather than the break radius.

To do so, we calculated the mean density of stars in the shells of radii 2 kpc 4 kpc as a function of time for different simulations. The results are shown in the upper panel of Figure 7. In the lower panel we plot the radius at which the density of stars drops below Mkpc as a function of time. In both plots there is a clear systematic difference between the measurements for different initial disk orientations: the stars are stripped more effectively on prograde orbits as demonstrated by the lower outer densities and smaller radii where the density drops to a fixed value. Although the values for the measurements were chosen in an arbitrary way, we expect the results to be similar if these parameters are slightly modified. They clearly confirm that the amount of tidal stripping depends very strongly on the initial inclination of the dwarf’s disk.

## 4. Comparison with semi-analytic predictions

To describe the encounter between a disky dwarf galaxy and the Milky Way we use the impulse approximation as discussed in Binney & Tremaine (1987) and D’Onghia et al. (2010). According to this approximation the th component of the acceleration of each dwarf’s star with respect to the center of mass is given by

 ˙vk=−∑j∂2ψ∂xk∂xj∣∣∣x=0xj, (1)

where is the Milky Way’s (i.e. perturber’s) potential and are Cartesian coordinates. We integrate equation (1) over a finite time period to obtain velocity increments that can be compared with increments measured from simulations

 Δvk=Δt∫0˙vkdt. (2)

We work in a frame centered on the center of the dwarf’s mass, the dwarf’s disk lies in plane and position of the Milky Way at the pericenter is . In this frame the trajectory of each star originating from the dwarf is

 x=[rcos(Ωt+ϕ0),rsin(Ωt+ϕ0),0], (3)

where is the radius, is the angular velocity and is the initial azimuthal angle of the star. The trajectory of Milky Way is given by

 X=Vt+X0, (4)

where and are constant vectors, fitted to mimic the perturber’s trajectory from simulations as a straight line during a given time period. Note that and are different for each of our simulations as they depend on the inclination . This dependence can be found by rotating the trajectory for the prograde case with matrix defined by Euler angles, in order to obtain trajectories in other cases. One of the Euler angles is the inclination and other two depend on the perturber’s orbit.

To approximate the gravitational potential of the Milky Way we sum the potential from its stars and the dark matter halo. The first part is represented as a point-mass potential, while the second is given by the NFW profile

 ψ=−GMD|x−X|−gGMHln(1+c|x−X|/rv)|x−X|, (5)

where is the mass of the Milky Way disk, is the virial mass of its halo, is the virial radius, is the concentration parameter and (see Łokas & Mamon 2001). Substituting equations (1), (3), (4) and (5) into (2) we obtain formulae which can be numerically integrated to get velocity increments.

In order to calculate the velocity increments described above we estimate the angular velocity in the dwarf as . We calculate integrals of equation (2) over a small period of time, so that the assumptions concerning the trajectories are valid. In our simulations the outputs were saved every 0.05 Gyr and we choose to integrate over this time to compare our predictions with velocity increments occurring in the simulations when the tidal force is the strongest, i.e. between the output preceding the first pericenter and the one as close as possible to this pericenter. We check how the distribution of increments changes with the distance from the center of the dwarf. The results for all simulations are summarized in Figure 8.

The red points in Figure 8 represent the values predicted by the impulse approximation, while the blue points correspond to the values measured from the simulations. At radii smaller than 2 kpc the semi-analytical predictions do not reproduce the distribution of velocity increments because at these radii the velocity changes are dominated by the gravitational potential of the dwarf which was not included in the predictions. However, for radii greater than 2 kpc the agreement between the full simulations and semi-analytic predictions for this short time period is very good. In particular, in the upper left panel of Figure 8 for simulation I0 we find two very well-defined branches corresponding to stars on different sides of the dwarf. Branches from simulations are not exactly symmetric with respect to zero while the branches from theoretical predictions are. The difference is due to the fact that the analytic predictions only take into account the lowest order terms.

For simulations I270 the increments are almost identical as for I90, as expected due to symmetry of the two configurations with respect to the orbital plane. In some of the panels in Figure 8 the velocity increments for radii larger than 2 kpc are approximately zero. However, the negligible values are consistently obtained both from the simulations and the semi-analytic calculations.

We further illustrate these results in Figure 9 where we plot the surface density maps of the stellar component and the velocity vectors for a random sample of a hundred stars at radii larger than 2 kpc. The plots show the dwarfs at the first pericenter passage, i.e. after they have been affected by a tidal impulse from the Milky Way for the first time. The coordinate system is as defined above, with coordinates in the plane of the dwarf’s disk.

The comparison of the upper left panel of Figure 9, corresponding to simulation I0, to the other three confirms that for the prograde encounter the effect of the tidal force in the strongest: the dwarf galaxy disk is already strongly distorted toward a bar-like shape and two tidal arms are formed. There are significant increments of velocity along the axis. Small increments of velocities are also visible along the axis in the right panels of Figure 9 corresponding to simulations I90 and I270. The effect of the tidal force is weakest for the retrograde case (lower left panel of Figure 9) where the dwarf’s initial disk is affected very little.

As discussed by D’Onghia et al. (2010) strongest tidal interactions occur when the intrinsic angular velocities of stars in the dwarf’s disk are comparable to the angular velocity of the satellite on its orbit. This condition can be written as

 Ωdisk≃Ωorb (6)

and we can define the resonance parameter

 α=|Ω|Ωorb, (7)

that should be of the order of unity for the strongest, resonant response.

To demonstrate the resonant nature of the tidal effects in our simulations we measured for each simulation output and found radii from the center of the dwarf at which . The time dependence of this radius for simulation I0 is shown in Figure 10. We can see that the variability of this radius reflects the varying orbital velocity of the dwarf resulting in rather large values at apocenters and much smaller ones at pericenters. The slow variation over time scales much larger than the orbital period is caused by the mass loss and decreasing . Whenever the dwarf reaches the pericenter of its orbit around the Milky Way, this characteristic radius drops below 2 kpc, and even down to 1 kpc at the later pericenters. Note that these radii are of the order of 2-3 half-light radii of the stellar component which means that a significant fraction of stars is affected. Comparing with Figure 8 showing the velocities at the first pericenter we confirm that this radius (equal to 1.7 kpc at this time) is exactly where the tidal effects start to prevail over the dwarf’s potential.

## 5. Conclusions

In this work we extended previous studies of the efficiency of the tidal stirring mechanism to include the dependence on the initial inclination of the dwarf galaxy disk with respect to its orbit around the Milky Way. Our simulation setups involved a dwarf galaxy placed on a typical, eccentric orbit around a Milky Way-like host and its evolution was followed for 10 Gyr. We considered four configurations, an exactly prograde, an exactly retrograde and two intermediate orientations of the disk. We found the efficiency of tidal stirring to be very strongly dependent on this inclination.

The effect of the tidal interaction turns out to be the strongest for the exactly prograde orientation of the dwarf’s disk (I0). In this case the disk transforms into a strong bar () at the first pericenter passage and a remnant bar is retained until the end of the evolution. Although the bar becomes weaker with time, the shape of the stellar component is consistently prolate, but tending to a spherical toward the end of the evolution. This morphological transformation is mirrored in the kinematics by the gradual decrease of the dwarf’s rotation velocity. Already at the first pericenter passage the rotation drops significantly and the velocity dispersion increases, mostly in the radial direction due to the formation of the bar. During subsequent pericenter passages the rotation is further decreased down to rather low values at the end of the simulation where . Thus the streaming motions are almost completely replaced by random motions of the stars. Interestingly, the radial velocity dispersion dominates the whole time, which manifests itself in the anisotropy parameter close to even at the end.

In the two cases of perpendicular orientations of the dwarf’s disk with respect to the orbit (I90 and I270) the evolution is to some extent similar to the I0 case. In these configurations the bar also forms after the first pericenter passage, but more slowly (over a time scale of about 1 Gyr corresponding to half the orbital period) and is significantly weaker (), also in the subsequent evolution. The overall shape of the stellar component can be characterized more as triaxial than decidedly prolate. The transition from the streaming to random motions of the stars also happens less efficiently with only slightly below unity at the final outputs. The anisotropy parameter is close to zero due to the contribution of the still significant rotation.

In the exactly retrograde case (I180) no strong evolution is present: the dwarf’s stellar component does not form a bar and remains disky. The only signatures of tidal evolution in this case are the mass loss (similar as in other cases, mostly in dark matter), small decrease of the rotation velocity, slight evolution of the anisotropy parameter from negative toward isotropic and a non-negligible thickening of the disk. The difference between this and the other cases is also visible in the stellar density profiles which are less affected and do not show any clear transition from the bound component to the tidal tails.

We have interpreted these changes in the context of the resonant stripping mechanism recently discussed by D’Onghia et al. (2009, 2010). In particular, we calculated the velocity increments the dwarf’s stars should experience at the first pericenter and compared them with the direct measurements from simulations. We find a very good agreement between the two, confirming the interpretation that the evolution we see in the full -body treatment is indeed due to the orientation of the dwarf’s disk. The resonant nature of the phenomenon is further confirmed by the behavior of the ratio between the angular velocity of the stars in the dwarf and the angular velocity of its orbital motion. This ratio turns out to be of the order of unity only near the pericenters and this is indeed when the tidal effects are the strongest.

The results presented here suggest that the most important mechanism underlying the tidal evolution of disky dwarfs orbiting a bigger galaxy is indeed of resonant nature. We propose to refer to the processes of morphological and dynamical evolution of the dwarfs we described as ‘resonant stirring’ in analogy to the ‘resonant stripping’ mechanism found by D’Onghia et al. (2009, 2010) to increase the mass loss in similar configurations. As discussed by D’Onghia et al. (2010), the resonance is broad, hence the name ‘quasi-resonant stirring’ would be more appropriate. In physical terms, this resonance can be traced to the fact that the stars with remain for an extended period of time on the line joining the dwarf galaxy and the perturber. For these stars the tidal force (which is strongest along this line) has the longest time to operate which results in the largest velocity increments and the largest stirring. In the context of the tidal radius calculations, the difference between the prograde and retrograde cases comes from the change of sign of the Coriolis force. The relation between the two approaches remains to be investigated and we plan to address this issue in our future work.

## Acknowledgments

This work was supported in part by PL-Grid Infrastructure, the Polish National Science Centre under grant 2013/10/A/ST9/00023, by US National Science Foundation Grant No. PHYS-1066293 and the hospitality of the Aspen Center for Physics. ED gratefully acknowledges the support of the Alfred P. Sloan Foundation and of the NSF Grant No. AST-1211258 and ATP-NASA Grant No. NNX14AP53G. MS benefited from the summer student program of the Copernicus Center. We thank L. Widrow for providing procedures to generate -body models of galaxies for initial conditions.

### Footnotes

1. affiliation: Nicolaus Copernicus Astronomical Center, Bartycka 18, 00-716 Warsaw, Poland
2. affiliation: Nicolaus Copernicus Astronomical Center, Bartycka 18, 00-716 Warsaw, Poland
3. affiliation: Warsaw University Observatory, Al. Ujazdowskie 4, 00-478 Warsaw, Poland
4. affiliation: Nicolaus Copernicus Astronomical Center, Bartycka 18, 00-716 Warsaw, Poland
5. affiliation: Department of Astronomy, University of Wisconsin-Madison, Madison, WI 53706, USA
6. affiliation: Alfred P. Sloan Fellow

### References

1. Binney, J., & Tremaine, S. 1987, Galactic Dynamics (Princeton, NJ: Princeton Univ. Press)
2. D’Onghia, E., Besla, G., Cox, T. J., & Hernquist, L. 2009, Nature, 460, 605
3. D’Onghia, E., Vogelsberger, M., Faucher-Giguere, C. A., & Hernquist, L. 2010, ApJ, 725, 353
4. Hénon, M. 1970, A&A, 9, 24
5. Holmberg, E. 1941, ApJ, 94, 385
6. Kazantzidis, S., Łokas, E. L., Callegari, S., Mayer, L., & Moustakas, L. A. 2011, ApJ, 726, 98
7. Keenan, D. W., & Innanen, K. A. 1975, AJ, 80, 290
8. Klimentowski, J., Łokas, E. L., Kazantzidis, S., Mayer, L., & Mamon, G. A. 2009, MNRAS, 397, 2015
9. Kozlov, N. N., Syunyaev, R. A., & Éneev, T. M. 1972, SPhD, 17, 413
10. Łokas, E. L., & Mamon, G. A. 2001, MNRAS, 321, 155
11. Łokas, E. L., Kazantzidis, S., & Mayer, L. 2011, ApJ, 739, 46
12. Łokas, E. L., Majewski, S. R., Kazantzidis, S., et al. 2012, ApJ, 751, 61
13. Łokas, E. L., Gajda, G., & Kazantzidis, S. 2013, MNRAS, 433, 878
14. Łokas, E. L., & Semczuk, M. 2014, in Proc. XXXVI Polish Astronomical Society Meeting, ed. A. Różańska, & M. Bejger (Warsaw, Poland: PAS), 53
15. Łokas, E. L., Athanassoula, E., Debattista, V. P., et al. 2014a, MNRAS, 445, 1339
16. Łokas, E. L., Ebrová, I., del Pino, A., & Semczuk, M. 2014b, MNRAS, 445, L6
17. Mayer, L., Governato F., Colpi, M., et al. 2001, ApJ, 559, 754
18. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 (NFW)
19. Pawlowski, M. S., & Kroupa, P. 2014, ApJ, 790, 74
20. Read, J. I., Wilkinson, M. I., Evans, N. W., Gilmore, G., & Kleyna J. T. 2006, MNRAS, 366, 429
21. Springel, V. 2005, MNRAS, 364, 1105
22. Springel, V., Yoshida, N., & White, S. D. M. 2001, New Astronomy, 6, 79
23. Toomre, A., & Toomre, J. 1972, ApJ, 178, 623
24. Widrow, L. M., & Dubinski, J. 2005, ApJ, 631, 838
25. Widrow, L. M., Pym, B., & Dubinski, J. 2008, ApJ, 679, 1239
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minumum 40 characters