Next: 6.4 Evaluation of Simulation
Up: 6.3 Molecular Dynamics Simulation
Previous: 6.3.2 Continuous Potentials
6.3.3 Beyond Basic Molecular Dynamics
- orientation dependent potentials
- ionic or multipolar interactions
- polymer chains (geometrically constrained molecules)
- nonequilibrium dynamics (e.g. laminary flow); see
[EVANS 86], [HOOVER 99]
- thermostatted dynamics (Nosé-Hoover)
Geometrical constraints / SHAKE method:
Internal geometrical constraints, e.g. rigid bond lengths
(see table 6.2):
- Stiff harmonic bonds - uneconomical
- SHAKE method by Ryckaert et al. [RYCKAERT 77] - better
SHAKE: Consider the smallest non-trivial Kramers chain
consisting of three atoms connected by massless rigid bonds of lengths
. Constraint equations:
with
etc.
Lagrange constraint forces: keep the bond lengths constant.
The constraint force on atom is parallel to
.
Atom is subject to two constraint forces along
and
. Atom is kept in line by
a force along
:
where
are the physical accelerations due to Lennard-Jones
or other pair potentials.
Procedure:
- Let the positions
be given at time .
Integrate the equations of motion for one time step with all
set to zero, i. e. without considering the constraint forces;
denote the resulting preliminary positions (at time ) as
. These will not yet
fulfill the constraint equations; instead, the values of
and
will have some nonzero values
, .
- Now we make the correction ansatz
with undetermined , requiring that the corrected positions
fulfill the constraint equations:
with
; the terms
are quadratic in .
Instead of solving these two quadratic equations for the unknowns
we ignore, for the time being, the small quadratic terms.
The remaining linear equations are solved iteratively,
meaning that this system
of linear equation is solved to arrive at an improved estimate
for which is again inserted in 6.4-6.6
leading to a new set of linearized equations etc., until the absolute
values of are negligible; generally, this will occur after a very
few iterations.
Solving the linearized equations involves a matrix inversion.
To avoid this we introduce one more simplification:
In passing through the chain from one end to the other
we consider only one constraint per atom:
- First the bond
is repaired by displacing and .
- By repairing the next bond
we disrupt the first bond
again; this is accepted in view of further iterations.
- By going through the chain several times we reduce both the error
introduced by neglecting the quadratic terms and
the error due to considering only one constraint at a time.
In our case the procedure is:
insert this in 6.4-6.6 and iterate until
are negligible.
The technique is easily generalized to longer chains.
Molecules and robots:
Robot arms made up of several successive links and joints bear
some resemblance to chain molecules.
A standard problem of robotics, the inverse kinematic problem,
may therefore be solved a la simulation.[KASTENMEIER 86]
The principle of the so-called constrained dynamics/stochastic
optimization technique is as follows:
- Define a required ``world trajectory''for the
final element in the chain.
- If the robot is redundant, meaning
that it has more degrees of freedom (links and joints) than necessary,
the problem is underdetermined: different combinations of
movements by the individual joints produce identical paths of the
final member.
- This redundancy may be exploited to fulfill additional requirements,
such as an overall minimum of angular accelerations (i. e. mechanical wear)
in the joints, avoidance of obstacles, etc.
- The last element is now moved, one time step at a time, along its
requested trajectory, and SHAKE is invoked to have the preceding joints
follow.
- The additional requirements are combined into a cost function
which is minimized, at each time step, using a stochastic search
(simulated annealing or simple random search).
Thermostats / Nosé-Hoover method:
In a nonequilibrium process work must be invested which systematically
heats up the system. In lab experiments a thermostat takes care of this;
in simulations the problem of thermostating is non-trivial.
The temperature of a molecular dynamics sample is not an input parameter
to be manipulated at will; rather, it is a quantity to be ``measured''
in terms of an average of the kinetic energy of the particles,
( dimension).
Suggestions as to how one can maintain a desired
temperature in a dynamical simulation range from repeatedly rescaling
all velocities (``brute force thermostat'') to adding a suitable
stochastic force acting on the molecules.
Such recipes are unphysical and introduce an artificial trait of
irreversibility and/or indeterminacy into the microscopic dynamics.
A interesting deterministic method of thermostating a simulation
sample was introduced by Shuichi Nosé and William Hoover.
Under very mild conditions the following augmented equations of motion
will lead to a correct sampling of the canonical phase space at a given
temperature :
The coupling parameter describes the inertia of the thermostat.
The generalized viscosity is temporally varying and may assume
negative values as well.
For systems of many degrees of freedom (particles) this is sufficient
to produce a pseudo-canonical sequence of states.
For low-dimensional systems it may be necessary to use two NH thermostats
in tandem [MARTYNA 92].
Next: 6.4 Evaluation of Simulation
Up: 6.3 Molecular Dynamics Simulation
Previous: 6.3.2 Continuous Potentials
Franz J. Vesely Oct 2005
See also: "Computational Physics - An Introduction," Kluwer-Plenum 2001