[gmx-developers] suggestion for PME for 3.3.1 for free_energy=yes

[David Mobley]

Dear all,

I wanted to write and suggest that the default PME parameters be set
differently when free_energy=yes in the 3.3.1 release. In particular,
I think the default pme_order needs to be 6 with a fourier spacing of
0.10 or better. ewald_rtol probably isn't as important but I would
suggest 1e-06.

Currently, in 3.3, the default pme_order=4 with ewald_rtol=1e-05 and
fourierspacing=0.12. Michael Shirts and I have done some testing (he
can comment on his testing separately if he likes) and found that with
pme_order=4, results for free energy calculations involving
modifications of electrostatics can be substantially wrong.

My particular testing has mostly been done with rcoulomb=1.0 or 0.9,
although Michael has looked more at varying this, as well. Regardless,
I find that results for computed free energies involving changes in
electrostatics can be wrong by 1-4 kcal/mol with the default pme
parameters as compared to their values with higher pme_order and lower
fourierspacing. For my own work I think I've ended up settling on
pme_order=6 and fourier_spacing=0.10 with ewald_rtol=1e-6, as my
results don't seem to change substantially if I further reduce
fourier_spacing or ewald_rtol, or even increase rcoulomb. Again,
though, I want to emphasize that the difference in free energies going
from the default PME parameters to the ones I suggest here is 1-4
kcal/mol (I tested this on some polar molecules including phenol; the
difference is smaller for less polar molecules), which is roughly
comparable to the total hydration free energy for such molecules.

If for some reason you decide it's not a good idea to change the
default PME parameters (at least for free_energy=yes), I would
strongly encourage adding some suitable warnings if free_energy=yes
with these parameters, and probably also in the documentation.

I can provide exact values for the computed free energies, etc. as
needed. It is possible of course that the *best* parameters for these
things are system-dependent, but it seems pretty clear that
pme_order=4 with fourier_spacing=0.12 is NOT adequate for free energy
calculations involving modification of charges, at least when the
charges are reasonably polar.

Thanks,
David Mobley
UCSF