Error in the invariant measure of numerical discretization schemes for canonical sampling of molecular dynamics
View/ Open
source_files.zip (13.90Mb)
Date
28/11/2013Author
Matthews, Charles
Metadata
Abstract
Molecular dynamics (MD) computations aim to simulate materials at the atomic level
by approximating molecular interactions classically, relying on the Born-Oppenheimer
approximation and semi-empirical potential energy functions as an alternative to solving
the difficult time-dependent Schrodinger equation. An approximate solution is
obtained by discretization in time, with an appropriate algorithm used to advance the
state of the system between successive timesteps. Modern MD simulations simulate
complex systems with as many as a trillion individual atoms in three spatial dimensions.
Many applications use MD to compute ensemble averages of molecular systems at
constant temperature. Langevin dynamics approximates the effects of weakly coupling
an external energy reservoir to a system of interest, by adding the stochastic Ornstein-
Uhlenbeck process to the system momenta, where the resulting trajectories are ergodic
with respect to the canonical (Boltzmann-Gibbs) distribution. By solving the resulting
stochastic differential equations (SDEs), we can compute trajectories that sample the
accessible states of a system at a constant temperature by evolving the dynamics in
time. The complexity of the classical potential energy function requires the use of
efficient discretization schemes to evolve the dynamics.In this thesis we provide a systematic evaluation of splitting-based methods for
the integration of Langevin dynamics. We focus on the weak properties of methods
for confiurational sampling in MD, given as the accuracy of averages computed via
numerical discretization. Our emphasis is on the application of discretization algorithms
to high performance computing (HPC) simulations of a wide variety of phenomena,
where configurational sampling is the goal.
Our first contribution is to give a framework for the analysis of stochastic splitting
methods in the spirit of backward error analysis, which provides, in certain cases,
explicit formulae required to correct the errors in observed averages. A second contribution
of this thesis is the investigation of the performance of schemes in the overdamped
limit of Langevin dynamics (Brownian or Smoluchowski dynamics), showing the inconsistency
of some numerical schemes in this limit. A new method is given that is
second-order accurate (in law) but requires only one force evaluation per timestep.
Finally we compare the performance of our derived schemes against those in common
use in MD codes, by comparing the observed errors introduced by each algorithm
when sampling a solvated alanine dipeptide molecule, based on our implementation of
the schemes in state-of-the-art molecular simulation software. One scheme is found to
give exceptional results for the computed averages of functions purely of position.