Bug #2215

The units of "rinv" and "sh_ewald" are inconsistent.

Added by Chong-Li Zhao about 3 years ago. Updated almost 3 years ago.

Target version:
Affected version - extra info:
Affected version:


The following expression in "gromacs/gmxlib/nonbonded/nb_generic.cpp":
rinvcorr = (fr->coulomb_modifier == eintmodPOTSHIFT) ? rinv-fr->ic->sh_ewald : rinv;

Unit of "rinv" is inverse of nanometer;
"sh_ewald" is just a number without unit in that
"ic->sh_ewald = gmx_erfc(ic->ewaldcoeff_q*ic->rcoulomb);"
is defined in

Associated revisions

Revision 3f970c00 (diff)
Added by Mark Abraham almost 3 years ago

Fix value of Ewald shift

In all the short-ranged kernel flavours, sh_ewald is subtracted from
rinv, which have inconsistent dimensions. Fortunately, rcutoff is
often close to 1, and the inconsistent shifts often cancel in
practice, and energy differences computed on neighbour lists of the
same size will have the error cancel. The difference doesn't even
show up in the regressiontests, but would if we had a unit test
of a single interaction.

Fixes #2215

Change-Id: Ia2ea57f3bd9d521879783b207353d9d6f4ccb4a8


#1 Updated by Mark Abraham about 3 years ago

There seems to be a minor, but embarrassing, issue here. sh_ewald is used as if it was computed as erfc(beta * cutoff)/ cutoff but computed only as erfc(beta * cutoff). Fortunately, cutoff is approximately one, and the errors from opposite charges approximately cancel, and for a given cutoff energy differences will have even more fortuitous cancelling.

I think we should subtract the cutoff distance from sh_ewald in forcerec.cpp. And plan to use a name that is more descriptive of the quantity, ie that it is a potential shift (albeit in 1/r form).

This issue was introduced along with the Verlet scheme, and now affects the group scheme also, in every Ewald kernel that I have seen.

#2 Updated by Gerrit Code Review Bot about 3 years ago

Gerrit received a related patchset '1' for Issue #2215.
Uploader: Mark Abraham ()
Change-Id: gromacs~release-5-1~Ia2ea57f3bd9d521879783b207353d9d6f4ccb4a8
Gerrit URL:

#3 Updated by Mark Abraham about 3 years ago

  • Target version set to 5.1.5
  • Affected version changed from 2016.3 to 5.1.5

#4 Updated by Berk Hess about 3 years ago

  • Category set to documentation
  • Status changed from New to Fix uploaded

Note that this bug only affects the reported energies, not the forces or the sampling.
And, as Mark noted, the effects on the energy are very small, so this is relatively harmless (but still embarrassing).

#5 Updated by Mark Abraham almost 3 years ago

  • Status changed from Fix uploaded to Resolved

#6 Updated by Mark Abraham almost 3 years ago

  • Status changed from Resolved to Closed

Also available in: Atom PDF