Projected Hessian Learning: Fast Curvature Supervision for Accurate Machine-Learning Interatomic Potentials
Authors: Austin Rodriguez, Justin S. Smith, Sakib Matin, Nicholas Lubbers, Kipton Barros, Jose L. Mendoza-Cortes
Organizations: Theoretical Division, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, United States. · Department of Chemical Engineering & Materials Science, Michigan State University, East Lansing, Michigan 48824, United States. · NVIDIA Corp., 2788 San Tomas Expy, Santa Clara, California 95051, United States. · Computer and Artificial Intelligence Division, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, United States. · Center for Nonlinear Studies, Los Alamos National Laboratory, Los Alamos, New Mexico 87545, United States. · Department of Physics and Astronomy, Michigan State University, East Lansing, MI, 48823, United States.
The Hessian matrix (second derivatives) encodes far richer local curvature of the potential energy surface than energies and forces alone. However, training machine-learning interatomic potentials (MLIPs) with full Hessians is often impractical because explicitly forming and storing Hessian matrices scales quadratically in cost and memory. We introduce Projected Hessian Learning (PHL), a scalable second-order training framework that injects curvature information using only Hessian-vector products (HVPs). Rather than constructing the Hessian, PHL projects curvature along stochastic probe directions and uses an unbiased stochastic trace-based loss with favorable system-size scaling, enabling curvature-informed training without quadratic memory growth. We benchmark PHL on a chemically diverse dataset of reactants, products, transition states, intrinsic reaction coordinates, and normal-mode sampled geometries computed at omegaB97XD/6-31G(d). We compare energy-force training (E-F), two HVP-based schemes (E-F-HVP with one-hot or randomized probes), and full energy-force-Hessian training (E-F-H). With randomized probes per minibatch, both HVP schemes match full-Hessian training in energy, force, and Hessian accuracy while delivering >24x epoch speedups for the small molecular systems studied. In a fixed-probe regime with one HVP per molecule, randomized projections consistently outperform one-column probing, especially for far-from-equilibrium geometries. Overall, PHL replaces explicit Hessian supervision with force-complexity curvature training, retaining most second-order accuracy gains while scaling to larger, more complex molecular systems.
Figures & tables
Figure 1: Conceptual comparison of curvature supervision strategies used to train MLIPs. Left: Full Hessian training explicitly uses all (3N)×(3N) second-derivative elements of M systems. Middle: One-column HVP probing samples a single Hessian column via a canonical basis vector. Right: Projected Hessian Learning (PHL) (this work). Curvature is enforced through stochastic projections using random probing vectors to form Hessian-vector products, yielding random weighted combinations of Hessian columns (Hutchinson-style probing). Both HVP-based strategies avoid explicit Hessian construction; PHL aggregates information across multiple curvature directions in expectation.
Figure 2: Validation loss curves for energy, force, and Hessian predictions under the fixed-vector approach. Each panel shows the RMSE of the validation set as a function of training epoch for energy and force training (E–F, red), energy, force, and Hessian–vector product training using the one-column method (E–F–HVP One-Column, orange), energy, force, and HVP training using the Hutchinson estimator (E–F–HVP PHL, purple), and full energy–force–Hessian training (E–F–H, blue).
Figure 3: Validation RMSE comparison across the Test Set, IRC dataset, and NMS dataset for models trained with different levels of information using (a) randomized probe vectors (resampled each minibatch) and (b) fixed probe vectors (one probe per system). Results are shown for energy and force training (E–F); energy, force, and Hessian–vector products (HVP) using the one-column method (E–F–HVP One-Column); energy, force, and HVP using the Hutchinson estimator (E–F–HVP PHL); and energy, force, and full Hessian training (E–F–H). From top to bottom, panels report RMSE for energies, forces, and Hessians. Bars represent mean RMSE values, with error bars indicating variability across ensembles of five independently trained models. Incorporating Hessian information improves accuracy relative to E-F across all datasets; under randomized probing, one-hot and PHL achieve statistically indistinguishable performance for the small molecular systems studied here, whereas under fixed probes PHL provides systematically lower errors, most notably for extrapolative NMS geometries.
Property
Dataset
Mean Difference (kcal/mol)
p -value
Significance
Energy
Test
-0.0039
0.950
–
Energy
IRC
-0.0775
0.348
–
Energy
NMS
-0.2320
0.342
–
Force
Test
0.0134
0.832
–
Force
IRC
-0.0198
0.732
–
Force
NMS
-0.1444
0.195
–
Table 1: Two-sample t -test results comparing RMSE values from the one-column and PHL estimators under the randomized-vector training approach. The mean difference is defined as (one-column – PHL). Negative values indicate lower RMSE for the one-column method. Asterisks denote statistically significant differences ( p<0.05 ).
Property
Dataset
Mean Difference (kcal/mol)
p -value
Significance
Energy
Test
0.2370
0.063
–
Energy
IRC
0.1303
0.399
–
Energy
NMS
0.8989
0.006
Force
Test
0.1268
0.097
–
Force
IRC
0.2076
0.013
Force
NMS
0.8086
0.006
Table 2: Two-sample t -test results comparing RMSE values from the one-column and PHL estimators under the fixed-vector training approach. The mean difference is defined as (one-column – PHL). Negative values indicate lower RMSE for the one-column method. Asterisks denote statistically significant differences ( p<0.05 ).
Figure 4: Execution time per training epoch for different methods: energy and force training (E–F); energy, force, and Hessian–vector products (HVP) using the one-column method (E–F–HVP One-Column); energy, force, and HVP using the PHL method (E–F–HVP PHL); and energy, force, and full Hessian training (E–F–H). Bars show mean times with error bars indicating variability across epochs. Both one-column and PHL estimators achieve more than a 24-fold speedup compared to full Hessian training, while E–F training is the fastest but lacks Hessian information.
Figure 5: Scaling of CPU time with system size for density functional theory (DFT) calculations using Gaussian16 on a semi-log plot (a) and a log-log plot (b). Shown are the computational costs for evaluating total energies (green squares), forces (red triangles), and full Hessians (blue diamonds), and the estimated cost for evaluating Hessian–vector products (HVP, purple line). While force and energy costs scale comparably, the cost of full Hessian evaluations grows much more steeply with the number of atoms. In contrast, HVP evaluations provide a significantly cheaper alternative, with scaling closer to that of two force calculations, highlighting their efficiency for incorporating second-derivative information.
Appendix figures & tables6 assets
Supplementary material from the paper’s appendix.
Appendix
Figure S1 : Bland–Altman analysis of energy RMSE differences for randomized probes per minibatch. Differences are defined as ΔRMSE=RMSEone-column−RMSEPHL across the Test, IRC, and NMS datasets. Each point represents a paired set of trained models. The green line denotes the mean difference and the shaded region its 95% CI; black lines indicate limits of agreement. The 95% CI includes zero for all datasets, indicating no statistically significant difference in energy accuracy between probing strategies when probes are randomized each minibatch.
Figure S2 : Bland–Altman analysis of force RMSE differences for randomized probes per minibatch, with ΔRMSE defined as in Fig. S1 . The 95% CI includes zero for all datasets, indicating no statistically significant difference in force accuracy between one-column and PHL methods when probes are randomized each minibatch.
Figure S3 : Bland–Altman analysis of Hessian RMSE differences for randomized probes per minibatch, with ΔRMSE defined as in Fig. S1 . The 95% CI includes zero for all datasets, indicating no statistically significant difference in Hessian accuracy between probing strategies in the randomized-probe regime.
Figure S4 : Bland–Altman analysis of energy RMSE differences for fixed probes per system. Differences are defined as ΔRMSE=RMSEone-column−RMSEPHL across the Test, IRC, and NMS datasets. The mean difference is not significant for Test and IRC (95% CI includes zero), but is significantly positive for NMS (95% CI excludes zero), indicating better energy extrapolation for Hutchinson-based PHL in the data-limited regime.
Figure S5 : Bland–Altman analysis of force RMSE differences for fixed probes per system, with ΔRMSE defined as in Fig. S4 . The 95% CI excludes zero for IRC and NMS, demonstrating that the one-column estimator yields significantly higher force errors than the Hutchinson-based PHL estimator on intermediate (IRC) and extrapolative (NMS) geometries.
Figure S6 : Bland–Altman analysis of Hessian RMSE differences for fixed probes per system, with ΔRMSE defined as in Fig. S4 . The 95% CI excludes zero for Test, IRC, and NMS, indicating that the one-column estimator produces significantly higher Hessian errors than the Hutchinson-based PHL estimator across all datasets in the fixed-probe regime, with the largest discrepancies observed on NMS.
Local curvature of potential energy surfaces is critical for predicting certain experimental observables of molecules and materials from first principles, yet it remains far beyond reach for complex systems. In this work, we introduce a Hessian-informed Machine Learning Interatomic Potential (Hi-MLIP) that captures such curvature reliably, thereby enabling accurate analysis of associated thermodynamic and kinetic phenomena. To make Hessian supervision practically viable, we develop a highly efficient training protocol, termed Hessian INformed Training (HINT), achieving two to four orders of magnitude reduction for the requirement of expensive Hessian labels. HINT integrates critical techniques, including Hessian pre-training, configuration sampling, curriculum learning and stochastic projection Hessian loss. Enabled by HINT, Hi-MLIP significantly improves transition-state search and brings Gibbs free-energy predictions close to chemical accuracy especially in data-scarce regimes. Our framework also enables accurate treatment of strongly anharmonic hydrides, reproducing phonon renormalization and superconducting critical temperatures in close agreement with experiment while bypassing the computational bottleneck of anharmonic calculations. These results establish a practical route to enhancing curvature awareness of machine learning interatomic potentials, bridging simulation and experimental observables across a wide range of systems.
Bangchen Yin, Jian Ouyang, Zhen Fan +7
ByteDance Seed · Tsinghua University · Peking University
While machine-learning interatomic potentials (MLIPs) have successfully learned potential energy surfaces (PES) and atomic forces, many practical applications, such as vibrational analysis and transition state search, rely heavily on the PES Hessian. Yet standard MLIPs are trained on energy and forces only, and using the Hessian as an additional training objective requires backpropagation whose memory cost grows quadratically with the number of atoms. To address these limitations, we propose two Hessian-derived data augmentation schemes: isotropic Gaussian displacement (UniAug) and normal mode-weighted displacement (ModeAug). Both methods utilize simple Taylor expansions, generating physically grounded molecular structures with highly accurate energy and force labels. By preserving the underlying MLIP architecture and training procedure, this approach contributes to improving energy and force predictions in various benchmarks. Comprehensive evaluations across non-equilibrium and equilibrium datasets demonstrate that our approach enhances model accuracy where reference forces are large and training data are limited.
Fundamental tasks in computational chemistry, from transition state search to vibrational analysis, rely on molecular Hessians, which are the second derivatives of the potential energy. Yet, Hessians are computationally expensive to calculate and scale poorly with system size, with both quantum mechanical methods and neural networks. In this work, we demonstrate that Hessians can be predicted directly from a deep learning model, without relying on automatic differentiation or finite differences. We observe that one can construct SE(3)-equivariant, symmetric Hessians from irreducible representations (irrep) features up to degree l=2 computed during message passing in graph neural networks. This makes HIP Hessians one to two orders of magnitude faster, more accurate, more memory efficient, easier to train, and enables more favorable scaling with system size. We validate our predictions across a wide range of downstream tasks, demonstrating consistently superior performance for transition state search, accelerated geometry optimization, zero-point energy corrections, and vibrational analysis benchmarks. We open-source the HIP codebase and model weights to enable further development of the direct prediction of Hessians at https://github.com/BurgerAndreas/hip
Andreas Burger, Luca Thiede, Nikolaj Rønne +5
University of Toronto · Vector Institute for Artificial Intelligence · Technical University of Denmark +4