[gmx-developers] milestoning

Mark Abraham mark.j.abraham at gmail.com
Fri Aug 2 16:40:46 CEST 2013

On Wed, Jul 31, 2013 at 7:16 PM, Sergio Decherchi
<Sergio.Decherchi at iit.it> wrote:
> Dear Gromacs Developers,
> I would like to allow Gromacs to do Markovian Milestoning because it seems a
> very
> effective algorithm for the estimation of both free energy and kinetics. To
> achieve
> this goal I have to confine the simulation inside a Voronoi cell.
> To do this confinment, in litterature it is suggested (see Venturoli and
> Vanden Eijden "Markovian Milestoning") that one can, upon Voronoi cell
> violation, invert velocities and reset the position to the previous
> integration step such that the system can remain inside the Voronoi cell.
> I decided to interface to Gromacs via plumed2 such that collective variables
> can be
> used to compute the metric that rules the Voronoi cell.
> Following plumed2 philosophy I got velocities just before the do_force call
> in md.c.
> Assuming a velocity verlet integrator which are the velocities that I get in
> this point of the code?
> I mean, which time step?

With VV, the positions and velocities are in sync at all times, by
construction. Except during the update stage, of course, but obviously
that is not happening right before do_force().

> Do you suggest a specific strategy to invert velocities and set the previous
> positions?

There's no elegant way to do that. GROMACS is currently designed
around the MD assumption that the next configuration will always
evolve from the former. One would have to checkpoint the entire state
in memory to be able to "roll back" like you might do with MC.

> Gromacs expects an atom in a domain decomposition cell; after the reset of
> the position to the previous step, this cell may have changed thus leading
> to an inconsitency. I got this kind of problems.
> How can I set the positions to the previous integration step and at the same
> time assuring
> that Gromacs is aware of this change of positions from the domain
> decomposition stand point?

I gather the "fix" doesn't need recalculating forces. The simplest
seam to use would be right before the update of positions and
velocities. Construct backups of x and v at that point, call the
update routines, check for whatever violations might arise, and fix
them (either by changing the result of the update or changing the
velocities that go into the update and then calling the update again).
AFAIK, DD does not balance load between neighbour-search steps, and
since the above intercept of the update phase completes before that
even starts.

> This problem can be avoided by using particle decomposition, but it is
> probably suboptimal.
> I have done many trials and I have other issues.
> Often lincs algorithm fails; when I compute
> the 'reflection' on the Voronoi facet I use subset of atoms because the
> metric is defined in a
> subset of atoms. To get the maximal stability should I change velocities and
> reset positions
> of the whole system? Still lincs should be enabled or disabled at all?

Not sure if your problems are intrinsic to the method or your attempts
to implement it. LINCS will be fine if you've fed something physically
reasonable into the update stage, but if you're abusing a bond by
moving its atoms in awkward directions, LINCS will need help (at


> Thanks,
> Cheers,
> Sergio Decherchi
> --
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-developers-request at gromacs.org.

More information about the gromacs.org_gmx-developers mailing list