Monte Carlo simulations of globular cluster evolution. I. Method and test calculations

We present a new parallel supercomputer implementation of the Monte Carlo method for simulating the dynamical evolution of globular star clusters. Our method is based on a modified version of Hénon's Monte Carlo algorithm for solving the Fokker-Planck equation. Our code allows us to follow the evolution of a cluster containing up to 5 × 10⁵ stars to core collapse in ≲40 hours of computing time. In this paper we present the results of test calculations for clusters with equal-mass stars, starting from both Plummer and King model initial conditions. We consider isolated as well as tidally truncated clusters. Our results are compared to those obtained from approximate, self-similar analytic solutions, from direct numerical integrations of the Fokker-Planck equation, and from direct N-body integrations performed on a GRAPE-4 special-purpose computer with N = 16384. In all cases we find excellent agreement with other methods, establishing our new code as a robust tool for the numerical study of globular cluster dynamics using a realistic number of stars.