Move the surface-chemistry Fortran codegen upstream into pyrometheus

[sbryngelson]

Follow-up to #1821, which added heterogeneous surface chemistry. Not a defect in that PR — a structural observation about where its code ended up.

The asymmetry

For the gas phase, MFC's toolchain generates Fortran in one line (toolchain/mfc/run/input.py:516):

thermochem_code = pyro.FortranCodeGenerator().generate("m_thermochem", sol, ...)

For the surface phase, #1821 added ~400 lines of hand-written Fortran emitted by string concatenation (generate_surface_thermochem, input.py:141-560). It produces get_surface_net_production_rates and get_surface_reaction_heat_flux — the direct analogues of pyrometheus's own get_net_production_rates.

That is a code generator for Cantera mechanisms living in a CFD solver's build toolchain.

Why this matters beyond tidiness

nongas_h_rt_expressions (input.py:240-264) takes a Cantera species' coeffs, splits it at tmid into low/high ranges, and evaluates the NASA 7-coefficient enthalpy polynomial — reimplementing what pyrometheus already generates as get_species_enthalpies_rt, only for bulk/surface species, because pyrometheus does not know those exist.

Two bugs already came out of that duplicated surface: a hasattr(rate, "pre_exponential_factor") guard that silently accepted sticking-coefficient rates and emitted them as plain Arrhenius (wrong by orders of magnitude — a sticking coefficient is a dimensionless probability, not a prefactor), and silently-dropped coverage dependencies. Both are fixed, but they are the failure mode of maintaining a second, partial chemistry code generator.

What would move upstream

• ct.Interface loading and phase resolution
• Interface-Arrhenius rate expressions with explicit reaction orders
• Sticking-coefficient conversion: gamma/(1 - gamma/2) * sqrt(RT/(2*pi*W)) / site_density^n
• Blowers-Masel rates (activation energy shifted by reaction enthalpy at runtime)
• Coverage-dependent rates
• NASA thermo for bulk/surface species

I checked pyrometheus (installed version): no surface support — ct.Interface, InterfaceKinetics, coverage and sticking have zero hits outside an unrelated flamelets/ module. So this is a genuine upstream feature, not a port.

It also belongs there on the merits: sticking-coefficient conversion and Blowers-Masel shifts are properties of a Cantera mechanism, not of MFC. Any pyrometheus consumer doing catalysis or ablation wants them.

Payoff

toolchain/mfc/run/input.py loses ~400 lines and goes back to delegating both phases. The surface generator gets pyrometheus's own test suite instead of MFC's integration tests as its only coverage.

Sequencing: this should land upstream before the MFC-side splits (#TBD), since it deletes the most code.

---

Tracked in #1893.

[sbryngelson]

The upstream half of this is now open: pyrometheus/pyrometheus#118.

It adds generate_surface to all three backends, so a ct.Interface generates code the same way a ct.Solution does:

FortranCodeGenerator.generate("thermochem", sol)                  # gas, unchanged
FortranCodeGenerator.generate_surface("surface_thermochem", interface,
                                      gas_module_name="thermochem")

The surface module holds only the surface mechanism and refers to the gas one, which is the shape MFC already consumes — so generate_surface_thermochem in toolchain/mfc/run/input.py can be replaced by a call rather than rewritten.

It covers everything the hand-written generator does and the parts it rejects: interface-Arrhenius, sticking coefficients with the kinetic-theory conversion and Motz-Wise correction, coverage-dependent rates, reversible surface reactions, bulk phases, explicit reaction orders, and ConstantCp thermo. Verified against Cantera to ≤1.5e-14 on two mechanisms — Cantera's ptcombust and this repo's own carbon mechanism — with the Fortran and C++ numbers taken from compiled binaries.

Worth noting for whoever does the MFC-side swap: the two bugs found in the hand-written generator during the #1821 review — sticking coefficients silently emitted as Arrhenius, and coverage dependence silently dropped — are both handled upstream rather than rejected, so the swap removes code and lifts restrictions.

Adding this repo's carbon mechanism as a second fixture also turned up an unrelated pyrometheus bug, fixed in pyrometheus/pyrometheus#120: generating from a mechanism with no reactions produced code that raised IndexError when called.

Blocked on #118 landing. Tracked in #1893.