# Velocity Verlet Algorithm

The goal of this section is to introduce the velocity Verlet algorithm for molecular dynamics (MD) simulations and to show how to implement it to compute time-correlation functions (e.g., ) and transport coefficients (e.g., the diffusion coefficient ).

Consider the task of computing the diffusion coefficient according to Eq. (496). The quantity of interest is, therefore, the equilibrium ensemble average of the velocity autocorrelation function

 (502)

The computation of  thus requires sampling initial conditions according to the ensemble distribution described by  and for each initial condition compute the value of the particle velocity  at time .

The velocity Verlet algorithm (W. C. Swope, H. C. Andersen, P. H. Berens, and K. R. Wilson, J. Chem. Phys. 76, 637 (1982)) computes the particle velocity  and position  at time t+ as follows:

 (503)

where  is a small time increment,  is the particle mass and  is the total force acting on the particle at time t. Given the initial conditions  and  one can compute  and  simply by applying Eqs.(503) successively n times, with .

Note that by implementing the algorithm one generates a phase space trajectory (i.e., a sequence of snapshots'' for the particle coordinates and velocities at all intermediate times  (with j=1,2,...,n)).

Molecular dynamics simulations thus provide the sequence of microscopic configurations through which the model system passes in time. Such detailed microscopic information allows one to compute the result of a measurement of an observable (i.e., an ensemble average) according to the time average introduced by Eq. (111) (i.e., simply by averaging the value of the observable throught the whole manifold of microscopic configurations generated during the time of the measurement). Therefore, another way of computing the ensemble average introduced by Eq. (502) is

 (504)

where  is the time of the measurement of the diffusion constant , a time that is much larger than the relaxation time of the velocity autocorrelation function.

Exercise:

Compute the velocity autocorrelation function for a fluid of argon atoms using the program developed for computing the radial distribution function. (Hint: substitute the Metropolis Monte Carlo algorithm by the velocity Verlet algorithm).