LINCS: A Linear Constraint Solver for Molecular Simulations

1
LINCS A Linear Constraint Solver for Molecular
Simulations

• Berk Hess, Hemk Bekker, Herman J.C.Berendsen,
  Johannes G.E.M.Fraaije
• Journal of Computational Chemistry, 1997

Ankur Dhanik
2
Outline

• Introduction to Molecular Dynamics
• Problem description
• Some solutions
• LINCS
• Results

3
Introduction to Molecular Dynamics

• A classical molecular simulations method
• Successive configurations of the system are
  generated by integrating Newtons laws of motion
• Integration algorithms should strive to reduce
  computation and conserve energy
• Various integration algorithms
• Standard Verlet algorithm
• Leap-frog algorithm

4
Introduction to Molecular Dynamics

• Standard Verlet Algorithm
• r(t?t) 2r(t) - r (t-?t) ?t2a(t)
• v(t) r(t?t) - r(t-?t)/2?t
• Velocities are not directly generated, and are
  one time step behind
• Leap-frog algorithm
• r(t?t) r(t) ?tv(t?t/2)
• v(t?t/2) v(t-?t/2) ?ta(t)
• v(t) v(t?t/2) v(t-?t/2)/2
• Velocities are half time step behind

5
Introduction to Molecular Dynamics

• Choosing the time step
• One order of magnitude smaller than the shortest
  motion (bond vibrations)
• Severe restriction as these high frequency
  motions have minimal effect on the overall
  behavior of the system
• Constrained dynamics

6
Introduction to Molecular Dynamics
System Type of motion present Suggested times step(s)
Atoms translation 10-14
Rigid molecules Translation and rotation 5 X 10-15
Rigid molecules, rigid bonds Translation, rotation, torsion 2 X 10-15
Rigid molecules, flexible bonds Translation, rotation, torsion, vibration 10-15or 5 X 10-16
The different types of motion present in various
systems together with suggested time steps
7
Introduction to Molecular Dynamics

• In constrained dynamics bonds and angles are
  forced to adopt specific values throughout a
  simulation
• Constraints are categorized as holonomic and
  non-holonomic.
• Suppose a particle on the surface of a sphere
• r2 a2 0 Holonomic
• r2 a2 gt 0 Non-holonomic
• In a constrained system
• Particles are not independent
• Magnitude of constraint forces are unknown

8
Introduction to Molecular Dynamics

• Integration algorithms
• Standard Verlet Algorithm
• Leap-frog algorithm
• Choosing time step
• Constrained dynamics
• Holonomic and non-holonomic constraints

9
Outline

• Introduction to Molecular Dynamics
• Problem description
• Some solutions
• LINCS
• Results

10
Problem Description

• Design an algorithm for solving constrained
  molecular dynamics, the constraints being
  holonomic
• The algorithm should strive for following
  features
• Numerical stability
• Energy conservation
• Computational efficiency

11
Some solutions

• Reset coupled constraints after an unconstrained
  update
• Non-linear problem
• SETTLE
• Solves analytically
• Very fast, but unsuited for large molecules
• SHAKE
• Iterative method
• Sequentially all the bonds are set to the correct
  length
• Simple and numerically stable
• No solutions may be found when displacements are
  large, difficult to parallelize

12
Some solutions

• EEM
• The second derivatives of constraint equations
  are set to zero
• All the constraints are dealt with simultaneously
• Corrections are necessary to achieve accuracy and
  stability

13
LINCS

• Does an unconstrained update
• Sets the projection of the new bonds onto the old
  directions of the bonds to the prescribed lengths
• Similar to EEM, with some practical differences
• Implements
• Efficient solver for the matrix equation
• A velocity correction that prevents rotational
  lengthening
• A length correction that improves accuracy and
  stability

14
LINCS
15
LINCS
Newtons equation of motion
(I-TB) is projection matrix which sets the
constrained coordinates to zero. T transforms
motions in the constrained coordinates into
motions in cartesian coordinates
16
LINCS
17
LINCS

• Numerical errors can accumulate which leads to
  drifts

Velocity correction
Position correction

• Drawback the projection of the new bonds onto
  the old directions rather than the new bonds are
  set to prescribed lengths

18
LINCS

• Correction for rotational lengthening
• To correct rotation of bond i, the projection of
  the bond on the old direction is set to

• The corrected positions are

19
LINCS
Constrained new position is given by
20
LINCS

• The first power of An gives the coupling effects
  of neighboring bonds
• The second power gives the coupling effect over a
  distance of two bonds
• The inversion through a series expansion makes
  parallelization easy
• In one timestep, the bonds influence each other
  when they are separated by fewer bonds than
  highest order in expansion

21
LINCS

• Parallelization
• Consider a linear-bond constrained molecule of
  100 atoms to be simulated on a dual processor
  computer
• Uses rotation correction and an expansion to the
  second power of An
• Because order of expansion is two, bonds
  influence each other over a distance of 6
• Update of position and call of LINCS algorithm
  must be done for atom 1-56 and 44-100 on
  processors 1 and 2 respectively
• 1-50 update from processor 1 and 50-100 from
  processor 2

22
Results
23
Results

• Solves constrained molecular dynamics
• Numerically stable
• Conserves energy
• Three to four times faster than SHAKE
• Can be easily parallelized