[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,
>>
>> 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.)
>>
> 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.
>>
>> 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.
>>
>
>
> --
> 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.
More information about the gromacs.org_gmx-developers
mailing list