[gmx-developers] Gromacs 3.3.4 PME energies keep changing with longer cutoffs

[Gabriel Rocklin]

Hi all,

I'm seeing some confusing behavior with Gromacs 3.3.4 and PME
electrostatics.  My potential energies seem to not to converge a single
"correct" value as I make the PME parameters more and more exact.

I took 10 snapshots from a single precision, constant volume simulation I'm
doing of a protein in water and saved them as .gro files with -ndec = 8.  I
then used mdrun -rerun to re-evaluate their energies using a number of
different PME parameters.  I found that the energies were not sensitive to
ewald_rtol or fourierspacing (below certain values), but they were very
sensitive to rlist and rcoulomb even out to very long distances.

For example, using pme_order=6, spacing=0.1, and tolerance=0.00000001, as I
increase rlist and rcoulomb together from 0.7 to 2.5nm, I get these
potential energies for my 10 snapshots:

0.7 -656927.375 -657225.4375 -656475.0 -656451.875 -657604.8125 -656090.6875
-657553.25 -657241.8125 -657747.9375 -656863.5625
0.9 -654605.375 -654916.1875 -654163.9375 -654132.125 -655300.875 -653771.5
-655239.875 -654936.875 -655420.75 -654539.25
1.0 -654553.5625 -654862.0625 -654108.9375 -654078.0 -655246.875
-653718.4375 -655182.6875 -654884.3125 -655365.8125 -654486.375
1.1 -654536.875 -654842.5 -654094.3125 -654061.25 -655227.9375 -653701.4375
-655166.5 -654864.3125 -655349.625 -654468.75
1.2 -654533.5625 -654842.625 -654087.1875 -654054.5 -655224.625 -653695.125
-655164.25 -654863.5625 -655340.5 -654463.125
1.3 -654531.9375 -654837.0625 -654087.5625 -654054.5625 -655223.6875
-653694.125 -655161.6875 -654857.5 -655345.0625 -654462.5
1.5 -654532.375 -654839.8125 -654086.1875 -654057.25 -655223.25 -653693.3125
-655164.4375 -654864.1875 -655340.1875 -654461.5625
1.8 -654527.125 -654840.3125 -654086.9375 -654054.0625 -655221.25
-653696.6875 -655155.75 -654860.0625 -655343.5 -654465.8125
2.0 -654532.5625 -654843.125 -654083.75 -654053.6875 -655223.375 -653695.375
-655158.8125 -654860.1875 -655342.3125 -654464.4375
2.5 -654535.875 -654838.75 -654089.875 -654055.0625 -655220.25 -653692.75
-655163.4375 -654857.375 -655337.125 -654461.375

So not only do the numbers keep changing, but the relative energies change
as well, so that for example going from 2.0 to 2.5 nm, the standard
deviation of the change in energy across the 10 snapshots is 4 kJ.

The system is a dodecahedral box with a volume of about 500 nm^3.

Is this normal behavior?  Should these values converge to specific energies?
 And how far out would you normally have to set rcoulomb/rlist to see this?

Thanks,
Gabriel

Gabriel Rocklin
UCSF Shoichet/Dill Labs
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-developers/attachments/20100216/995db0fd/attachment.html>