Accessing Numerical Energy Hessians with Graph Neural Network Potentials and Their Application in Heterogeneous Catalysis
Brook Wander, Joseph Musielewicz, Raffaele Cheula and John R. Kitchin, “Accessing Numerical Energy Hessians with Graph Neural Network Potentials and Their Application in Heterogeneous Catalysis”, The Journal of Physical Chemistry C 129(7), 3510-3521 (2025). https://doi.org/10.1021/acs.jpcc.4c07477
This summary was written by a language model from the paper and has not been checked line by line. Rely on the published paper for details.
Note: the extracted text available for this summary is the arXiv preprint (2410.01650v2) with its Supplementary Information, not the typeset journal version; numbers below are from that text.
Question
The second derivative of the potential energy (the Hessian) gives vibrational frequencies, and from them zero-point energies and vibrational entropies for Gibbs free energies, and it also drives Hessian-based saddle point optimizers. Computing it by DFT finite differences is costly. The paper asks whether graph neural network potentials pretrained on Open Catalyst 2020 (OC20) forces and energies, which were never trained on curvature, can supply usable Hessians for adsorbates on catalyst surfaces; whether cheap potentials make entropy treatments beyond the harmonic approximation practical; what the easiest route to Gibbs corrections is; and whether model Hessians help transition state optimization (Wander et al. 2025).
Approach
The authors built a DFT (VASP, RPBE, OC20 settings) dataset of 4,599 finite-difference Hessians: 1,812 relaxed systems from the OC20 validation set, 2,151 from the OC20 training set (under 100 atoms), and 636 transition states from the OC20NEB set, using two displacements of 0.0075 Å per direction. Five pretrained models (SchNet, PaiNN, GemNet-T, and EquiformerV2 with 31M and 153M parameters) produced Hessians by the same finite-difference scheme; finite differences were preferred over automatic differentiation, which was slower and no more accurate against VASP. Accuracy was measured by MAE of sorted frequencies, of the Gibbs correction at 300 K (computed with ASE HarmonicThermo, discarding imaginary modes), and of the largest-imaginary or smallest-real (LI/SR) frequency that identifies a saddle point. The 153M EquiformerV2 was also fine-tuned on energies and forces of the small-displacement single points from training-set minima (learning rate 5e-6, best at 12 epochs). For 100 validation systems, complete potential energy sampling (CPES) with the 31M EquiformerV2 mapped the adsorbate’s 2D energy surface on a 0.5 Å grid of constrained relaxations to obtain a translational partition function and entropy, checked with DFT single points. Finally, the Sella optimizer, with a model Hessian recomputed every step, refined ML-NEB transition states from the CatTSunami workflow on OC20NEB. Code and data are in the gibby repository.
Key findings
The best pretrained model (EquiformerV2 153M) reproduced DFT frequencies with an MAE of 58 cm-1 on relaxed systems and 88 cm-1 on transition states, with a small systematic overestimate of frequencies; model ranking followed OC20 force and energy accuracy. With a constant bias correction the Gibbs correction MAE at 300 K was 0.042 eV for minima and 0.064 eV for transition states, beating the common practice of computing a correction for each adsorbate on one surface and reusing it elsewhere (“adopt-a-value”, 0.058 eV) and slightly beating a lookup table of mean per-adsorbate frequencies (0.045 eV), which the SI tabulates but which cannot cover arbitrary adsorbates or transition states. Fine-tuning lowered the all-frequency MAE on minima to about 29 cm-1 and the Gibbs MAE to 0.031 eV (0.068 eV for transition states), but it worsened the LI/SR frequency on transition states (from about 103 to 202 cm-1) because training only near minima biased the model against imaginary modes.
CPES showed that neglecting translational freedom misses on average 0.16 eV (maximum 0.27 eV) of the Gibbs correction at 300 K, with 94% of sampled systems above 0.1 eV. On 17 systems checked by DFT single points, the ML translational term had an MAE of 0.023 eV, though a few spurious low ML energies caused outliers, and constraint failures or desorption degraded some grids. In transition state search, Sella refinement raised the converged fraction from 80.2% (CatTSunami baseline) to 93.0%, a 65% cut in unconverged cases, at a small cost in energy MAE (0.069 versus 0.060 eV); keeping only cases where both NEB and Sella converged gave 0.041 eV MAE and 95.2% success but 70.6% coverage.
Significance
The work shows that force-trained catalysis potentials carry useful curvature information for free, which lets thermochemical corrections and Hessian-driven saddle searches be computed per system rather than borrowed across surfaces. It also flags translational entropy as a sizable, often ignored term that cheap potentials can now estimate, and gives a practical caution that fine-tuning data must resemble the target task (saddle points) to help it.