[gmx-developers] Re: water optimization?

[Berk Hess]

Michael Shirts wrote:

>Hi, Erik-
>
>I'm working with David Mobley and John Chodera at UCSF to port my free
>energy changes that I put in 3.1.4 into 3.3.  Here are some of the
>things we've been thinking about -- it would be great to hear your
>thoughts (and I'm cc'ing this to the developer's list to get
>additional thoughts). Are there any suggestions or thoughts on these
>subjects as we work on getting them in?
>
>Here are our thoughts:
>
>1) A toggle for decoupling versus mutating interactions of perterbed atoms.
>
>2) Being able to print out the energy differences to arbitrary numbers
>of different lambda values from the equilibrium sampled value, for use
>in FEP, or BAR/WHAM
>
>3) Eventually, a new g_bar (or perhaps modification of g_wham) to
>handle these calculations with the standard gromacs tools.
>
>4) Evaluating the free energy only every N steps (it becomes expensive
>when multiple PME calls are required)
>
>5) A good tutorial for free energy calculations, with sample files, etc.
>
>David, John, am I missing anything?  Anyone, any thoughts?
>
>  
>
I don't understand how you want to do the decoupling exactly.

I have implemented free-energy for PME in 3.3.
Currently we interpolate the interactions, not the charges,
so we always have to do two PME calculations per step.
However, the quite expensive setup of the charge interpolation
functions only needs to be done once. With the new 64-bit loops
on an AMD64 FE wiht PME runs about a factor of 0.8 slower than a normal run.
I think this penalty is so small that it is not worth implementing
point 4 (which would also require implementing charge interpolation).

For decoupling you want a decoupled state where the decoupled
molecule has plain Coulomb interactions without cut-off.
This is not possible with the normal neighborlists, as there is
no guarantee that the molecule is smaller than the cut-off length.
I think the simplest way to accomplish decoupling is with pair
interactions.

My proposition is:

Add a pairs interaction type 2 without parameters that uses
the normal LJ parameters and a FudgeQQ of 1.
Add pairs type 3 and 4 that are equivalent to 1 and 2,
except that they always use the A state parameters.

Add a decouple option to the mdp file where the user can
select a moleculetype.
For this moleculetype grompp should then change all pairs
type 1 to type 3 and for all non-excluded atom pairs add
exclusions and pairs type 4.

Berk.