[gmx-developers] Re: Various integrator questions

[Michel Cuendet]

Hi

Michael wrote :

>... the overhead of keeping
>Berendsen is very low, and it can easily be included at this time. 
>Erik Lindhal pointed out in a separate e-mail that it is very stable
>and quick in relaxing to near-equilibrium ensembles.
>
In my opinion it is mandatory to keep the berendsen thermostat as an 
option. Nose-Hoover is good if you want to reproduce the right ensemble, 
and if you have enough sampling time. But NH is not appropriate for 
equilibration purposes. When the initial temperature is far from the 
target temperature, NH can give nasty oscillations, and due to it's 
second order structure will require much more time to stabilize than the 
berendsen thermostat (first order).

By the way, it would be nice (and not too difficult) to implement 
Nose-Hoover Chains. First, only NHC, and not NH, mathematically gives 
the canonical ensemble in the case where the center of mass is fixed 
(the physicist would say : if the system has enough d.f. this doesn't 
matter). But more importantly, the chains prevent the thermostat 
variable to enter a resonant mode with periodic oscillations, and allows 
for faster equilibration.

>>I'll put in a plug again for Velocity Verlet, which avoids this whole
>>>issue, as the velocities are determined at the same time as the
>>>coordinates.  I'll start working on implementing my local version on
>>>the 3.3 codebase, so that people can check it out.
>>>
>>>Best,
>>>Michael
>>>
>>    
>>
>
>Although by itself correct, there is another issue here.
>
>Although Velocity Verlet avoid kinetic energy and pressure complications,
>it introduces complications with the constraint algorithms as the velocities
>need to be constrained as well, whereas there are automatically correct 
>with leap-frog.
>I guess there is no conceptual problem and one can always just constrain
>the full-step velocities. But this does involve coding a velocity 
>version for all
>constraint algorithms. (I have already implemented LINCS code for this).
>
>Berk.
>

The integrators I was advertizing for are precisely generalizations of 
VV, based on the Trotter expansion of the Liouvilean. They have been 
explicitely derived also in the case of constraints. Martyna et al. were 
even so much concerned that these integrators should be implemented in 
MD packages, that they gave the full fortran code for them in the annex 
of their paper !! Certainly a good source of inspiration for the 
implementation, not only on the theoretical level.

Well, I just warmly recommend the three papers below...

Martyna et al. "Explicit reversible integrators for extended systems dynamics", Molecular Phys. 87(5) : 1117, 1996.

Ciccotti and Kalibaeva, "Molecular Dynamics of Complex Systems: ...", Lecture Notes in Physics 640 : 150, 2004.

Kalibaeva et al. "Constant pressure - Constant Temperature Molecular Dynamics: a correct constrained ensemble using the molecular virial", Mol. Phys. 101(6) : 765, 2003

Cheers,

Michel