Generating random numbers for long molecular dynamics simulations

Question

I've written a basic 2D Langevin dynamics simulator in C++, for a particle in a potential, solving the equation:

$$M\ddot{X} = - \nabla U(X) - \gamma M \dot{X} + \sqrt{2 \gamma k_B T M} R(t)$$

This relies on a source of random numbers for $R(t)$, and I've been using std::random (in particular, a Mersenne Twister RNG).

My current approach is to generate the random numbers all in one go at the beginning of the simulation, rather than at each time step, as it was (for me at least) faster. This is fine for runs up to $10^{8}$ time steps, but eventually I run out of RAM to store the random numbers.

I'd like to at least have the ability to run longer simulations (I'm investigating diffusion behaviour of particles in simple potentials), and would like to know what potential pitfalls to look out for when either

• splitting the simulation into segments, or
• breaking up the calculation of the random numbers into smaller chunks (say do a block of numbers ever $10^{6}$ steps.

One thing I noticed was the risk of repeating sequences - see "Vulnerability in Popular Molecular Dynamics Packages Concerning Langevin and Andersen Dynamics"

Answer 1

I suspect the reason why generating the random numbers on the fly is slower for you is due to the rather large state of the Mersenne Twister.

Switching to something like the PCG or XorShift+ random number generator would have several advantages for you:

1. Higher quality of randomness (Mersenne Twister fails several tests for randomness)
2. Smaller state, so generating the random number on the fly should be faster than storing in memory, thus solving your problem.
3. General speed - both PCG and XorShift+ are around 2x faster than the Mersenne Twister.

As for repeating sequences, this should not happen if you initialize the generator properly, using std::seed_seq. See the PCG C++ documentation for instance.