Machine-learning interatomic potentials with MACE¶

ACEsuit MACE is a pre-trained E(3)-equivariant graph-neural-network interatomic potential. It predicts energy, forces, and stress on a DFT-fitted potential-energy surface — orders of magnitude faster than an SCF, at close-to-DFT accuracy within its trained domain. vibe-qc exposes it as a first-class method="mace", plus direct drivers for periodic systems.

MACE is a learned potential, not a quantum-chemistry method: it solves no Schrödinger equation and produces no wavefunction. vibe-qc drives MACE’s pre-trained forward pass and attributes it as an external model — the energy is MACE’s, not a vibe-qc ab-initio total. This tutorial walks from a single molecule through geometry optimization, periodic crystals, molecules adsorbed on oxide and metal surfaces, and a comparison with vibe-qc’s own semi-empirical methods.

What you’ll learn: run MACE on molecules and crystals; optimize geometries and relax unit cells; compute adsorption energies of water and ammonia on MgO, CaO, Al₂O₃, and metal surfaces; choose a model and respect its license; and when to reach for MACE versus ab-initio or semi-empirical methods.

Installation¶

MACE is an optional dependency (a heavy PyTorch + e3nn stack), behind the [mace] extra, on Python ≤ 3.13:

pip install 'vibe-qc[mace]'

On Python 3.14 the runtime is import-gated with a clear message (MACE’s matscipy dependency has no 3.14 wheel yet). On macOS, vibe-qc auto-sets KMP_DUPLICATE_LIB_OK=TRUE (its native core and PyTorch each link an OpenMP runtime) — verified not to change results. Foundation-model weights download on first use into ~/.cache/mace/.

1. A first calculation — an organic molecule¶

method="mace" plugs into run_job exactly like an SCF method, minus the basis set (MACE needs none):

from pathlib import Path
from vibeqc import Molecule, run_job

# ethanol, C2H5OH
mol = Molecule.from_xyz("ethanol.xyz")
run_job(mol, method="mace", output="ethanol_mace")

ethanol_mace.out carries the energy, an MACE provenance block (model, license, training/theory, the energy-scale caveat), and the references; ethanol_mace.bibtex / .references cite the MACE method paper + the foundation-model paper. There is no .molden — MACE has no orbitals.

Energy scale. A MACE energy lives on a model-specific reference scale (each model subtracts its own per-element atomic energies), not a vibe-qc total and not comparable across models. For H₂O, MACE-MPA-0 gives ≈ −0.51 Ha; the organic MACE-OFF23 gives ≈ −76.5 Ha (near its ωB97M total). Use MACE energies for relative quantities — geometry optimization, reaction energies at fixed composition + model, adsorption energies, MD — never as absolute totals.

2. Geometry optimization¶

Add optimize=True — ASE’s BFGS drives the MACE calculator directly:

run_job(mol, method="mace", optimize=True, fmax=0.02, output="ethanol_opt")

You get ethanol_opt.traj (one frame per step) as usual. MACE’s geometries are close to experiment — relaxing water and ammonia:

molecule	bond	MACE-MPA-0	experiment
H₂O	O–H / H–O–H	0.970 Å / 103.9°	0.958 Å / 104.5°
NH₃	N–H / H–N–H	1.015 Å / 107.4°	1.012 Å / 106.7°

3. Choosing a model — and its license¶

The default is the MIT-licensed MACE-MPA-0 (materials, 89 elements). Pass MLIPOptions to select another model, device, or dtype:

from vibeqc.mlip import MLIPOptions

run_job(mol, method="mace",
        mlip_options=MLIPOptions(model="off23-medium",
                                 accept_academic_license=True),
        output="ethanol_off23")

model key	family	domain	license
`medium-mpa-0` (default), `small`, `medium`, `large`	MACE-MP / MPA	materials (89 elts)	MIT
`off23-small`, `off23-medium`, `off23-large`	MACE-OFF23	organic (H,C,N,O,F,P,S,Cl,Br,I)	ASL

The organic MACE-OFF23 weights are under the Academic Software License — academic, non-commercial use only. vibe-qc raises a PermissionError if you select an ASL model without acknowledging it (accept_academic_license=True or VIBEQC_ACCEPT_ASL=1). MIT models are never gated. Picking an organic model for an out-of-domain element (say a metal) raises a clear error pointing you at a materials model.

4. Periodic crystals¶

Periodic systems use direct drivers (MACE is not wired into run_periodic_job, which is SCF-only) — following the periodic semi-empirical precedent. run_periodic_mace returns energy, forces, and stress; optimize_periodic_mace_cell does a variable-cell relaxation:

import numpy as np
from ase.build import bulk
from vibeqc.ase_periodic import atoms_to_periodic_system
from vibeqc.mlip.mace import run_periodic_mace, optimize_periodic_mace_cell

system = atoms_to_periodic_system(bulk("Si", "diamond", a=5.43))

res = run_periodic_mace(system)
print(res.energy)                       # Ha (model reference scale)
print(np.diag(res.stress()) * 29421.0)  # GPa  (Ha/bohr³ → GPa)

relaxed = optimize_periodic_mace_cell(system, fmax=0.01)

For bulk silicon this relaxes the (slightly compressed) cell — the stress drops from −1.85 GPa to ≈ 0, the energy lowers, and the lattice constant moves from 5.43 Å to MACE’s ≈ 5.47 Å equilibrium. Forces are zero by the diamond symmetry, as they should be.

5. Molecules on surfaces — adsorption energies¶

This is where a fast, broadly-trained materials model shines: thousands of single-points per second let you screen adsorption across many surfaces. The recipe is the standard one — adsorption energy

E_ads = E(slab + molecule) − E(slab) − E(molecule)

all evaluated with the same MACE model so the per-element reference energies cancel. Here is a complete, runnable example for water on MgO(100):

import numpy as np
from ase.build import bulk, surface, molecule, add_adsorbate
from ase.optimize import BFGS
from ase.constraints import FixAtoms
from vibeqc.mlip.mace import mace_calculator

# one MACE-MPA-0 calculator, reused for every structure
calc = mace_calculator()

def energy(atoms, relax=True):
    a = atoms.copy(); a.calc = calc
    if relax:
        z = a.positions[:, 2]
        a.set_constraint(FixAtoms(mask=z < z.min() + 1.2))   # freeze the bottom
        BFGS(a, logfile=None).run(fmax=0.15, steps=80)
    return a.get_potential_energy()                          # eV

# gas-phase water in a box
h2o = molecule("H2O"); h2o.set_cell([14, 14, 14]); h2o.center(); h2o.pbc = True
e_gas = energy(h2o)

# MgO(100) slab (3 layers, 2×2), relaxed
slab = surface(bulk("MgO", "rocksalt", a=4.212), (1, 0, 0), 3, vacuum=8.0) * (2, 2, 1)
e_slab = energy(slab)

# water on top of a surface Mg, relaxed
ads = slab.copy()
mg = max((i for i, z in enumerate(ads.numbers) if z == 12), key=lambda i: ads.positions[i, 2])
add_adsorbate(ads, molecule("H2O"), height=2.1, position=ads.positions[mg, :2], mol_index=0)
e_total = energy(ads)

print(f"E_ads(H2O/MgO) = {e_total - e_slab - e_gas:.3f} eV")   # ≈ -0.81 eV

Running the same recipe for water and ammonia across a range of oxide and metal surfaces (MACE-MPA-0, relaxed, bottom layer fixed):

surface	atoms	E_ads(H₂O) / eV	E_ads(NH₃) / eV
MgO(100)	24	−0.81	−1.16
CaO(100)	24	−3.18 †	−1.35
Al₂O₃(0001)	60	−1.50	−1.64
TiO₂(110)	36	−0.66	−1.08
Pt(111)	27	−0.37	−1.42
Cu(111)	27	−0.11	−0.45

The trends are physically sensible: ammonia (a stronger Lewis base) binds more strongly than water on every surface except CaO. The standout is † water on CaO(100): it binds far more strongly (≈ −3.2 eV) than on MgO because CaO is much more basic — in the relaxed structure one O–H points down and forms a very short hydrogen bond (≈ 1.5 Å) to a surface oxygen, the onset of the dissociative surface-hydroxyl chemistry CaO is known for. (This is a stable, reproducible MACE result — ≈ −3.16 eV across starting heights of 1.8–2.5 Å — not a stray local minimum.) The less basic MgO binds water moderately (≈ 0.8 eV); the more reactive Al₂O₃ binds both strongly; on the noble-ish Cu(111) water barely physisorbs (≈ 0.1 eV) while it adsorbs more strongly on Pt(111). A fixed-charge force field would miss the MgO < CaO basicity trend that MACE reproduces here. A handful of these single-points take seconds — the same screen at DFT would be hours to days.

These are illustrative numbers from a quick relaxation (loose fmax, thin slabs, a single adsorption site). For publication you would converge slab thickness, k-points-equivalent supercell size, the adsorption site, and fmax — but the workflow is exactly the one above. On strongly basic surfaces like CaO, check whether the water stayed molecular or dissociated (inspect the relaxed O–H distances) — the two are different physical states with different energies. The full script is in examples/mlip/04_surface_adsorption.py.

6. MACE versus vibe-qc’s semi-empirical methods¶

vibe-qc ships its own semi-empirical platform (DFTB, GFN2-xTB, PM6, OMx — see semiempirical). Both are fast, non-ab-initio energy engines, but they sit in different places:

	MACE (this tutorial)	semi-empirical (vibe-qc’s own)
nature	learned potential (pre-trained GNN)	physics-based, vibe-qc’s own code
accuracy (in domain)	close to the reference DFT	lower; qualitative–semi-quantitative
electronic structure	none (no orbitals/charges/gap)	yes — charges, orbitals, HOMO–LUMO
charge / spin states	ignored (neutral potential)	represented
elements	trained set (10 organic / 89 materials)	any parameterised element
transferability	in-distribution only	broad (physics-based)
typical use	accurate energies/forces/MD, screening	electronic properties, any element

Concretely, MACE reproduces the water and ammonia geometries above to within ~0.01 Å and ~1° of experiment — materially closer than typical semi-empirical methods (the MACE-OFF23 paper reports MACE-OFF beating both GFN2-xTB and ANI-1ccx on organic benchmarks). The trade-off is that semi-empirical methods give you a wavefunction (Mulliken charges, a gap, frontier orbitals) and work for any parameterised element and charge state, where MACE gives none of those and only within its training domain. A detailed, numbers-backed comparison lives in semiempirical_mlip_comparison.

Rule of thumb: reach for MACE when you want DFT-quality energies/forces/geometries/MD on systems inside its training domain and don’t need electronic structure; reach for semi-empirical (or ab-initio) when you need charges/orbitals/gaps, unusual elements or charge states, or a physically transferable model.

7. Visualization¶

A MACE run writes the same structural output as any vibe-qc job, so geometry and (with optimize=True) the relaxation trajectory visualize out of the box: moltui ethanol_opt.xyz (terminal), or pass output_qvf=True and open the .qvf archive in vibe-view for the 3-D structure + trajectory animation + citations. There is no orbital/density/band visualization — there is no wavefunction to plot. See mlip § Visualization.

8. Limitations + good practice¶

No electronic structure — no charges, orbitals, gaps, spectra.
Trained-domain only — accuracy degrades out of distribution; MACE cannot tell you when you’ve left its domain. Sanity-check against DFT for unusual chemistry.
Charge / spin ignored — vibe-qc warns if you pass a non-neutral or open-shell system; the energy is the neutral value regardless.
Energy scale — model-specific reference, never an absolute total or comparable across models.
Attribution — every run cites the MACE papers; cite them in published work (CLAUDE.md § 10: vibe-qc drives MACE, it does not claim MACE’s energy as its own).