Berk,
I have added the check.
The problem was that Gromacs did not truncate the potential
at rcoulomb with PME, only rlist was used. However rcoulomb
was used in the determination of beta.
To avoid that people thought that rcoulomb had an effect
on the cut-off, I have added the check.
OK, so, in other words, rcoulomb now only affects beta with PME; rlist
is used for the electrostatic cutoff.
Is there a straightforward way this could be changed in future
versions so that the cutoff used for electrostatics could be different
from rlist? In particular, I am concerned that this means that one
can't vary the vdw cutoff independently of the electrostatics cutoff
beyond some limit. That is, rlist needs to be >rvdw with shift/switch
vdw potentials:
WARNING 1 [file equil_constp0.7.mdp, line unknown]:
For energy conservation with switch/shift potentials, rlist should be 0.1
to 0.3 nm larger than rcoulomb/rvdw.
Suppose I want to run with a 1.4 nm rvdw -- then I now have to use
rlist=1.5 nm and rcoulomb=1.5 nm? That doesn't seem good. One would
like to be able to change rvdw independently of the electrostatics --
for example, to explore the effects of long range van der Waals
cutoffs. In fact, Michael Shirts and I and others had been trying to
do exactly that in GROMACS 3.1.4 and 3.3 (adjust rvdw and rlist
independently of rcoulomb in binding free energy calculations, in
order to see how strongly the choice of vdW cutoff affects computed
binding free energies). If changing rlist also changes the real space
cutoff for electrostatics, things are substantially more complicated
and it seems unclear how to even address this issue.
Maybe a solution would be to NOT use rlist for beta and for the PME
potential, and instead to use rcoulomb or some extra variable... Then
one could change the neighbor list (rlist) and rvdw without having to
tweak the electrostatic interactions...
What is the optimal cut-off scheme is a different issue.
Indeed one would always want the force to go smoothly
at the cut-off, or before the cut-off in case one has charge groups
or nstlist>1.
However for PME one can not have 'exact' electrostatics while
the particle-particle force is zero at the cut-off since the reciprocal
space requires an error function contribution in real space.
Therefore the real space interaction has infinite range, but decays
very rapidly.
In PME the real space interaction is erfc(beta r)/r, which in Gromacs
is applied to all atom pairs in the neighborlist. The cut-off error
made here is very small, since ewald_rtol=1e-5 (I also often use 1e-6).
The alternative would be to somehow make the direct space interaction
go exactly to zero before the cut-off. This would lead to larger
errors, as one then misses more of the electrostatic interactions.
The only advantage would be that the integration could be more
reversible (assuming that all other algorithms are well reversible).
Given the accuracies of all other algorithms I would say it does not
make sense to remove some electrostatic interactions at the cut-off
when they are already 1e-5 or 1-e6 smaller than the Coulomb
interaction at that distance.
So, if I understand properly, you're saying that the notion of an
electrostatics "cutoff" doesn't really make sense with PME, so my
logic about having the neighborlist extend past the "cutoff" doesn't
really apply. Rather, one just wants to ensure that the "cutoff" is
sufficiently large that electrostatic interactions are fairly small
and thus any inaccuracies due to neighbor list issues will be
basically negligible...?
Thanks,
David
Berk.
_______________________________________________
gmx-users mailing list gmx-users@gromacs.org
http://www.gromacs.org/mailman/listinfo/gmx-users
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
_______________________________________________
gmx-users mailing list gmx-users@gromacs.org
http://www.gromacs.org/mailman/listinfo/gmx-users
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
