1 Introduction

Collisionless N-body simulations are heavily used to study the dynamical evolution of galaxies, or systems of galaxies, and have so far been remarkably successful in producing many interesting results. Yet proper care has to be taken to eliminate possible sources of numerical errors which could, if present, lead to erroneous results. One of the possible sources of errors stems from the fact that the number of particles in a simulation is several orders of magnitude less than the number of stars in a typical galaxy, or, in other words, that the graininess in a computer realisation is much higher than that of the galactic system it is meant to represent. This could lead to errors since a particle moving through an N-body representation of a given continuous system representing a galaxy is deflected from the orbit it would have had in the corresponding smooth continuous medium, due to two-body encounters. This effect is known as two-body relaxation and the characteristic time linked to it as two-body relaxation time (hereafter $T_{\rm relax})$.

The relaxation time will obviously increase with the number of particles in the configuration, and tend to infinity as the number of particles tends to infinity, in which limit the evolution of the system will follow the collisionless Boltzmann equation. Thus the relaxation time of simulations will be much shorter than that of galaxies, which is much longer than the age of the universe. Relaxation leads to a loss of memory of the initial conditions and an evolution of the system towards a state of higher entropy. It is thus necessary to have good estimates of the relaxation times of N-body simulations, since we can trust their results only for times considerably shorter than that.

In the early times of N-body simulations, when the number of particles used was of the order of a few hundred, the authors by necessity gave estimates of relaxation times in order to enhance the credibility of their results. Unfortunately in most cases only simple analytical estimates were used and the corresponding relaxation times were found to be comfortably, although perhaps unrealistically, high. As computers became faster, the number of particles used in simulations was increased. Authors using several tens or hundreds of thousands particles deemed it unnecessary to include such simple estimates of the relaxation times, since it was well known that the simple analytical estimates would give reassuringly high relaxation times. Nevertheless, it is not clear whether the simple analytical estimates are in all cases sufficiently near the true values. This could well be doubted since the simple analytical estimates rely on a number of approximations, which are not in all cases valid.

Since different N-body methods may lead to different relaxation rates, it is of interest to discuss relaxation times when introducing a new method. It could thus have been feared that in a tree code (Barnes & Hut 1986) the relaxation time would not be determined by the number of particles, but by the number of nodes, which would then act as "super-particles''. Since the number of nodes is always much smaller than the number of particles, this would entail considerably shorter relaxation times than direct summation with the same number of particles, and thus constitute a major disadvantage of the tree code. This fear was put to rest by Hernquist (1987) who showed that the relaxation for tree code calculations does not differ greatly from that obtained by direct summation provided the tolerance parameter is less than 1.2. Similarly Hernquist & Barnes (1990) compare relaxation rates in direct summation, tree and spherical harmonic N-body codes, while Weinberg (1996) introduces a modification of the orthogonal function potential solver that minimises relaxation.

Several methods have been used to measure two-body relaxation. Standish & Aksnes (1969), Lecar & Cruz-Conzález (1972) and Hernquist (1987) have measured the angular deflection of test particles moving in a configuration of N field particles. Although this method has the disadvantage of not including collective effects, it has the advantage that all the parameters can be changed independently of each other, and that the results are easy to interpret. Theis (1998) performed semi-analytical calculations, assuming a homogeneous medium and also ignoring cumulative effects. The most widely used approach is to monitor the energy conservation of individual particles in systems in which, had it not been for the individual encounters, the individual energies would have been conserved. This method, which includes collective effects, is well suited for testing relaxation rates introduced by different codes, but can only be used with systems in equilibrium. It has been used e.g. by Hernquist & Barnes (1990), Hernquist & Ostriker (1992), Huang et al. (1993) and Weinberg (1996). Theuns (1996) measured the diffusion coefficients as a function of the energy in a direct summation N-body simulation by studying the properties of the random walk in energy space for particles of given energy and found very good agreement with theoretically calculated diffusion coefficients. Finally a number of studies (e.g. Farouki & Salpeter 1982, 1994; Smith 1992) rely on a measurement of the mass segregation, i.e. on the fact that, due to two-body encounters, high mass particles lose energy and spiral towards the center, while light ones gain energy and move to larger radii. Thus the configuration of the high mass particles contracts, while that of the light particles expands, and from the rate at which this happens we can calculate the two-body relaxation time.

In this paper we will calculate the relaxation times in a large number of cases, using the first of the methods mentioned above, i.e. by measuring deflection angles of individual trajectories of test particles in a configuration of rigid field particles. We will cover a much larger part of the parameter space than was done so far, and we will also extend to larger number of particles. All the calculations presented in this paper were made on the Marseille GRAPE-3, GRAPE-4 and GRAPE-5 systems. The Marseille GRAPE-3 systems have been described by Athanassoula et al. (1998), while a general description of the GRAPE-4 systems and their PCI interfaces has been given by Makino et al. (1997) and Kawai et al. (1997) respectively. A description of the GRAPE-5 board can be found in Kawai et al. (2000). Opting for a GRAPE system restricts us to a single type of softening, the standard Plummer softening, but has the big advantage of allowing us to make a very large number of trials, covering well the relatively large parameter space. Theis (1998) compared the relaxation rates obtained with the standard Plummer softening to those given by a spline (Hernquist & Katz 1989) and showed that the differences between the two are only of the order of 20-40%.

This paper is organised as follows: in Sect. 2 we briefly summarise the simple analytical estimates of the relaxation time. In Sect. 3 we describe the numerical methods used in this paper and discuss the validity of their approximations. Here we also introduce the mass models which will be used throughout this paper. The values of the parameters to be used, and in particular the values of the softening, are derived and discussed in Sect. 4. In Sect. 5 we give results for the relaxation time. We specifically discuss the effect of number of particles, of the velocity and of the softening, and compare results obtained with GRAPE-3 and GRAPE-4. We also give a prediction for the relaxation time in an N-body simulation. We summarise in Sect. 6.