Source-linked AI summary
Learning consistent molecular mechanics force fields from first principles
Berkay Günes, Leif Seute, Jigyasa Nigam, Frauke Gräter
TL;DR
Classical force fields are efficient for large-scale simulations, but empirical nonbonded parameters limit their adaptability, even in hybrid machine-learning approaches. grappa-fullFF learns bonded and nonbonded parameters together from ab initio data with physics-inspired regularization. It reports accurate electric response, geometry optimization, and conformational sampling without externally assigned nonbonded parameters.
Problem
Classical force fields support large-scale simulations through computational efficiency, but empirically assigned nonbonded parameters remain a limitation in hybrid ML/MM approaches.
Source: §1 · pp. 1–2Method
grappa-fullFF jointly learns bonded and nonbonded force-field parameters from ab initio data, using electrostatic-potential and volume-based physical regularization.
Source: §1 · p. 2, §2.2 · p. 4Results
grappa-fullFF recovers electrostatic potentials near MBIS accuracy, retains grappa-bonded geometry-optimization performance, and reproduces alanine-dipeptide free-energy surfaces about twice as accurately as espaloma-0.3.2.
Source: §1 · p. 2Takeaways & Limitations
Abstract
from arXiv · showhide
Classical force fields (FFs) remain the workhorse for large-scale simulations even as machine-learned interatomic potentials (MLIPs) approach ab initio accuracy. They decompose total configuration energies into simple effective interactions whose parameters are traditionally assigned based on atom or bond types, enabling efficient simulations but also limiting their ability to adapt across configurations. Recent machine learning approaches have improved the accuracy and transferability of bonded parameters in these FFs by inferring them as functions of local atomic environments, but still rely on empirical nonbonded parameters for practical simulations. In this work, we introduce a unified approach, \texttt{grappa-fullFF}, which learns both bonded and nonbonded parameters \emph{consistently} and simultaneously from ab initio reference data. By incorporating physically inspired regularization via supervision of the electrostatic potential and an architecture that facilitates charge equilibration, our model recovers accurate electric response properties, achieves state-of-the-art accuracy on geometry optimization benchmarks, and reproduces the conformational sampling of both classical and existing machine-learned FFs, without relying on externally assigned nonbonded parameters.
1 Introduction
Classical MM remains efficient for large-scale simulations, but empirically assigned parameters limit accuracy and transferability, and hybrid ML/MM approaches still rely on empirical nonbonded parameters. grappa-fullFF learns bonded and nonbonded parameters jointly from ab initio data and reports accurate electric response, geometry optimization, and alanine-dipeptide sampling.
1 Introduction
- Classical MM remains widely used for large-scale, long-timescale simulations because it is exceptionally efficient, while MLIPs remain computationally prohibitive at that scale.Source: §1 · p. 1
- Hybrid ML/MM methods improve bonded parameters but still rely on empirical nonbonded parameters or separate learned and classical descriptions, limiting generalization across configurations.Source: §1 · p. 2
- grappa-fullFF recovers QM electrostatic potentials near MBIS accuracy, retains grappa-bonded geometry performance, and is about twice as accurate as espaloma-0.3.2 for alanine-dipeptide free-energy surfaces.The alanine-dipeptide comparison is relative to AMBER99SB-ILDN.Source: §1 · p. 2
2 Methodology
The methodology addresses ambiguity in decomposing total energies into force-field components by jointly predicting bonded and nonbonded parameters and regularizing the nonbonded terms with physical information. Charge fluxes conserve total charge, while volume ratios provide a physically motivated basis for LJ parameters.
2 Methodology
- Energy degeneracy means different bonded and nonbonded interactions can yield similar total energies and forces, leaving learned parameters potentially nonphysical and unreliable across environments.The method addresses this ambiguity with additional physical constraints.Source: §2 · p. 2
2.1 Model Architecture
- The model extends grappa-bonded with charge-flux and volume-ratio heads for nonbonded parameters, while leaving bonded parameterization unchanged.Source: §2.1 · p. 3
- Antisymmetric charge fluxes transfer charge between bonded atoms while conserving total charge at every refinement step.Later flux predictions also condition on the current charge state and its local change.Source: §2.1 · p. 3
- Per-atom volume ratios provide a physically motivated coupling of the predicted LJ σ and ϵ parameters.The model predicts positive volume ratios from latent atomic representations.Source: §2.1 · p. 3
2.2 Model training through differentiable MM and physics-informed regularization
- A differentiable MM functional trains parameters on QM energies and forces, while ESP and volume-ratio losses regularize only the charge-flux and volume heads.ESP matching constrains charges without depending on a chosen charge-partitioning scheme.Source: §2.2 · pp. 3–4
3 Results
The results evaluate grappa-fullFF on electrostatic properties, optimized geometries, and explicit-solvent conformational thermodynamics. It approaches MBIS electrostatics, retains the bonded-only model’s geometric accuracy, and reproduces its conformational sampling.
Accuracy of the Learned Electrostatics
- grappa-fullFF’s predicted charges nearly match MBIS accuracy on QM electrostatic properties, while outperforming AM1-BCC and espaloma-charge on ESP and dipole-moment errors.Charges are inferred through ESP supervision rather than trained against reference atomic charges.Source: §3 · p. 4
Structural Accuracy in Geometry Optimization
- In a standard geometry-optimization benchmark, starting from unrelaxed conformations, espaloma-0.3.2, grappa-bonded, and grappa-fullFF outperform GAFF-2.11 and OpenFF-2.1.0 across geometry RMSE, TFD, and relative-energy errors.Source: §3 · p. 5, Table 2 · p. 5
Conformational Thermodynamics in Explicit Solvent
- Both Grappa models achieve JS divergences ≲0.08 nats against AMBER99SB-ILDN for capped alanine dipeptide conformational sampling, with only ∆JS = 0.013 nats between them.The comparison uses 500 ns explicit-solvent MD trajectories and 5° bins for the two backbone dihedral angles.Source: §3 · p. 5, Figure 2 · p. 5
4 Discussion and Conclusions
grappa-fullFF learns bonded and nonbonded parameters from ab initio energies and forces, using physics-inspired regularization to address energy degeneracy. It retains geometric accuracy and reproduces explicit-solvent conformational thermodynamics, though training used isolated gas-phase monomers.
4 Discussion and Conclusions
- Physics-inspired constraints on charge-conserving bond fluxes through ESP and LJ parameters through TS-vdW volume ratios resolve ambiguity among interaction terms.This regularization addresses energy degeneracy in decomposing the force-field energy.Source: §4 · p. 5
- The learned charges and LJ parameters emerge without optimization against reference parameter values, and polar hydrogens receive reduced LJ σ values for hydrogen bonding.The resulting nonbonded parameters offer a cheaper alternative to electron-density-based MBIS calculations.Source: §4 · p. 6
- grappa-fullFF matches grappa-bonded’s geometry-optimization performance and reproduces explicit-solvent conformational thermodynamics despite training exclusively on isolated gas-phase monomers.The authors identify explicit-solvent training and further refinement for larger, more delocalized molecules as next steps.Source: §4 · p. 6
A.1 Qualitative Analysis of Electrostatic Potential Surfaces
The ESP visualization compares predicted and QM electrostatic potentials on a representative molecule’s 1.4× van der Waals surface. grappa-fullFF closely reproduces the QM surface, whereas AM1-BCC and espaloma-0.3.2 show larger deviations.
A.1 Qualitative Analysis of Electrostatic Potential Surfaces
- The 1.4× vdW surface approximates intermolecular contact distances, making its ESP a view of electrostatics relevant to nonbonded interactions.The figure shows a representative polar molecule.Source: §A.1 · p. 11, Figure 3 · p. 11
- grappa-fullFF reproduces the QM ESP surface’s spatial structure and magnitude, comparable to MBIS-derived charges.The reference ESP is derived from DFT.Source: §A.1 · p. 11, Figure 3 · p. 11
- AM1-BCC and espaloma-0.3.2 show larger deviations from QM in several surface regions, sometimes over- or under-polarizing.The authors interpret the visualizations as showing that ESP supervision constrains predicted charges to reproduce electrostatics from the ab initio electron density.Source: §A.1 · p. 11, Figure 3 · p. 11
A.2 Cumulative Error Distributions for Geometry Optimization
Across the OpenFF Industry Benchmark geometry-optimization metrics, grappa-fullFF’s error distributions nearly match grappa-bonded’s. The authors attribute this preservation of bonded accuracy to physical regularization of the additional nonbonded parameters.
A.2 Cumulative Error Distributions for Geometry Optimization
- Physical regularization constrains the learnable nonbonded degrees of freedom so they do not substantially alter bonded interactions despite energy degeneracy.This explains why joint parameter learning preserves the bonded-only model’s accuracy.Source: §A.2 · p. 11, Figure 4 · p. 12
- The CDFs compare geometry RMSE, TFD, and relative conformer-energy errors across OpenFF Industry Benchmark models.The distributions for grappa-fullFF and grappa-bonded are nearly indistinguishable.Source: §A.2 · p. 11, Figure 4 · p. 12
- Machine-learned force fields outperform GAFF-2.11 and OpenFF-2.1.0 on geometric fidelity, while grappa-fullFF matches its bonded-only predecessor.The reported metrics are geometry RMSE, TFD, and relative conformer energies.Source: §A.2 · p. 11, Figure 4 · p. 12
A.3 Quantitative Divergence of Ramachandran Free-Energy Surfaces
The section quantifies how closely force fields reproduce AMBER99SB-ILDN’s Ramachandran distribution using divergence metrics from explicit-solvent simulations. Both Grappa models are close to the reference, and grappa-fullFF is substantially closer than espaloma-0.3.2.
A.3 Quantitative Divergence of Ramachandran Free-Energy Surfaces
- The divergences compare each force field’s Ramachandran free-energy surface with AMBER99SB-ILDN using symmetric KL and JS metrics.The probability distributions are computed from 500 ns trajectories with 5° × 5° binning.Source: §A.3 · p. 12, Table 3 · p. 12
- 0.0820 nats: grappa-fullFF’s JS divergence is approximately twice as close to the AMBER99SB-ILDN distribution as espaloma-0.3.2.The comparison uses 5° × 5° histograms of (ϕ, ψ) from 500 ns of unbiased TIP3P MD; smaller divergence means closer to AMBER99SB-ILDN.Source: §A.3 · p. 12, Table 3 · p. 12
- Both grappa-bonded and grappa-fullFF closely reproduce the AMBER distribution.The reported JS divergence difference between the Grappa models is negligible, approximately 0.013 nats.Source: §A.3 · p. 12
B Ablation of the physics-informed regularization
This ablation compares electrostatic-potential regularization with direct supervision of predicted charges against MBIS reference charges. ESP regularization alone gives the lowest ESP errors across the three partitions, while charge supervision alone recovers charge magnitude but raises ESP errors.
B Ablation of the physics-informed regularization
- Without either charge regularization, ESP RMSE increases by 35–56% across all three partitions, and the charge-magnitude ratio is 0.74.The ratio is defined for mean absolute charge magnitudes on 30 PubChem molecules.Source: B Ablation of the physics-informed… · p. 13
- Charge supervision alone (λq = 1000, λESP = 0) recovers a charge-magnitude ratio of approximately 1 but increases ESP errors.This objective supervises predicted charges against MBIS reference charges through an MSE loss.Source: B Ablation of the physics-informed… · pp. 12–13
- ESP regularization alone gives the lowest ESP errors across all three dataset partitions, although the charge-magnitude ratio is 0.78.The authors conclude that ESP regularization can improve electrostatics and electric response properties without exactly recovering MBIS charges.Source: B Ablation of the physics-informed… · p. 13
- The ablation varies the ESP and MBIS charge-loss weights while holding the energy, force, and volume-ratio loss weights fixed at 2, 1, and 100.The final row of Table 4 is the configuration used to train grappa-fullFF.Source: B Ablation of the physics-informed… · p. 12, Table 4 · p. 13
C Implicit Learning of Polar Hydrogen Radii
grappa-fullFF learns environment-dependent Lennard-Jones radii for polar hydrogens from ab initio energies and force gradients. Its learned radii depart from TS-vdW references in a way that supports hydrogen bonding without unphysical Pauli repulsion.
C Implicit Learning of Polar Hydrogen Radii
- Regularizing LJ parameters toward TS-vdW volume ratios while training on ab initio energies and forces yields environment-dependent radii for hydrogen bonding.The learned values avoid the unphysical Pauli repulsion associated with larger polar-hydrogen radii.Source: C Implicit Learning of Polar Hydrogen… · p. 13, Figure 5 · p. 13
- grappa-fullFF learns systematically lower Lennard-Jones σ values for hydrogen than the relatively uniform TS-vdW reference values.The TS-vdW radii are approximately 2.3–2.5, whereas the learned values spread from approximately 1.1–2.5.Source: Figure 5 · p. 13
- Training on ab initio energy and force gradients alone teaches the model the downward spread in hydrogen σ values.The figure passage states that no reference radii are used to teach this correction.Source: Figure 5 · p. 13
D Molecular dynamics protocol
The MD protocol compares four force fields on the same solvated capped alanine dipeptide system. Simulations run for 500 ns in explicit TIP3P solvent, with the learned model parameterizing the solute and conventional parameters used for water.
D Molecular dynamics protocol
- All four force fields simulate the same capped alanine dipeptide in a 2.4 nm cubic box with 418 TIP3P waters, ions at 0.15 M ionic strength, and 1,278 total atoms.The system topology and starting coordinates are generated once and reused across the force fields.Source: D Molecular dynamics protocol · pp. 13–14
- grappa-fullFF parameterizes the solute entirely with the neural network while TIP3P water retains conventional parameters.The setup tests stability of the learned nonbonded parameters at the solute–solvent interface.Source: D Molecular dynamics protocol · p. 14
- After 100 ps NVT and 1 ns NPT equilibration, each system runs for 500 ns of NPT dynamics at 300 K and 1 atm with a 2 fs timestep.Long-range electrostatics use particle-mesh Ewald; frames are saved every 1 ps.Source: D Molecular dynamics protocol · p. 14