a self-starting scheme and therefore much easier to implement than the multistep ABM scheme of the same order. Its local truncation error has much smaller coefficient, which allows a larger timestep for the same accuracy. On general-purpose computers, this advantage is partly cancelled by the additional CPU time to calculate the first time derivative, but on GRAPE-type machines this is an significant advantage ([@Makino1992], [@Makino1993]).

The output are taken at intervals of a fixed time, which is a certain fraction of the crossing time. We recorded central density, core radius, number of particles in the core, radii of Lagrangian shells, velocity dispersion within Lagrangian shells, binaries and their binding energies. The core parameters are calculated following Casertano and Hut (1985), with modifications described in McMillan (1990).

The accuracy parameter of the time integration is adjusted so that the energy error between two outputs is smaller than a certain prescribed value. The value we used is $1\times 10^{-5}\sim 1\times 10^{-6}$ depending on $N$. We required higher accuracy for larger $N$, since the duration of the simulation is longer. When the energy error is very large, the program automatically reads the output at the previous checkpoint and restarts with a reduced accuracy 