by Professor Niels Grønbech Jensen, Department of Mathematics, Department of Mechanical and Aerospace Engineering, UC Davis, USA.
Numerical simulations of atomic and molecular ensembles by Molecular Dynamics always involve discretization of time, and as the time step is increased the discrete-time behavior becomes increasingly different from that of the anticipated continuous-time dynamics. This creates a dilemma for any simulation of a dynamical system: use a small time-step, resulting in dynamics that resemble the desired continuous-time behavior at the expense of computational efficiency; or use a large time step that makes the simulation finish sooner at the expense of accuracy. Adding to this dilemma is a much overlooked fundamental problem with discrete-time dynamics; namely that the so-called velocity variable does, in fact, not produce the conjugate momentum variable to the simulated position in discrete time.
Thus, the commonly used Verlet algorithm in Molecular Simulations, is not only producing errors that depend on the square of the simulated time step, it is further inducing those errors differently in velocity and position, such that the self-consistency between kinetic and configurational measures is broken. Effects arising from this problem include that kinetic energy cannot be used as a precise indicator of configurational sampling temperature, and this results in a number of potentially severe issues if the time step is challenged. For example, isothermal ensembles that are simulated by methods, in which the measured kinetic temperature is used to control the thermodynamics, will inevitably lead to erroneous configurational sampling. It is therefore essential to understand the features of different algorithms such that optimal properties can be chosen for a given set of problems and objectives.
Our aim is here to investigate and improve the simulation techniques for systems in thermal equilibrium. We briefly review our stochastic Størmer-Verlet algorithm for the evolution of Langevin equations in a manner that preserves proper configurational sampling (diffusion and Boltzmann distribution) in discrete time. The resulting method, which is as simple as conventional Verlet schemes, has been numerically tested on both low-dimensional nonlinear systems as well as more complex molecular ensembles with many degrees of freedom. In light of the fundamental artifacts introduced by discrete time, we provide a simple intuitive picture of the unique benefits of our algorithm. We then discuss the solution to the “velocity problem”, and we suggest a constructive approach to achieving simultaneous measures of high quality for both configurational and kinetic sampling in discrete time. Time permitting, we finally introduce a companion algorithm for controlling pressure in molecular ensembles; i.e., a barostat for so-called NPT simulations in Molecular Dynamics.