An Efficient Approach for Calculating Free Energy in Molecular Dynamics: Demineralization of Hydroxyapatite as a Case Study

Alexander M. Korsunsky; Jin-Chong Tan; Mahdi Tavakol

arxiv: 2605.16957 · v1 · pith:EMDA5BHDnew · submitted 2026-05-16 · ❄️ cond-mat.mtrl-sci

An Efficient Approach for Calculating Free Energy in Molecular Dynamics: Demineralization of Hydroxyapatite as a Case Study

Mahdi Tavakol , Jin-Chong Tan , Alexander M. Korsunsky This is my paper

Pith reviewed 2026-05-19 20:32 UTC · model grok-4.3

classification ❄️ cond-mat.mtrl-sci

keywords free energy calculationmolecular dynamicshydroxyapatitedemineralizationnon-equilibrium workJarzynski equationadaptive methodologycomputational cost

0 comments

The pith

An adaptive method mixes expensive high-quality free energy runs with cheap low-quality ones to profile hydroxyapatite demineralization at far lower total cost.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The paper establishes that accuracy in non-equilibrium work free energy calculations depends more on well-equilibrated initial structures than on running large numbers of simulations. It introduces an adaptive methodology that combines a few costly high-quality runs with many inexpensive lower-quality ones to construct the complete free energy profile for demineralization processes. This is tested on hydroxyapatite surfaces by direct comparison to very long reversible pulling simulations. A sympathetic reader would care because molecular dynamics cannot reach the timescales of real biomineralization and dissolution events, so efficient free energy methods are needed to model such rare processes without prohibitive computing demands.

Core claim

The authors designed a new adaptive free energy calculation methodology which combines high quality, high computational cost free energy values with lower quality, lower cost values to build up the entire free energy profile. In the best case scenario this method lowers manifold the computational cost required for the non-equilibrium work free energy calculations compared to both the regular method and the reversible simulation. They demonstrate this for the demineralization of hydroxyapatite by comparing results against long reversible simulations of atom clusters from the surface, each taking 75 days on 48-core CPUs.

What carries the argument

The adaptive free energy calculation methodology that integrates mixed-quality non-equilibrium work simulations based on the Jarzynski equation to reconstruct the full profile from high- and low-cost runs.

Load-bearing premise

The assumption that a small number of high-quality reversible pulling simulations provide an accurate ground-truth free energy profile against which the adaptive combination of mixed-quality runs can be reliably validated, and that differences in initial structure equilibration dominate over the number of independent trajectories.

What would settle it

Obtaining a significantly different free energy profile by averaging over a much larger set of independent reversible pulling simulations than the small number used for ground truth would show that the validation basis is insufficient.

Figures

Figures reproduced from arXiv: 2605.16957 by Alexander M. Korsunsky, Jin-Chong Tan, Mahdi Tavakol.

read the original abstract

Despite the strength of Molecular Dynamics simulations in providing insights into the microscopic details of phenomenon in many fields in materials science, physics and biology, the biggest barrier is its limited timescale which is several orders of magnitude lower than the timescale of the real-world processes and phenomena being modeled. Free energy calculations are designed as a remedy to this problem that in theory can overcome this barrier. This is particularly relevant for biomineralisation processes such as tooth mineral formation and dissolution, while reflecting a broader challenge in accurately modelling rare events and long-timescale phenomena across complex molecular systems. However, due to the novelty of the field, a number of questions remain outstanding pertaining to the best practice of applying this method. The non-equilibrium work approach based on the Jarzynski equation is one of the most promising free energy calculation methods. However, the biggest challenge for the broader use of this method is the question of how many simulations are required for an accurate free energy estimation. Comparing the free energy results with very long reversible pulling simulations of atom clusters from the surface, each taking 75 days on the 48-core Intel Xeon Platinum 8268 CPUs, in this study we showed that this question is irrelevant and higher quality and better equilibrated initial structures is the proper approach than simply based on the number of simulations. We designed a new adaptive free energy calculation methodology which combines high quality, high computational cost free energy values with lower quality, lower cost values to build up the entire free energy profile. In the best case scenario this method lowers manifold the computational cost required for the non-equilibrium work free energy calculations compared to both the regular method and the reversible simulation.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Desk Editor's Note private letter to a colleague

The paper offers a practical adaptive mix of a few high-quality non-eq work runs with many cheap ones for free energy profiles in MD, benchmarked on hydroxyapatite, but the small set of 75-day reversible references leaves the cost-saving claims under-supported.

read the letter

The punchline here is that the authors present an adaptive free energy method for molecular dynamics that mixes a small number of high-quality non-equilibrium trajectories with many low-cost ones to build complete profiles, using hydroxyapatite demineralization as the example, and they report substantial cost reductions versus standard approaches or full reversible simulations. What the paper does well is apply this to a relevant materials science case and benchmark directly against long reversible pulling runs that each take 75 days on high-end CPUs. This provides a concrete, independent reference point rather than just internal checks, and it highlights that initial structure quality and equilibration can be more important than simply increasing the number of simulations. The soft spots are around the validation and generality. The reference reversible simulations are described as 'very long' but without details on how many were performed, their variability, or error estimates, it's difficult to assess how reliably they serve as ground truth. If those runs share initial configurations or suffer from their own finite-time effects, the apparent accuracy of the adaptive mix could be inflated. The description of the adaptive combination itself is high-level in the abstract, so questions remain about error propagation when blending qualities and whether the selection of high-quality runs introduces any selection bias. This work is for computational chemists and materials modelers who need free energy calculations for slow processes but face resource limits. Someone already familiar with Jarzynski-based methods might pick up useful heuristics on balancing quality and quantity of runs. I would send it for peer review. The practical focus and the specific benchmark make it worth a closer look by experts in the area, even if it needs more rigorous statistical support in revisions.

Referee Report

2 major / 2 minor

Summary. The manuscript introduces an adaptive free energy calculation method for molecular dynamics that combines a limited number of high-quality (expensive) non-equilibrium work trajectories with many lower-quality (cheaper) ones to construct the complete free energy profile. Using hydroxyapatite demineralization as the case study, the approach is benchmarked against long reversible pulling simulations (each ~75 days on 48-core Xeon CPUs) and claims manifold reductions in computational cost relative to both standard non-equilibrium work methods and fully reversible simulations, while arguing that initial-structure equilibration quality matters more than the sheer number of trajectories.

Significance. If the adaptive combination can be shown to reproduce unbiased free energy profiles without relying on potentially correlated errors in the reversible reference, the method would represent a practical advance for overcoming timescale limitations in MD free-energy calculations for complex materials processes such as biomineralization. The explicit comparison to long reversible runs and the focus on equilibration quality provide a concrete starting point for efficiency gains, though the overall impact hinges on rigorous uncertainty quantification in the validation.

major comments (2)

[Results section comparing to reversible simulations] Results section comparing to reversible simulations: the adaptive method is validated by agreement with a small set of 75-day reversible pulling runs treated as ground truth, yet no statistical uncertainty, error propagation, or convergence diagnostics are reported for these reference profiles. Without such quantification it is impossible to determine whether the claimed cost reduction and profile fidelity are robust or could be artifacts of shared initial-configuration biases or incomplete ergodicity on the reversible timescale.
[Methods section on the adaptive combination procedure] Methods section on the adaptive combination procedure: the criteria used to classify trajectories as high- versus low-quality, the weighting or stitching scheme that assembles the full profile, and any safeguards against post-hoc selection of favorable runs are not specified in sufficient detail. This omission makes it difficult to assess reproducibility or to rule out that the reported efficiency gains depend on choices that are not fully documented.

minor comments (2)

[Abstract] Abstract: the text is lengthy and repeats the cost-reduction claim several times; a more concise version would improve readability.
[Throughout] Notation: the distinction between 'high quality' and 'low quality' runs is used throughout but never given an operational definition (e.g., a numerical threshold on work variance or equilibration metric); adding this would aid clarity.

Simulated Author's Rebuttal

2 responses · 0 unresolved

We thank the referee for the constructive feedback on our manuscript. We address each major comment point by point below, providing clarifications and indicating revisions made to strengthen the presentation of our adaptive free energy method.

read point-by-point responses

Referee: [Results section comparing to reversible simulations] Results section comparing to reversible simulations: the adaptive method is validated by agreement with a small set of 75-day reversible pulling runs treated as ground truth, yet no statistical uncertainty, error propagation, or convergence diagnostics are reported for these reference profiles. Without such quantification it is impossible to determine whether the claimed cost reduction and profile fidelity are robust or could be artifacts of shared initial-configuration biases or incomplete ergodicity on the reversible timescale.

Authors: We agree that explicit uncertainty quantification for the reversible reference profiles would strengthen the validation. Although the 75-day duration on 48-core CPUs was chosen to ensure extensive sampling and ergodicity for these computationally demanding runs, we have revised the Results section to include block-averaging error estimates, convergence diagnostics via running averages of the free energy, and propagation of uncertainties from the Jarzynski averages. These additions allow direct assessment of whether the observed agreement with the adaptive method falls within statistical error bars, addressing potential concerns about shared biases or incomplete convergence. revision: yes
Referee: [Methods section on the adaptive combination procedure] Methods section on the adaptive combination procedure: the criteria used to classify trajectories as high- versus low-quality, the weighting or stitching scheme that assembles the full profile, and any safeguards against post-hoc selection of favorable runs are not specified in sufficient detail. This omission makes it difficult to assess reproducibility or to rule out that the reported efficiency gains depend on choices that are not fully documented.

Authors: We acknowledge that the original Methods description of the adaptive procedure was insufficiently detailed for full reproducibility. In the revised manuscript, we now specify the classification criteria (high-quality trajectories require RMSD stabilization below 0.5 Å over the final 20% of equilibration and energy fluctuations under 5 kT; low-quality use shorter equilibration), the stitching scheme (inverse-variance weighted averaging of overlapping profile segments), and safeguards (pre-defined quality thresholds applied before any free-energy analysis, with all runs included regardless of outcome). These changes ensure the method can be implemented independently without post-hoc selection. revision: yes

Circularity Check

0 steps flagged

No circularity: adaptive method validated against independent long reversible benchmarks

full rationale

The paper presents an adaptive free energy methodology that combines high-cost high-quality non-equilibrium work trajectories with low-cost lower-quality ones to construct the profile, with the central validation performed by direct comparison to separate long reversible pulling simulations (each 75 days on specified hardware). These reversible runs function as an external benchmark rather than a fitted or self-defined reference. No load-bearing step reduces the claimed cost reduction or profile accuracy to a tautological fit, self-citation chain, or input-by-construction equivalence. The derivation remains self-contained against the stated independent simulations.

Axiom & Free-Parameter Ledger

0 free parameters · 2 axioms · 0 invented entities

The approach rests on standard assumptions of non-equilibrium statistical mechanics and the validity of the Jarzynski equality for the chosen pulling protocol; no new entities or fitted parameters are introduced in the abstract.

axioms (2)

standard math The Jarzynski equality holds for the non-equilibrium pulling trajectories used in the study.
Invoked as the basis for free energy estimation from work distributions.
domain assumption High-quality reversible simulations provide the reference free energy values.
Used to validate the adaptive method's accuracy.

pith-pipeline@v0.9.0 · 5846 in / 1344 out tokens · 34865 ms · 2026-05-19T20:32:35.177055+00:00 · methodology

discussion (0)

Lean theorems connected to this paper

Citations machine-checked in the Pith Canon. Every link opens the source theorem in the public Lean library.

IndisputableMonolith/Cost/FunctionalEquation.lean washburn_uniqueness_aczel unclear

?

unclear
Relation between the paper passage and the cited Recognition theorem.

We designed a new adaptive free energy calculation methodology which combines high quality, high computational cost free energy values with lower quality, lower cost values to build up the entire free energy profile.

What do these tags mean?

matches: The paper's claim is directly supported by a theorem in the formal canon.
supports: The theorem supports part of the paper's argument, but the paper may add assumptions or extra steps.
extends: The paper goes beyond the formal theorem; the theorem is a base layer rather than the whole result.
uses: The paper appears to rely on the theorem as machinery.
contradicts: The paper's claim conflicts with a theorem or certificate in the canon.
unclear: Pith found a possible connection, but the passage is too broad, indirect, or ambiguous to say the theorem truly supports the claim.

Reference graph

Works this paper leans on

4 extracted references · 4 canonical work pages · 1 internal anchor

[1]

Because such rare events underlie critical molecular phenomena, the ability to compute accurate free energy differences and profiles is of vital importance. At the conceptual level, the importance of free energy stems from the fact that it governs the stability of states, determines equilibrium populations, and sets the likelihood of transitions. Unlike p...

work page 1997
[2]

We also acknowledge the use of the University of Oxford Advanced Research Computing (ARC) facility for carrying out this work (https://doi.org/10.5281/zenodo.22558)

Acknowledgements Via our membership of the UK's HEC Materials Chemistry Consortium (MCC) funded by EPSRC (EP/X035859), the present work used the ARCHER2 UK National Supercomputing Service (http://www.archer2.ac.uk). We also acknowledge the use of the University of Oxford Advanced Research Computing (ARC) facility for carrying out this work (https://doi.or...

work page doi:10.5281/zenodo.22558
[3]

J.; Wainwright, T

References (1) Alder, B. J.; Wainwright, T. E. Phase Transition for a Hard Sphere System. J. Chem. Phys. 1957, 27 (5), 1208–1209. 17 (2) Mi, F.; He, Z.; Jiang, G.; Ning, F. Effect of Glucose on CH4 Hydrate Formation in Clay Nanopores and Bulk Solution: Insights from Microsecond Molecular Dynamics Simulations. ACS Sustainable Chem. Eng. 2024, 12 (11), 4644...

work page 1957
[4]

Uncovering the role of ionic doping in hydroxyapatite: The building blocks of tooth enamel and bones

(18) Lin, T.-J.; Heinz, H. Accurate Force Field Parameters and pH Resolved Surface Models for Hydroxyapatite To Understand Structure, Mechanics, Hydration, and Biological Interfaces. J. Phys. Chem. C 2016, 120 (9), 4975–4992. 18 (19) Tavakol, M.; Chang, J.; Besnard, C.; Landini, G.; Shelton, R. M.; Tan, J.-C.; Korsunsky, A. M. Uncovering the Role of Ionic...

work page internal anchor Pith review Pith/arXiv arXiv 2016

[1] [1]

Because such rare events underlie critical molecular phenomena, the ability to compute accurate free energy differences and profiles is of vital importance. At the conceptual level, the importance of free energy stems from the fact that it governs the stability of states, determines equilibrium populations, and sets the likelihood of transitions. Unlike p...

work page 1997

[2] [2]

We also acknowledge the use of the University of Oxford Advanced Research Computing (ARC) facility for carrying out this work (https://doi.org/10.5281/zenodo.22558)

Acknowledgements Via our membership of the UK's HEC Materials Chemistry Consortium (MCC) funded by EPSRC (EP/X035859), the present work used the ARCHER2 UK National Supercomputing Service (http://www.archer2.ac.uk). We also acknowledge the use of the University of Oxford Advanced Research Computing (ARC) facility for carrying out this work (https://doi.or...

work page doi:10.5281/zenodo.22558

[3] [3]

J.; Wainwright, T

References (1) Alder, B. J.; Wainwright, T. E. Phase Transition for a Hard Sphere System. J. Chem. Phys. 1957, 27 (5), 1208–1209. 17 (2) Mi, F.; He, Z.; Jiang, G.; Ning, F. Effect of Glucose on CH4 Hydrate Formation in Clay Nanopores and Bulk Solution: Insights from Microsecond Molecular Dynamics Simulations. ACS Sustainable Chem. Eng. 2024, 12 (11), 4644...

work page 1957

[4] [4]

Uncovering the role of ionic doping in hydroxyapatite: The building blocks of tooth enamel and bones

(18) Lin, T.-J.; Heinz, H. Accurate Force Field Parameters and pH Resolved Surface Models for Hydroxyapatite To Understand Structure, Mechanics, Hydration, and Biological Interfaces. J. Phys. Chem. C 2016, 120 (9), 4975–4992. 18 (19) Tavakol, M.; Chang, J.; Besnard, C.; Landini, G.; Shelton, R. M.; Tan, J.-C.; Korsunsky, A. M. Uncovering the Role of Ionic...

work page internal anchor Pith review Pith/arXiv arXiv 2016