Dear Grace,
> Does anyone know how the density output in a NPT simulation (constant
> pressure temperature) is calculated?
>
> I used the formula:
>
> (density in N/V)*10*weight/6.022 = density in g/mL where
> weight is the
> molecular weight (eg. weight_water=18)
>
> This does not seem to give me the correct density that is outputted!
In the sander code the density is stored in ener(42) and is set on line 555
of runmd.f to:
ener(42) = tmass/(0.602204*vol)
Where tmass is the total mass set in rdparm to the sum of the atomic masses
from the prmtop file. Which are in grams per mol - hence H2O = 1.008 + 1.008
+ 16.000 = 18.016 (note mass of H is not 1 if you have a lot of hydrogens
you could introduce quite a large error by assuming water is mass 18.)
The volume, I believe is simply the volume of the periodic box. So if the
periodic box volume is in A^3 = 1 x 10^-30 m^3 = 1x10^-27 L = 1x10^-24 cm^3
- hence the reason avogadro's number is specified in the denominator as
0.602204. Hence we have (g*mol^-1) / (mol-1 * A^3 * conversion from A^3 to
cm-3)
Hence the density is expressed in g/cm-3.
So try re-working you calculation into g/cm-3 (using the masses in the
prmtop files) and see if you get the same answer. Note, the volume of your
periodic box will change during the simulation so you may want to either
calculate it every ntwx steps (and extract the box coordinates from the
mdcrd file). Or add a print statement to the code that prints out the box
volume for you.
All the best
Ross
/\
\/
|\oss Walker
| Department of Molecular Biology TPC15 |
| The Scripps Research Institute |
| Tel:- +1 858 784 8889 | EMail:- ross.rosswalker.co.uk |
|
http://www.rosswalker.co.uk/ | PGP Key available on request |
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu
Received on Thu Sep 02 2004 - 01:15:10 PDT
