[gmx-users] Is the AMBER or CHARMM Force field running in GROMACS

[Berk Hess]

Hi,

I have found the problem.
There is a subtle rounding problem is single precision.
In the settle code there are two variables which are used
to determine the center of mass of a water molecule:
wO = mO/(mO+2*mH)
wH = mH/(mO+2*mH)
where mO and mH are the mass of the water Oxygen and Hydrogen.
wO+2*wH should be exactly one, but numerically it can differ one or
two bits from one. This results in a consistent velocity shift
for all water molecules. In single precision this can result in a
significant center of mass motion and pressure deviations.
With the current code the Gromacs choice of masses gives exactly
one for the mentioned sum, but the Amber tip3p masses do not.

A simple solution is to change in src/mdlib/csettle.c the line:
  static real wo,wh,wohh;
to
  static double wo,wh,wohh;
fsettle.f can be corrected analogously.

Berk.

_________________________________________________________________
Talk with your online friends with MSN Messenger http://messenger.msn.nl/