The galaxy models used in these experiments have a central bulge, an exponential disk, and a spherical dark halo. Disks are actually represented in two different ways; as flat, rotating structures corresponding to the usual notion of a disk, or as spherical, isotropic distributions with the same cumulative mass profiles as an exponential disk. When the latter disk representation is used, distribution functions for all components may be calculated exactly using Abel integrals.
For ease of comparison with earlier work, we initially adopted galaxy models similar to those used by Velazquez & White (1996), with spherical bulges following Hernquist (1990) profiles, exponential & isothermal disks, and non-singular isothermal halos (Hernquist 1993). But a problem came to light in trying to realize the halo using an isotropic distribution function; the r−1 profile of the inner bulge creates a potential with a nonzero gradient as r → 0, and in such a potential well a halo with a constant-density core can only be realized with an anisotropic distribution function. We therefore replaced the isothermal halo profile with an NFW profile; the latter has the same profile as the bulge as r → 0, so it can be realized with an isotropic distribution function.
For numerical convenience, the bulge profile was tapered off smoothly
at large radii:
The NFW halo has a logarithmically diverging mass, so it must
be tapered at large radii; we used the form adopted by Springel &
The disk's density profile depends on cylindrical radius R and
Parameters for the galaxy model are listed below. For each component, M is the mass, a is the scale radius, and b is the taper radius.
|Component||Integrated Mass||Scale Radius||Taper Radius|
|bulge||Mb = 0.3125||ab = 0.15||bb = 8|
|disk||Md = 1||Rd = 1|
|halo||Mh = 14||ah = 2.5||bh = 24.9485|
These quantities are reported in arbitrary units with G = 1. Any such system of units is internally self-consistent, and the results may be rescaled as needed. To scale this model to the Milky Way, for example, we may equate the disk's scale length Rd and mass Md to the scale length and mass of the MW's disk, Rd = 3.5 kpc and Md = 6 ×1010 M\odot.
Fig. 1 shows density and mass profiles for each component of the galaxy model. Here the disk is represented by a spherical configuration which has the same mass profile M(r). At all radii the spherically-averaged disk density and mass profiles fall below the corresponding halo profiles. This is not the case in the V&W model - there the constant density of the halo core allows the disk to dominate (slightly) at intermediate radii.
Initial data are generated from isotropic distribution functions fb(E), fs(E), and fh(E), which represent the bulge, `spherical disk', and halo, respectively. We calculated these distribution functions using an Abel equation which relates component c's density profile ρc(r) and distribution function fc(E) to the potential generated by all components, Φ(r). To allow for the finite resolution of the N-body force calculation, we smoothed the total density profile ρ(r) with a Plummer kernel before evaluating the potential function Φ(r). This smoothing procedure employs a semi-empirical fitting function with a free parameter κ; the tests below were run to establish the optimal value of this parameter.
The N-body realizations used in these tests had Nb = 20480 bulge particles, Ns = 65536 `spherical disk' particles, and Nh = 229376 halo particles. Plummer smoothing with a length scale of ε = 0.05 length units was adopted. Forces were calculated with a modified tree-code. Bodies were advanced with a time-centered leap-frog; for these initial tests a time-step ∆t = 1/16 was used.
|κ = 1.15||κ = 1.25||κ = 1.30||κ = 1.35|
Fig. 2 and the associated animations shows how bodies evolve in binding energy Ei during a fairly short run of 8 time units. As the animations make clear, bodies are scattered in binding energy by stochastic events and by large-scale potential fluctuations; the former give rise to a diffusive process, while the latter create coherent patterns of motion. These coherent motions are very nearly absent in the κ = 1.35 model, indicating that this model is closest to equilibrium.
Initial data for the disk is generated using a routine which sets up a self-gravitating disk in approximate equilibrium with the external potential due to the bulge and halo. Initial data for the bulge and halo are generated using the same Abel integral technique described above; to calculate the required distribution functions, we assume that the total potential is spherical. Of course, the actual potential is aspherical because of the disk, but this does not seem a significant source of error in the present experiments.
Last modified: September 24, 2002