Time length molecular dynamics

Question

What limits time length of simulation in molecular dynamics to femto and in coarse graining to nano or micro second?

Answer 1

There is a practical and a "theoretical" limit to why you can't just run longer with molecular dynamics simulations.

First, the practical limit. If you are looking at an atom and you want to accurately resolve something, you need to have a sufficient number of temporal integration points per time scale of interest. Since we are looking at atoms and molecules, everything we might want to know is related to the motion and vibration of the atoms in the lattice. So let's just take a really simple example of simulating a single Cesium-133 atom. This is the atom used in atomic clocks and so we know it oscillates at roughly $1\text{e}9~\text{Hz}$. So to resolve 1 period of motion for the atom, we might need, say, 100 temporal points per period. We're looking at something like a timestep of roughly 10 picoseconds. To run "longer" physical times, you need a whole lot of time steps to get anywhere more practical. In problems like this example, it's easy to just run longer. But if it took you 6 months to get your 100 billion atom system to reach a point where you can get some statistical information, another 6 months to get double the physical time may not be something you can afford. Many of the simulations I compare against are run at tenths or hundredths of femtoseconds in order to be stable and accurate, so just running longer is ridiculously expensive.

Now the "theoretical" limit. I put it in quotes because it's really due to the implementation in hardware that limits it. You should know (and if you don't, start digging into it) that the Explicit Euler scheme is bad. It's simple, sure, but it's bad -- it imposes a rather harsh time step restriction, but most importantly it is unconditionally unstable. It doesn't conserve energy and it will blow up eventually. And being first order, it's not really accurate.

An alternative is a fully implicit method. These are unconditionally stable, allowing for arbitrarily large time steps, but in practice you can't do that. You need to resolve the physical processes of interest so you still need to select reasonable time scales. So before when we took 100 temporal steps per period, that is probably due to stability constaints -- maybe we can only take 10 steps and still get a good answer if our method is stable. The downside is fully implicit methods typically require a global matrix inversion. These don't scale well in parallel as processor counts get large. So they tend to be avoided.

What is used instead are semi-implicit methods, or symplectic methods. The most common in MD is the Verlet integration (either position or velocity depending on the needs). Looking at Verlet specifically, it's a second-order accurate method and allows much larger time steps compared to a fully explicit scheme, although smaller than a fully implicit scheme. Most importantly, it is energy conserving in theory.

But here's the problem -- although energy conserving in theory, no scheme is actually energy conserving in practice. Finite precision in hardware means errors accumulate over time. So if you are running at 0.01 femtoseconds and you ran 1 billion time steps, you probably accumulated some numerical error due just to finite precision. At best this means you answer drifted and is wrong, but at worst your simulation will be unstable and blow up. Verlet has a tendency to do the latter and it will be unstable for large numbers of timesteps. So you cannot just keep running longer -- you are either inaccurate or unstable.

As an aside, you also run into more numerical precision errors -- trying to add something like $1e\text{-}17$ to $1e\text{-}9$ could lead to underflows. This can be fixed by using higher precisions numbers (quad instead of double) but then it just doubled memory space and probably compute time, so it will take many times longer in CPU time to get to the same physical time. And now we're back at the practical problem.

When you coarse-grain your simulation, you can take larger time steps because you no longer are interested in the smallest physical time scales. Rather than looking at the motion of individual atoms relative to one another, you look at the motion of individual molecules relative to one another. In length and time scales, this may be orders of magnitude larger. So now you can run at nano- or picosecond time steps instead of femtoseconds. When you ran your original simulation 1 billion steps and only got to 1 nanosecond, now you can run for 10000 time steps to get to 1 nanosecond. Fewer time steps means less accumulated numerical error; alternatively, if you ran the same number of time steps until it was unstable, you've gone way further in time in the coarse simulation.

All of these issues (and more) are not confined to MD. The same is true in any numerical simulation.