204 Molecular dynamics simulations
check whether their separation is larger thanrcut-off. In the same paper in which he
introduced the midpoint integration algorithm into MD, Verlet [8] proposed keeping
a list of particle pairs whose separation lies within some maximum distancermaxand
updating this list at intervals of a fixed number of steps – this number lies typically
between 10 and 20. The radiusrmaxis taken larger thanrcut-offand must be chosen
such that between two table updates it is unlikely for a pair not in the list to come
closer thanrcut-off. If both distances are chosen carefully, the accuracy can remain
very high and the increase in efficiency is of the order of a factor of 10 (the typical
relative accuracy in macroscopic quantities in a MD simulation is of order 10−^4 ).
There exists another method for keeping track of which pairs are within a certain
distance of each other: thelinked-cell method. In this method, the system is divided
up into (rectangular) cells. Each cell is characterized by its integer coordinates
IX,IY,IZ in the grid of cells. The cell size is chosen larger than the interaction range
which is about the size ofrmax>rcut-offin the Verlet method. If we wanted a list
of particles for each cell, we could simply restrict the interactions to particle pairs
in the same, or in neighbouring cells. However, as particles will leave and enter the
cells, the bookkeeping of these lists becomes a bit cumbersome. This bookkeeping
can however be done very efficient by using a list of particle indices. The procedure
is reminiscent of the use of pointers in a linked list. We need two ingredients: we
must have a routine which generates a sort of table containing information about
which particle is in what cell, and we need to organise the force calculation such
that it uses this information.
To be specific, let us assume that there areM×M×Mcells. The particles are
numbered 1 throughN, so each particle has a definite index. We use an integer array
called ‘Header’ which is of sizeM×M×M: Header(IX,IY,IZ) tells us thehighest
particle index to be found in cell IX,IY,IZ. We also introduce an integer array ‘Link’
which is of sizeN. The arrays Header and Link are filled in the following code:
dimension header(M,M,M), link(N)
Set Header (IX,IY,IZ) to 0
Set Link(I) to 0
FOR I = 1,N DO
IX = int(M*x(I)/L)+1
IY = int(M*y(I)/L)+1
IZ = int(M*z(I)/L)+1
link(i) = header(IX,IY,IZ)
header(IX,IY,IZ) = I
END FORNow, Header contains the highest index present in all cells. Furthermore, for particle
I, Link(I) is another particlein the same cell. To find all particles in cell IX, IY, IZ,