Linear Analysis of Stochastic Verlet-Type Integrators for Langevin Equations

We provide an analytical framework for analyzing the quality of stochastic Verlet-type integrators for simulating the Langevin equation. Focusing only on basic objective measures, we consider the ability of an integrator to correctly simulate two characteristic configurational quantities of transport, a) diffusion on a flat surface and b) drift on a tilted planar surface, as well as c) statistical sampling of a harmonic potential. For any stochastic Verlet-type integrator expressed in its configurational form, we develop closed form expressions to directly assess these three most basic quantities as a function of the applied time step. The applicability of the analysis is exemplified through twelve representative integrators developed over the past five decades, and algorithm performance is conveniently visualized through the three characteristic measures for each integrator. The GJ set of integrators stands out as the only option for correctly simulating diffusion, drift, and Boltzmann distribution in linear systems, and we therefore suggest that this general method is the one best suited for high quality thermodynamic simulations of nonlinear and complex systems, including for relatively high time steps compared to simulations with other integrators.

💡 Research Summary

This paper develops a rigorous analytical framework for evaluating stochastic Verlet‑type integrators used to simulate Langevin dynamics. The authors focus on three elementary configurational observables that any correct thermostat must reproduce: (i) the Einstein diffusion constant D = kBT/α for a particle on a flat potential, (ii) the drift velocity vd = f/α for motion on a tilted planar potential, and (iii) the configurational temperature (or variance) of a harmonic oscillator, κ⟨r²⟩ = kBT. These three relations constitute the benchmark for linear systems.

The authors then express any stochastic Verlet‑type integrator in a unified configurational form:

rₙ₊₁ = 2c₁rₙ − c₂rₙ₋₁ + Δt²/m (c₃fₙ + c₄fₙ₋₁) + Δt²/m (c₅βₙ⁻ + c₆βₙ⁺),

where the coefficients cᵢ(γΔt) and the noise correlation ζ(γΔt) are the only free parameters; γ = α/m is the reduced friction. The stochastic increments βₙ⁻, βₙ⁺ are constructed from a weight function ψ(s) on the interval