244 Molecular dynamics simulations
In a careful treatment of the Ewald technique, the sum is carried out formally
by taking a large volume of some particular shape (e.g. a sphere) containing the
system replicas and then this shape is increased. The reason for this is that the
sum over the interactions is conditionally convergent, i.e. it depends on the order
in which the various contributions are taken into account. This is explained by the
fact that the system replicas all have a dipole moment and will hence build up a
surface charge at the boundary of the huge volume. The most natural boundary
condition (the one which is arrived at in more pedestrian derivations) is consistent
with the sphere being embedded in a perfectly conducting medium. For the sphere
embedded in a dielectric, a correction must be included [56]. It is important to be
aware of this when calculating (say) dielectric properties of a charged system.
8.7.2 Efficient evaluation of forces and potentialsAs a result of the long range of the forces, all interacting pairs must be taken
into account in the calculation of forces or potentials. The straightforward imple-
mentation, considered in the previous sections of this chapter, also called the
particle–particle method(PP) because all pairs are considered explicitly, would
requireO(N^2 )steps, but it turns out possible to reduce this to a more favourable
scaling. We shall briefly sketch two other methods, and then consider a third one,
thetree methodin greater detail.
In theparticle-mesh(PM) method, a (usually cubic) grid in space is defined. A
mass (or charge) distributionρis then defined by assigning part of each particle’s
mass (or charge) to its four neighbouring grid points according to some suitable
scheme. The potential can then be found by solving Poisson’s equation on the grid
∇D^2 U(r)=− 4 πρ(r) (8.137)where∇D^2 is the finite-difference version of the Laplace operator. Using fast Fourier
transforms (seeAppendix A9), this calculation can be carried out in a number of
steps proportional toMlogMwhereMis the number of grid points. Knowing
the potential, the force at any position can be found by taking the finite difference
gradient of the potential, after a suitable correction for the self-energy resulting
from the inclusion of the interaction of a particle with itself in this procedure.
This method obviously becomes less accurate for pairs of particles with a small
separation, as in that case the Coulomb/gravitation potential is not sufficiently
smooth to be represented accurately on the grid. Therefore it is sensible to treat
particles within some small range (for example a range comparable to the grid
constant) by the PP method. This can be done by splitting the force into a smooth
long range (LR) and a short range (SR) part:
F=FLR+FSR, (8.138)