Previously, the effects of the shocks were calculated in a linear regime, assuming that the potential of the cluster does not change during or after the perturbation. We also consider whether an anisotropy of stellar velocities develops as a result of preferential loss of stars on radial orbits. Finally, we investigate how effects of the shocks are diminished, because of conservation of adiabatic invariants of the stars in the center of the cluster. The gravitational potential and density of the system are expanded in a series of basis functions, and some first terms of the expansion are retained. The code is much faster than direct N-body integration and allows simulations with 1,000,000 particles or more.
Computation time scales linearly with the number of particles.
The code is also almost perfectly parallelizable, since particles effectively interact with the mean field and not with each other. We developed an efficient MPI version of the code for the IBM SP, as well as a simpler version on the shared memory architecture for SGI Origin 2000. We simulate a globular cluster with 1,000,000 stars, randomly drawn to satisfy the initial distribution function. In the initial equilibrium model, the outer radius of the system is seven times the core radius. A subset of 10,000 stars is illustrated in Figure 1.