Ajith Peter Posters

1
A Penalty Function Method for Constrained
Molecular Dynamics Simulation
Ajith Gunaratne and Zhijun Wu Department of
Mathematics, Program on Bioinformatics and
Computational Biology, Iowa State University
Constrained Dynamics MD simulation can be used to
study various dynamic properties of proteins, but
a long sequence of iterations has to be carried
out even for small protein motions due to the
small time step (10-15 sec) required. The bonding
forces are among those causing fast protein
vibrations that require small time steps to
integrate, and they may be replaced by a set of
bond length constraints, to increase the step
size and hence the simulation speed. Lagrange
multiplier methods have been developed for
constrained dynamics simulation. The multipliers
are determined in every step to satisfy the
constraints through the solution of a nonlinear
system of equations.
Penalty Function Method
Results
Principle of Least Action
The penalty function method (Penalty) is
implemented in CHARMM and applied to protein
BPTI. The simulation results are analyzed and
compared with Verlet and Shake.
We consider the least action problem with
constraints for a given physical system as a
constrained optimization problem and define a
quadratic penalty function for the problem. The
equations of motion then become
Principle of Variation
Shake
Cooper et al., U. Southampton, UK
Euler / Lagrange Equations
The simulation was carried out on SGI
workstations. Heating and equilibrium were
conducted to bring Verlet, Shake, and Penalty
simulations to the same point. The following
simulations were recorded for analysis. The
results showed that the penalty function method
had high correlations with Shake and outperformed
Verlet.
Equations of Motion
m
F
Folding of protein chicken villin headpiece
(HP-36) simulated by Duan and Kollman (Science
282, 1998). The simulation ran for 100 days on
256 processor Cray T3E.
Concluding Remarks We have implemented a penalty
function method in CHARMM and tested it on
protein BPTI. The dynamic properties obtained
were highly correlated with Shake and showed
clear advantages over Verlet. The method is
equivalent to an unconstrained method with the
bonding forces gradually strengthened as the
simulation proceeds. It promises faster
convergence to the limit of the system. The
penalty parameter has to be selected
appropriately (not altering the system too much).
The ultimate effect appeals further
investigations with longer simulations. Contact
ajith_at_iastate.edu, zhijun_at_iastate.edu
Theorem Every limit point x of the set of global
minimizers xk of the penalty functions pk
corresponding to the penalty parameters µk, as
µk?8, is a global minimizer of the original
objective function.
Theorem Every limit point x of the set of KKT
points xk of the penalty functions pk, µk?8,
that satisfies the constraint qualification is a
KKT point of the original objective function.
Velocity Auto Correlation (51 CYS C? )
Equations of Motion
The RMS fluctuations of backbone and non-backbone
atoms and the velocity auto correlations of Ca
atoms calculated by the penalty function method
agreed well with those by Shake.
The simulation with the penalty function method
can be done by using a conventional unconstrained
solver such as Verlet, only with the penalty
parameter increased in an appropriate manner as
the simulation proceeds.
Verlet Algorithm
Note that we can scale the constraints with their
force constants Kj when forming the penalty
terms. The resulting force function can then be
viewed as a smooth continuation of the original
force field.
Shake Algorithm
Backbone RMS Fluctuations
Non-backbone RMS Fluctuations
Acknowledgements We would like to thank the
Department of Mathematics, Iowa State University,
for providing research support for Ajith
Gunaratne.
aijk
The penalty function method is easy to implement
and costs less than a Lagrange multiplier method,
which requires the solution of a nonlinear system
of equations in every time step.
xl
ajkl
aijkl
aklm
xm
xk
ajklm