On-the-Fly MLFFs: Adaptive MD Force Fields

• On-the-fly MLFFs are adaptive interatomic models that update predictions during molecular dynamics using active learning and uncertainty quantification.
• They employ Bayesian inference with local descriptor regression to selectively trigger first-principles calculations, drastically reducing computational costs.
• These models are applied to phase transitions, extreme conditions, and disordered systems, though challenges remain in uncertainty calibration and long-range interaction capture.

On-the-fly machine-learned force fields (MLFFs) are interatomic models that are generated or refined during molecular dynamics (MD) rather than being fully prepared in advance. In the canonical form, the current MLFF predicts energies, forces, and often stresses together with an uncertainty estimate; a new first-principles calculation is then executed only when the current configuration is judged insufficiently represented by the existing training data. This closes an adaptive loop between simulation, error assessment, and reference labeling, and contrasts with the conventional train-then-deploy workflow of pre-fitted or offline-trained force fields (Jinnouchi et al., 2019, Unke et al., 2020, Smidstrup et al., 15 Sep 2025).

1. Conceptual scope and defining characteristics

The defining feature of an on-the-fly MLFF is not merely that an ML potential is evaluated during MD, but that the potential is updated during the simulation itself through active acquisition of new first-principles labels. In the review literature, this appears as an iterative loop: train a preliminary ML-FF on a small initial reference set, run MD, collect additional conformations whenever model predictions become unreliable according to an uncertainty criterion, perform new reference calculations, retrain, and repeat until no more unreliable regions are discovered (Unke et al., 2020). In the VASP implementation for melting-point calculations, this logic is operationalized as automatic force-field generation on the basis of Bayesian inference during MD, with first-principles calculations executed only when new configurations outside the already sampled dataset appear (Jinnouchi et al., 2019).

A persistent source of confusion is the use of “on-the-fly” to describe any MLFF evaluated inside an MD loop. Several papers in the broader MLFF literature explicitly do not satisfy the adaptive-learning criterion. QuantumATK, for example, “leverages machine-learned force-fields to simulate thermal properties and to generate realistic amorphous geometries in large-scale systems,” but the paper describes customized fitted MLFFs and pre-fitted general-purpose MLFFs in an offline train-then-use paradigm rather than a live active-learning loop (Smidstrup et al., 15 Sep 2025). ALIGNN-FF is a unified pretrained graph neural network force field spanning 89 elements, but it does not implement online retraining, uncertainty-based query selection, or autonomous dataset growth during simulation (Choudhary et al., 2022). Likewise, ASTEROID is a multi-stage data-cost-aware training framework for MLFFs that is explicitly framed as complementary to on-the-fly pipelines rather than a complete online loop (Bukharin et al., 2023).

A compact classification of the literature discussed here is therefore useful:

Category	Representative papers	Characterization
Genuine on-the-fly / adaptive MLFFs	(Jinnouchi et al., 2019, Liu et al., 2021, Kumar et al., 2023, Kumar et al., 2024, Timmerman et al., 2024, Kumar et al., 2024, Liu et al., 2020, Eggestad et al., 29 Oct 2025, Srivastava et al., 23 Jan 2026)	MD-time uncertainty assessment triggers new first-principles labeling and model update
Offline or adjacent MLFF methods	(Smidstrup et al., 15 Sep 2025, Choudhary et al., 2022, Bukharin et al., 2023, Alwis et al., 4 Mar 2026, Yin et al., 2024, Fonseca et al., 2023)	Pre-fitted deployment, transfer-learning, runtime profiling, ensemble fusion, or post hoc validation

This distinction matters methodologically. A paper may be highly relevant to on-the-fly MLFFs while not itself presenting adaptive retraining. The GPU study of ANI2x and TorchMD-Net, for instance, analyzes the cost of running pretrained MLFFs inside MD time stepping and is therefore an enabling deployment study rather than a paper on online model construction (Alwis et al., 4 Mar 2026). The same applies to validation tools such as FFAST, which address where MLFFs fail without implementing live acquisition or retraining (Fonseca et al., 2023).

2. Statistical and mathematical foundations

Most on-the-fly MLFFs in this corpus adopt a local energy decomposition. In the VASP melting-point framework, the total energy is written as

$U=\sum_{i=1}^{N_a} U_i,$

with each local contribution represented through kernel regression on local descriptors (Jinnouchi et al., 2019). The zirconium phase-transition work uses the same extensive structure, $U=\sum_i U_i$, with $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$ and a kernel expansion over selected basis environments (Liu et al., 2020). In the warm dense matter literature, the same idea is written elementwise as

$$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$

with SOAP descriptors and a polynomial kernel (Kumar et al., 2024). The Hugoniot framework uses the analogous expression for the electronic free energy $\mathcal{F}^{\rm MLFF}$ (Kumar et al., 2024).

The practical reason for this local form is that it yields energies, forces, and stresses from a single differentiable model. In the melting-point implementation, the target vector for each structure contains one energy-per-atom equation, $3N_a$ force components, and six independent stress components, and prediction is written compactly as

$\mathbf{y}^{\alpha}=\bm{\phi}^{\alpha}\mathbf{w},$

with the global training problem assembled as $\mathbf{y}=\mathbf{\Phi}\mathbf{w}$ (Jinnouchi et al., 2019). The zirconium paper uses the same linear-in-parameters structure and emphasizes that force rows are derivatives with respect to atomic coordinates and stress rows are derivatives with respect to unit-cell coordinates (Liu et al., 2020).

The adaptive element of these frameworks rests on Bayesian linear regression. In the VASP on-the-fly method, the posterior mean and covariance of the coefficients are

$$\bar{\mathbf w}=\frac{1}{\sigma_v^2}\mathbf{\Sigma}\mathbf{\Phi}^{\mathrm T}\mathbf Y,\qquad \mathbf{\Sigma}^{-1}=\frac{1}{\sigma_w^2}\mathbf I+\frac{1}{\sigma_v^2}\mathbf{\Phi}^{\mathrm T}\mathbf{\Phi},$$

and the predictive posterior for a new structure is

$$p(\mathbf y\mid \mathbf Y)=\mathcal N(\bm{\phi}\bar{\mathbf w},\bm{\sigma}), \qquad \bm{\sigma}=\sigma_v^2\mathbf I+\bm{\phi}^{\mathrm T}\mathbf{\Sigma}\bm{\phi}.$$

The diagonal elements of $U=\sum_i U_i$0 are used as the Bayesian errors in energies, forces, and stresses (Jinnouchi et al., 2019). Closely related BLR expressions are used in the WDM and Hugoniot frameworks, where the predictive uncertainty is written as the diagonal of a posterior covariance expression involving $U=\sum_i U_i$1, $U=\sum_i U_i$2, the kernel matrix, and the design matrix (Kumar et al., 2024, Kumar et al., 2024).

A major extension of this formalism is $U=\sum_i U_i$3-learning. In the orbital-free-to-Kohn–Sham framework, the learned target is not the full Kohn–Sham potential energy surface but the correction

$U=\sum_i U_i$4

with runtime predictions assembled as

$U=\sum_i U_i$5

This residual-learning design is explicitly motivated by the claim that the OF-to-KS correction is smoother and more local than the full KS target (Kumar et al., 2023). A related beyond-DFT strategy appears in zirconia, where an on-the-fly DFT-level MLFF is corrected by a second MLFF trained on $U=\sum_i U_i$6, $U=\sum_i U_i$7, and $U=\sum_i U_i$8 between RPA and DFT (Liu et al., 2021).

Numerical linear algebra can materially affect the final model. The zirconium study shows that solving the overdetermined least-squares problem by SVD and pseudoinversion of the design matrix improves the MLFF relative to the usual inversion of the squared matrix in regularized Bayesian regression, because the design matrix is strongly ill-conditioned and the normal equation squares the condition number (Liu et al., 2020). This is a methodological point about on-the-fly MLFFs themselves, not only about zirconium.

3. Uncertainty, triggering rules, and adaptive data acquisition

The distinctive operational question in on-the-fly MLFFs is when to call the reference electronic-structure method. In the VASP melting-point framework, a new first-principles calculation is performed when either the maximum Bayesian force uncertainty or the maximum spilling factor exceeds its threshold, subject to a short history-based rule that avoids oversampling nearby configurations (Jinnouchi et al., 2019). The threshold for the Bayesian error is itself adaptive: it is determined from recent post-refinement uncertainty values rather than fixed once and for all (Jinnouchi et al., 2019). In the orbital-free $U=\sum_i U_i$9-MLFF framework, new KS calculations are performed whenever the uncertainty in the force prediction exceeds a threshold $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$0, and $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$1 is set adaptively from the first post-training MLFF step (Kumar et al., 2023). The WDM and Hugoniot implementations in SPARC use the same logic, phrased as a force-uncertainty threshold $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$2 that triggers a new SQ-DFT or KS-DFT label when exceeded (Kumar et al., 2024, Kumar et al., 2024).

Not all adaptive workflows use the same representation, and chemical complexity can become the dominant bottleneck. In SPARC’s on-the-fly Bayesian MLFF framework, SOAP-based featurization causes the descriptor dimensionality to grow with the number of chemical species, which in turn increases the number of DFT calls required during active learning. Replacing SOAP with normalized Gaussian multipole features changes that scaling: for bulk systems with up to six elements, the number of DFT calls remains approximately independent of the number of chemical species, with alloy training counts of 63, 55, 70, 49, and 52 Kohn–Sham calls for PtAg, PtAgAu, PtAgAuIr, PtAgAuIrPd, and PtAgAuIrPdRh using GMP, versus 124, 147, 145, 153, and 170 with SOAP (Timmerman et al., 2024). This is one of the clearest examples in the literature of an on-the-fly MLFF bottleneck being shifted from the regression scheme to the descriptor layer.

At the same time, the calibration of the uncertainty model is a recognized weak point. The GMP paper explicitly shows that the Bayesian uncertainty estimate can be poorly calibrated and, in some systems, not even correlated with the true error; Pearson correlations between maximum predicted uncertainty and actual RMS force error span roughly $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$3 to $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$4 for GMP and from near zero to about one for SOAP (Timmerman et al., 2024). It also identifies a failure mode in which a dynamically updated threshold becomes too permissive after a large error, so later DFT checks are no longer triggered even as the trajectory enters unphysical regions (Timmerman et al., 2024). This is not a minor implementation detail: it is a central controversy in the design of uncertainty-driven on-the-fly learning.

The acquisition rule may also be embedded in a broader data-management strategy. In the WDM and Hugoniot frameworks, only atoms whose Bayesian force errors exceed a threshold are first added to the training set, and then CUR decomposition is used to downsample the selected environments (Kumar et al., 2024, Kumar et al., 2024). In the RPA zirconia workflow, a second active-learning pass collects 1275 small-cell candidate structures, after which a leverage-score CUR algorithm compresses them to 168 representative configurations for expensive high-level labeling (Liu et al., 2021). This suggests a general pattern: on-the-fly MLFFs increasingly combine adaptive acquisition with post-acquisition compression.

4. Representative scientific applications

The earliest explicit large-scale demonstration in this set is melting-point prediction. In VASP, on-the-fly MLFF generation was applied to Al, Si, Ge, Sn, and MgO, with more than 99% of first-principles calculations bypassed during force-field generation and production MD accelerated by 2000–5000 times per MD step relative to first-principles MD (Jinnouchi et al., 2019). The generated force fields remained compact—typically fewer than 500 structure datasets and fewer than 1000 local reference configurations—and thermodynamic perturbation theory showed that the MLFF melting points quantitatively reproduce first-principles melting points (Jinnouchi et al., 2019).

Finite-temperature solid-state phase transitions have become a major proving ground. For zirconium, on-the-fly MLFFs reproduce the first-order displacive $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$5-$U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$6 transition, including an abrupt jump of the volume and a cooperative displacement of atoms, with a predicted transition temperature of about 1049 K after free-energy correction, compared with the experimental 1136 K (Liu et al., 2020). For zirconia, on-the-fly DFT-level learning plus $U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$7-learning to RPA accuracy yields a single production-ready MLFF capable of reproducing high-level energies, forces, and stresses; the tetragonal–monoclinic transition temperature moves from 1492 K with MLFF-SCAN to 1415 K with MLFF-RPA$U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$8, close to the experimental 1400 K (Liu et al., 2021). In ferroelectric oxides, VASP’s on-the-fly MLFF implementation generated force fields starting only from the 0 K ground-state structure and qualitatively recovered the phase sequences of BaTiO$U_i=F[\rho_i^{(2)},\rho_i^{(3)}]$9, PbTiO$$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$0, LiNbO$$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$1, and BiFeO$$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$2, while also revealing order-disorder and displacive features of local cation displacements; however, the transition temperatures remained sensitive to the exchange-correlation functional, with PBEsol judged more robust than LDA and $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$3SCAN (Eggestad et al., 29 Oct 2025).

Extreme conditions constitute another major domain. In warm dense carbon, the SQ-MLFF framework reproduces diffusion and viscosity from direct SQ-DFT while reducing 500,000-step trajectories to only 306–110 SQ calls depending on temperature and achieving speedups of roughly $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$4 (Kumar et al., 2024). The Hugoniot extension adds an internal-energy model to the on-the-fly free-energy MLFF formalism, accelerates Kohn–Sham Hugoniot calculations by up to two orders of magnitude, and computes Hugoniots for 14 materials in the FPEOS database between 10 kK and 2 MK, with less than 1.4% of MD steps requiring KS labeling in the carbon validation (Kumar et al., 2024).

Disordered soft matter and polymers show that the approach is not confined to crystalline inorganic systems. For polymer glass transitions, VASP’s Bayesian on-the-fly framework was used to generate MLFFs from only about 1000 AIMD-sampled configurations per polymer during 200 ps high-temperature NPT runs on roughly 200–250 atom cells; these models were then deployed in 4600–8000 atom cooling simulations over 23 temperature points, leading to approximately six orders of magnitude cost reduction relative to direct AIMD and $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$5 predictions across twelve polymers in good agreement with experiment (Srivastava et al., 23 Jan 2026). This suggests that on-the-fly MLFFs are now credible for long-time disordered-matter simulations, provided the training trajectory is thermodynamically relevant.

The same theme appears in materials for quantum devices, though there the MLFF role is not on-the-fly. QuantumATK positions MLFFs as enabling realistic amorphous Al/AlO$$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$6/Al junction structures, thermal-property simulations, and phonon/electron-phonon calculations relevant to superconducting qubits and sensors, but the paper does not describe adaptive retraining during production simulation (Smidstrup et al., 15 Sep 2025). This is a reminder that the scientific success of MLFFs in a domain does not imply the use of on-the-fly methodology.

5. Performance, hardware, and deployment constraints

On-the-fly MLFF research sits at the intersection of statistical learning and high-performance MD. Descriptor choice directly affects the number of expensive reference calls. In the multicomponent alloy setting, normalized Gaussian multipole descriptors not only keep DFT calls roughly flat with element count but also maintain compact design matrices: in the alloy benchmarks, GMP used 39–55 columns versus 203–276 for SOAP (Timmerman et al., 2024). This reduction is operationally important because the total cost of optimized on-the-fly implementations is dominated by the number of first-principles calls (Timmerman et al., 2024).

Once the force field exists, runtime execution becomes its own bottleneck. The GPU workload study of ANI2x and TorchMD-Net shows that MLFF-based MD is not “just inference”: each MD step requires both forward and backward computations because forces are energy gradients (Alwis et al., 4 Mar 2026). For ANI2x, radial AEV construction scales approximately as $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$7, angular AEV construction as $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$8, and the descriptor has 1008 components (Alwis et al., 4 Mar 2026). For the TorchMD-Net equivariant transformer, per-layer complexity scales approximately as $$\mathcal{E}^{\rm MLFF}=\sum_E \sum_{i=1}^{N_A^E}\sum_{t=1}^{N_T^E}\tilde w_t^E\,k(\mathbf{x}_i^E,\tilde{\mathbf{x}}_t^E),$$9 (Alwis et al., 4 Mar 2026). In both cases, the dominant performance problems are irregular memory access, scatter/gather operations, poor cache reuse, and many short-lived kernels. The paper therefore treats runtime MLFF execution as a distinct workload class rather than a simple extension of classical MD or standard neural inference (Alwis et al., 4 Mar 2026).

This systems perspective matters for on-the-fly MLFFs even when the paper itself does not implement adaptive learning. If uncertainty-triggered DFT calls become infrequent, then the critical path shifts toward descriptor construction, neural or kernel inference, and force backpropagation. The GPU paper shows that kernel fusion can radically improve this regime: NNPOps reduces ANI AEV evaluation from roughly 60–65 CUDA kernel launches per AEV to four fused kernels, yielding nearly 100× faster AEV evaluation and about 17×–20× speedup including energy/force computation (Alwis et al., 4 Mar 2026). A plausible implication is that future on-the-fly MLFF frameworks will need co-design across acquisition logic, descriptor engineering, and low-level MD implementation, not merely better uncertainty criteria.

Current deployment stacks also impose practical limits. The same GPU study notes that current ANI2x integrations into OpenMM and GROMACS lack domain decomposition support for ML potentials, preventing scalable multi-GPU MLFF MD in the same way as classical FFs (Alwis et al., 4 Mar 2026). This is not an algorithmic limitation of on-the-fly learning per se, but it constrains the system sizes and trajectory lengths at which adaptive MLFFs can be exploited efficiently.

6. Validation, limitations, and adjacent research directions

A recurring lesson is that average error metrics are inadequate for judging practical safety. FFAST was developed precisely because global MAE and RMSE can hide the localized force failures that destabilize MD or corrupt observables (Fonseca et al., 2023). In stachyose, atom-projected errors showed that carbons and oxygens near glycosidic bonds were substantially worse predicted than the hydrogen-dominated average suggested; in DHA, errors rose as the molecule folded, especially near the carboxylic group (Fonseca et al., 2023). For on-the-fly MLFFs, this means that active-learning success cannot be assessed only by aggregate force RMSE. Local motifs and trajectory-specific outliers matter more.

Long-range physics remains a structural limitation of many local on-the-fly models. The ferroelectric phase-transition paper uses an 8 Å cutoff for local environments and explicitly states that the learned force field does not directly include electrostatic interactions beyond that radius. The authors connect this locality assumption to the systematic underestimation of transition temperatures in ferroics, where long-range Coulomb interactions are essential (Eggestad et al., 29 Oct 2025). The same paper also highlights finite-size effects and incomplete ergodicity over short MD windows as additional sources of bias (Eggestad et al., 29 Oct 2025). These are not easily fixed by adding more labels; they reflect representational and sampling limits.

Not every strategy that reduces high-fidelity label demand is itself on-the-fly. ASTEROID lowers data cost by pretraining on cheap inaccurate or unlabeled data, then fine-tuning on a small accurate set. It improves sample efficiency and downstream MD robustness, but it does not include uncertainty estimation, streaming updates, adaptive frame selection, or online retraining during MD (Bukharin et al., 2023). Its relevance is therefore as an offline or semi-offline complement to on-the-fly campaigns, not a substitute for adaptive acquisition. Similarly, EL-MLFFs stack several pre-trained MLFFs with a GNN meta-model and improve force RMSE substantially on methane and methanol/Cu(100), but the method does not implement active learning, online retraining, or a conservative energy model, and the paper does not demonstrate stable long-time MD (Yin et al., 2024).

A final misconception concerns “world models” or universal pretrained MLFFs. ALIGNN-FF demonstrates that one can train a single graph-based force field across 89 elements and obtain usable energies, forces, relaxation behavior, and short MD for structurally diverse materials (Choudhary et al., 2022). This suggests a possible role for universal models as priors or warm starts in adaptive workflows. But the paper itself does not provide uncertainty quantification or online correction, and it therefore leaves open the core on-the-fly question of when to trust the model outside its training manifold (Choudhary et al., 2022).

Taken together, these studies show that on-the-fly MLFFs are best understood not as a single model family but as a construction strategy: a symmetry-respecting, differentiable interatomic model is coupled to MD, monitored for unreliability, selectively corrected by first-principles calls, and iteratively refined until the intended portion of configuration space is covered (Unke et al., 2020). The most mature implementations already support phase transitions, warm dense matter, shock Hugoniots, multicomponent alloys, and polymer glass transitions, but the central open problems remain uncertainty calibration, locality versus long-range physics, stable large-scale deployment, and validation beyond global average errors (Timmerman et al., 2024, Fonseca et al., 2023).