# [gmx-users] Re: Wrong heat capacity (water)

Andreas Kring akring at fys.ku.dk
Thu Apr 17 12:20:11 CEST 2008

```Continuation of thread
http://www.gromacs.org/index.php?option=com_wrapper&Itemid=84
(answered by David van der Spoel)

>> Hello all,
>>
>> I ran a MD simulation of TIP4P water in the NVT ensemble (976
molecules,
>> density=0.98 g/cm3, T=300K (Nose-Hoover), SETTLE and PME enabled). This
>> works out fine, i.e. I get a correct radial distribution function
>> (g_OO), but the heat capacity (C_V) is wrong. I calculated it using the
>> command:
>>
>> \$  g_energy -f ener.edr -ndf 6 -nmol 976
>>
>> and get the value below (which sould have been approximately 75
J/(mol*K)).
>>
>> Is the command above incorrect?'

>Not necessarily. It would be good if you check the following things:
>
>- is the simulation in equilibrium?
>- how large are the fluctuations for the equilibrium part?
>- is the density correct?
>
>You can compute the cV from the fluctuations manually as well. Does that
>match the value computed by the program?
>

Thanks for you answer. Yes, you were right. The simulation was not in
equilibrium. I've run a new simulation which I have made sure is in
equilibrium. I have also checked that the density is correct, as you
suggest above, using the g_density command.

The numbers for the new simulation can be seen below (same simulation
conditions as stated in the top of this mail), and the answer is still
not right? With the command

\$ g_energy -f ener.edr -ndf 6 -nmol 976

I get Cv = 305.413 J/(mol*K), while a manual calculation gives

Cv = 1/(k_B*T*T)*Var(E) = 85.6 J/(mol*K), where Var(E) was found using
the RMSD numbers remembering that the unit of these are kJ/mole. The
table value is Cv equal to approximately 75 J/(mole*K).

Any ideas?

Kind regards
Andreas

Energy                      Average       RMSD     Fluct.      Drift
Tot-Drift
-------------------------------------------------------------------------------
Potential                  -41.3012   0.219318          0    8.76027
87.6114
Kinetic En.                 7.48311   0.135439          0    4.46649
44.6693
Total Energy               -33.8181   0.256267          0    13.2266
132.28
Temperature                 300.002    5.42982    5.40392   0.183467
1.83485
Pressure (bar)             -269.723    384.283    384.133    3.71038
37.1075
Heat Capacity Cv:      305.413 J/mol K (factor = 0.639443)

```