[gmx-developers] details about Heat capacity caculation

[Haiqing Zhao]

Hi, David, 

thanks for reply.

I checked it. Version 4.5.5 really makes a difference from 4.5.3 which I used.
In 4.5.5, Cv works exactly well by  nmol*(RMSD(E))^2/RT^2. (Even though I dont know here why need to multiply the # of molecular.)

But still I'm sure the Cp is not from simple Enthalpy Fluctuation,  as you see from the example below--167296 J/mol, and there is several orders mistake.

Please help me check it.
Thanks.

e.g.

Enthalpy                   -17216.8        4.2    167.296    22.4922  (kJ/mol)

You may want to use the -driftcorr flag in order to correct
for spurious drift in the graphs. Note that this is not
a substitute for proper equilibration and sampling!

Temperature dependent fluctuation properties at T = 300.076.

Heat capacities obtained from fluctuations do *not* include
quantum corrections. If you want to get a more accurate estimate
please use the g_dos program.

WARNING: Please verify that your simulations are converged and perform
a block-averaging error analysis (not implemented in g_energy yet)
Volume                                   = 1.8548e-05 m^3/mol
Enthalpy                                 =   -37.0255 kJ/mol
Coefficient of Thermal Expansion Alpha_P =  0.0129106 (1/K)
Isothermal Compressibility Kappa         = 3.00555e-10 (J/m^3)
Adiabatic bulk modulus                   = 3.32718e+09 (m^3/J)
Heat capacity at constant pressure Cp    =    80.3937 J/mol K



On Nov 1, 2011, at 2:05 PM, David van der Spoel wrote:

> On 2011-11-01 16:54, Haiqing Zhao wrote:
>> 
>> Dear Gmxers,
>> 
>> I have one question about how Gromacs calculate the heat capacity. I dont think it is only Enthalpy fluctuation÷(kT^2) in NPT ensemble. And also I checked it in NVT, it is not simply Total energy(/internal energy) fluctuation÷(kT^2).
>> 
> Which gmx version? It's more or less correct in 4.5.4 but you really need qm correctiions. Check g_dos for that
> 
> 
> 
> 
>> One data for NPT is: ( -nmol   465  -nconstr    3 )
>> 
>> Energy                      Average   Err.Est.       RMSD  Tot-Drift
>> -------------------------------------------------------------------------------
>> Potential                  -44.6187      0.007   0.293318  0.0281861  (kJ/mol)
>> Kinetic En.                 7.59226     0.0032    0.19886  0.0201858  (kJ/mol)
>> Total Energy               -37.0264     0.0091   0.359776  0.0483719  (kJ/mol)
>> Temperature                 300.076       0.13    7.85971   0.797822  (K)
>> Pressure                   0.663728       0.16    447.392  -0.388627  (bar)
>> Volume                      14.3219     0.0093   0.133543 -0.0667073  (nm^3)
>> pV                         0.438734    9.5e-05 0.00136276 -0.000679764  (kJ/mol)
>> Enthalpy                   -17216.8        4.2    167.296    22.4922  (kJ/mol)
>> 
>> Temperature dependent fluctuation properties at T = 300.076. #constr/mol = 3
>> Isothermal Compressibility: 3.00555e-05 /bar
>> Adiabatic bulk modulus:        33271.8  bar
>> Heat capacity at constant pressure Cp:    67.9219 J/mol K
>> Thermal expansion coefficient alphaP: 9.66233e-05 1/K
>> 
>> And data for NVT, ( -nmol   7631  -nconstr    3 )
>> 
>> Energy                      Average   Err.Est.       RMSD  Tot-Drift
>> -------------------------------------------------------------------------------
>> Total Energy               -33.2264     0.0003  0.0900307 -0.000219237  (kJ/mol)
>> Temperature                 300.008     0.0063    1.95324 -0.0148441  (K)
>> 
>> Temperature dependent fluctuation properties at T = 300.008. #constr/mol = 3
>> Heat capacity at constant volume Cv:    70.1813 J/mol K
>> 
>> (BTW, I found that (1)if you choose Enthalpy as output, g_energy will always give Cp, no matter what ensemble it is. (2)And only you choose tot energy and temp as output, no other terms, g_energy will give you Cv. If not, you cannot get Cv by g_energy.)
>> 
>> I appreciate if anybody could help me to check it.
>> Thanks!!
>> 
>> ---------------------------------------------
>> Haiqing Zhao
>> PH.D. candidate
>> in Computational Biophysics
>> Michigan Technological University
>> http://www.phy.mtu.edu/~haizhao
>> 
> 
> 
> -- 
> David van der Spoel, Ph.D., Professor of Biology
> Dept. of Cell & Molec. Biol., Uppsala University.
> Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205.
> spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
> -- 
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-developers-request at gromacs.org.