VMD-L Mailing List

From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Tue Jun 28 2011 - 18:08:31 CDT

On Tue, Jun 28, 2011 at 5:05 PM, <gaborekt_at_duq.edu> wrote:
> Hi Axel,
>
> 1.  Yes, the systems I have been calculating RDFs for are all periodic.

ok.

> 2.  The basic differences I've seen so far are not large, however they may
> be significant for the analysis I plan on doing.  Some of the peak heights
> and well depths are slightly different in VMD than in CHARMM (only
> different by ~ 0.8 or so), however the most noticeable difference I've
> seen so far deals with the convergence of the RDF.

this difference may be due to choice of histogram bins and
which distance "r" they represent. this paper should have
the details on the choices the code in VMD is making.
http://dx.doi.org/10.1016/j.jcp.2011.01.048

> One thing I've observed is that the RDF plateaus at different values in
> each program.... for example, when I calculate a specific RDF out to 20
> angstroms, in VMD it plateaus at about 0.97 while in CHARMM it plateaus at
> approximately 1.06.... this worries me because I'm planning on doing a

with the code in VMD you should be careful with distances
beyond half the box diameter. periodic replica are not considered,
so your statistical sample gets smaller up to sqrt(1/2) * box.
for any reasonable statistics, the code in VMD should go exactly to 1.0.
please make sure that you have the current version of VMD when
testing this. there has been the occasional obscure and rare bug
that has surfaced over the last years that may matter in your case.

> Kirkwood-Buff analysis where the convergence of the RDF will have a huge
> effect on my results.
>
> Also, I noticed that in CHARMM the RDF tails off once I reach the edge of
> my water box... for example, it will plateau at 1.06 up til 18 angstroms
> (my box is only 45 A X 38 A X 37 A) and then the RDF will dip down towards
> zero the rest of the way.... meanwhile, for the same RDF in VMD it remains
> plateaued at 0.97 all the way up to 20 angstroms.... I'm assuming this
> just has something to do with the way each program deals with the periodic
> cell.

yes, the code in VMD allows you to go "into the edges", i.e. it correctly
normalizes, if your sample is not a complete spherical shell anymore,
but you have to factor in the reduced statistical sample size, _and_
that for calculations with ewald summation (or equivalent), the potential
due to the reciprocal space contribution is supposed to have the largest
errors the further you get away from an individual charge.

> I'm also wondering if the sloppiness you mentioned with the normalization
> constant may be showing here since my systems are rather small (~
> 4000-5000 atoms).

yes. sloppy normalization will show the more, the smaller the system is.
it will also show the more the larger your histogram bins are. however,
particularly the latter is a balancing act between accurate representation
of the g(r) and good statistics.

> Thanks for your help so far,

no problem. please keep me posted, if you notice anything odd.
getting the g(r) code in VMD as close to perfection as possible
has become a matter of pride for me. since it is something that
i have been working and improving on (in increments) for over
15 years now (since i was an undergrad and programmed my
first (rather sloppy, but less so than that my colleague!) g(r) code. ;-)

cheers,
     axel.

> Tim
> 

-- 
Dr. Axel Kohlmeyer
akohlmey_at_gmail.com  http://goo.gl/1wk0
Institute for Computational Molecular Science
Temple University, Philadelphia PA, USA.