[gmx-users] free energy with TIP4P bug fixed

[chris.neale at utoronto.ca]

I just realized that my previous post was incorrect on its final test.  
The correct conclusion is:

I have done some tests that, for a single system, empirically indicate
that the coulomb(SR) discrepancy occurs only with use of the free
energy code and only when tip4p is in the same energygrps as the atoms  
that are perturbed. The discrepancy is NOT observed when non-perturbed  
atoms from the perturbed molecule are included with tip4p in a single  
energygrps selection.

energygrps = SOL_DPC_DPN_&_!C1_C2_O7      DPN_&_C1_C2_O7
-3.22653e+05

In my haste, I had inadvertently constructed energygrps .ndx entries  
utilizing a non-existant "C7", when I really intended to include "O7"

My Appologies,
Chris.


-- original message (Note: author acknowledges incorrect statements below) --


I have done some tests that, for a single system, empirically indicate
that the coulomb(SR) discrepancy occurs only with use of the free
energy code and only when tip4p is in the same energygrps as any part
of the molecule which is perturbed, even if the common atoms are not
perturbed themselves.

This system is composed of:
SOL = 7180 tip4p
DPC = 49 detergents
DPN = perturbed detergent in which the charges are turned off.

general free energy mdp options in these runs:
init_lambda         = 0.00
sc_alpha            = 0.0
nsteps              = 0

Different mdp options and resulting Coulomb(SR) values are shown below:

## Reproduce the discrepancy:
free_energy=no; energygrps=System
-3.22653e+05
free_energy=yes; energygrps=System
-3.23064e+05

## Discrepancy is resolved with separate energygrps for tip4p,
perturbed molecule, and rest of system:
free_energy=no; energygrps=SOL     DPC       DPN
-3.22653e+05
free_energy=yes; energygrps=SOL    DPC       DPN
-3.22653e+05

## Discrepancy is resolved when energygrps separate tip4p from the
perturbed molecule, with the remainder of the system being flexible to
join with either tip4p or the perturbed molecule:
free_energy=yes; energygrps=SOL         DPN_DPC
-3.22657e+05
free_energy=yes; energygrps=SOL_DPC     DPN
-3.22653e+05

## Reproduce the discrepancy when energygrps treat tip4p together with
the perturbed molecule, even when the remainder of the system is
separated:
free_energy=yes; energygrps=SOL_DPN       DPC
-3.23061e+05

## Reproduce the discrepancy when only 3 atoms of DPN are perturbed
while those 3 perturbed atoms are in a separate energygrps section
from tip4p, while tip4p is in the same energygrps section as the
non-perturbed atoms from DPN (note that the perturbed .itp is
different in this case in order to (i) perturb only 3 atoms, and (ii)
to ensure that perturbed atoms are not in the same charge group as any
non-perturbed atoms, which is required to pass grompp due to the
energygrps selection):
free_energy=yes; energygrps=energygrps = SOL_DPC_DPN_&_!C1_C2_C7
DPN_&_C1_C2_C7
-3.22886e+05

Again, this is not exhaustive in that it is only a single system, so I
can not be sure that there are no other issues, but at least the
discrepancies noted here probably hold for all systems.

Chris.