[gmx-users] "Coul. Recip." in Logfile for Classical Ewald

[Barnett, James W.]

I'm doing a classical Ewald simulation with GROMACS (obviously I should 
be using something other than Ewald in my simulations like PME, but I am 
just trying to understand the concept and write some analysis code.). In 
the logfile does "Coul. Recip." include the self correction (or any 
other kind of correction)? I'm having a little trouble seeing where the 
specific number comes from and am trying to understand. Adding up the 
kspace contribution (from GROMACS) and self contribution (my own 
calculation) doesn't seem to get me the number from the logfile. I know 
I'm missing something.

As an example I ran 33 SPC waters in a 1 x 1 x 1 box with cutoffs at 0.4 
(this is just a small quick simulation to get some numbers for this 
purpose). Again this is classical Ewald, fourierspacing is the default 
0.12, and ewald-rtol is the default 1e-5.

I added the following line to "src/mdlib/ewald.c" just before "return 
energy;" in the "do_ewald" function and recompiled :

    printf("%5.6f\n",energy);

This prints out the k-space contribution to the electrostatic energy at 
each step. As an example, at step 0 this it printed out 1008.472290 
(kJ/mol).

After getting this I calculated the self contribution as follows. The 
logfile says that 1/beta (the Gaussian width) is 0.128065 nm, so beta = 
7.808535 for this simulation. There are 33 water oxygens with a charge 
of -0.82 and 66 water hydrogens with a charge of 0.41. The following 
equation for the self term comes from section 4.9.1 of the manual (f = 
138.935485):

-f * beta / sqrt(pi) * [ (33)(-0.82)^2 + (66)(0.41)^2 ]  = -20372.33 kJ/mol

The logfile states for step 0 that Coulomb (Recip.) is -1.23120e+03 
kJ/mol, but this does not equal the total of the k-space and self 
correction terms:

1008.5 + -20372.3 = -19363.8

I appreciate anybody's help with understanding what I'm missing here in 
understanding how GROMACS calculates the self correction or if Coul. 
Recip. is something else.

Thanks,
Wes Barnett