On the destruction and over-merging of dark halos in dissipationless N-body simulations

Ben Moore, Neal Katz & George Lake
Department of Astronomy, FM-20,
University of Washington, Seattle, WA 98195, USA.

Abstract:

1 INTRODUCTION

Understanding how structure evolves in the universe is a fundamental problem in cosmology. Cold dark matter dominated universes have been extensively investigated with some success at reproducing observations from galactic to cluster scales (e.g. Davis et al 1985, White et al 1987). Structure formation within these models proceeds in a bottom up fashion, with small over-densities in the mass distribution collapsing first and subsequently merging hierarchically to form larger objects. Galaxies form from gas that has cooled within the deep potential wells provided by the dark matter halos and the accretion of gas rich sub-clumps that are part of the merging hierarchy (White & Rees 1978). The dominant force that drives structure formation is gravity, hence numerical simulations have proved a useful tool to study the non-linear growth of structure formation over a wide range of dynamical scales.

Several authors have attempted to study the formation of single rich galaxy clusters within a cosmological context in order to obtain very high spatial and mass resolution (Carlberg & Dubinski 1991, Summers 1993, Carlberg 1994, Evrard et al 1994). Their aim was to resolve individual halos in clusters at the present epoch. However, even with mass and length resolution and kpc respectively, no halos were found to survive within a distance of at least half the virial radius from the cluster center. The so called `over-merging' problem, or the inability to resolve sub-structure within dense environments, was first noted by White et al (1987) and Frenk et al (1988). This problem is most serious when comparing cosmological models with the observations. In order to create artificial galaxy catalogues it is customary to identify galaxies with over-dense peaks in the simulated mass distribution; i.e. peaks with 1-d velocity dispersion km/s could be identified with galaxies (the Milky Way is close to an spiral galaxy). However, the cluster sized objects that actually form lack significant substructure or smaller peaks, leading investigators to artificially populate dense regions in the dark matter distribution with `galaxies' (Davis et al 1985, Gelb & Bertschinger 1994).

Clusters in the real universe contain many galaxies that have retained their dark halos (at least within their optical radii). When the over-merging problem was first noticed, it was suggested that cooling gas would increase the central density of the halos, helping them stay intact within the cluster. Preliminary results from smoothed particle hydrodynamic simulations show that galaxy like clumps of gas and dark matter do survive to the present epoch (Katz et al 1992, Katz & White 1993, Summers 1993, Evrard et al 1994, Navarro et al 1995). These techniques are potentially very valuable, but are presently limited by their low resolution. Summers, Davis, and Evrard (1995) discuss different methods for identifying galaxy tracers within dissipationless simulations. They argue that the dynamics of the galaxy tracers identified as Carlberg (1994), is different from the dynamics of the clumps that form after a gaseous component is included.

In this letter we discuss the mechanisms that cause dark halos to dissolve in dense environments within a Hubble time. Artificial numerical effects include two-body relaxation of particles within the softened dark halos and the heating of halos by artificially massive N-body particles. Two physical heating effects will always be present, even in the limit of infinite numerical resolution: tidal heating of halos on eccentric orbits and impulsive heating from rapid encounters between halos. We shall use numerical simulations to follow the evolution of halos within a cluster environment in order to isolate these effects and test analytic estimates of the dissolution time-scales. Our aim is to determine what erases substructure in present simulations and to determine whether or not future simulations at higher resolution will be able to resolve halos within dense environments.

2 DISRUPTION MECHANISMS

2.1 NUMERICAL EFFECTS

2.1a Relaxation

A galaxy halo with internal velocity dispersion at a distance from a cluster with velocity dispersion , will have a tidally limited mass

Here we have assumed isothermal cluster and halo potentials so that the tidal radius of a halo is simply . Typical particle masses in cosmological simulations are between and , therefore the number of particles, , within an halo at kpc is in the range 50-500. The evaporation time-scale from two-body encounters is of order , where the median relaxation time-scale for a halo with half mass radius can be written (Spitzer & Hart 1971), i.e. for isothermal potentials

Clearly, is larger than a Hubble time for moderately large .

Evaporation is accelerated in the presence of a tidal field (c.f. Chernoff & Weinberg 1990). To test the evaporation time-scale of a halo orbiting within a cluster, we constructed equilibrium halos with truncated isothermal density profiles. Our aim is simulate both the evaporation rate from halos that form in current simulations and those which might form given much better length resolution. Current cosmological simulations have force resolutions limited by gravitational softening kpc, and produce halos with correspondingly soft, low density cores with a size comparable to the force resolution. Within the next few years a resolution of kpc will be obtained and the structure of individual halos should be better resolved.

To test the ideas presented above, we evolved a series of model halos constructed with a total mass of within a outer radius of 30 kpc. Such a halo might surround an galaxy 300 kpc from the center of a rich cluster. The halos were placed on circular orbits at 600 kpc in an isothermal cluster potential with a 1-d velocity dispersion of 1000 km/s. This distance was chosen so that the halos's outer radius was well inside the tidal radius imposed by the cluster, in order to minimize the effects of tidal heating. Figure 1 shows the fraction of mass remaining within a fixed distance of 30 kpc from the center of each halo over a period of 100 Gyrs. We evolved the particle distributions using the TREECODE (Barnes & Hut 1986, Hernquist 1987) with 20,000 time-steps. The number of particles, , spline force softening, (Hernquist & Katz 1989), and core radius, , of each model are also indicated in Figure 1.

  
Figure 1: Mass loss rates from dark halos owing to evaporation. The halos are constructed using the indicated model parameters for the core radius, , and the softening length , and particle number . Each halo was placed on a circular orbit at 600 kpc from the center of an isothermal potential for 100 Gyrs.

We find that the initial rate of mass loss is consistent with the standard evaporation formula above. However, the rate of mass loss is not linear with time, but rapidly slows down as the physical size of the halo approaches the softening or core size. Hence, halos with large softening and as few as 20 particles can survive many Hubble times. The rate of evaporation increases as the softening, hence core size, is reduced, i.e. individual encounters between particles can transfer more energy at the same impact parameter. Halos with the same softening will evaporate faster if the physical core is large since the binding energy is correspondingly lower. From these tests we conclude that relaxation effects are not important at driving mass loss from halos within current simulations.

2.1b Particle - halo heating

Carlberg (1994) explained the halo disruption in his particle cluster simulation as a consequence of heating by massive N-body particles, particles much more massive than any viable dark matter candidate. We can make an analytic estimate for using the impulsive approximation, since the typical encounter time-scale is smaller than the galaxy's internal time-scale, i.e. , where is the impact parameter. This calculation was originally performed for the disruption of open clusters by giant molecular clouds (Spitzer 1958). As a cluster particle passes by or through a halo it tidally distorts the system and increases its internal kinetic energy. Following Binney & Tremaine (1987), we equate the disruption time-scale as the time for the halos binding energy, , to change by order of itself owing to impulsive energy input:

where is the perturbers mass and is the perturbers half mass radius, which we equate to the gravitational softening length. Here, we have assumed that the galaxy's mean square radius is similar to the half mass radius . Therefore, we derive

by substituting the tidal radius of the halo for the half mass radius and noting that the number density of perturbers within an isothermal potential, .

The analogous analytic estimate for globular cluster dissolution by black holes in the halo of the Milky Way was tested using N-body simulations by Moore (1993) and found to be accurate within a factor of about two. We therefore conclude that present dissipationless cluster simulations have obtained sufficient resolution to avoid this problem. For example, Carlberg's (1994) cluster simulations have and kpc, and Summers et al (1995) have and kpc, yielding 's of about 65 and 23 Gyrs, respectively.

2.2 PHYSICAL EFFECTS

2.2a Tidal heating

At a certain distance from the halo, particles will escape and become bound to the cluster potential. Our truncated softened isothermal halos will have somewhat smaller tidal radii than those given by the above formula for true isothermal potentials but are more realistic halo potentials. We test the analytic formula by placing model halos on circular orbits at different radii within a cluster and following their evolution with the TREECODE as in Section 2.1a. We define the tidal destruction radius to be the distance from the cluster center at which the halos lose 50% of their mass over 5 Gyrs. (Typical halos completely disrupt at this cluster-centric distance over a Hubble time.)

We evolve several model halos for 10 Gyrs with 2000 time-steps, varying the force softening and core radius. The initial number of particles is 10,000 so as to minimize relaxation effects. Our model halo with a core radius kpc survives to kpc. Halos with kpc and kpc will survive to within about kpc. Halos with kpc can survive at kpc.

Halos that survive several orbits have tidal radii close to the value given by equation (1) using set equal to their pericentric distance. Hence, in the limit of infinite resolution, singular isothermal halos will be tidal stripped to a radius where their outer density is approximately equal to the density of the cluster at pericenter.

Halos become very unstable when the tidally imposed radius approaches a value smaller than two or three times the core radius. Furthermore, halos with a physical core will be more unstable with a larger softening. Hence, in the absence of all other disruption mechanisms, tidal disruption from the cluster potential is sufficient to erase all halos within 500 kpc for a Plummer force softening of 10 kpc. Note that a Plummer softening equivalent to our spline softening would be 30% smaller. Hence, simulations with 20 kpc Plummer softening will dissolve halos at about 1 Mpc within a Coma sized cluster! Any additional source of heating will drive this survival radius even further outwards.

Weinberg (1994) has stressed the importance of resonant orbit coupling between stellar dynamical systems orbiting within a tidal field. The orbital time of a halo at half the clusters virial radius is of order years, a time-scale similar to the internal orbital period within an halo at its half mass radius. Halos on eccentric orbits within a cluster will lose mass as a consequence of heating by the tidal field, even though the `tidal shock' varies quite slowly with time. We test the importance of this effect numerically by placing the above model halos in eccentric orbits within the same cluster potential. Typically over a Hubble time, between 10 and 20 percent of the mass is lost owing to the slow tidal shocks over a period of about 10 Gyrs. Halos with larger cores lose more mass than more tightly bound systems, but tidal shocking plays a relatively minor role in erasing halos.

2.2b Halo-halo heating

An additional heating term arises from impulsive encounters between halos within the cluster. The rate of energy input via impulsive heating is , where is the fraction of the cluster mass contained within the perturbers. Therefore, the relative importance of heating from cluster particles versus individual halos is . Halos are typically between 10 and 100 times as massive as the particles in current simulations so we expect that the heating rate from halo-halo encounters will be about an order of magnitude more important than that from particle-halo encounters.

Performing the same calculation for estimating , but replacing the softened particle perturbers with halos we find

Assuming tidally truncated isothermal halos (c.f. Section 2.1a), varies linearly with (at the cluster virial radius ) and we can write equation (5) as

The scaling of this formula shows that is independent of position within the cluster. However, note that this derivation breaks down when the perturber mass is much larger than the halo mass and when the impact parameter is so large that use of the impulse approximation breaks down, i.e. when the encounter occurs beyond about 75 kpc from the halo. Hence, at large distances from the cluster center, very few encounters occur within this impact parameter and the heating rate falls to zero. Furthermore, in deriving this time-scale we have ignored the effect of the tidal field in estimating the binding energy of the halo. As we demonstrate below, halos closer to the center have lower binding energies and hence will be easier to disrupt. Thus, we would expect halos to survive at the edge of clusters and to be disrupted with increasing ease towards the cluster center.

We test this analytic estimate by constructing a cluster of galaxies similar to the Coma cluster. Within the virial radius (), galaxies are drawn from a Schechter luminosity function with and ( and ), including all galaxies brighter than . Each galaxy is an isothermal potential that is tidally limited using equation (1) with equal to its pericentric distance. Velocity dispersions are assigned using the Faber-Jackson relation (1976); our minimum luminosity corresponds to a halo with 1-d velocity dispersion km/s. For km/s the cluster mass within the virial radius of 1.5 Mpc is , hence for a mass to light ratio of 250, the total cluster luminosity within this radius is . This normalization yields approximately 950 galaxies brighter than our minimum luminosity and 30 brighter than . At 500 kpc from the cluster center, about 25% of the total cluster mass remains associated with the cluster galaxies.

Our model halos are placed on circular orbits at 300 kpc and 600 kpc and feel the potential field of all the cluster members and the analytic potential of the remaining cluster background. The 950 perturbing cluster galaxies orbit on eccentric orbits within the cluster. Since halos lose up to half of their mass over a Hubble time, we multiply the initial perturbers mass by 0.75. Since the impulsive energy input is proportional to the square of the perturber mass this is a conservative means of accounting for the gradual mass loss from the other cluster halos. The typical encounter velocity is equal to km/s requiring 10,000 time-steps over a Hubble time to model the impulsive encounters between galaxies.

Our treatment of the cluster environment contains all of the important dynamical effects. Its major shortcoming is that we assume that the cluster is in place at the beginning of the simulation. Figure 2 shows the results from several numerical calculations. We find that the survival of a given halo depends sensitively on the ratio of the core radius to the tidal radius. As the core radius approaches the tidal radius the halos disrupt completely in a relatively short time-scale. Our estimate of the halo-halo disruption time-scale from equation (5) is in relatively good agreement with the numerical results for halos with small core radii.

  
Figure 2: Mass loss rates from dark halos owing to halo-halo heating. For each model, circular orbits in different directions were followed to demonstrate the stochastic nature of halo-halo heating. The dotted curves show models with kpc in orbit at 300 kpc. The solid curves show models with kpc at 300 kpc. The short dashed curves show models with kpc at 600 kpc. The long dashed curves show models with kpc at 600 kpc.

The curves show a stochastic evolution that leads to a dissolution time-scale of halos with fixed physical properties that varies by a factor . The heating from a single encounter depends upon the square of the perturber mass, hence given a Schechter luminosity function, most of the heating is due to galaxies with luminosities . This explains why the evolution is fairly chaotic since the number of massive perturbers is relatively small. Typically, the evolution of a halo is driven by between 5 and 10 large impulsive encounters.

3 DISCUSSION AND CONCLUSIONS

Current dissipationless N-body simulations have achieved a resolution sufficient to resolve large halos within cluster environments. However, little substructure and no galaxy halos were found in these cluster simulations. We find that two body evaporation is unimportant for halos with more than 30 particles. A second artificial numerical effect, particle-halo heating, does not pose a problem within numerical simulations as long as the particle mass is kept below . However, present cosmological simulations have softening lengths of order 5 - 20 kpc leading to halos with large low-density cores. As the limiting tidal radius approaches the halo core radius the dissolution of the halo occurs very rapidly. Hence, the over-merging problem within large over-dense regions is due to the large force softening combined with tidal heating from the mean field of the cluster and encounters with other dissolving halos.

In the limit of infinite numerical resolution, the survival of individual halos within dense environments depends critically on their inner structure. If halos form with singular isothermal density profiles, then they can always be resolved at some level. Alternatively, including a gaseous component that dissipates energy can increase the central density of halos, effectively decreasing their disruption time-scale under halo-halo collisions. Recent results on the structure of dark matter halos give conflicting results (Moore 1995). The highest resolution simulations of Carlberg (1994), Warren et al (1992) and Crone et al (1994) give very steep inner density profiles. Furthermore, after correcting for force softening, Crone et al show that profiles are singular and fall steeper than on all scales. These results disagree with those of Katz & White (1993) and Navarro et al (1994), who find that the density profiles of halos fall with radius as over a large central region. If the former is true then it will be possible to identify halos within dissipationless simulations, although the halos themselves would be inconsistent with some observations of galaxy rotation curves (Moore 1994). If the simulated halos have shallow inner density profiles as some observations indicate, or as the latter authors find, then the over-merging problem can only be resolved by including a gaseous component.

If dissipationless dark matter forms singular isothermal halos, then in order to resolve galaxies in a cluster environment, force softenings less than a kpc must be adopted. In this case the computational cost becomes enormous if sufficient numbers of particles are used to resolve high density regions and to minimize relaxation effects. Including a dissipational component effectively increases the density of galaxy sized halos with far fewer particles and with only a relatively small increase in the computations run time. Therefore, apart from a consistency check, it does not make sense to continue increasing the resolution of dark matter only calculations which are designed for the sole purpose of identifying dark halos within dense environments.

Acknowledgments

We thank Ray Carlberg for interesting discussions stimulated by his million particle cluster simulation and the referee, Frank Summers, for a thorough reading and constructive comments. This work was funded by NASA through the LTSA and HPCC/ESS programs.

References

Barnes J. & Hut P. 1986, Nature, 324, 446.

Binney, J. and Tremaine, S. 1987, Galactic Dynamics, Princeton University Press, page 442, ed. J. Ostriker.

Carlberg R. & Dubinski J. 1991, Ap.J., 369, 13.

Carlberg R. 1994, Ap.J., 433, 468.

Chernoff D.F. & Weinberg M.D. 1990, Ap.J., 351, 121.

Crone M.M., Evrard A. & Richstone D. 1994, Ap.J., 434, 402.

Davis M., Efstathiou G., Frenk C.S. & White S.D.M. 1985, Ap.J., 292, 371.

Dubinski J. & Carlberg R. 1991, Ap.J., 378, 496.

Evrard A.E., Summers FJ & Davis M. 1994, Ap.J., 422, 11.

Faber S.M. & Jackson R.E. 1976, Ap.J., 204, 668.

Frenk C.S., White S.D.M., Davis M. & Efstathiou G. 1988, Ap.J., 327, 507.

Gelb J.M. & Bertschinger E. 1994, Ap.J., 436, 467.

Hernquist L. 1987, Ap.J.Supp., 64, 715.

Hernquist L. & Katz N. 1989, Ap.J.Supp., 70, 419.

Katz N., Hernquist L. & Weinberg D. 1992, 399, L109.

Katz N. & White S.D.M. 1993, Ap.J., 412, 455.

Moore B. 1993, Ap.J.Lett., 413, L96.

Moore B. 1994, Nature, 370, 629.

Moore B. 1995, Conf. Proc. 5th Annual Maryland Astrophysics conference, `Dark Matter'.

Navarro J., Frenk, C.S. & White S.D.M. 1995, M.N.R.A.S. , in press.

Spitzer L. 1958, Ap.J., 127, 17.

Spitzer L. & Hart M.H. 1971, Ap.J., 164, 399.

Summers F.J. 1993, Ph.D. thesis, Berkeley.

Summers F.J., Davis M., & Evrard A.E. 1995, Ap.J., submitted.

Warren S.W., Quinn P.J., Salmon J.K. & Zurek H.W. 1992, Ap.J., 399, 405.

Weinberg M.D. 1994, A.J., in press.

Warren S.W., Quinn P.J., Salmon J.K. & Zurek H.W. 1992, Ap.J., 399, 405.

White S.D.M., Davis M., Efstathiou G. & Frenk C.S. 1987, Nature, 330, 451.

White S.D.M. & Rees M.J. 1978, M.N.R.A.S., 183 341.