[gmx-developers] Problems with artificial minima when using soft-core potentials

[Bert de Groot]

Dear Alan,

thank you for bringing this up.

Just to put a few things straight:
-the first thing we tried when running into softcore instabilities was to adopt 
the Beutler approach with separate alpha values for the LJ and Coulombic 
interactions. This indeed helps to a certain extent (as does Berk's fix for 
hydrogens). However, both approaches still lead to additional minima at 
intermediate lambda values for certain atom types, which cannot be avoided by 
any choice of separate alpha values for LJ and Coulomb interactions (or at least 
not for any combination we tried). We explore this in Fig S4 of our paper. The 
artificial minima particularly occur when using fast-growth methods and 
frequently involved . Unfortunately the effect can easily remain unnoticed and 
therefore lead to inaccurate free energy estimates. This issue is solved in our 
softcore approach.
It may be that the systematic approach by Steinbrecher and Case mentioned by 
Michael Shirts yields suitable alpha values for the original Beutler softcore 
also in the cases that we found to be problematic. This would be interesting to 
test.

-you mention that the potential has a discontinuity at r=0 in our approach, 
which may lead to potentially unstable dynamics. This is not the case. Only the 
force is discontinuous (just like the hardcore LJ and Coulomb forces are 
discontinuous at r=0). This does not lead to unstable dynamics. In addition, 
also the discontinuous force is not an issue, simply because r=0 is never 
reached. This holds in the hardcore case as well as in our case. Now one could 
argue that with a softcore potential r=0 might be more easily reached than in 
the hardcore case, because of a lower enthalpic cost. This is indeed true, but 
as we explore in section S4 of our paper, even in single precision the 
probability to hit r=0 is negligible for all practical purposes, not only for 
enthalpic but particularly for entropic reasons.

Our softcore implementation is currently obtainable as patch for gromacs 4.0.7, 
and will soon (when we're done testing) be available for 4.5 and commited to git 
for future gromacs releases (4.6 or master).

best regards,

Bert

______________________________________
Prof. Bert de Groot, Ph.D.

Max Planck Institute for Biophysical Chemistry
Computational biomolecular dynamics group
Am Fassberg 11
37077 Goettingen, Germany

tel: +49-551-2012308, fax: +49-551-2012302
http://www.mpibpc.mpg.de/groups/de_groot

On 09/03/2012 05:08 PM, Alan E. Mark wrote:
> Dear Berk,
>              This has been coming up in various circumstances for a long time.
> Bert de Groot recently however published a paper in JCTC highlighting the
> problem in GROMACS and proposing a complex solution which in fact makes the
> potential have a discontinuity at r=0. This also makes the dynamics potentially
> unstable. Stable dynamics was exactly why the soft-core was introduced in this
> form in the first place.  He claims there is a problem with our soft-core when
> in reality the problem is in the form implemented in GROMACS which is a pity.  I
> would also note that  by being able to scale the alpha parameters for the LJ and
> electrostatics separately (i.e. using a larger shifting distance for the
> electrostatics) avoids any need to turn off the electrostatics first. Turning
> the electrostatics off completely first can add the the range of non-physical
> states that need to be sampled making the calculations slower and more complex
> than necessary.
>
> Alpesh and I are writing a response to Bert's paper and a few other factors. I
> will send you a copy when we have a good draft.
>
> Alan
>
>