SPARTA

Steric Path Alignment for Reactive Transformation Assessment

SPARTA generates reactant-product conformer pairs and iteratively aligns them using an Image-Dependent Pair Potential (IDPP) nudged-elastic-band (NEB) path optimization.

Path generation and alignment are purely geometric and need no electronic-structure program. Evaluating electronic energies along the aligned path gives a barrier estimate that acts as a fast proxy for kinetic feasibility, at a small fraction of the cost of locating a transition state. SPARTA will use any calculator you give it; see Electronic energies.

Installation

pip install sparta-chem

The import name is sparta:

from sparta import Reaction

Python API

Quick start

from sparta import Reaction, XTBCalculation, run_sparta

rxn = Reaction(
    reactant_smiles="[C:1]([H:2])([H:3])([H:4])[H:5].[F:6][F:7]",
    product_smiles="[C:1]([H:2])([H:3])([H:4])[F:6].[H:5][F:7]",
)

rxn = run_sparta(rxn, calculator=XTBCalculation())
print(rxn.barrier_kcal)          # SPARTA barrier, kcal/mol

run_sparta is the whole pipeline in one call. It scans stretch factors from 1.0 to 2.0, aligns one endpoint pair at each, evaluates electronic energies along the path, and keeps the factor giving the lowest barrier.

Both SMILES must be fully atom-mapped, including hydrogens, and the two mappings must cover the same set of map numbers, that correspondence is what defines which bonds change. Hydrogens have to be written as explicit mapped atoms ([C:1]([H:2])...) rather than left implicit ([CH4:1]): SPARTA parses with removeHs=False, so an implicit hydrogen is not an atom and will be missing from the generated geometry.

Geometric scoring only

Without a calculator nothing electronic is computed, and the result carries the IDPP score alone. This is the cheaper, purely geometric mode:

from sparta import Reaction, generate_initial_conformers, run_iterative_alignment

rxn = Reaction(reactant_smiles=..., product_smiles=..., stretch_factors=[1.3, 1.5, 1.8])
rxn = generate_initial_conformers(rxn)
rxn = run_iterative_alignment(rxn)

print(rxn.aligned_results[0]["score"])     # IDPP score, not an energy

Each function returns the reaction, so these also chain: Reaction(...).generate_initial_conformers().run_iterative_alignment().

Accessing results

# aligned_results is sorted by score (ascending, lower is better)
best = rxn.aligned_results[0]

print(best["score"])            # float: max interior IDPP energy
print(best["reactant_coords"])  # np.ndarray (n_atoms, 3)
print(best["product_coords"])   # np.ndarray (n_atoms, 3)
print(best["images"])           # list of ase.Atoms, full NEB path
print(best["energies"])         # list of float, per-image IDPP energies
print(best["seed"])             # int
print(best["stretch_factor"])   # float
print(best["conformer_idx"])    # int: index into rxn.top_k_conformers

# present only when a calculator was supplied
print(best["barrier_kcal"])         # SPARTA barrier, kcal/mol
print(best["electronic_energies"])  # list of float, per interior image
print(best["peak_image_idx"])       # index of the highest-energy image
print(rxn.barrier_kcal)             # same as aligned_results[0]["barrier_kcal"]
print(rxn.reference_energy)         # relaxed reactant reference, Hartree

Two different quantities live here and they are not interchangeable:

aligned_results is sorted by barrier_kcal when a calculator was given, and by score otherwise, so aligned_results[0] is always the best result by whichever metric was available.

The reactant reference relaxes each connected fragment separately and sums the energies. Optimizing a multi-fragment reactant as one unit would fold an arbitrary inter-fragment interaction into every barrier.

All intermediate conformers are also available:

rxn.all_conformers      # every generated conformer pair, sorted by initial NEB score
rxn.top_k_conformers    # the top-k subset passed to iterative alignment

Electronic energies

The IDPP score is geometric. To obtain an electronic barrier, evaluate single-point energies on the images of the aligned path with a calculator of your choice.

SPARTA ships an XTBCalculation wrapper for convenience, which shells out to the xtb executable. xtb is not a Python dependency and is not installed by pip, install it separately and make sure it is on PATH, or pass an explicit path via xtb_cmd. Nothing else in SPARTA requires it.

from sparta.calculation import XTBCalculation

calc = XTBCalculation(charge=0, multiplicity=1)       # or xtb_cmd="/path/to/xtb"
# XTBCalculation also takes `timeout` (seconds, default 120)
energies = [calc.energy(img) for img in best["images"]]

ASECalculation wraps any ASE calculator, which covers most machine-learned potentials:

from sparta.calculation import ASECalculation

calc = ASECalculation(ase_calculator)                 # also takes `fmax`, `max_steps`

Any object exposing the same energy(atoms) interface works in their place, so the path code never needs to know what is underneath.

Reaction parameters

generate_initial_conformers parameters

run_iterative_alignment parameters

Choosing a stretch factor automatically

A single stretch factor is not optimal for every reaction. find_optimal_stretch_factor scans upward until the score stops improving and leaves the best result on the reaction object:

from sparta import find_optimal_stretch_factor

rxn = find_optimal_stretch_factor(rxn, start=1.0, step=0.1, max_stretch=3.0)

It also accepts every alignment parameter of run_iterative_alignment.

Getting the barrier for an aligned reaction

If a reaction has already been aligned, sparta_barrier returns its barrier in kcal/mol, computing the electronic energies only if they are not already there:

from sparta import sparta_barrier

barrier = sparta_barrier(rxn, calculator=XTBCalculation())

Reproducing the published results

The defaults of run_sparta are the published settings. Stated explicitly:

from sparta import Reaction, XTBCalculation, run_sparta

rxn = Reaction(
    reactant_smiles=...,
    product_smiles=...,
    n_conformers=30,          # per stretch factor
    n_seeds=1,
)

rxn = run_sparta(
    rxn,
    calculator=XTBCalculation(),
    stretch_factors=None,     # scan 1.0 to 2.0 in steps of 0.1
    patience=3,               # stop after 3 factors without improvement
)

run_sparta sets top_k=1, so one endpoint pair per stretch factor carries through to alignment. The reported barriers used g-xTB; XTBCalculation invokes whichever xtb build is on PATH, and GFN2 and g-xTB give different numbers for the same reaction.

To evaluate a fixed set of factors instead of scanning, pass them: run_sparta(rxn, calc, stretch_factors=[1.0, 1.1, ..., 2.0]).

Command line

Installing the package provides a sparta command that processes a CSV of reactions and writes a results CSV.

Input file format

A CSV with at minimum these two columns:

reactant_smiles,product_smiles
[C:1]([H:2])([H:3])([H:4])[H:5].[F:6][F:7],[C:1]([H:2])([H:3])([H:4])[F:6].[H:5][F:7]
[C:1]([H:3])([H:4])=[C:2]([H:5])[H:6].[H:7][Cl:8],[C:1]([H:3])([H:4])([H:7])[C:2]([H:5])([H:6])[Cl:8]

As above, every atom carries a map number and hydrogens are explicit.

Additional columns are preserved in the output.

Basic usage

# Defaults
sparta reactions.csv

# Specify output file
sparta reactions.csv -o results.csv

# Override hyperparameters via flags
sparta reactions.csv \
    --n-conformers 10 --n-seeds 3 \
    --stretch-factors 1.3 1.5 1.8 \
    --top-k 5 --idpp-objective peak

# Load hyperparameters from a JSON config
sparta reactions.csv --config config.json

# Compute SPARTA barriers rather than geometric scores alone
sparta reactions.csv --calculator xtb

# Save the best aligned conformer pair as XYZ files
sparta reactions.csv --save-conformers --conformer-dir ./conformers

Config file format

CLI flags override config file values. All keys are optional.

{
    "n_conformers": 10,
    "n_seeds": 3,
    "stretch_factors": [1.3, 1.5, 1.8],
    "top_k": 5,
    "n_images": 15,
    "neb_fmax": 0.001,
    "neb_steps": 200,
    "propagation_threshold": 0.0001,
    "propagation_step_size": 2.0,
    "max_iterations": 40,
    "objective": "peak",
    "d_power": 4
}

Output

The output CSV contains all original columns plus:

All CLI flags

positional arguments:
  input                       CSV file with reactant_smiles and product_smiles columns

options:
  -o, --output FILE           Output CSV (default: <input>_results.csv)
  -c, --config FILE           JSON config file for hyperparameters

  --n-conformers INT          Conformers per (seed, stretch_factor) [default: 10]
  --n-seeds INT               Number of random seeds [default: 1]
  --stretch-factors SF [...]  Stretch factor(s) [default: 1.5]
  --top-k INT                 Top-k conformers to align [default: 5]
  --n-images INT              NEB images [default: 15]

  --neb-fmax FLOAT            NEB force convergence [default: 0.001]
  --neb-steps INT             Max NEB steps per run [default: 200]
  --propagation-threshold F   Endpoint force convergence [default: 0.0001]
  --propagation-step-size F   Initial endpoint step size Å [default: 2.0]
  --max-iterations INT        Max NEB+propagation iterations [default: 40]
  --d-power INT               IDPP distance exponent [default: 4]

  --save-conformers           Save best XYZ conformer pair per reaction
  --conformer-dir DIR         Directory for XYZ files [default: conformers]
  --verbose                   Print per-conformer optimization progress

License

MIT, see LICENSE.