Equations

Force calculation

The force for each particle is calculated with a direct N-Body scheme, where all forces are calculated accurately from each object to each other object. The equation for the force is (boldface parameters are Cartesian vectors, see Wikipedia for more details.)

The acceleration follows from ${\bf a}_i = {\bf F}_i/m_i$.

Integration

The equations of motion are integrated via the Verlet-leapfrog
predictor corrector integrator, see Wikipedia for more details.

The equation to update the position:

The equation to update the velocity:

Units

The units of the integration are so called Heggie-Units, in which case M=G=1. Here

As a consequence, the total energy of the system is Etot = 0.25. A consise description of these units can be found here.

The Algorithm

The gravitational N-body problem
read t=0 snapshot from input
calculate energy dt= 0.001
initialize tend = 1
WHILE t < tend DO
Calculate ai(ri)
Advance positions ri → ri+1
Calculate ai+1(ri+1)
Advance velocities vi → vi+1
t += dt
nsteps++
IF (nsteps%10 = 0)
calculate energy
print diagnostics (energy, energy error)
ENDIF
END WHILE
dump final snapshot