Learning consistent molecular mechanics force fields from first principles

Berkay Günes; Leif Seute; Jigyasa Nigam; Frauke Gräter

physics.chem-phOct 6, 2026

Learning consistent molecular mechanics force fields from first principles

Authors: Berkay Günes, Leif Seute, Jigyasa Nigam, Frauke Gräter

Organizations: Max Planck Institute for Polymer Research Mainz, Germany · Heidelberg University, Germany · Heidelberg Institute for Theoretical Studies, Germany · Initiative for Computational Catalysis, Flatiron Institute, United States

Abstract

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.

Figures & tables

Figure 1: Grappa architecture and training routine. The shared graph-attention encoder feeds the atom embeddings in the individual parameter prediction heads. All parameters are trained simultaneously on QM energies and forces. The nonbonded parameter heads are additionally regularized by the electrostatic potential (ESP) and the ratio of the effective and free atomic volume, $v_{i}$ .

	PubChem		Dipeptide		DES Monomers
Charge Model	ESP	Dipole	ESP	Dipole	ESP	Dipole
MBIS 1 1 1 Require DFT calculations	7.157(54)	$1.190(26)$	8.47(37)	$1.724(79)$	$6.22(25)$	$0.622(50)$
AM1-BCC 1 1 1 Require DFT calculations	$10.414(87)$	$1.156(25)$	$12.32(54)$	$1.88(10)$	$8.50(40)$	$0.597(47)$
espaloma-charge	$35.7(10)$	$2.084(40)$	$30.3(11)$	$2.174(71)$	$18.2(17)$	$0.758(66)$
grappa-fullFF	$8.90(58)$	1.048(24)	$12.7(13)$	1.584(81)	6.12(26)	0.561(46)

Table 1: Comparison of ESP (

10^{-3}e/\text{\AA}

) and dipole moment (

e\text{\AA}

) RMSE across charge models and datasets. Parenthetical values give the standard deviation on the last digits.

Force Field	RMSE med (Å)	TFD med	ddE MAE	ddE RMSE
gaff-2.11	0.4997	0.0476	2.341	3.892
openff-2.1.0	0.3688	0.0355	2.036	3.495
espaloma-0.3.2	0.2725	0.0241	1.639	3.064
grappa-bonded	0.2571	0.0232	1.608	2.773
grappa-fullFF	0.2661	0.0221	1.589	2.973

Table 2: Comparison of force field performance metrics: median of the Root Mean Square Error (RMSE med) in geometries between MM-optimized and QM-optimized conformers, Torsion Fingerprint Deviation (TFD), and the Mean Absolute Error (MAE) and RMSE of the conformer-averaged relative energy (ddE).

Figure 2: Ramachandran FES of capped Alanine dipeptide from

500

ns explicit-solvent (TIP3P) MD.

Appendix figures & tables5 assets

Supplementary material from the paper’s appendix.

Appendix

Figure 3: ESP mapped onto the $1.4\times$ vdW surface of a representative molecule. grappa-fullFF closely recovers the ESP surface derived from DFT with high accuracy, matching the computationally more expensive partial charges derived from the MBIS density partitioning. Standard empirical methods such as AM1-BCC and machine-learned charge models (espaloma-0.3.2) show noticeably higher deviations from the QM reference.

Figure 4: Cumulative distribution functions (CDFs) for the geometry Root Mean Square Error (RMSE), Torsion Fingerprint Deviation (TFD), and relative conformer energies (|ddE|) evaluated on the OpenFF Industry Benchmark Season 1 v1.1. The machine-learned force fields (espaloma-0.3.2, grappa-bonded , grappa-fullFF ) consistently outperform classical force fields (GAFF-2.11, OpenFF-2.1.0) on geometric fidelity. grappa-fullFF matches the structural accuracy of its bonded-only predecessor, confirming that the joint prediction of nonbonded parameters does not degrade bonded performance.

Force field	Symmetric KL [nats]	Jensen–Shannon [nats]
AMBER99SB-ILDN (reference)	0.000	0.0000
espaloma-0.3.2	1.752	0.1596
grappa-bonded	0.432	0.0691
grappa-fullFF	0.502	0.0820

Appendix

Table 3: Divergence of each force field’s Ramachandran free-energy surface from AMBER99SB-ILDN, computed from a

5^{\circ}\!\times\!5^{\circ}

histogram of

(\phi,\psi)

over

500

ns of unbiased TIP3P MD. Smaller is closer to AMBER99SB-ILDN.

Config	$\lambda_{q}$	$\lambda_{\text{ESP}}$	ESP ${}_{\text{PubChem}}$	ESP ${}_{\text{Dipep.}}$	ESP ${}_{\text{DES}}$	$\dfrac{\langle q_{i}\rangle}{\langle q_{i}^{\text{MBIS}}\rangle}$	ddE MAE
No regularization	0	0	13.84	17.35	9.10	0.74	—
$\lambda_{q}\neq 0,\lambda_{\text{ESP}}=0$	1000	0	11.80	16.52	8.71	0.99	1.72
$\lambda_{q},\lambda_{\text{ESP}}\neq 0$	1000	100	11.07	20.24	7.55	$\approx\!1.00$	—
$\lambda_{q}=0,\lambda_{\text{ESP}}\neq 0$	0	1000	8.86	12.84	6.13	0.78	1.58
ESP in units of $10^{-3}e/\text{\AA}$

Appendix

Table 4: Ablation of the two nonbonded regularization terms.

\lambda_{q}

weighs the loss against per-atom MBIS charges as defined in Eq. ( 7 ),

\lambda_{\text{ESP}}

weighs the ESP loss defined in Eq. ( 6 ). The energy, force, and volume-ratio loss weights are held fixed at

(\lambda_{E}=2,\lambda_{F}=1,\lambda_{v}\!=\!100)

\langle q_{i}\rangle/\langle q_{i}^{\text{MBIS}}\rangle

is the ratio of mean (absolute) charge magnitudes on 30 PubChem molecules. ddE MAE is from the geometry-optimization benchmark of Table 2 .

Figure 5: Comparison of predicted Lennard-Jones $\sigma$ parameters against the TS-vdW reference values for the three OMol25 SPICE subsets. The TS-vdW references assign relatively uniform radii to all hydrogen atoms ( $\sigma\approx 2.3-2.5\AA$ ), while grappa-fullFF learned a systematic downward spread ( $\sigma\approx 1.1-2.5\AA$ ). This correction is learned purely by training the model on ab initio energy and force gradients, which allows grappa-fullFF to accurately model hydrogen bonds without unphysical Pauli repulsion.

Explore similar work

Implicit Machine Learning Force Fields Accelerate Molecular Dynamics Simulations
Jul 31, 2026Johannes Maeß, Leon Werner, J. Thorben Frank +5Molecular Dynamics
BIFOLD, Berlin, Germany. · Machine Learning Group, TU Berlin, Berlin, Germany. · Department of Artificial Intelligence, Korea University, Seoul, South Korea. +1
Physics-Aligned Electronic Ground-State Learning Improves Generalization
Oct 7, 2026Eike S. Eberhard, Xaver Kainz, Viktor Kotsev +2
Technical University of Munich (TUM) · Munich Data Science Institute (MDSI) · Munich Center for Machine Learning (MCML) +1
EquiFiLM: Charge-Conditioned Equivariant Force Fields via Feature-wise Linear Modulation
Jul 6, 2026Samuel Sahel-Schackis, Ken-ichi Nomura, Aiichiro Nakano +2Quantum Chemistry Wireless Foundation Models
Department of Physics, Stanford University · Collaboratory for Advanced Computing and Simulation, University of Southern California · Linac Coherent Light Source, SLAC National Accelerator Laboratory +1

Jul 31, 2026, cs.LGJ/K move · Enter open · S save

Implicit Machine Learning Force Fields Accelerate Molecular Dynamics Simulations

Johannes Maeß, Leon Werner, J. Thorben Frank, Winfried Ripken+4

BIFOLD, Berlin, Germany. · Machine Learning Group, TU Berlin, Berlin, Germany. · Department of Artificial Intelligence, Korea University, Seoul, South Korea.+1

Molecular Dynamics

Open paper AskPDF arXiv

We introduce implicit machine learning force fields (I-MLFFs), which replace explicit stacks of neural network layers with self-consistent fixed-point equations. In molecular simulations, this formulation enables intermediate representations to be reused across successive timesteps, thereby warm-starting force evaluation. The resulting models effectively combine the computational footprint of a shallow, single-layer MLFF with the representational capacity and accuracy of a deep neural network. Our approach unlocks architecture-agnostic efficiency gains that are inaccessible when force prediction and trajectory integration are considered separately. We demonstrate this across three major classes of graph neural networks: invariant, equivariant Cartesian tensor, and SO(3)-equivariant spherical-tensor architectures. Each yields a two- to five-fold reduction in compute and memory footprint. Crucially, these gains are achieved while retaining full atomistic resolution and the original integration timestep, avoiding spatial or temporal coarse graining. Our contribution therefore advances the scaling frontier of quantum-mechanically faithful molecular simulation, enabling longer trajectories and larger atomistic systems within fixed GPU memory and compute budgets, and thereby opening access to new insights across biomolecular and material systems.

Abstract

Explore similar work

Implicit Machine Learning Force Fields Accelerate Molecular Dynamics Simulations

Physics-Aligned Electronic Ground-State Learning Improves Generalization

EquiFiLM: Charge-Conditioned Equivariant Force Fields via Feature-wise Linear Modulation