AMBER Archive (2007)Subject: AMBER: Still puzzled about this virial computation
From: David Cerutti (dcerutti_at_mccammon.ucsd.edu) 
Date: Thu Feb 01 2007 - 02:00:10 CST
 
 
 
 
Hello,
 
     Thanks to Rob and Ross for their helpful inputs.  I'm still puzzled a 
 
bit about some of the computations, though.
 
 I am using the following in my run file, but it appears to be ineffective 
 
at removing the various long-range corrections (for which I do not yet 
 
have an equivalent in my code) from SANDER:
 
 &ewald
 
   use_pme=0,
 
   eedmeth=4,
 
   vdwmeth=1,
 
&end
 
 1.) I can't seem to get my (real-space, no 1:4 interactions, no long-range 
 
approximations, strict cutoff, periodic) van-der Waals energies to agree 
 
with SANDER even when I use vdwmeth=0.  This seems to be the case no 
 
matter how many particles I use--2, 216, 1728.  It's very bad with only 
 
two particles--a difference of around 30%.  I get better and better 
 
agreement as I use a larger and larger cutoff with a bigger and bigger 
 
system--presumably because I'm minimizing the effect of any long-range 
 
correction, eventually obtaining that 1% discrepancy in the case of a 
 
1728-particle system with an 18.0A cutoff.
 
 2.) If I work hard to get the periodicity out of the system (by leaving 
 
the configuration as is and making the box much larger so I can use a 
 
200.0A cutoff, which is essentially doing an all:all calculation of a 
 
non-periodic system), I get very good agreement (only 0.03% discrepancy, 
 
which can be attributed to slight imprecision in my energy computer).
 
 3.) In fact... I find that setting vdwmeth=0 has absolutely no effect on 
 
the computed energy, even with a small 9.0A cutoff--it's the same as if I 
 
do not specify a vdwmeth parameter, or explicitly set it to 1.  Could this 
 
be a bug?
 
 4.) It seems that I get a consistent relationship between my computed 
 
virial and the SANDER value, particularly in the case of that 
 
"non-periodic" system--my computed virial is two times the SANDER value to 
 
six significant digits.  The system I'm working with is a collection of 
 
monatomic particles, so I know there's no constraints contributing to the 
 
virial.  I'm not sure how I can get energies that essentially match and 
 
yet have a virial that has some factor gone awry--which is why I really 
 
need someone to explain the formula that AMBER uses when it computes the 
 
virial.
 
 5.) Like my ability to get the van-der Waals terms to agree under certain 
 
conditions, I also get good agreement with the electrostatics under some 
 
conditions--two water molecules in a big box or 1728 waters in a huge 
 
box--essentially non-periodic systems.  If I try to run a system where 
 
periodicity would matter (one with a net charge or a neutral system 
 
without a vast gulf separating the various pieces), I run into trouble.
 
 6.) I don't know about constraints contributing to the virial, so to test 
 
whether my code agrees with AMBER when electrostatics get into the 
 
picture, I used a gas of "sodium" particles (I set the charges to 0.1e and 
 
turned off the Lennard-Jones interactions to focus on the electrostatics 
 
this time).  When doing so, I finally get a hint as to what the formula 
 
for computing that virial is: the virial comes out as negative 0.5 times 
 
the energy; I would have expected it to be negative 1.0 times the energy, 
 
so that's some explanation for where the factor of two in my computation 
 
of the virial for a gas of Lennard-Jones particles came from.  This 
 
relationship between the SANDER-computed energy and virial remains if I 
 
change the cutoff (which, as I would expect in a system without Ewald 
 
terms, has a significant effect on the total energy computed as well as 
 
the virial--but they vary in that 2:1 proportion).
 
 7.) Why does the total energy seem to be greatly affected by some sort of 
 
long-range electrostatic correction, even when I've set eedmeth=4 and 
 
use_pme=0 (see point 5)?
 
 There seem to be conundrums around every corner.  Can anyone help 
 
straighten this out?
 
 Thanks!
 
 Dave
 
-----------------------------------------------------------------------
 
The AMBER Mail Reflector
 
To post, send mail to amber_at_scripps.edu
 
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu
 
 
  
 |