Saturday, March 16, 2013

Calculation of Molecular Entropies Using Temperature Integration

Kim Sharp Journal of Chemical Theory and Computations 2013, 9, 1164
Contributed by +Jan Jensen

Extracting accurate entropy changes from molecular dynamics (MD) simulations is notoriously difficult.  The entropy is usually computed using the Boltzmann-Gibbs-Planck expression$$S=-k \int p(\mathbf{q}) \ln p(\mathbf{q}) d\mathbf{q} $$where $\mathbf{q}$ is the postion vector and $p(\mathbf{q})$ the corresponding probability distribution function, which converges slowly for large systems.  In this paper Kim Sharp makes a convincing argument for using the Clausius-van't Hoff (CvH) expression$$S_{CvH}=\frac{\Delta E}{T_{CvH}}=\frac{\Delta E \left( \ln T_2 - \ln T_1\right)}{T_2-T_1}$$There are two issues which must be considered when using this equation.  One issue is that it assumes the heat capacity is constant between $T_1$ and $T_2$ placing a limit on how large the difference in $T$ can be.  The other issue is that the entropy does not converge to 0 as $T$ approaches 0 for simulations using classical energy functions such as the harmonic oscillator.

To address the latter problem Sharp argues that at sufficiently low $T$ classical systems will behave harmonically so that the entropy can be estimated by either the quasiharmonic method or using normal mode frequencies.  

The approach is tested in two ways. CvH conformational entropy changes on going from 280 K to 320 K for five ligands are shown to be close to corresponding values computed using probability density functions, while converging an order of magnitude faster. Furthermore, water vaporization  entropies at 300 K computed using the CvH approach were within 1 cal/mol/K of experiment and conventional Monte Carlo simulations.  This result can be obtained by cooling in as little as three steps to about 50 K and applying the quasiharmonic method for lower temperatures.

In addition to a more rapid convergence compared to existing methods the method also obviates any need to define suitable coordinates for the probability density functions, which can present a problem for some applications such as QM/MM based MD.  Further, MD simulations are often started by heating from an energy minimized structure and this procedure can be conveniently used to extract the entropy of the system.

Thanks to +Casper Steinmann for alerting me to this paper

 Creative Commons License
This work is licensed under a Creative Commons Attribution 3.0 Unported License