Statistical Mechanics is the basis of molecular level calculations of properties and averages. Nowadays free energy calculations are fairly "turn-the-crank", which is not necessarily a good thing. It is easy to perform a free energy calculation and get a number.

In the molecular dynamics package GROMACS, for instance, alchemical free energy calculations are only valid if done in the NVT ensemble. This is because GROMACS samples potential energies, thus one can calculate (or press a button on pymbar or another piece of software)

begin{equation}

Delta A = -kT ln langle exp(-beta Delta mathcal{U}) rangle

end{equation}

Where $Delta A$ is the difference in the Helmholtz free energy between two states, $beta = frac{1}{kT}$, $mathcal{U}$ is the difference in sampled potential energy between two states, and $langle rangle$ is an ensemble average.

If you run an NPT simulation you must instead sample

begin{equation}

Delta G = -kT ln frac{langle V exp(-beta Delta mathcal{U}) rangle}{langle V rangle}

end{equation}

While it is tempting to say that the volumes cancel so you only need the ensemble average of the NVT Boltzmann factor, this is not at all correct. GROMACS is aware that it is not correctly sampling the NPT and has tucked away a note in the user manual to the effect that for water at 298.15 K the error is small i.e., less than a kilojoule.

- a kilojoule can make a large difference to an equilibrium constant
- Chemistry does not happen at only 298.15! if you increase the temperature or density of the solvent (not all of us do water) the error quickly gets out of hand – this is particularly true for ionic liquids.

Is this assumption (possibly error in other codes) present in other codes? Are there any MD codes that properly sample the NPT ensemble for free energy?

While Godzilla has shown entirely correctly how one *should* sample the free energy in the NPT, neither his equation, or mine, which is equivalent for $Delta G$, are implemented in GROMACS, and possibly others. I am looking for if any MD programs do implement the proper free energy sampling. They all properly do NPT sampling (I assume), but sampling the free energy is a step above and beyond.

Matter Modeling Asked by B. Kelly on August 19, 2021

1 AnswersIt is straightforward to show that in a typical $NPT$ setting the Zwanzig equation still only depends on the energy difference and not on the volume (here I define $H$ to be the Hamiltonian of each system, respectively and $x$ to represent all variables over phase space):

$$frac{Z_{B,NPT}}{Z_{A,NPT}} = int int e^{-beta Delta H(x)} frac{e^{-beta H_A(x)}}{Z_{A,NPT}} e^{-beta PV} dxdV = left<e^{-beta Delta H}right>_{A,NPT}$$

so I am not sure where you get this equation from $-$ it seems wrong to me. Of course, you still get $NPT$ contribution in this way, since the exponential averages have to be performed in the isothermal-isobaric ensemble, and these will generally be different to averages obtained in $NVT$ (although in many practical purposes, virtually the same, especially in dense systems that are not near any critical points).

I am not sure which part of the GROMACS manual you are referring to, but I suspect you are talking about the Berendsen barostat. The Berendsen barostat is known to be not rigorous but very good for equilibration (again, for many systems it is practically the same as a rigorous barostat and usually not the limiting error). For a production run you should use the Parrinello-Rahman barostat, which is considered to be rigorous. However, it is deterministic, which can lead to some practical problems regarding non-physical oscillations. This is the "rigorous" barostat implemented in GROMACS and to my knowledge this way of calculating free energies in $NPT$ is completely valid from a theoretical standpoint.

Another rigorous barostat which is also stochastic is the Monte Carlo barostat. Unfortunately, it is only implemented in AMBER and OpenMM and not in GROMACS. As far as I know, NAMD comes with a Nosé-Hoover Langevin piston barostat which is also rigorous. You can see that pretty much every MD engine comes with a Berendsen barostat and a single rigorous barostat and they all use different rigorous barostats.

As a final remark, if you are planning to run simulations at high temperatures, I'd be *much* more worried about the validity of your model, presumably a force field, at these extreme conditions. Also, while it is true that 1 kJ/mol can be considered significant in some cases, there are so many errors when calculating free energies, including force field accuracy, insufficient sampling, estimator convergence (e.g. Zwanzig), or even the fact that you restrict your system to be a tiny box, that I would be much less worried about the impact of the barostat, or indeed the choice of thermodynamic ensemble on my results in most conceivable applications.

Answered by Godzilla on August 19, 2021

1 Asked on January 5, 2022 by quantumx

computational chemistry coupled cluster gaussian quantum chemistry

1 Asked on January 5, 2022

ab initio calculations density functional theory functional many body perturbation theory quantum chemistry

1 Asked on January 3, 2022

density functional theory mathematical modeling quantum chemistry spectroscopy td dft

3 Asked on January 1, 2022 by anyon

1 Asked on January 1, 2022 by bnd

1 Asked on December 22, 2021

lammps molecular dynamics molecular interactions molecular modeling software

2 Asked on December 22, 2021

1 Asked on December 20, 2021

electronic structure fcidump integrals one topic per answer software

1 Asked on December 20, 2021 by c-alexander

density functional theory from phys se quantum chemistry spectroscopy td dft

1 Asked on December 20, 2021

1 Asked on December 18, 2021 by kfc

1 Asked on December 16, 2021

density functional theory from phys se mathematical modeling statistical thermodynamics variational principle

1 Asked on December 16, 2021 by luphys

alloys density functional theory from phys se relativistic software

1 Asked on December 16, 2021

1 Asked on December 13, 2021 by roman-korol

atomic physics full configuration interaction high precision molecular modeling periodic system

1 Asked on December 13, 2021

band structure density functional theory from phys se quantum espresso

9 Asked on December 11, 2021

high performance computing matter modeling specific technology one topic per answer open source software

4 Asked on December 9, 2021 by georgios-karaiskakis

computational chemistry density functional theory functional history

2 Asked on December 9, 2021 by koroma

band structure lattice dynamics perovskites phonons symmetry

Get help from others!

Recent Answers

- eric_kernfeld on How to test consistency of responses?
- kjetil b halvorsen on How to test consistency of responses?
- Philipp on How do i draw a ray in unity
- Justin Markwell on Unity app crashes when using unmodified custom Android manifest (didn’t find class “UnityPlayerActivity”)
- DMGregory on MouseLook Script “Pops” back to the last value when the script is enabled after being disabled or destroyed

Recent Questions

- MouseLook Script “Pops” back to the last value when the script is enabled after being disabled or destroyed
- Unity app crashes when using unmodified custom Android manifest (didn’t find class “UnityPlayerActivity”)
- How do i draw a ray in unity
- How to test consistency of responses?
- How can I understand these variograms?

© 2022 AnswerBun.com. All rights reserved.