Module 5:

Molecular Dynamics Simulations

MD simulation trajectories of the trimeric SARS-CoV-2 spike protein and FDA approved or
investigational drug molecules binding to a construct of the SARS-CoV-2 trimeric spike protein.
 Why Not Quantum Mechanics?
• Modeling the motion of a complex molecule by solving the wave functions of
the various subatomic particles would be accurate…

− 2 2
  + U ( x, y, z ) ( x, y, z ) = E ( x, y, z )
2m

• But it would also be very hard to program and take more computing power
than anyone has!
 Classical or Molecular Mechanics

• Instead of using Quantum mechanics, we can use classical Newtonian mechanics

to model our system.

• Calculate the structure and energy of molecules based on nuclear motions.

• electrons are not considered explicitly

• will find optimum distribution once position of nuclei are known

• Born-Oppenheimer approximation of Schrödinger equation

➢ nuclei are heavier and move slower than electrons

➢ electrons move fast enough to adjust to any nuclei movement

 Molecular Mechanics

Molecular mechanics treats a molecule as a collection of

weights connected with springs, where the weights represent the
nuclei and the springs represent the bonds.

• Collection of atom types (to define the atoms in

a molecule), parameters (for bond lengths, bond
angles, etc.) and equations (to calculate the
energy of a molecule)
 MM: Atom Types

In MM, a given element may have several atom types.

• sp3-Hybridized carbons have a tetrahedral bonding geometry
• aromatic carbons have a trigonal bonding geometry.
• C-C bond in the ethyl group differs from a C-C bond in the phenyl ring
• C-C bond between the phenyl ring and the ethyl group differs from all other C-C bonds in ethylbenzene.
 Force
Bonded Terms V ( R ) = Ebonded + E non−bonded Non-bonded Terms
Bond Stretch van der Waals Interaction

Bond Angle

1 1 2
E str = k str ,ij (r ij − r 0)2 Ebend = k b , ijk (θijk − θ0)
2 2
Dihedral

Electrostatic Interaction
Improper
q1q2
Dihedral
F=
E d = k d, ijkl (1+cos (nΦ− Φs)) Eid =
1
k id ,ijkl (Φijkl − Φ0)2
4 0 r 2

2
 MM: Force Field
Force Field used to calculate the energy and geometry of a molecule.

• Total energy of a molecule is divided into several parts called force potentials, or potential energy
equations.
• Force potentials are calculated independently, and summed to give the total energy of the molecule.
• These equations define the potential energy surface of a molecule.
Cars cruise down an expressway at 25 m/s. Engineers are designing an off-ramp in an
interchange with a deceleration of −2.0 m/s2 that lasts 3.0 s.
• What velocity will cars have at the end of the off-ramp?
• What minimum length should the ramp have?

v = v0 + at
v = (25 m/s) + (−2.0 m/s2)(3.0 s)
v = 19 m/s

Δs = v0t + ½at2
Δs = (25 m/s)(3.0 s) + ½(−2.0 m/s2)(3.0 s)2
Δs = 66 m
“Dynamics”
 Equation of Motion

• Newton’s second law or the equation of

motion, F=ma, where F is the force exerted on
the particle, m is its mass and a is its
acceleration.

• Knowing the force on each atom, one can

determine the acceleration of each.

• Integration of the equations of motion then

yields a trajectory that describes the positions,
velocities and accelerations of the particles as
they vary with time.

• From this trajectory, the average values of

properties can be determined.
 Integration Algorithms

• The potential energy is a function of the atomic positions (3N) of all the atoms in
the system.

• Due to the complicated nature of this function, there is no analytical solution to the
equations of motion; they must be solved numerically.

• Numerous numerical algorithms have been developed for integrating the equations
of motion, such as
Verlet algorithm
Velocity Verlet
Leap-frog algorithm
Three Atoms
Progress in time
Atoms move forward in time
 The idea
Molecular Dynamics
Simulations
c what atoms do in real life, assuming a
•
n potential energy
Build realistic function
atomistic model of the system
e energy
• function
Simulate allowsofusyour
the behavior to system
calculate
over the
time force
using specific conditions
perienced by any pressure,
(temperature, atom given
volume,theetc.)
positions of the
er atoms
•
wton’s MD
lawssimulations modelthose
tell us how the motion of some
forces group ofthe
will affect particles by solving the classical
equations
tions of of motion.
the atoms
Energy (U)

Position Position 4
MD Systems
 Molecular Dynamics

MD is a computational microscope. Understanding how things work on a

detailed molecular level ranging from electronic structures to long-time phase
behavior of molecules
 MD algorithm

Set Initial conditions (positions, velocities)

Calculate forces

Solve equations of motion

Move atoms

Increase time

Repeat as long as you need

Calculate the desired physical quantities

 Time Step dt

• Time step dt should be as large as possible to still get accurate trajectories (on the time
scale needed) and conserve of energy.

• In general, dt should be ≈ 0.01 x the fastest behavior of your system

Choosing the timestep: for numerical stability and accuracy in conserving

the energy, one typically needs to pick a time step that is at least one order
of magnitude lower than the fastest time scale in the system.

Atoms oscillate about once every 10-14 s in a solid 

MD time steps are ≈ 10-15 s in simulations of solids.
 Geometry and boundary conditions
• Most simulations probe the structural and
thermodynamic properties of a system of a few
thousand particles.
• For such small systems, the choice of boundary
conditions can strongly affect the measured
properties.

• If nothing special is done, atoms

near the boundary would have less
neighbours than atoms inside.

• This causes surface effects in the

simulation to be much more
important than they are in the real
system.
 Periodic boundary conditions

Consider box of size L,

repeat box infinitely many Basic box
times in all directions
copies
a particle which goes out from the simulation box by one side is reintroduced in the box by the
opposite side.
 Solvation
• Biological activity is the result of interactions between molecules and occurs at the
interfaces between molecules (protein-protein, protein-DNA, protein-solvent, DNA-
solvent, etc).

• Many biological processes occur in aqueous solution

• Solvation effects play a crucial role in determining molecular conformation, electronic

properties, binding energies, etc

• How to model solvation?

• explicit treatment: solvent molecules are added to the molecular system

• implicit treatment: solvent is modeled as a continuum dielectric
 MD – Shortcomings
• Quality of the forcefield

• Size and Time – atomistic simulations can be performed only for systems of
a few tenths of angstroms on the length scale and for a few nanoseconds on
the time scale

• Conformational freedom of the molecule – the number of possible

conformations a molecule can adopt is enormous, growing exponentially with
the number or rotatable bonds.

• Only applicable to systems that have been parameterized

• Connectivity of atoms cannot change during dynamics – no chemical

reactions.
 MD
Programs:

AMBER
CHARMM
VASP (DFT)
XMD
CPMD(DFT) (Car-Parrinello MD)
NAMD
Web-servers for running MD Simulations
• [Link]
• [Link]
• [Link]
• [Link]