Fix LJPME Lorentz-Berthelot combination rules

Christian found a bug in the way we evaluate
arithmetic combination rules for LJ-PME reciprocal
space. With this change, both group and verlet kernels
produce virtually identical total energies for both
types of reciprocal-space combination rules. It does
not affect the (default) geometric LJPME case.

Change-Id: I4bb0de1df9588dce3ca94e39a370a0b994b25b53

diff --git a/src/gromacs/mdlib/forcerec.c b/src/gromacs/mdlib/forcerec.c
index 3b2d69ca048b06ff393bba9d74c54d8412632532..2ea1e549c84130320b034463525f6ceeb419c4be 100644 (file)
--- a/src/gromacs/mdlib/forcerec.c
+++ b/src/gromacs/mdlib/forcerec.c
@@ -218,8 +218,8 @@ static real *mk_nbfp_combination_rule(const gmx_ffparams_t *idef, int comb_rule)
                 sigmaj = pow(c12j / c6j, 1.0/6.0);
                 epsi   = c6i * c6i / c12i;
                 epsj   = c6j * c6j / c12j;
-                c6     = epsi * epsj * pow(0.5*(sigmai+sigmaj), 6);
-                c12    = epsi * epsj * pow(0.5*(sigmai+sigmaj), 12);
+                c6     = sqrt(epsi * epsj) * pow(0.5*(sigmai+sigmaj), 6);
+                c12    = sqrt(epsi * epsj) * pow(0.5*(sigmai+sigmaj), 12);
             }
             C6(nbfp, atnr, i, j)   = c6*6.0;
             C12(nbfp, atnr, i, j)  = c12*12.0;