From: "David Mobley" <[EMAIL PROTECTED]>
Reply-To: Discussion list for GROMACS users <gmx-users@gromacs.org>
To: "Discussion list for GROMACS users" <gmx-users@gromacs.org>
Subject: [gmx-users] confused about rcoulomb<=rlist
Date: Wed, 12 Jul 2006 07:58:23 -0700
Dear all,
I desperately need some clarification about how GROMACS handles
coulomb interactions with rcoulomb and rlist, and why it makes sense
to do things this way.
In particular, there was a thread on the revisions list
(http://www.gromacs.org/pipermail/gmx-developers/2006-April/001581.html)
a while back explaining that in 3.3.1, grompp checks to ensure that
rcoulomb=rlist for PME/Ewald/PPM. I was very surprised by this, as I'd
always been using rlist>rcoulomb, which I thought was best. Can
someone explain why one would want rcoulomb=rlist? I don't quite
understand the explanation on the revisions/developers list.
My confusion is exacerbated by the fact that the GROMACS 3.3 manual
includes the following statement:
"The neighbor search cut-off rlist should be 0.1 to 0.3 nm larger than rvdw
to accommodate for the size of charge groups and diffusion between neighbor
list updates."
This is exactly the same logic I was applying to rcoulomb when I
concluded that I wanted to make sure rlist was larger than rcoulomb --
basically, that I didn't want things outside rlist diffusing to within
rcoulomb between neighbor list updates, and hence not being counted in
the real-space interactions when they should be. Is there some other
concern that outweighs this or means that this point is irrelevant?
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.
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.
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
