Molecular Dynamics Group, University of Queensland |

Molecular Dynamics Group, University of Groningen |

MD

Tutorial

Although normally represented as static structures, molecules such as lysozyme are in fact dynamic. Most experimental properties, for example, measure a time average or an ensemble average over the range of possible configurations the molecule can adopt. One way to investigate the range of accessible configurations is to simulate the motions or dynamics of a molecule numerically. This can be done by computing a trajectory, a series of molecular configurations as a function of time, by the simultaneous integration of Newton's equations of motion

(eqn. 1)

and

(eqn. 2)

for all atoms (i = 1, 2 ,...,N) of the molecular system. The atomic
coordinates, r, and the velocity, v, of atom, i, with mass,
m_{i}, thus become functions of time. The force F_{i}
exerted on atom i by the other atoms in the system is given by the
negative gradient of the potential energy function V which in turn
depends on the coordinates of all N atoms in the system:

(eqn. 3)

For small time steps δt, eqn. (2) can be approximated by

(eqn. 4)

and eqn. (1) likewise by

(eqn. 5)

Eqns (4) and (5) form the so-called leap-frog scheme for
integrating Newton's equations of motion. Typically a time step of 1
to 10 fs is used for molecular systems. Thus a 100 ps (10^{-10} seconds)
molecular dynamics simulation involves 10^{5} to
10^{4} integration steps. Even using the fastest computers
only very rapid molecular processes can be simulated at an atomic
level. As with any aspect of modelling, the accuracy of the predicted
dynamics will depend on the validity of the underlying assumptions of
the model. In this case the model is essentially defined by the force
field that is used. For this exercise we will be using the GROMOS96
empirical force field.

This technique is commonly referred to as Molecular Dynamics (MD). A detailed explanation of the concepts behind this numerical technique can be found at http://www.fisica.uniud.it/~ercolessi/md/md/, a web site developed by Furio Ercolessi. A printable version can be downloaded from his page here.

Factors that govern the outcome of MD simulations are:

- choice of the degrees of freedom
- force field parameters
- treatment of non-bonded interactions
- solvation effects
- boundary conditions
- treatment of temperature and pressure
- integration time step
- starting configuration

As with any aspect of modeling, the accuracy of the predicted dynamics will depend on the validity of the underlying assumptions of the model. In this case this is essentially defined by the model for the intermolecular interactions (or potential energy) used. That model is a mathematical function (force field) that describes how the value for the potential energy depends on the spatial arrangement of all the atoms.

The following is designed to acquaint you with the general features of the molecular dynamics software package Gromacs. This will be done by using lysozyme as an example. Gromacs is a widely used molecular dynamics simulation package developed at the University of Groningen. Information on Gromacs can be found at http://www.gromacs.org/.

To run a simulation several things are needed:

- a file containing the coordinates for all atoms
- information on the interactions (bond angles, charges, Van der Waals)
- parameters to control the simulation.

The .pdb or .gro file contains the coordinates for all atoms and is the input structure file for MD simulation. The interactions are listed in the topology (.top) file and the input parameters are put into a .mdp file. To get an idea of the different file types processed by Gromacs, follow this link.

The exercise falls apart in four sections, corresponding to the actual steps in an MD simulation.

- Conversion of the pdb structure file to a Gromacs structure file, with the simultaneous generation of a descriptive topology file.
- Energy minimization of the structure to release strain.
- Running a full simulations.
- Analyzing results.

©2005 T.A.Wassenaar