8.4 Integration methods: symplectic integrators 211
than the typical collision time (time of free flight), you should find
〈x^2 〉∝t^2 , (8.24)
and this crosses over to diffusive behaviour
〈x^2 〉=Dt, (8.25)
withDthe diffusion constant. In the solid phase, the diffusion constant is 0. In the
gas phase, the diffusive behaviour sets in at later times than in the fluid.
If the program works properly, keeping a Verlet neighbour list as discussed
in the previous section can be implemented. Verlet [8] usedrcut-off=2.5σand
rmax =3.3σ. A more detailed analysis of the increase in efficiency for various
values ofrmaxwithrcut-off=2.5σshows thatrmax=3.0σwith the neighbour list
updated once every 25 integration steps is indeed most efficient [ 2 , 15 ].
programming exercise
Implement the neighbourlist in your program and check whether the results
remain essentially the same. Determine the increase in efficiency.
8.4 Integration methods: symplectic integrators
There exist many algorithms for integrating ordinary differential equations, and a
few of these are described inAppendix A. In this section, we consider the particular
case of numerically integrating the equations of motion for a dynamical system
described by a time-independent Hamiltonian, of which the classical many-particle
system at constant energy is an example. Throughout this section we consider the
equation of motion for a single particle in one dimension. The discussion is easily
generalised to more particles in more dimensions.
The Verlet algorithm is the most popular algorithm for molecular dynamics and
we shall consider it in more detail in the next subsection. Before doing so, we
describe a few criteria which were formulated by Berendsen and van Gunsteren
[16]for integration methods for molecular dynamics. First of all,accuracyis an
important criterion: it tells us to which power of the time step the numerical traject-
ory will deviate from the exact one after one integration step (see alsoAppendix A).
Note that the prefactor of this may diverge if the algorithm is unstable (e.g. close to a
singularity of the trajectory). The accuracy is the criterion that is usually considered
in numerical analysis in connection with integration methods.
Two further criteria are related to the behaviour of the energy and other con-
served quantities of a mechanical system which are related to symmetries of the
interactions. Along the exact trajectory, energy is conserved as a result of the time-
translation invariance of the Hamiltonian, but the energy of the numerical trajectory
will deviate from the initial value and this deviation can be characterised by itsdrift,