[gmx-users] ewald_rtol and convergence

[Michael Shirts]

Hi, all-

So, I'm wondering what the effect of ewald_rtol is.  My understanding
is that it should just be shifting the calculation between the real
space and the reciprocal space summations.  But the total
electrostatic is changing non-negligibly as well.  My system is a box
of waters with a protein -- total charge is zero.  PME Cutoff is 0.9,
fourier_spacing is 0.02 (so well converged in all cases), PME order 6.
Octahedral box, (5.78643,5.45550,4.72461) box vector, all the same
configuration.

Ewald_rtol                   total                  Short range       
 Long range         1/beta (nm)
1e-02                  -53195.725708   -46946.948640     -6248.777068 
     0.494129
1e-03                  -53432.220675   -45014.961538     -8417.259138 
     0.386805
1e-04                  -53652.773007   -43093.541701   -10559.231306  
    0.327146
1e-05                  -53872.777698   -41119.255332   -12753.522366  
    0.288146
1e-06                  -54092.976378   -39070.390423   -15022.585955  
    0.260198

I did some comparison to straight Ewald, as well.

Ewald_rtol                   total                  Short range       
 Long range
1e-02                   -52899.579997   -46946.948640    -5952.631357
1e-03                   -53001.090875   -45014.961538    -7986.129337
1e-05                   -53182.993892   -41119.255332   -12063.738560

We see that the short range part is exactly the same in PME and
straight Ewald, as it should, though the long range part is not the
same, even when the PME part is well converged (i.e., I reduced the
fourier_spacing until any changes were sub 0.002 kJ/mol, which is one
part in 25 million).

I'm unsure why the differences.  It there coding issue with the
octahedral box? Is there some effect with the dielectic at infinity
(i'm using epsilon_surface = 0).  What is the 'true' electrostatic
energy?
I tried different cutoffs -- presumably, as we take a low rtol,
pushing the calculation further into the direct space, and choose a
longer cutoff, the correction will get more accurate.

PME: (order 6, ewald_rtol = 1e-02)
cutoff                     total                  Short range        
Long range             1/beta

0.9                  -53195.725708   -46946.948640     -6248.777068   
   0.494129
1.8                  -52875.688467   -49910.307947     -2965.380520   
   0.988258
2.7                  -52798.032604   -50836.767614     -1961.264991   
   1.48239
3.6                  -52764.124490   -51296.887618     -1467.236872   
   1.97652
4.5                  -52745.340510   -51572.825196     -1172.515314   
   2.47065

Ewald: (fourier_nx,ny,nz = 25)

0.9                  -52899.579997   -46946.948640    -5952.631357    
    0.494129
1.8                  -52785.799065   -49910.307947    -2875.491118    
    0.988258
2.7                  -52761.626727   -50836.767614    -1924.859113    
    1.48239
3.6                  -52747.927053   -51296.887618    -1451.039435    
    1.97652
4.5                  -52738.223387   -51572.825196    -1165.398191    
    2.47065

So it seems that -52730 is closer to the "right" answer.  Why isn't
this more possible to with other combinations of terms?  Am I missing
something?  Is this level of agreement (both for Ewald, and between
PME and Ewald) an unreasonable expectation? :)

Thanks,
Michael