Lum Nforbi wrote:
Hi Andrew Paluch,
You had asked me to supply the formulas that I attempted using to
calculate Cv. The formula from J. M. Haile is labelled isometric heat
capacity and given by Cv = Nk/(N-NT*(3(N/2)-1)<Ek*^-1>) and the one in
Allen and Tildesley is quite complicated to write out (has some special
characters in it).
I also found in one of the volumes of gmx-users where someone had
used the formula
Cv = 1/(k_B*T*T)*Var(E) and instead of having the expected 75 J/Kmol, he
got 86 J/Kmol. He also mentioned that Var(E) can be obtained from RMSD
numbers. I know that RMSD numbers are displayed from g_energy, but how
do you get Var(E) from these numbers?
The RMSD is the square root of the variance.
The RMSD depends strongly on the temperature coupling. Berendsen
coupling damps the fluctuations and hence gives too low Cp/Cv.
Thanks,
Lum
Date: Wed, 9 Dec 2009 20:45:57 -0500
From: Andrew Paluch <[email protected] <mailto:[email protected]>>
Subject: Re: [gmx-users] Problems with calculating Cv and Cp
To: Discussion list for GROMACS users <[email protected]
<mailto:[email protected]>>
Message-ID: <[email protected]
<mailto:[email protected]>>
Content-Type: text/plain; charset="us-ascii"
Lum,
You'll have better luck if you perform a few simulations at different
temperatures at the same pressure and mole number, and then
numerically differentiate the resulting enthalpy. I believe that the
default heat capacity that Gromacs will print out is only valid for
NVE simulations, and I am not sure of the formulas that you are
referring to from Allen and Tildesley & J.M. Haile. Also, don't
expect your calculated heat capacity with TIP3P to agree the
experimental value. You should under-predict the experimental value.
Lastly, depending on what you are interested in, you may have better
luck with a water model other than TIP3P.
Hope this helps,
Andrew
_______________________________________________
_______________________________________________
Andrew Paluch
Department of Chemical and Biomolecular Engineering
University of Notre Dame du Lac
[email protected] <mailto:[email protected]>
_______________________________________________
_______________________________________________
On Dec 9, 2009, at 8:18 PM, Lum Nforbi wrote:
> Dear all,
>
> I have run an 8 ns NPT simulation of 2000 molecules of TIP3P
> water and I have a very low Cv value of 12.4748 J/mol K (factor =
> 0.000164481). The result is below. I ignored this value and have
> tried using formulas for Cv that I found in the two books: Allen
> and Tildesley & J. M. Haile but I can't come out with the right
> answer.
> Has anyone ever calculated Cv or Cp for water manually from
> scratch and gotten the right answer? If so, please, could you give
> me the details of what you did?
>
> Statistics over 4000001 steps [ 0.0000 thru 8000.0005 ps ], 10 data
> sets
> All averages are exact over 4000001 steps
>
> Energy Average RMSD Fluct.
> Drift Tot-Drift
> ------------------------------
> -------------------------------------------------
> Potential -79498.8 285.331 285.331
> 0.000110205 0.881641
> Kinetic En. 14960.5 191.868 191.854
> -0.0010135 -8.10799
> Total Energy -64538.3 353.238 353.232
> -0.000903288 -7.22631
> Temperature 299.962 3.84702 3.84673
> -2.0321e-05 -0.162568
> Pressure (bar) 0.945534 194.613 194.613
> 0.000164281 1.31425
> Box-X 3.94528 0.00514348 0.00514348
> 0 0
> Box-Y 3.94528 0.00514348 0.00514348
> 0 0
> Box-Z 3.94528 0.00514348 0.00514348
> 0 0
> Volume 61.4097 0.240325 0.24032
> 6.42593e-07 0.00514075
> Density (SI) 974.3 3.80629 3.80622
> -1.04747e-05 -0.0837974
> Heat Capacity Cv: 12.4748 J/mol K (factor = 0.000164481)
> Isothermal Compressibility: 2.27095e-05 /bar
> Adiabatic bulk modulus: 44034.4 bar
>
>
> Thank you,
>
> Lum
>
> --
--
David van der Spoel, Ph.D., Professor of Biology
Molec. Biophys. group, Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone: +46184714205. Fax: +4618511755.
[email protected] [email protected] http://folding.bmc.uu.se
--
gmx-users mailing list [email protected]
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
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
