Size of the system limited by storage on computer and speed of programs.
A. Nonbond cutoffs
-The force/energy loop is almost always the rate-limiting step O(
).
-The nonbonding terms are the most time-consuming.
To save computational time, apply a spherical cutoff for nonbonding terms:
V(r)=0 for
Advantages
-a cutoff can reduce the number of neighbors explicitly considered by a factor of approximately
(where R is radius of entire system)
-
should be large enough so that the spherical cutoff should be a small perturbation
-for LJ atoms, a typical cutoff radius is
,
which is just 1.6% of the well depth
-both the van der Waals and electrostatic interactions are significant up to 15
or more
Disadvantages
-calculated properties are for truncated potential and may be different from properties of non-truncated potential
must test for each new system
-sudden cutoff of nonbonded interactions at
discontinuities in energy and forces
1. Shifted and shifted-force potentials
Energy is discontinuous at
(i.e. when a pair of molecules crosses
this boundary the total energy will not be conserved)
Avoid this by shifting the potential function V(r) by an amount
- does not affect the forces so does not affect equations of motion
-does affect the total energy, and its contribution varies depending
on the particular configuration
Force is discontinuous at
(i.e. the force is zero for
), which can cause numerical instability
Avoid this by using shifted-force potential to make derivative zero at cutoff distance
2. Switching functions
Multiply V(r) by a switching function
-smoothly turns off interactions over a range of distances
to
-
decreases from 1 to 0 between
and
Example:
Brooks et al., J. Comp Chem. 4, 187-217 (1983)
This brings the potential V and the force
down to zero in a smooth way between
and
3. Neighbor lists
-Buffer region is created between
and
-A neighbor list is kept for all pairs of atoms closer than
-Neighbor list is not updated every step
-Only those pairs of atoms in neighbor list are considered during calculation
of nonbonded interactions (saves time in calculating distances
between pairs of atoms to test if closer than
)
-Neighbor list is updated often enough to ensure that no atom outside the
buffer
region can become closer than
before the neighbor list is updated
Example: update the neighbor list whenever any atom moves more than 1/2 the buffer width
-Buffer width and velocities determine time between updates
4. Charge groups and switching atoms
London-van der Waals: often reasonable to truncate at 8-10
Electrostatic: charge-charge interactions go as 1/r
very long-ranged
-Since most molecules are composed of neutral fragments with
dipoles and quadrupoles, the leading term in the electrostatic
interaction between molecules or parts of molecules is dipole-dipole,
which goes as
(still long-ranged but falls off much faster)
-If cutoffs are applied on atom-by-atom basis, they could split dipoles,
which would artificially introduce a large charge-charge interaction;
Avoid splitting dipoles by applying cutoffs over charge groups
-charge group: small group of atoms near each other that have a net charge of zero
-one atom in each group is designated the switching atom, and cutoffs
are applied for distances between switching atoms in different groups
-size of charge group must be smaller than cutoff distance;
Example:
For a typical charge group of size 1-3
,
cutoffs larger than 7 or 8
are reasonable
B. Boundaries
Problem: want to simulate bulk liquid but can only include a relatively small number of atoms
large fraction of molecules would lie on the surface and thus experience
different forces than those in bulk liquid
1. Periodic boundary conditions
-Replicate cubic box to form an infinite 3-dimensional lattice
-As a molecule moves in the original box, all its periodic images move in exactly the same way
-If a molecule leaves the central box, an image enters through the opposite face
-Only need to store coordinates of central box
-Assumes that the properties of a small, infinitely periodic system are the same as those of a macroscopic system
must check by performing simulation with different box lengths L
Minimum image model:
-each molecule can interact only with the molecule or molecular image closest to it
-must truncate potentials so that cutoff distances are less than L/2
(i.e.
)
-maximum of N(N-1)/2 pairwise-additive interaction terms
Explicit image model:
-if
must generate periodic images of molecules in central box (ghost moleculeses) to as great a distance as necessary for cutoff (i.e. include all ghost molecules that are within
of a real particle)
-ghost molecules shadow their respective particles in central box
-the molecules that are ghost molecules can change during the simulation
2. Stochastic forces at spherical boundary
-System is represented by a sphere, and the boundary potential mimics
the potential that a solvent or lattice would produce
-Outer shell (
) of the sphere is designated as a buffer region, and molecules in this region obey a stochastic (Langevin) equation of motion rather than the standard Newtonian equation of motion
-To prevent evaporation either
1. include an additional shell (
) of frozen molecules
2. include an effective boundary potential
C. Long-range forces
Long-range forces: potentials that fall off no faster than
(i.e. charge-charge (
); charge-dipole (
);
and dipole-dipole (
))
Problem:
conditionally convergent sums
different answer depending on method
of grouping terms
Disadvantage of minimum image method
-Cube will be electrically neutral, but periodic structure will be imposed on the liquid
(i.e. similary charged ions will occupy positions in opposite corners of the cube)
1. Ewald sums
-Technique for summing interaction between an ion and all its periodic images
-Convert one conditionally convergent sum into two absolutely convergent sums
-Surround the sphere with a conducting material in order to cancel net dipole
Physical description:
-Surround each ion with a Gaussian charge distribution of opposite sign and
equal magnitude to screen the interactions so that they are now
short-ranged (and the sum of interactions is absolutely convergent)
-Add a cancelling charge distribution so that the overall potential is identical to the original one
-Sum this cancelling distribution in reciprocal space so that it is absolutely convergent
Mathematical description:
-Multiply lattice sum by a convergence function
to make the sum absolutely
convergent
-Add a term equal to the product of the lattice sum and
to preserve equality
-Fourier transform this second term so that it is also absolutely convergent
General lattice sum:
where the first summation is over all simple cubic lattice points
and the second two sums are over the N molecules in the central box. The prime on the first sum indicates that i=j
is omitted for
.
Ewald sum:
is the complementary error function:
first term: real space sum
second term: reciprocal space sum
last term: cancels the self term that is included in the reciprocal space sum
-Sums in first and second terms are absolutely convergent
-
is chosen so that both sums converge appropriately
2. Reaction field method
-The particles within a cutoff sphere are treated directly, and the charges
outside the cavity are treated as if they form a dielectric continuum
that produces a reaction field within the cavity
-Assumes that the interaction from molecules beyond a cutoff distance can be handled in an average way, using macroscopic electrostatics
|
|
|
|
|
|
|
|
|
|