Skip to content

Chemistry file formats — mechanisms & species database

aeolus1d reads chemistry from three on-disk artefacts. This page is the authoritative schema for each.

Which one do I want?

I want to… Use File
Define custom global / semi-global mass-action kinetics (gas, or gas + surface coverage) [pipe.*.species] reaction = "mechanism" → mechanism_file a TOML MechanismSpec (§1)
Import a published elementary surface mechanism [pipe.*.catalyst] kind = "mechanism" → mechanism_path a CHEMKIN .inp or Cantera .yaml file (§2)
Supply per-species thermodynamics (NASA-7 cp, molar mass, composition) the species thermo database one .npz per species (§3)

The two mechanism systems are independent — different TOML fields, different file formats, different solver seams. Don't confuse mechanism_file (§1) with mechanism_path (§2).

Background on how these plug into the solver: Theory §7 Strang-split chemistry, §3.5 Species transport, and the Component index — species & reactions.

1. Declarative mechanism spec (mechanism_file, TOML)

The lightweight format: a mass-action mechanism you write by hand. A pipe opts in through its species block:

[[pipe]]
id = "doc"
# … geometry, mesh …

[pipe.species]
name    = "CO"                 # primary species  = mechanism species[0]
extras  = ["O2", "CO2", "C3H6", "H2O"]   # = mechanism species[1:]
reaction = "mechanism"
mechanism_file = "mechanisms/co_hc_oxidation.toml"

A relative mechanism_file resolves against the case file's directory (the example above is a case in examples/). Since 2026-09-29; before, it resolved against the working directory. A path found only there still loads, with a FutureWarning naming both places, and that fallback will be removed. An absolute path is used as written.

Alignment rule (enforced at build time): the mechanism's species list must equal [primary, *extras] in the same order — otherwise the build raises a ValueError. This keeps the reaction columns aligned with the pipe's transported-species layout under every driver. The file is read with tomllib, so mechanism_file is TOML only (the CHEMKIN/Cantera path in §2 is the one that takes .inp/.yaml).

1.1 Gas-phase example

examples/mechanisms/co_hc_oxidation.toml — DOC-style CO + C₃H₆ oxidation:

species = ["CO", "O2", "CO2", "C3H6", "H2O"]

[[reactions]]
# CO + 1/2 O2 -> CO2
stoich = { CO = -1.0, O2 = -0.5, CO2 = 1.0 }
A  = 1.0e10
Ea = 9.0e4
dH = 2.83e5

[[reactions]]
# C3H6 + 9/2 O2 -> 3 CO2 + 3 H2O
stoich = { C3H6 = -1.0, O2 = -4.5, CO2 = 3.0, H2O = 3.0 }
A  = 2.5e9
Ea = 1.3e5
dH = 1.93e6

1.2 Schema (MechanismSpec)

Top-level keys (parsed by aeolus1d.closures.mechanism.MechanismSpec.from_dict):

key type required meaning
species list[str] yes gas species, canonical order; [0] is the pipe's primary
reactions array of tables yes one [[reactions]] block per reaction (below)
molar_mass table {species = kg/mol} no overrides; common species default from the builtin table, so usually omitted
surface_species list[str] no coverage variables θ_k ∈ [0,1] (see §1.3)
site_capacity table {θ = mol-sites/m³} no site density Ω_k per coverage
theta_init table {θ = value} no initial coverage, default 0

Each [[reactions]] table:

key type required meaning
stoich inline table yes stoichiometry; negative = consumed, positive = produced. Gas and surface (θ_k) species allowed
A float yes Arrhenius pre-exponential. Mass-action units — they depend on the reaction order (e.g. m³·mol⁻¹·s⁻¹ for a single first-order reactant)
Ea float yes activation energy [J/mol]
dH float no (def 0) heat of reaction [J per mol of extent]; positive = exothermic
orders inline table no override reaction orders; default is first order in each consumed species (occupied/concentration order for surfaces)
vacancy inline table no (1−θ_k) orders for surface sites (site-availability dependence)
A_rev, Ea_rev float no (def 0) reverse rate; A_rev = 0 ⇒ irreversible
orders_rev, vacancy_rev inline table no reverse-direction orders (default: first order in each produced gas species)

1.3 Surface coverage

Add surface_species to model adsorption/desorption. A coverage θ_k carries an occupied concentration Ω_k·θ_k and a vacant concentration Ω_k·(1−θ_k); a surface stoich coefficient is the coverage change per unit extent (dθ = extent / Ω). orders give the occupied (θ) order, vacancy the (1−θ) order.

examples/mechanisms/nh3_ads_des_surface.toml — Langmuir NH₃ storage (reproduces the bundled nh3_ads_des closure):

species = ["NH3"]
surface_species = ["theta_NH3"]

[site_capacity]
theta_NH3 = 100.0          # mol-sites / m^3 (typical Cu-zeolite)

[theta_init]
theta_NH3 = 0.0

[[reactions]]
# adsorption: NH3(gas) + free site -> NH3*   rate = k·c_NH3·(1-θ)·Ω
stoich  = { NH3 = -1.0, theta_NH3 = 1.0 }
orders  = { NH3 = 1.0 }
vacancy = { theta_NH3 = 1.0 }
A  = 1.0e3                  # m^3 / (mol-site . s); non-activated
Ea = 0.0
dH = 8.0e4                  # +exothermic on adsorption

[[reactions]]
# desorption: NH3* -> NH3(gas) + free site   rate = k·θ·Ω
stoich = { NH3 = 1.0, theta_NH3 = -1.0 }
orders = { theta_NH3 = 1.0 }
A  = 1.0e10
Ea = 8.0e4
dH = -8.0e4                 # desorption is endothermic

When surface_species is present the mechanism compiles to a surface kinetics object (py_reaction_surface) and each coverage is seeded into the pipe's surface_state from theta_init; a pure-gas spec compiles to py_reaction_multi. Working end-to-end case: examples/35_mechanism_pipe.toml (plot: examples/plots/plot_35_mechanism_pipe.py).

GIL-free integration (opt-in)

The compiled surface mechanism (py_reaction_surface) is the one reaction seam wired to the optional GIL-free chemistry backend: set [case] chem_solver = "numbalsoda" (requires the optional numbalsoda package) to integrate each cell's stiff mass-action ODE through numbalsoda's nogil LSODA in a prange instead of the batched SciPy LSODA. The per-cell right-hand side is the same kernel (aeolus1d.closures.mechanism._mech_cell_dstate_jit), so the numerics agree to integrator tolerance; any pipe it can't handle falls back to SciPy. Bespoke catalyst closures (§2, Brinkmeier/Voltz/SCR) are not on this path.

2. Elementary mechanism import (mechanism_path, CHEMKIN / Cantera)

For detailed mechanisms, the catalyst closure imports a CHEMKIN or Cantera file and builds surface kinetics over the brick:

[pipe.scr.catalyst]
kind = "mechanism"
mechanism_path = "mechanisms/scr_surface.inp"   # .inp | .yaml | .yml
T_solid_init = 300.0

[pipe.scr.catalyst.surface_init]
theta_NH3 = 0.0
theta_NO  = 0.0
field type meaning
mechanism_path str path to the mechanism file; extension selects the importer
surface_init table {θ = value} initial surface coverage per species
T_solid_init float (def 300.0) initial substrate temperature [K]

The extension dispatches the importer (aeolus1d/io/build.py → _attach_mechanism_catalyst):

  • .inp → aeolus1d.chemistry.importers.chemkin_mech.import_chemkin_mech — standard CHEMKIN-II/III ELEMENTS / SPECIES / REACTIONS sections, plus SITE blocks for surface species (CHEMKIN-Surface). Reaction lines carry A b Ea; the REACTIONS header unit token (CAL/MOLE, KCAL/MOLE, JOULES/MOLE, KJOULES/MOLE, KELVINS) sets Ea units (converted to J/mol on import). =/<=> are reversible, => irreversible; the non-standard <= is read as a right-to-left irreversible reaction (sides swapped on import). Rate coefficients accept a Fortran D/d exponent (1.0D13 → 1e13).
  • .yaml / .yml → aeolus1d.chemistry.importers.cantera_yaml.import_cantera_yaml_as_mechanism — a subset of Cantera ≥2.5 YAML: a units block, species with NASA-7 thermo, and a reactions list with rate-constant: {A, b, Ea}.

Both importers return the Mechanism dataclass (aeolus1d.chemistry.reaction.Mechanism: elements, species, reactions, surface_sites, …). Advanced kinetics (third-body / (+M) falloff, LOW/TROE/SRI, PLOG, DUPLICATE) are recognised where supported; Chebyshev rates and non-NASA-7 thermo are not. The importer modules in src/aeolus1d/chemistry/importers/ are the authoritative spec for the exact accepted subset — consult them before relying on an exotic feature.

3. Species thermodynamic database (.npz)

Per-species NASA-7 thermodynamics live as one NumPy .npz per species under src/aeolus1d/chemistry/db_data/. These supply molar mass, atomic composition, and cp(T)/h(T)/s(T) for the imported mechanisms in §2.

Bundled species (14): CH4 CO CO2 H H2 H2O H2O2 HO2 N2 NH3 NO O O2 OH (GRI-Mech 3.0 / Burcat-Goos thermo).

3.1 Per-file array schema

Each <name>.npz (verified against NH3.npz) contains:

array shape dtype meaning
nasa7_low (7,) float64 NASA-7 coefficients a₁…a₇, low-T range
nasa7_high (7,) float64 NASA-7 coefficients, high-T range
T_min () float64 low end of validity [K]
T_mid () float64 low/high polynomial switch [K]
T_max () float64 high end of validity [K]
M_w () float64 molar mass [g/mol] (= kg/kmol; e.g. NH₃ → 17.03056)
elements_keys (n,) <U str element symbols, e.g. ['N','H']
elements_counts (n,) int64 per-element atom counts, e.g. [1,3]
phase () <U str phase tag, e.g. 'G' (gas)

In memory this maps to the frozen SpeciesThermo dataclass (aeolus1d.chemistry.species.SpeciesThermo).

Pickle-free by policy

The .npz files store plain typed arrays only — no pickled Python objects — so they load with allow_pickle=False (a security gate from the 2026-06-10 audit). Keep any species DB you author in the same plain-typed form.

3.2 Resolving & adding species

from aeolus1d.chemistry.db import load_species, save_species

nh3 = load_species("NH3")                 # bundled db_data/
mine = load_species("C2H4", db_dir="my_db")  # user directory
save_species(species_obj, db_dir="my_db")    # write my_db/C2H4.npz

load_species(name, db_dir=None) reads <db_dir>/<name>.npz, falling back to the bundled db_data/ when db_dir is None. New species come from a thermo importer (aeolus1d.chemistry.importers.nasa7 / janaf) then save_species. To re-emit the bundled DB in the pickle-free format, run python scripts/regen_species_db.py [db_dir].

Burcat's database (~3700 species, fetched on demand)

The bundled seed covers ten common species. When a mechanism names something else, Burcat & Ruscic's Third Millennium database is the usual source — it carries NASA-7 polynomials for roughly 3700 species, with Active Thermochemical Tables updates.

It is not bundled. Burcat's terms are "free of charge on condition that proper quotation is given" and "the sale of this information is strictly forbidden and it cannot be listed, printed or copied in any commercial publication form". aeolus1d is MIT licensed, which grants commercial use, so shipping the file inside the wheel would hand users a bundle whose parts disagree about what they may do with it. It is fetched to a user-local cache instead — the same treatment local_cases/ gets for licensed vendor data.

aeolus1d thermo-db fetch     # download once (~2.3 MB) into the cache
aeolus1d thermo-db info      # cache location, species count, citation

The cache defaults to ~/.local/share/aeolus1d/burcat (honouring XDG_DATA_HOME), overridable with AEOLUS_BURCAT_DIR. It sits deliberately outside the installed package so licensed data never lands where git add can reach it.

from aeolus1d.chemistry.burcat_db import ensure_burcat_db, load_burcat_species

ensure_burcat_db()                     # no-op once cached
c2h4 = load_burcat_species("C2H4")     # SpeciesThermo

Species are keyed on the identifier alone (N2), with Burcat's full label (N2 REF ELEMENT) registered as a second key — some species differ only by that qualifier, e.g. CH4 RRHO versus CH4 ANHARMONIC.

Parsing uses aeolus1d.chemistry.importers.burcat, not the CHEMKIN nasa7 importer. The file is NASA-7-like but differs on the first line of every record — float element counts (N 2.), the molecular weight where T_common belongs (so the split temperature is the database's 1000 K convention), a fifth number on line 4, CRLF endings, and no THERMO/END wrapper. The CHEMKIN parser recovers exactly one species from it.

If you publish results that use it, cite:

A. Burcat and B. Ruscic, "Third Millennium Ideal Gas and Condensed Phase Thermochemical Database for Combustion with Updates from Active Thermochemical Tables", ANL-05/20 / TAE 960 (2005).

See also