Biomedical Engineering Reference
In-Depth Information
Subtracting and rearranging we get
v
A
t
v
A
t
t
2
t
2
+
=
−
+
a
A
(
t
)
t
+···
(10.8)
I will switch back and forth between (for example)
v
and d
r
/d
t
in order to try and improve
the readability of the equations. Also, I could have written
v
A
(
t
) or (
v
A
)
t
to mean the
instantaneous velocity of particleAat time
t
. The acceleration
a
is calculated from the force.
Using the same procedure for the Taylor expansion of
r
A
at the time point
t
+
t
/2 we get
v
A
t
t
t
2
+
=
+
+
+···
r
A
(
t
t
)
r
A
(
t
)
(10.9)
Equations (10.8) and (10.9) form the so-called
leapfrog algorithm
, which is reputed to be
one of the most accurate and stable techniques for use in molecular dynamics.
A suitable time increment
t
for molecular dynamics is a femtosecond (10
−
15
s). In the
leapfrog scheme, the velocities are first calculated at time
t
+
t
/2. These are used to
+
calculate the positions of the particles at time
t
t
and so on. In this way the velocities
leap over the positions and then the positions leap over the velocities.
10.4.2 Verlet Algorithm
If instead we start from the Taylor expansion of
r
A
(
t
) we have
d
r
A
d
t
d
2
r
A
d
t
2
1
2
(
t
)
2
r
A
(
t
+
t
)
=
r
A
(
t
)
+
t
+
+···
t
t
(10.10)
d
r
A
d
t
d
2
r
A
d
t
2
1
2
(
t
)
2
r
A
(
t
−
t
)
=
r
A
(
t
)
−
t
+
+···
t
t
which gives (assuming that third-order and higher terms are negligible)
d
2
r
A
d
t
2
(
t
)
2
r
A
(
t
+
t
)
=
2
r
A
(
t
)
−
r
A
(
t
−
t
)
+
(10.11)
t
This is known as the
Verlet algorithm
. The acceleration is obtained from the force exper-
ienced by atom A at time
t
. The velocity does not appear in the expression, but it can be
obtained from the finite difference formula
r
A
(
t
+
t
)
−
r
A
(
t
−
t
)
v
A
(
t
)
=
(10.12)
2
t
TheVerlet algorithm uses positions and accelerations at time
t
and the position at time
t
−
t
to calculate a new position at time
t
t
. All these have to be stored at every iteration.
A variant is the
velocity Verlet
algorithm, which requires only the storage of positions,
velocities and accelerations that all correspond to the same time step. It takes the form
+
d
r
A
d
t
d
2
r
A
d
t
2
1
2
(
t
)
2
r
A
(
t
+
t
)
=
r
A
(
t
)
+
t
+
t
t
(10.13)
d
r
A
d
t
t
+
d
2
r
A
d
t
2
d
2
r
A
d
t
2
t
+
t
t
1
2
v
A
(
t
+
t
)
=
t
+
There are many other algorithms in the literature, each with their own strengths and
weaknesses.
