Locally Trained SOAP-GAP Model

• The paper introduces a machine-learned interatomic potential that uses SOAP descriptors and Gaussian process regression to predict energies, forces, and virials with quantum accuracy.
• The methodology leverages active learning and local filtering of datasets to enhance accuracy and transferability within specific configurational regions.
• Practical applications include large-scale molecular dynamics and Monte Carlo simulations, balancing computational efficiency with high-fidelity results.

A Locally Trained SOAP-GAP (Smooth Overlap of Atomic Positions - Gaussian Approximation Potential) model is a machine-learned interatomic potential constructed by fitting Gaussian process regression to local atomic environments, defined via the SOAP descriptor. This model targets specific configurational or spatial regions (such as defect cores, surfaces, or alloy concentrations) to achieve enhanced accuracy and transferability within the selected local domain. The model enables quantum-accurate prediction of energies, forces, and virials for atomistic simulations, with practical deployment in large-scale molecular dynamics, Monte Carlo sampling, and related techniques.

1. Theoretical Foundations

A SOAP-GAP model is predicated on the principle that the total potential energy of an atomic system can be decomposed as a sum over local atomic contributions, each represented as a feature vector derived from the local environment:

$E_\text{tot} = \sum_i E_i \,,$

where each $E_i$ is modeled as a function of the local atomic neighborhood using the SOAP representation and a kernel-based regression within the Gaussian process (GP) framework.

Neighbor Density and Basis Expansion

The local environment for atom $i$ is encoded via a smeared neighbor density

$$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$

where $z_j$ is the chemical identity, $\sigma$ is the atomic Gaussian width, and $f_\text{cut}$ is a smooth radial cutoff. This density is projected onto a set of orthogonal radial basis functions $g_n(r)$ and spherical harmonics $Y_{lm}(\hat{r})$:

$$\rho^\alpha(\mathbf{r}) = \sum_{nlm} c_{nlm}^\alpha\, g_n(|r|)\, Y_{lm}(\hat{r})$$

SOAP Power Spectrum

Permutationally, translationally, and rotationally invariant features are constructed from the expansion coefficients as the SOAP power spectrum:

$E_i$0

All $E_i$1 are concatenated into a single descriptor vector $E_i$2 for each atomic environment.

SOAP-GAP Kernel

Similarity between environments is measured by a normalized dot-product polynomial kernel:

$E_i$3

with typical exponents $E_i$4. Normalization ensures $E_i$5.

Gaussian Process Regression

Each local atomic energy is expanded over $E_i$6 "sparse" (inducing) environments $E_i$7 as

$E_i$8

with weights $E_i$9 found by solving the regularized linear system:

$i$0

where $i$1, $i$2 contains target observables (energies, forces, virials), and $i$3 encapsulates regularization parameters for different observables.

2. Data Preparation and Active Learning Workflow

Locally trained SOAP-GAP models require an appropriately curated and labeled dataset. The principal steps in database construction, as implemented for Ag-Pd alloys (Rosenbrock et al., 2019), are:

1. Initial Dataset Generation
   • Enumerate fcc- and bcc-based supercells up to 4 atoms for multiple compositions.
2. Preliminary MTP Fit
   • Fit a polynomial MTP to DFT-computed structures.
3. Active Learning Iteration
   • Relax each structure; high extrapolation grades trigger DFT calculation and inclusion in the training set.
4. Expansion
   • Extend to larger cells (up to 12 atoms), yielding 10,850 structures and 774 unique DFT-relaxed configurations.
5. High-Precision DFT Evaluation
   • All 774 training configurations re-evaluated at tight DFT settings (k-point density $i$4 Å$i$5, EDIFF=$i$6).
6. Locality Filtering
   • Restrict training references by "config_type" (e.g., "surf" for surface, or select regions in alloys), if fitting a local model (Klawohn et al., 2023).

No off-lattice or liquid data is included during training for solid-phase models.

3. Descriptor, Sparsification, and Kernel Construction

Descriptor construction and sparsification are the computational bottlenecks for SOAP-GAP models. The procedure is as follows:

• SOAP Descriptor Evaluation
  • For each atom in each frame, evaluate $i$7 within a cutoff $i$8 (e.g., $i$9–$$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$0 Å).
  • Descriptor dimensionality scales as $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$1, where $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$2 is the number of elements.
• Compression Techniques
  • Apply embedding (e.g., $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$3, $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$4) to map chemical and radial indices to lower $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$5-dimensional representation.
  • Use reduced forms such as $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$6, $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$7 schemes to achieve linear scaling with $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$8.
• Sparse Point Selection
  • Choose $$\rho_i^\alpha(\mathbf{r}) = \sum_{j \in \text{nbrs}} \delta_{z_j, \alpha} \exp\left( -\frac{\left|\mathbf{r} - \mathbf{r}_{ij}\right|^2}{2\sigma^2} \right) f_\text{cut}(|\mathbf{r}_{ij}|)$$9 representative environments (e.g., $z_j$0=500–2000) via CUR, $z_j$1-means, random, or uniform sampling.
  • Remove near-duplicates according to a jitter threshold.
• Kernel Assembly
  • Build kernel blocks $z_j$2 ($z_j$3) and $z_j$4 ($z_j$5).
  • Compute derivatives for force and stress learning.

4. Model Fitting, Regularization, and Validation

The fitting procedure specifically addresses the high-dimensional kernel regression while ensuring stability and generalization:

• Regularization
  • Distinct regularization ($z_j$6) per observable: typical settings are $z_j$7 eV/atom, $z_j$8 eV/Å, $z_j$9 eV/atom.
  • Additional small jitter ($\sigma$0) added to kernel diagonals.
• Parallel Fit and Hyperparameter Optimization
  • Distributed descriptor computation and linear algebra across $\sigma$1 MPI ranks. Use ScaLAPACK QR for scalable solves (Klawohn et al., 2023).
  • Hyperparameters (cutoff, basis resolution, kernel exponent, regularization strength) tuned by cross-validation on held-out frames (10–20% test split).
• Model Validation
  • Reported benchmarks (Rosenbrock et al., 2019):
  • On liquid MD test set (6000 snapshots): energy RMSE $\sigma$2 meV/atom, force RMSE $\sigma$3 meV/Å, virial RMSE $\sigma$4 meV/Å$\sigma$5.
  • Phonon spectrum integrated RMSE: $\sigma$6 THz, train/test parity.
  • Transition path "swap" test: SOAP-GAP correctly predicts physical barriers without spurious minima.
• Posterior Variance Diagnostics
  • The prediction variance, $\sigma$7, identifies out-of-sample areas where model uncertainty is large.

5. Implementation and Practical Usage

A locally trained SOAP-GAP model is typically implemented using the QUIP/GAP suite (Klawohn et al., 2023), with a workflow as follows:

• Configuration File Example:

$f_\text{cut}$9

• Fitting Invocation:

$g_n(r)$0

When focusing on a local region, ensure that configurations are tagged accordingly (e.g., config_type="surf").

• Hardware and Scaling:
  • Descriptor memory scales as $\sigma$8, with kernel memory as $\sigma$9.
  • Fitting $f_\text{cut}$0 with $f_\text{cut}$1 is feasible using MPI and descriptor compression.
  • Training is typically completed within hours on $f_\text{cut}$2 CPU cores.
• Post-Fitting Deployment:
  • Validated XML models are compatible with ASE, LAMMPS, and other atomistic simulation environments.

6. Computational Cost and Comparative Performance

SOAP-GAP offers DFT-like accuracy but at significantly reduced computational cost relative to ab initio methods:

Model	Typical RMSE (energy)	Relative Evaluation Speed	Suitable for
SOAP-GAP	$f_\text{cut}$3 meV/atom	10³–10⁴× faster than DFT	Off-lattice modeling
MTP	$f_\text{cut}$4–$f_\text{cut}$5 meV/atom	$f_\text{cut}$6 faster than GAP	Phase diagrams, large sampling

Evaluating a slice with $f_\text{cut}$7 billion potential calls is feasible with MTP but remains prohibitive with SOAP-GAP. GAP evaluations scale as $f_\text{cut}$8 and are several times slower per atom than MTP, but orders of magnitude faster than on-the-fly DFT (Rosenbrock et al., 2019).

7. Limitations and Domain of Applicability

Locally trained SOAP-GAP models, while excelling in local accuracy and transferability within their training domain, are limited in their extrapolation ability outside of that space. The exclusion of off-lattice and liquid configurations in training restricts the model’s generalization to high-temperature or disordered phases. Model quality is strongly dependent on choice and diversity of training data, regularization strength, and descriptor completeness. The computational requirements for fitting SOAP-GAP scale steeply with descriptor and sparse set size, motivating compression and parallelism techniques introduced in recent GAP frameworks (Klawohn et al., 2023).

A plausible implication is that, for applications requiring extensive sampling (e.g., nested sampling, large-scale Monte Carlo), alternative polynomial potentials such as MTP offer a better cost-precision trade-off, while SOAP-GAP is preferred for simulations where off-lattice or local accuracy is paramount, and model uncertainty needs quantification.