# [gmx-developers] details about Heat capacity caculation

Haiqing Zhao haizhao at mtu.edu
Tue Nov 1 21:00:11 CET 2011

```Hi,

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)
pV                              0.438734    9.5e-05 0.00136276 -0.000679764  (kJ/mol)
Enthalpy                   -17216.8        4.2    167.296    22.4922  (kJ/mol)

Yes, I did not doubt Cp=Var(Enthalpy)/k T^2.
Actually my question is only about Enthalpy. As we see in the example, the total energy= -37.0264KJ/mol, pv = 0.4387 Kj/mol, How can enthalpy=  -17216.8 KJ/mol?

sorry for the time.

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

> On 2011-11-01 20:15, Haiqing Zhao wrote:
>>
>> Hi, David,
>>
>>
>> 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.)
>>
> This is supposed to work for liquids, hence you need to normalize.
> How many molecules do you have, 465?
> Heat capacity of 80 resembles water. I'm pretty sure it works correctly, since we just submitted a paper about it and I have looked into this in detail. We used:
>
> kB T^2 c_P = Var (Enthalpy)
>
> which comes straight from Allen and Tildesley.
>
>> 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.
>>
>> 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.
>>
>> 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.
>>
>
>
> --
> 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.

```