Skip to content

Phase log

Detailed per-phase technical history of every milestone aeolus1d has shipped, lifted from the original README to keep the README itself focused on capabilities and on-ramp. Read top-to-bottom for the canonical chronology, or jump to the table of contents the heading anchors provide.

For a higher-level summary by milestone, see the Changelog. For what's planned next, see the Roadmap.


Status

Milestones 1–6 shipped + the Phase 25 engine fuel chain (25k → 25q + 25o), Phase 26 case-level species metadata, Phase 26.5 component protocols, audio probe / WAV export, parallel case runner. Plus Milestone 7 (largely shipped) through 2026-05-06: alternative- fuel kinetics, full aftertreatment depth (DPF active regen + controller, SCRF, two-stage SCR, urea injector), VGT 3-D map + boost controller, Livengood-Wu knock + closed-loop retard + TOML support, FMI export scaffolding, cycle-resolved post- processing, PINN closure adapter. 1380+ passing tests (933 solver + 455 GUI), 1 known xfail (Helmholtz resonance — see milestone 3 below), ruff check clean.

Recent polish (2026-05-14) — * Phase 7c-doc-S — DOC sulfur-poisoning model. New aeolus1d.closures.doc_sulfur module with two pieces: DOCSulfurOxidation (a stateful drop-in replacement for the base make_doc_oxidation closure with a multiplicative deactivation factor f_S = 1 − alpha_S · theta_S on both HC and CO branches) and DOCSulfurController (an on_step callable that integrates the lumped sulfate balance against mean DOC gas temperature using Arrhenius storage (low Ea — proceeds at exhaust T) and reduction (high Ea — meaningful only above ~870 K). Wired through TOML as species.reaction = "doc_oxidation_S" plus [[controller]] kind = "doc_sulfur"; the validator pins the controller-to-closure binding. Bundled examples/61_doc_sulfur_poisoning.toml + examples/plots/plot_61_doc_sulfur.py exercise a storage + active-regen cycle: 1 s of high-S exhaust at 550 K takes theta_S from 0 → 0.72 (50 % activity loss); a follow-up 1 s at 1100 K with clean fuel takes it back from 0.72 → 0.02 (full re-activation). Closes the DOC sulfur Phase 7c roadmap item. * Phase 7c-surge — compressor surge / choke regime diagnostics. Compressor now exposes last_regime (one of "normal" / "idle" / "backflow" / "surge" / "choke") plus normalised last_surge_margin / last_choke_margin distances to the top / bottom edges of the map's PR grid (positive when inside the envelope, ≤ 0 at or past the edge). Replaces the silent "wall both ends" off-design fall-through with informative state — controllers and post-processors can read the flag to detect when an operating point has hit a surge-line or saturated at choke. The wall-close behaviour is unchanged (existing callers don't regress); the new fields are pure diagnostics. Closes the surge/choke Phase 7c roadmap item per docs/roadmap.md. * Phase twin-turbo — sequential / two-stage turbocharger topology template. New bundled examples/60_twin_sequential_turbo.toml + examples/plots/plot_60_twin_sequential_turbo.py that exercises HP and LP turbo-shaft pairs plumbed in series (intake_atm → LP_compressor → intake_mid → HP_compressor → intake_manifold; exhaust cascades the other way through HP_turbine → LP_turbine). Two independent shafts with shaft_HP carrying smaller inertia for transient response. The kernels (Turbine / Compressor / Shaft) already supported N turbomachines + N shafts before this phase; the work was the bundled topology template + a regression test (tests/integration/test_twin_sequential_turbo.py) that pins the cascade direction as a load-bearing contract. Switching to a parallel-twin pattern (two compressors feeding the same manifold) is the same TOML pattern with different pipe wiring — left as a user-extension exercise. Closes the twin/sequential-turbo Phase 7c roadmap item. * Phase 7c-wastegate — WastegateController. Velocity-form PI controller in aeolus1d.control that drives Turbine.bypass_open (the parallel-bypass fraction the Turbine has carried since Phase 7c) to maintain a target manifold pressure. Sign-inverted vs BoostController — when boost overshoots the target the controller OPENS the wastegate to bleed exhaust around the turbine; when underboost it closes. [[controller]] kind = "wastegate" exposes the controller in TOML with turbine_id / manifold_pipe_id / p_target / Kp / Ki / bypass_min / bypass_max fields, validated by _validate_controller_refs. Bundled examples/59_wastegate_boost_ctrl.toml + examples/plots/plot_59_wastegate.py exercise an overboost scenario (manifold seeded at 2.4 bar, target 1.8 bar) — the wastegate opens from 0.10 toward bypass_max as the integral term accumulates. Closes the last open Phase 7c item per docs/roadmap.md. * Phase Perf-7 — packed-kernel parallel-vs-serial dispatcher. Empirical full-app benches showed the existing @njit(parallel=True) packed kernels with default NUMBA_NUM_THREADS = cpu_count (16 on the audit laptop) oversubscribe the small-P workload at typical aeolus1d network sizes. Engine-case measurement at default 16 threads vs 1 thread: case 50 (P=6) 15.5 s vs 7.0 s; case 52 (P=2) 1.85 s vs 0.91 s — parallel was ~2 × WORSE. Added serial twins _packed_residual_fo_ser / _packed_min_dt_ser in aeolus1d/solver/packed_pipes.py and a Python dispatcher that picks parallel vs serial based on a fast-pipe-count threshold (_PARALLEL_THRESHOLD_P, default 12). Wall-clock on the engine cases: 50 16.78 s → 6.24 s (2.69 ×), 52 2.07 s → 0.95 s (2.18 ×); case 54 unchanged (uses advance_network_implicit, different hot path). Override via aeolus1d.solver.packed_pipes.set_parallel_threshold(P) (set to 0 to always use parallel — legacy behaviour — or to a very large number to always use serial). The SIMD-rewrite alternative was tried and reverted: bit-exact regression (fastmath FMA reordering) and a wall-clock regression (16-33 ×) on the engine cases. * Phase Perf-6 — chemistry advection JIT. Two new @njit kernels in aeolus1d/solver/species_advection_njit.py (species_advection_primary_jit, species_advection_extras_jit) replace the per-cell Python loop in network_timeloop._species_rate. Reaction callable + species-inflow BC overlay + DPF filter-trap stay in the Python wrapper (they need Python-side dict / callback access). Triggered by Phase Perf-history-2 cProfile that flagged _species_rate as 56 % of case 41 wall (1.19 s of 2.13 s) — same Python-loop- around-JIT-kernels pattern Phase Perf-2 fixed for the closure-Jacobian dispatcher. Apples-to-apples wall-clock on 41_dpf_regen.toml (default t_end = 0.06): 4.63 s → 1.74 s (2.66 ×). _species_rate self-time dropped from 1.19 s → 0.27 s; the new JIT kernel adds ~0.11 s, so end-to-end the chemistry hot path is ~3 × faster. Bit-exact algebraically (loop body unchanged, just moved into numba).

Recent polish (2026-05-07) — * Phase I-8-ilu0 — ILU(0) JFNK preconditioner. New NewtonConfig.jfnk_precondition_ilu0 (and [case] jfnk_precondition_ilu0 in TOML) reuses _block_thomas_solve on the chord-aged (L_J, D_J, Up_J) as M^{-1} — for block-tridiag J the ILU(0) factorisation is the block-Thomas factorisation, so the precond apply collapses to one Thomas pass. GMRES typically converges in 1–2 Krylov iters. Long-pipe bench at N=500: ILU(0) JFNK 0.398 s, pure JFNK 0.423 s, direct 0.227 s — ILU(0) is the cleanest matrix-free path (1.75–2.24 × slower than direct at all N tested). Mutually exclusive with jfnk_precondition (block-Jacobi); ILU(0) takes priority if both are set, pinned by test_jfnk_ilu0_precedence_over_block_jacobi. * Long-pipe bench + JFNK preconditioner reversal. New scripts/bench_long_pipe.py measures direct vs JFNK (preconditioned + pure-matrix-free) at user-controllable N. At N = 80: direct 0.009 s; jfnk_pure 0.024 s (2.6 ×); jfnk_precond 2.5 s (272 ×). At N = 500: direct 0.235 s; jfnk_pure 0.440 s (1.9 ×); jfnk_precond 121.9 s (519 ×). The block-Jacobi preconditioner that landed in the I-8-precond MVP is the wrong shape for 1D hyperbolic systems — it discards the off-diagonal acoustic + convective coupling that dominates J's spectral content, and GMRES iter count goes from 4 (no precond) to 42 (with precond) on the same matvec. Default flipped to jfnk_precondition = False; machinery kept opt-in for users with less-hyperbolic workloads. scripts/check_perf_history_refs.py + matching pytest pin the contract that every commit hash referenced in docs/perf-history.md resolves in the repo. * Phase I-8 (MVP) — opt-in Jacobian-Free Newton-Krylov. New NewtonConfig.use_jfnk flag (default off) replaces the direct block-tridiag-Thomas solve with GMRES (scipy) using FD-approximated J·v products: J·v ≈ (F(U + ε·v) − F(U)) / ε with Brown-Saad scaling. jfnk_precondition (default on when JFNK is on) builds the chord-aged diagonal-block D_J and uses it as a right-side block-Jacobi preconditioner via N independent 3×3 inverses (_block_jacobi_apply_jit); turn off for memory-bound pure-matrix-free runs. Slower than the direct path on the canonical case 54 — the direct block-tridiag-Thomas solve is the right cost shape at small N — so default stays off; JFNK is in the codebase mainly for users who push to networks where direct Jacobian assembly runs out of memory. [case] use_jfnk / jfnk_tol / jfnk_max_iters / jfnk_precondition surface in TOML. * Phase 27c soak test. New @pytest.mark.slow tests/integration/test_parallel_pipes_soak.py runs examples/54_adaptive_picard_turbo.toml end-to-end with both serial and parallel_pipes=True, asserts state bit-equivalence (rtol=1e-9), and pins repeat-run determinism. The soak measurement (full t_end = 50 ms): serial 7.2 s, parallel 18.8 s — parallel is ~2.6× slower on this benchmark, even stronger than the 2× the short cProfile run showed. Confirms parallel_pipes = False is the right default; the path stays for genuinely large networks where per-pipe JIT work clearly outweighs thread coordination overhead. * Phase Perf-history doc. docs/perf-history.md collects all seven Perf phases (Perf-1 → Perf-5 + 27c + 27c-cells) with the audit-driven before/after numbers and the empirical-vs-speculative pattern (JIT'd BCs sped up serial more than parallel; prange at N=80 was net-negative). Cumulative summary: case 40 0.323 → 0.020 s (94 % — 16× speedup); case 54 39.66 → 2.67 s (93 % — 15× speedup). * Phase Perf-4 — JIT'd BC dispatchers + chord knobs on TOML. make_stagnation_inlet_bc and make_pressure_outlet_bc delegate to @njit kernels (_stagnation_inlet_apply_jit, _pressure_outlet_apply_jit) so the BC apply releases the GIL — supports the Phase 27c parallel-pipes path and also speeds up serial. New [case] chord_max_age and chord_rate_threshold fields surface the Phase Perf-3 chord knobs to TOML cases. Wall-clock impact: * example 40 (single-pipe implicit): 0.048 → 0.020 s (~58 %) * example 54 (network, serial): 8.24 → 2.67 s (~68 %) Cumulative since the audit started: * example 40: 0.323 → 0.020 s (94 % — 16 × speedup) * example 54: 39.66 → 2.67 s (93 % — 15 × speedup) Note: Phase 27c parallel-pipes is now slower than serial on this benchmark (5.10 s parallel vs 2.67 s serial post-Perf-4) because thread-coordination overhead dominates once the BC + JIT-kernel work is fully GIL-released. Stays opt-in for very large networks where per-pipe JIT work clearly outweighs thread setup; default off. * Phase 27c — opt-in parallel per-pipe Newton solves. NewtonConfig.parallel_pipes and [case] parallel_pipes in TOML run each pipe's Newton step on a thread inside the Picard inner loop. Pipes are independent during a Picard iter (junction state is frozen between _fill_component_ghosts calls); numba kernels release the GIL inside the JIT regions so the Jacobian assembly + block-Thomas + closure-Jac work scales sublinearly with pipe count. Failure path falls back to sequential dt-halving, preserving legacy semantics. Wall-clock impact: example 54 8.24 → 5.73 s (~30 %). With this, cumulative case 54 39.66 → 5.73 s (~86 %). Default off until more soak time on the network regression suite. * Phase Perf-3 — chord (lazy-Jacobian) reuse. Hold J = I − dt · ∂R/∂U fixed across 2–3 inner Newton iterations, refresh on rate-degradation (F_norm_new / F_norm_old > 0.5) or line-search rejection. Two new NewtonConfig knobs: chord_max_age (default 2) and chord_rate_threshold (default 0.5). Skips the JIT Jacobian kernel + closure-Jac augmentation + numpy broadcasts on chord-reused iterations. Wall-clock impact: example 40 0.110 → 0.048 s (~56 %); example 54 18.09 → 8.24 s (~54 %). * Phase Perf-1-toml — [gas.eos] kind = "tabulated". The kernel-level EOS_TABULATED dispatch is now reachable from TOML via a new EOSSpec arm with CoolProp + grid-bounds fields. Builds the (T, ρ) → (p, e, ∂p/∂T, ∂p/∂ρ) panel at network construction time (lazy CoolProp import). New bundled example examples/58_tabulated_eos_implicit.toml exercises the path. Closes the verification loop on Phase Perf-1's batched eos_a_and_derivs_jit win. * Phase Perf-2 — JIT-ified closure-Jacobian dispatch. The Python for i in range(pipe.N) loop in _augment_jacobian_with_closures (which the post-Perf-1 cProfile flagged at 25–37 % of single-pipe wall time) moved into a single @njit kernel _augment_closure_jac_analytic_jit. Five real-gas helper functions (_real_gas_T_app_chain, _real_gas_heat_jacobian{,_from_chain}, _real_gas_porous_sutherland_jacobian{,_from_chain}) were plain Python despite calling @njit kernels internally — decorated @njit(cache=True, fastmath=True, inline="always") so the JIT loop inlines them. FD fallback for Python-override closures stays in the wrapper. Wall-clock impact: example 40 0.241 → 0.112 s (~54 %); example 54 27.7 → 18.3 s (~34 %). Cumulative across Perf-1 + Perf-2-quick + Perf-2: ex 40 0.323 → 0.112 s (~65 %), ex 54 39.7 → 18.3 s (~54 %). * Phase Perf-2-quick — broadcast D_J = I − dt·D_blk. Per-Newton-iter Python loop for i: D_J[i] = I − dt·D_blk[i] in _newton_step replaced with a single numpy broadcast (D_J = -dt·D_blk; D_J[:, _DIAG3, _DIAG3] += 1.0). Bit-equivalent to the legacy expression. Wall-clock drop measured by scripts/profile_perf1.py: ~25% on examples/40_implicit_adaptive_pipe.toml (0.323 s → 0.241 s) and ~30% on the 4-pipe network case examples/54_adaptive_picard_turbo.toml (39.7 s → 27.7 s). _newton_step's Python self-time dropped 29 ms → 1 ms (single- pipe) and 15.5 s → 3.1 s (network). Phase Perf-2 (full JIT- ification of the Newton orchestrator + closure-Jacobian dispatch loop) is queued in BACKLOG.md. * Phase Perf-1 — EOS / solver hot-path optimization sweep. Six concrete wins on the per-Newton-iteration hot path: (1) tabulated ∂p/∂T, ∂p/∂ρ panels (Phase B23 data) now wired directly into the Jacobian dispatch — replaces a 4-call central-difference fallback in eos_dpdT_dpdrho_jit; (2) eos_dedT_dedrho_jit for tabulated EOS uses the closed-form bilinear gradient of the energy panel (no more 4-call FD); (3) batched _bilinear_two / _bilinear_grad_two helpers share the (T, ρ) cell-index work across the dpdT and dpdρ panels; (4) new fused eos_a_and_derivs_jit returns (a, ∂p/∂T, ∂p/∂ρ, ∂a/∂T, ∂a/∂ρ) in one dispatch — collapses three separate calls per face per Newton iter in _eos_primitive_gradients; (5) fused tabulated (p, a) recovery in primitives_from_U_eos_jit after Newton T inversion; (6) block-Thomas _inv3x3_into writes through a hoisted scratch buffer instead of allocating per-cell. Plus Picard-loop bookkeeping skip when picard_max_iter == 1 and Tikhonov regularization in chemistry-BS _solve_2x2 so near-singular rows stop NaN-ing out and forcing rejected steps. * Closure-Jacobian programme closed (Phases I-Y + B16 → B23). Phase I-Y wraps the pipe Newton + ghost refresh in a Picard outer iteration so Compressor / Turbine / LODI ghost coupling is no longer split across the step. Phases B16 → B23 brought every shipped closure × kind × viscosity × EOS combination onto an analytic Jacobian path: friction (DARCY_CONST, COLEBROOK; both EOSes), heat (CONSTANT_TWALL; ideal + real gas), porous (DARCY_FORCHHEIMER + KC; constant μ + Sutherland μ; ideal + real gas), closed-form ∂a/∂U for cubic EOS (RK / PR, Phase B20) and tabulated EOS via differentiated bilinear interpolant (Phase B23), plus a shared T_app chain (Phase B21) that halves the EOS work when both real-gas heat + Sutherland porous fire on the same cell. Bench (@pytest.mark.slow): ideal-gas analytic ≈ 7.5× faster than the FD fallback; real- gas analytic ≈ 3×; tabulated EOS sound-speed Jacobian ≈ 4.3×. View → Solver Internals in the GUI shows the live per-pipe routing. * Adaptive-Picard tuning recipe (f9aa257) — new examples/57_picard_stiffness_step.toml + plot_57_* showing the controller's full life cycle (cold-start saturates cap, ratchets to floor, bumps for stiffness step, re-settles). Two explicit tuning signals documented in docs/phase-iy-picard.md. * Phase B22 — Wiebe combustion + closed-loop knock-controller TOML support. New [cylinder.combustion] block, all six knock_* fields parsed/saved. Bundled 52_knock_event.toml exercises the chain. * Phase 51d — full multi-species transport across 0-D plena. Junction.m_species now tracks every species column on connected pipes, not just the fuel column; gas leaving a plenum carries the well-mixed composition through to the downstream pipe. * CLI Cruise M import — aeolus1d import avl_cruise_m project.proj out.toml now works (was previously GUI-only). * GUI MVC rewrite (8 phases) — MainWindow is now a controller; case + layout + history live on aeolus1d_gui.model.CaseModel, run-result state on RunModel. Views subscribe to model signals (caseReplaced, caseMutated, runLoaded, canUndoChanged, …) instead of being manually re-synced after every action. See docs/architecture.md.

What's left. docs/roadmap.md is the canonical backlog. Open items: wastegate component, twin / sequential turbo topology, HCCI / RCCI multi-zone combustion, DOC sulfur-poisoning model, multi-pipe prange master kernel (Phase 27c), full FMI 2.0 round-trip.

The speculative bucket (JAX-mirror solver for AD, AMR, 3D-Euler core, cycle-resolved post) is in docs/roadmap.md. Pick anything there — they're independent.

Milestone 1 (closed) — 1D compressible-Euler pipe network core

  • Finite-volume 1D Euler solver on a network of pipes
  • HLLC Riemann solver, MUSCL 2nd-order reconstruction, SSP-RK2 time stepping
  • Boundary conditions: stagnation inlet, static-pressure outlet, wall, non-reflecting, transmissive
  • Junctions / flow splits through a 0D plenum with plenum-ghost port coupling (mean-flow mass and energy balance exact)
  • Quasi-1D area variation and integer-dispatched friction + heat closures (DARCY_CONST, CONSTANT_TWALL)
  • TOML case files, CLI (aeolus1d run, aeolus1d bench, aeolus1d examples), HDF5 output, benchmark suite (Sod, C-D nozzle, T-junction)

Milestone 2 (closed) — network capabilities

  • Phase 7 — port-BC investigation (Helmholtz xfail stays). Explored port-cell + pressure-match replacements for the plenum-ghost; a forced- response sweep confirmed the plenum-ghost is a first-order low-pass rather than a resonator. Physical diagnosis + measured transfer function logged in the xfail reason. Real fix needs a compliant / sub-cycled characteristic port protocol — scheduled for milestone 3.
  • Phase 8 — MUSCL + quasi-1D. The 2nd-order reconstruction now runs with variable area + closure sources. C-D nozzle at 400 cells: 0.265% Mach error (1st-order 0.332%).
  • Phase 9 — transient (time-varying) BCs. bc/transient.py with a Clock and piecewise-linear p₀(t) / T₀(t) / p_back(t) tables. Integration test recovers the imposed driving frequency to < 2%.
  • Phase 10 — Lax + shock-reflection benches. Lax L1(ρ) = 1.79% at MUSCL/400; shock reflection within 0.02% of Rankine-Hugoniot.
  • Phase 11 — time-series HDF5 output. aeolus1d run … --trace-dt DT writes a /time_series/ group with per-pipe and per-junction traces.
  • Phase 12 — SSP-RK3 + cleanup. RK3 (Shu-Osher) wired into both time loops; unused BCKind.JUNCTION_PORT retired.
  • Phase 13 — Colebrook friction + porous-media + ML closure seam.
  • FRICTION_COLEBROOK via the Haaland explicit approximation with laminar (64/Re) fallback below Re = 2300.
  • New PorousKind.DARCY_FORCHHEIMER volumetric drag closure (s_mom = −(μ/K) u − C_F ρ |u| u) — foundation for catalysts, DPFs, mufflers, flow straighteners. Porous-plug bench matches theory to 0.00–0.03%.
  • Python closure plug-in seam: attach py_friction / py_heat / py_porous callables to a Pipe and the timeloop overlays them on top of the @njit kernels — ready for ML surrogates.

Milestone 3 (closed) — engine/powertrain-simulator scaffolding

Phases 14–19d shipped. Each phase has a matching example (TOML or plot script) under examples/. 181+ tests pass + 1 xfail (the 0-D plenum Helmholtz, an approximation-specific limitation — the quasi-2D path measures to 0.01 %). Highlights:

  • Phase 14 — Winterbone neck-pipe junction. Each port can own a short internal FV pipe (neck_cells > 0) of Rayleigh-corrected length, providing the slug inertia missing from the bare plenum- ghost port.
  • Phase 14b — Corberán CP+h N-branch junction. Algebraic constant-pressure + total-enthalpy closure (no chamber state), faster + more accurate than the plenum for runner manifolds.
  • Phase 14c — Bassett-style per-port pressure-loss coefficient K_loss on the Corberán junction.
  • Phase 14d — Extended Helmholtz Resonator (EHR). Per-port u_n ODE state coupled through ρ_V·L_eff·du_n/dt = σ·(p_pipe − p_V). Right primitive for perforated-liner CAA work.
  • Phase 15 — Two-phase catalyst. Per-cell T_solid(x, t) field + Nusselt heat coupling h·a_v·(T_s − T_g) between gas and solid.
  • Phase 15b — Species transport + py_reaction ML seam. Single- scalar Y advected by the gas; Arrhenius demo closure ships as proxy for production SCR/combustion kinetics.
  • Phase 16 — Quasi-2D Helmholtz. Spatially-resolved chamber reproduces f_H to 0.01 % with tanh_step area profile and regional initial conditions.
  • Phase 17 — py_junction_bc ML seam parallel to the per-cell closures. Convergence point with of_gui/ OpenFOAM-trained plenum surrogates.
  • Phase 18 — Strang-split IMEX chemistry scaffold. Outer T(Δt/2) ∘ C(Δt) ∘ T(Δt/2) loop, per-cell scipy.solve_ivp(LSODA) for the stiff (Y, E) ODE.
  • Phase 19 / 19b / 19c — Real-gas EOS. Mass-specific Redlich- Kwong + Peng-Robinson; opt-in Pipe.eos_kind; @njit HLLC + Newton- based T-from-energy inverse, ≈160× faster than the Python-mode reference at runtime.
  • Phase 19d — @njit Bader-Deuflhard semi-implicit Bulirsch-Stoer chemistry sub-stepper. Drop-in alternative to LSODA via chem_solver = "bs_jit". Sub-step sequence {2, 6, 10, 14, 22, 34, 50}, frozen 2×2 closed-form Jacobian, Cramer's-rule solves. ~12–14× faster than LSODA on stiff Arrhenius cases.

Milestone 4 (closed) — engine-side stack + 2nd-order numerics + optimization

  • Phase 20 — EOS + closure stack. Lift the _check_eos_pipe_supported guard and wire EOS dispatch through the quasi-1D kernel (friction / heat / porous / catalyst). Real-gas shock tubes can already run; the next step is real-gas catalyst beds and porous plugs.
  • Phase 21 — EOS-aware junction ghosts. Plenum and EHR junctions fill the ghost ρ·E via the pipe's EOS internal energy. (Corberán remains ideal-gas-only — its CP+h closure uses Riemann invariants derived from ideal-gas isentropic relations; EOS-Corberán is a research extension.)
  • Phase 21b — EOS-aware Volume. The 0-D plenum's (m, E) → (p, T) recovery now goes through the EOS Newton inverse for real-gas plena. Closes the on-disc consistency story for plenum + EHR networks.
  • Phase 22 — Multi-species (passive extras). Pipe.Y_extras is a (N, n_extra) array of additional species advected with the same HLLC mass flux as the primary Y.
  • Phase 22b — Multi-species kinetics. Pipe.py_reaction_multi is a callable on the full (Y, Y_extras) stack returning rate arrays for every species. make_arrhenius_scr_multi_demo ships a 2-species standard-SCR rate law (NH3 + NO consumed in 1:1 molar ratio).
  • Phase 23 — Characteristic (LODI / Poinsot-Lele) port BC. Ships Port(lodi=True) — standard Poinsot-Lele subsonic outflow that pins p_b = p_V and extrapolates ρ, u from the pipe's outgoing invariant. Confirmed by experiment that LODI alone does NOT recover the Helmholtz mode on a 0-D plenum (the plenum's missing acoustic compliance is fundamental); examples/06_helmholtz_quasi2d.toml remains the path for Helmholtz-accurate plena.
  • Phase 24 — MP5 5th-order monotonicity-preserving reconstruction (Suresh-Huynh 1997). Drop-in alternative to MUSCL — scheme = "mp5" in the case header. Cheaper than WENO-5, simpler to @njit, and about 50× more accurate than MUSCL on smooth flow at the same resolution while preserving monotonicity at shocks. Requires Ng=3 (handled automatically by the loader). Pure ideal-gas + uniform area in this phase; MP5 + closures is a follow-on.
  • Phase 24b — HLLC-LM low-Mach corrected HLLC (Thornber et al., Comp. Fluids 2008). Drop-in alternative to standard HLLC that rescales the face velocity reconstructions toward the average by φ = min(1, M_local / M_threshold) before passing them to the star-state computation. At Mach ≥ threshold (= 0.1 default) the flux is bit-exactly standard HLLC; at low Mach it preserves significantly more amplitude on acoustic features and idle-state flows. Opt-in via [case] hllc_low_mach = true (or Pipe.hllc_low_mach); first-order ideal-gas uniform-area path for now. [2026-09-13] The kernel restriction stated here held. What this entry did not say, because it was not known, is that the flag reached only the explicit network driver's residual — a single-pipe run accepted it and ignored it. Fixed 2026-09-13; see docs/changelog.md.
  • Phase 24c/d/e/f — MP5 closes the dispatch corners. Phase 24 shipped MP5 only on uniform-area, closure-free, ideal-gas pipes; this set of follow-ons brings MP5 to feature-parity with MUSCL. 24c adds pipe_residual_quasi1d_mp5 — variable area + friction / heat / porous closures in the residual loop, same MP5 face-state reconstruction. 24d adds the per-cell-γ companion for mixture-tracking pipes. 24e extends the network-timeloop MP5 dispatch to use those kernels (and catalyst comes for free since _catalyst_rates_into overlays its source on top of any scheme's residual). 24f lifts the ideal-gas restriction with pipe_residual_eos_mp5_jit / pipe_residual_quasi1d_eos_mp5_jit — same Suresh-Huynh stencil applied to W = (ρ, u, T) (closed-form W → U for cubic EOSs). After this set the GUI's pre-flight scheme-compatibility helper has nothing left to flag. Tests: bit-exact reduction back to the simpler kernel for each combination's identity case (uniform γ_pad → quasi1d_mp5, EOS_IDEAL → ideal-gas MP5, etc.) plus end-to-end runs that verify the dispatch isn't silently downgrading to MUSCL.
  • Phase 25 — Chalet-Chesse-Hetet valve effective-area BC. Ships make_valve_inlet_bc(clock, lift_t, lift_y, lift_max, p0, T0, …) and the kind = "valve_inlet" TOML BC. A periodic lift table modulates an underlying stagnation inlet — at lift = 0 the BC collapses to a wall (closed valve), at lift = lift_max it behaves as a full stagnation inlet, intermediate fractions scale the ghost velocity by f. The model is phenomenological — a rigorous Cd · A_curtain orifice formulation is a future extension.
  • Phase 25b — Rigorous orifice-flow valve BC (Winterbone-Pearson Ch. 7; Heywood §6.2). Replaces the Phase-25 phenomenological velocity-scaling with the standard isentropic-orifice mass-flow formula (subsonic + choked branches). kind = "valve_orifice" in the TOML BC with A_max, A_pipe, Cd parameters that set the curtain area and discharge coefficient. The orifice BC respects the throat-vs-pipe area restriction, giving ~3× lower peak mass flow than Phase-25 at the same lift schedule when Cd·A_max ≪ A_pipe — matching engine-textbook physics.
  • Phase 25c — Quasi-steady-map compressor BC. Ships Compressor and TurboMap in aeolus1d.components.turbo. Two pipe ends are coupled via a tabulated centrifugal-compressor map (Serrano-Payri / Heywood Ch. 6 convention): inlet stagnation state + outlet back-pressure feed a bilinear (PR, N_corr) → (ṁ_corr, η_is) lookup, which sets the per-step mass flow and outlet stagnation temperature. The component is duck-typed as a Junction (zero-rate dummy volume) so it plugs straight into the existing advance_network(junctions=[...]) argument with no timeloop changes.
  • Phase 25d — Turbine + shaft-inertia ODE coupling. Adds Turbine (symmetric to Compressor with sign-flipped T0_out / shaft-power conventions and a synthetic_turbine swallowing-capacity map factory) and Shaft — a 1-DOF rotational inertia integrating I·dω/dt = (Σ Ẇ_T − Σ Ẇ_C)/ω via the existing per-port-u_n SSP-RK2 / RK3 hook (same machinery the EHR neck- velocity ODE uses). Compressor and Turbine accept an optional shaft= field; when set, they read ω from the shaft and contribute signed power back. The shaft enforces an omega_min floor to keep the ODE integrable across zero-speed crossings. Demo: paired hot-exhaust turbine + cool-intake compressor on a shared shaft — ω equilibrates when turbine power matches compressor load.
  • Phase 25e — Learned (trained) port-state BC surrogate. Ships LearnedJunctionBC in aeolus1d.closures.learned_junction_bc: a small from-scratch numpy MLP (5 → H → H → 3 with tanh + Adam) that obeys the Phase-17 py_junction_bc callable signature and is trained by SGD on synthetic samples from any target callable (e.g. resistive_demo, plenum-ghost, or CFD-derived data). Once trained, it is a drop-in replacement for the hand-coded callable — the Junction port-bc dispatch never knows the difference. Pure numpy, deterministic given a seed, no third- party deep-learning dependency. Demo ships a 600-epoch fit (~100 ms) on the existing resistive_demo target — the trained network reproduces the hand-coded BC's full transient to within 0.4 %.
  • Phase 25f — Learned compressor/turbine map surrogate. Same numpy-MLP / Adam pattern as Phase 25e, lifted from port-state to the TurboMap.lookup interface. Ships LearnedTurboMap in aeolus1d.closures.learned_turbo_map — a 2 → H → H → 2 network trained on samples from any TurboMap (a synthetic factory or a manufacturer table). Because Compressor and Turbine consume their map field via pure duck typing (.lookup, .p_ref, .T_ref), the trained network is a drop-in replacement — no Compressor / Turbine change needed. Demo: 600 epochs (~60 ms) on the synthetic_centrifugal factory, steady-state Compressor ṁ deviates by 0.3 % from the hand-coded reference run.
  • Phase 25g — Reciprocating-cylinder volume BC. Ships Cylinder and CylinderPort in aeolus1d.components.cylinder — a 0-D variable-volume plenum driven by a slider-crank kinematic (bore / stroke / con-rod / clearance volume / crank speed) with ≥1 valve ports on crank-angle-windowed open / close schedules (binary 1 / 0 — continuous poppet-lift schedule is the natural Phase-25h follow-on). The piston p · dV/dt work is added as an explicit source on the volume's energy balance, so the cylinder pumps gas through its ports while compressing / expanding the trapped charge. Closed-valve compression matches the adiabatic pV^γ = const textbook relation bit-near- exactly. Demo: a 4-stroke single cylinder (80 mm × 90 mm, CR = 10, 1800 rpm) breathing through paired intake (0°-180°) and exhaust (540°-720°) valve windows — peak compression pressure 26 atm, mass swing 11.7×.
  • Phase 25h — Wiebe combustion on the cylinder. Adds WiebeCombustion (textbook Heywood §9.2.2 single-zone model Q̇(θ) = m_fuel · LHV · dx_b/dθ · ω) as an optional field on Cylinder. The piston pdV-work term is joined by a chemical heat-release source on the cylinder's energy ODE. Demo: 80 mm bore × 90 mm stroke, CR = 10 SI engine at 1800 rpm with θ_start = 355°, Δθ = 60°, m_fuel = 50 mg, iso-octane LHV — peak pressure 68 bar (vs 26 bar pure compression), peak T 4902 K, peak Q̇ 796 kW. Closes the reciprocating-engine arc — combined with the Phase-25c/d turbo
  • Phase-25g cylinder + Phase-25b orifice valve infrastructure this gives a complete IC-engine 1-D simulation.
  • Phase 25i — Continuous poppet-lift schedule on Cylinder. Per-port piecewise-linear lift(θ) table with A_eff = Cd · A_port · lift(θ) / lift_max so the cylinder-to- pipe flux scales smoothly through the cam profile.
  • Phase 25j — Two-zone Wiebe combustion diagnostic. Optional WiebeCombustion(two_zone=True) partitions the bulk volume.E_internal into burnt + unburnt zones via the cumulative Wiebe burn fraction x_b(θ). Burnt-zone temperature from the energy balance; unburnt-zone temperature from adiabatic compression of the trapped charge to the cylinder's bulk pressure. Single-γ approximation (Heywood §9.4).
  • Phase 25k — Per-zone γ tracking in the two-zone diagnostic. Optional WiebeCombustion.gamma_burnt / WiebeCombustion.gamma_unburnt overrides drive separate cv_b = R/(γ_b−1) and cv_u = R/(γ_u−1) in the energy partition: E = m_b·cv_b·T_b + m_u·cv_u·T_u. Adiabatic- compression formula for T_u reads γ_u instead of the cylinder's single γ. With both fields unset the diagnostic collapses bit-for-bit to Phase 25j; setting realistic values (γ_u ≈ 1.40 cold-charge, γ_b ≈ 1.27 hot products) drops the burnt-zone temperature ~15-25 % at mid-burn — the textbook correction the Phase-25j docstring flagged.
  • Phase 44 — Case-driven optimization driver. A thin Python API (aeolus1d.optimize.optimize_case) wraps scipy.optimize.minimize with two pieces of glue: path-based parameter binding into the Case dataclass tree (e.g. "pipes[0].init.p" or 'pipes["tube"].length', with id- keyed list lookup) and probe-based objectives that compare a traced pipe time-series against reference data (ProbeObjective + make_probe_match_objective + load_reference_csv, with l1 / l2 / rmse / max metrics). The full version adds an [optimize] TOML block + aeolus1d optimize <case.toml> [--out tuned.toml] CLI subcommand that prints the per-step (theta, J) trajectory and substitutes the optimum back into a sibling TOML, plus a worked examples/optimize/tune_initial_pressure.toml demo (recovers a 1.50e5 Pa truth from a 1.10e5 Pa guess in ~32 evaluations). GUI dock (OptimizeDock) tabbed alongside the other inspect docks: subprocess-based QProcessOptimizeRunner streams the CLI's stdout into structured step(n, J, theta) signals, inline-editable parameter / probe tables, Browse CSV picker that stores paths relative to the loaded TOML, Save button that rewrites the [optimize] block in place. 47 tests.
  • Phase 43 — Higher-order RK time integrators (RK4 + RK5). Adds classical 4-stage 4th-order Runge-Kutta and Cash-Karp 6-stage 5th-order RK to the existing SSP-RK2 / SSP-RK3 set on advance_single_pipe. Generic _butcher_step helper walks any explicit Butcher tableau row by row so future schemes (Verner 6(5), Tsitouras, etc.) can drop in without reworking the loop. RK4/RK5 are non-SSP — appropriate for smooth-flow benchmark studies, not for shock-rich runs where TVD matters; the shock-rich engine path stays on SSP-RK2/3 and advance_network raises a clear error if the user picks RK4/5 there. Validated on smooth-wave advection: the canonical "MP5 spatial 5th-order finally exposed" result (order ≈ 1.99 on RK2, 3.42 on RK3, 5.00 on RK4 / RK5) lands cleanly. 8 integration tests pin every (scheme, RK) combination's observed accuracy floor and verify RK4 / RK5 agree on a well-resolved smooth grid.
  • Phase 42 — Multi-pipe FEM eigensolver with junction DOFs. Extends aeolus1d.modes from single-pipe to a 1-D network. Each pipe contributes its existing local M / K matrices; pipe-end DOFs that participate in a junction collapse via union-find into one shared "junction DOF" (continuity of velocity potential). The case's [[bc]] list maps to FEM BCs (wall / transmissive → Neumann, pressure_outlet / stagnation_inlet / nonreflecting → Dirichlet); bc_overrides lets callers force a specific BC at any (pipe_id, end) pair without rewriting the case. Returns a NetworkModesResult with one frequencies array plus per-pipe sampled mode shapes (per_pipe_shapes / per_pipe_x) and per-junction shape values for plotting. 6 integration tests pin: two equal pipes joined at a zero-volume junction reproduce the analytical modes of the combined L_total pipe to ≤ 1 % rel-err on the first 3 physical modes; a 3-pipe Y-network has matching mode values at the junction DOF (continuity to 1e-12); BC mapping produces the correct quarter-wave fundamental on a wall-open pipe; bc_overrides flips a closed-open quarter-wave into an open-open half-wave at twice the frequency; per-pipe x arrays span [0, L_pipe]; empty case raises.
  • Phase 41 — Mode-projection initial conditions. New init.kind = "mode" for [pipe.init] blocks: the loader runs aeolus1d.modes.solve_modes on the pipe geometry to get the mode-N velocity-potential shape φ_n(x), samples it at the pipe's cell centres, and applies an isentropic over-pressure p(x) = p_ref · (1 + ε · φ_n(x)) (with matching ρ; u = 0). Closes the spectrum calibration loop between the FEM eigensolver and the time-domain Euler solver: a clean mode-N IC must produce an FFT peak at f_n to within the time-domain truncation error. Schema: ModeInit dataclass with mode_index (1-based, skipping any rigid-body mode at f ≈ 0), amplitude, p_ref, T_ref, bc_left, bc_right, n_elements, p_order. Mirrors the FEM module's BC vocabulary so the IC and the eigensolve agree. GUI save_case writes the full [pipe.init] table for round-trip. 5 integration tests pin: spatial-shape projection (mode 1 of a wall-open pipe peaks at the wall), invalid mode index raises, mode-1 IC produces an FFT peak within 5 % of f_1, mode-2 IC produces an FFT peak within 5 % of f_2, and a [pipe.init.kind = "mode"] TOML dict parses through _parse_init into a ModeInit dataclass.
  • Phase 25m — Pollutant emissions (Zeldovich NOx, Hiroyasu soot, HC slip). Three production-relevant emission models ship as py_reaction_multi callable factories in aeolus1d.closures.pollutants: make_zeldovich_nox for thermal NO formation in the post-flame zone (extended Zeldovich quasi-steady-state collapse to d[NO]/dt = 2·k_1f·[N₂]·[O] with [O] from the ½ O₂ ⇌ O partial equilibrium), make_hiroyasu_soot for Hiroyasu-Kadota two-step soot (formation rate A_f·ρ_fv·p^0.5·exp(−Ea_f/RT) minus oxidation rate A_o·ρ_s·X_O2·p/√T·exp(−Ea_o/RT)), and make_hc_slip for Lyford-Pike-style first-order HC relaxation (Arrhenius decay toward an equilibrium mass fraction; freezes at exhaust-pipe T, fast at combustion T). Heywood §11.3 textbook forms throughout; defaults are typical engine-research values. All three have zero heat release on output (pollutant rates are tiny vs. the main combustion). 15 integration tests pin the property checks the literature requires: NOx Arrhenius scaling with T (rate at 2500 K vs 2000 K within a 30-1000× window), zero rate at zero O₂, monotone increase with [N₂] / [O₂]; soot formation dominates in fuel-rich + low-O₂ zones, oxidation dominates in fuel-lean + high-O₂ zones; HC slip decays at high T, freezes at low T, fixed point at Y == Y_eq.
  • Phase 25l — Cylinder wall heat transfer (Woschni / Hohenberg). Adds HeatTransferModel (abstract base) and two concrete subclasses, WoschniHeatTransfer (Heywood eq. 12.30, h = 3.26 · B^-0.2 · p_kPa^0.8 · T^-0.55 · w^0.8 with w = C1 · S_p_mean, default C1 = 2.28 — full-cycle simplified) and HohenbergHeatTransfer (1979 SAE 790825 form, h = 130 · V^-0.06 · p_bar^0.8 · T^-0.4 · (S_p + 1.4)^0.8). Plugged into Cylinder as an optional heat_transfer field; when set, the gas-side heat-loss term Q̇_wall = h · A_wall(θ) · (T_wall − T_gas) joins p · dV/dt and the optional Wiebe heat release on the cylinder energy ODE. Instantaneous wall area tracks the piston position (head + piston + crank-angle-dependent liner). 10 integration tests cover: closed-form correlation values, sign of the heat flux for hot-gas-cold-wall and cold-gas-hot-wall scenarios, zero flux at ΔT = 0, peak-pressure drop on closed-valve compression vs adiabatic baseline, Woschni / Hohenberg agreement to within a factor of 2 on instantaneous h, and monotone wall-area sweep TDC → BDC. Adiabatic default preserves the Phase-25g bit-identity behaviour.

Milestone 5 (closed) — advanced multi-species chemistry + catalyst stack

Combustion- and catalysis-grade species thermo + reaction kinetics, with native database persistence and standard-format importers (CHEMKIN / Cantera YAML / JANAF). Lays the groundwork for replacing the Phase-22b py_reaction_multi callable with mechanism-driven kinetics where cp(T), γ(T, Y) etc. come from real species thermo rather than scalar pipe-level constants.

  • Phase 26 — Species thermo database. Ships aeolus1d.chemistry package: SpeciesThermo (NASA-7 polynomial evaluator with cp(T), h(T), s(T) in both molar and per-kg form, T-range clamping), importers.nasa7 (parses standard CHEMKIN 4-line thermo blocks — same format Cantera / GRI-Mech / NIST files use), and a native NPZ-per-species database layer (db.save_species / db.load_species / db.list_species). The package ships a built-in default DB of 10 common species (N2, O2, H2, H2O, CO, CO2, CH4, OH, NO, NH3) with NASA-7 polynomials sourced from GRI-Mech 3.0 + Burcat — every bundled species reproduces the textbook Cp_mass(298.15 K) within 1.6 %. 17 unit tests + plot demo of all bundled species' Cp/R(T), Cp_mass(T), and h(T) − h(298.15) curves over 200-3000 K.
  • Phase 26b — JANAF tabulated importer. Ships aeolus1d.chemistry.importers.janaf with import_janaf_table(path) (parses the standard 4-column T / Cp/R / H/RT / S/R tab- or whitespace-separated text file every JANAF / NIST entry exports) and fit_nasa7_to_table (5-coefficient linear least-squares fit on each T-range with the integration constants a6 / a7 set from a single reference (T_ref = 298.15) row). Returns a SpeciesThermo with the same cp(T) / h(T) / s(T) API the rest of the package consumes — so a JANAF-imported species drops straight into Mechanism / MechanismKinetics / mixture_gamma. Header lines (# name=… M_w=… elements=…) carry the metadata the table format itself doesn't. 8 integration tests cover the round-trip on a synthetic table, file-format parsing of metadata, error paths (T_mid out of range, too few points), and a textbook diatomic recovery to better than 0.5 %.
  • Phase 27 — CHEMKIN mechanism reader. Ships aeolus1d.chemistry.reaction (Reaction + Mechanism dataclasses, modified-Arrhenius k(T) = A·T^b·exp(-Ea/RT) evaluator, stoichiometry-matrix builder), importers.chemkin_mech (parser for ELEMENTS / SPECIES / REACTIONS sections — handles = / <=> / => direction tokens, integer or floating-point stoichiometric coefficients, CAL/MOLE / KCAL/MOLE / JOULES/MOLE / KELVINS global Ea unit conversions to SI), and mech_db (native JSON-per-mechanism persistence — JSON because reactions are heterogeneous nested structures NPZ doesn't fit). The package ships a bundled ammonia_scr_demo mechanism in aeolus1d/chemistry/mech_db/. 12 unit tests + plot demo (parse → save → load round-trip; k(T) curves + stoichiometry heatmap + activation-energy bar chart). Continuation-line extensions (LOW / TROE / PLOG / three-body / DUPLICATE / REV) shipped in Phase 27b.
  • Phase 27b — CHEMKIN parser extensions. Adds full parsing of the rate-form extensions every published combustion mechanism (GRI-Mech, Konnov, Burcat, …) uses:
    • three-body reactions (+ M token, with per-species H2/2.0/ N2/1.0/ efficiency lists),
    • Lindemann / Troe / SRI falloff ((+M) marker plus LOW / TROE / SRI continuation lines),
    • pressure-dependent rate tables (PLOG rows interpolated in log-pressure),
    • explicit reverse rates (REV),
    • DUPLICATE markers for multi-channel reactions. The parsed extensions round-trip through the native JSON format unchanged. Ships an h2_oxidation_demo mechanism in the bundled DB exercising every rate type. 11 unit tests cover the parser and JSON round-trip.
  • Phase 27c — Full kinetics evaluation of all rate forms. The Phase-27b parser already produced Reaction objects with low_pressure / troe / sri / plog / rev / is_third_body / is_duplicate populated; this phase closes the loop so MechanismKinetics actually evaluates every form. Already shipped before 27c started: three-body [M]_eff = Σ α_i [X_i], Lindemann + Troe falloff, log-pressure PLOG interpolation. 27c adds:
    • SRI broadening F = d · T^e · F_cent^X with F_cent = a·exp(-b/T) + exp(-T/c) and X = 1/(1 + (log Pr)²) (CHEMKIN manual §3.5.4);
    • REV explicit reverse Arrhenius — when r.rev is set, k_r uses the user-supplied coefficients instead of k_f / Kc (with [M]_eff scaling for three-body REV);
    • DUPLICATE multi-channel additive rates — verified that paired duplicate reactions produce per-species rates equal to the sum of their individual contributions, which is what CHEMKIN's DUPLICATE marker requires. Both the numpy reference path and the @njit kernel (per_cell_rates_njit_full) handle every form, with the JIT path matching the numpy path bit-near-equivalent (1e-12 rel) across SRI / REV / REV+third-body parametric tests. 6 numpy + 3 JIT integration tests cover the new branches.
  • Phase 28 — Mechanism → solver wiring. Ships aeolus1d.chemistry.kinetics.MechanismKinetics that consumes a Phase-26 species DB + Phase-27 Mechanism and computes per-cell reaction rates q_j = k_f(T) · Π_i [X_i]^ν^R_{ij} with reverse rates q_r = q_f / Kc(T) derived from species enthalpies / entropies (Kc = Kp · (P_atm / RT)^Δν, Kp = exp(ΔS/R − ΔH/RT)). Per-species mass-fraction rate dY_i/dt = M_w_i · Σ_j ν^net_{ij} q_j / ρ and volumetric heat release s_E = − Σ_j q_j · ΔH_j(T) slot into the existing Phase-22b py_reaction_multi seam unchanged via the make_py_reaction_multi(tracked_names) factory. CHEMKIN-cgs → SI A-factor conversion uses (1e-6)^(n−1) with reaction order n; out-of-mechanism (untracked) reactive species are rejected at factory time. 10 unit tests + plot demo (0-D batch reactor on the bundled ammonia mechanism: NH3 + O2 fully consumed in 0.6 s, heat release sign / stoichiometric balance both correct).
  • Phase 28b — JIT kinetics kernel + variable-γ diagnostic. Numba @njit rate evaluator (aeolus1d.chemistry.kinetics_njit) that mirrors the pure-numpy MechanismKinetics.per_cell_rates math on flat arrays — bit-near-equivalent (1e-12 relative) at ~15-95× speedup depending on mechanism size. Opt-in via MechanismKinetics(..., use_jit=True). Adds aeolus1d.chemistry.mixture with diagnostic helpers: mixture_cp(T, Y), mixture_cv(T, Y), mixture_gamma(T, Y), mixture_R(Y) — fresh diatomic charge gives γ ≈ 1.40, polyatomic combustion products at 1500 K give γ ≈ 1.29, an 8 % gap that maps directly into the pipe solver's energy → T recovery error (the solver still uses the scalar pipe.gamma; full HLLC wiring is the natural Phase 28d follow-on). 11 unit tests + plot demo with side-by-side JIT vs numpy benchmark and γ_eff(T, Y) trace through the heat-release transient.
  • Phase 28d — Variable-γ scalar wiring (cell-averaged mixture). Pipe.use_mixture_gamma opt-in: each pre-stage refresh recomputes self.gamma from the cell-averaged Y vector, so closures / junctions reading the scalar see a γ that tracks the composition.
  • Phase 28e — Variable-γ in the HLLC flux kernel (first-order). Pipe.use_per_cell_gamma opt-in: the residual now uses each cell's own γ for pressure recovery / sound-speed / star-state algebra. hllc_flux_var_gamma(UL, UR, γL, γR) is bit-exact to hllc_flux for γL = γR. Restricted to the simple first-order ideal-gas path; MUSCL / MP5 / EOS / quasi-1D variants land in Phase 28f.
  • Phase 28f — Variable-γ in the higher-order reconstruction tier. Completes the variable-γ stack: pipe_residual_muscl_var_gamma and pipe_residual_mp5_var_gamma consume a per-cell γ_pad through the primitive-variable reconstruction (U_to_W_var_gamma), the face-state recovery, and the hllc_flux_var_gamma interface flux. pipe_residual_quasi1d_muscl_var_gamma handles the varying-area / closure case. Single-pipe advance_single_pipe dispatches through the variable-γ kernels when both use_per_cell_gamma=True and scheme="muscl" / "mp5", with a strict bit-exact fallthrough to the scalar kernels when either flag is off. 8 integration tests cover MUSCL / MP5 / quasi-1D parity, residual divergence under a real γ-step (1.40 ↔ 1.27), and end-to-end timeloop dispatch.
  • Phase 29b — Tabulated real-gas EOS (CoolProp-friendly). New EOS_TABULATED = 3 kind that packs a (T, ρ) → (p, e, ∂p/∂T, ∂p/∂ρ) table into the existing length-N params flat-array contract — same dispatch shape as RK / PR, no kernel signature changes. aeolus1d.closures.eos_tabulated.make_table_from_callables builds the table from any Python p(T, ρ) / e(T, ρ) pair; make_table_from_coolprop(fluid="CO2", T_min=…, T_max=…, …) wraps CoolProp's PropsSI (lazy import — install via pip install -e ".[real-gas]"). The @njit lookup kernel (aeolus1d.closures.eos_tabulated_njit) does bilinear interpolation in (T, ρ); pressure / internal-energy / sound-speed / T(ρ, e) inverse are wired through the standard eos_pressure_jit family, so any kernel that already takes (eos_kind, eos_params) works with tabulated EOS unchanged. Default 64×64 grid is ~128 kB per pipe and gives ~0.5 % bilinear- interp accuracy on a moderate-pressure box. 12 unit tests + 1 optional CoolProp roundtrip cover packing / metadata, on-grid exactness, off-grid bilinear accuracy vs RK, sound-speed FD derivatives, T(ρ, e) inverse round-trip, edge clamping, and a CO2 air-box CoolProp build (skipped when CoolProp isn't installed).
  • Phase 32 — MP5 in the network timeloop. Adds 5th-order reconstruction to the network path. Per-pipe scheme = "mp5" (requires Ng >= 3); boundary_face_flux_mp5_jit mirrors the same Suresh-Huynh stencil at the boundary face so junction fluxes match. Ideal-gas + uniform-area only (matches the existing pipe_residual_mp5 coverage). 5 integration tests.
  • Phase 33 — Variable-γ MUSCL in the network timeloop. Closes the per-cell-γ + MUSCL gap left by Phase 30. When scheme="muscl" AND use_per_cell_gamma are both set, the network residual routes through pipe_residual_muscl_var_gamma and the port-flux helper uses boundary_face_flux_muscl_var_gamma_jit. Mass conservation validated to ≤ 1e-9 across a junction with species-dependent γ. 5 integration tests.
  • Phase 34 — Test-suite parallelisation. pytest-xdist was already in the [dev] extra; this phase adds a Makefile with make test (parallel, skips @pytest.mark.slow), make test-all (full sweep), and make test-fast (with --durations=10 for spotting slow tests). Walltime drops from ~6 min to ~3 min on a typical workstation. Marks the slowest bench (Helmholtz quasi-2D, 20-period FFT integration) as slow.
  • Phase 35 — Validation suite. Introduces tests/validation/ for published-data / closed-form comparisons. Ships test_sod_against_analytical.py: a from-scratch implementation of the Sod-tube Riemann solution (Toro Ch. 4 algorithm with Newton on the star-state pressure), then asserts L1 density error against the numerical run for both first-order (≤ 10 %) and MUSCL (≤ 4 %) plus the qualitative MUSCL-better-than-first- order property.
  • Phase 36 — Sensitivity utilities. aeolus1d.sensitivity ships fd_gradient(objective, theta, h, kind) (central + forward finite-difference gradients) and parameter_sweep (a Cartesian-product sweep that returns a labelled grid) — the pragmatic first step for closure-parameter calibration, design of experiments, and ML loss-landscape inspection. Full AD through the @njit kernels is left as a future direction (Numba
  • autograd is limited; a JAX-mirror solver is its own project).
  • Phase 31 — Real-gas MUSCL in the network timeloop. Extends Phase 30's MUSCL plumbing to the EOS dispatch (RK / PR / tabulated CoolProp). The kernels live in hllc_eos_njit.py next to the existing first-order EOS ones: pipe_residual_eos_muscl_jit (uniform-area) and pipe_residual_quasi1d_eos_muscl_jit (variable-area + closures). Reconstruction is on W = (ρ, u, T) — temperature is the natural primitive for cubic EOS because W → U collapses to one eos_internal_energy call (closed-form for ideal gas, polynomial for RK / PR), avoiding a Newton iteration per face. The junction port-flux helper boundary_face_flux_eos_muscl_jit mirrors the same boundary reconstruction so pipe-side and volume-side fluxes remain identical. For eos_kind == EOS_IDEAL the kernel reduces to the regular MUSCL kernel to floating-point precision; for PR / RK / tabulated, the network now runs 2nd-order on refrigerant / cryo / CO2 pipework. 6 integration tests cover ideal-gas parity, PR-vs-first-order divergence, mass conservation across a CO2 plenum (≤ 1e-10 drift), quasi-1D dispatch, and RK dispatch.
  • Phase 30 — MUSCL in the network timeloop. The single-pipe advance_single_pipe already supported MUSCL / MP5; the network advance_network was first-order only because of pipe-vs-volume port-flux consistency concerns. Phase 30 adds a per-pipe Pipe.scheme flag ("first_order" default, "muscl" opt-in) and plumbs MUSCL through both the pipe residual and the volume's port-flux helper. The new _pipe_port_flux_muscl performs the same primitive-variable reconstruction at the boundary face that the network MUSCL kernel uses inside the pipe — so the pipe-side flux at a junction equals the volume-side flux exactly, and total mass across pipes + volumes is conserved to numerical precision (verified to drift ≤ 1e-10 relative on a wall-bounded two-pipe + plenum case). Quasi-1D pipes route through pipe_residual_quasi1d_muscl when scheme="muscl"; real-gas + per-cell-γ MUSCL combinations stay on first-order (their port-flux reconstruction is a follow-on). 8 integration tests cover the dispatch, mass conservation, and the MUSCL port-flux helper's collapse-to-first-order under uniform state.
  • Phase 37 — Cantera YAML mechanism import. New importer at aeolus1d.chemistry.importers.cantera_yaml; reads species (NASA-7 thermo, composition → M_w from atomic-mass table) and reactions (elementary Arrhenius, three-body, falloff with Troe, PLOG) into the same SpeciesThermo / Reaction dataclasses the existing CHEMKIN parser produces, so downstream MechanismKinetics works unchanged. Honours the per-mechanism units: block. [chemistry-yaml] extra (pyyaml). 12 unit tests + bundled H₂/O₂/H₂O/N₂ fixture.
  • Phase 38 — CGNS result writer. Optional CGNS export alongside the existing HDF5 (aeolus1d run --cgns OUT.cgns). One Structured Zone per pipe, single-cell zones for junctions, optional BaseIterativeData for the time-series trace. Opens natively in ParaView / Tecplot / EnSight / GT-SUITE Post. [interop] extra (pycgns) — lazy import; the writer raises a clear ImportError pointing to the extra when missing. 5 integration tests (skipped in installs without pycgns).
  • Phase 39 — Commercial 1D-tool case import / export. Best-effort structural translation between aeolus1d TOML and GT-Power / AVL Boost / Ricardo WAVE XML formats. A neutral Schematic intermediate (aeolus1d.io.interop.schematic) captures pipe geometry + junction connectivity + BC kinds every tool understands; per-tool modules (gt_power.py / avl_boost.py / wave.py) handle the XML element-name differences. Anything proprietary (turbo maps, fuel injection, valve cam tables) is dropped at import with a structured InteropReport listing what was lost. Symmetric: aeolus1d → tool drops aeolus1d-specific fields (real-gas EOS, MUSCL/MP5, LODI, chemistry, Wiebe combustion, learned closures) and reports them. CLI: aeolus1d import <tool> <in> <out.toml>, aeolus1d export <tool> <case.toml> <out>. 10 integration tests + a synthetic GT-Power XML fixture.
  • Phase 40 — FMI 2.0 FMU export. Wrap an aeolus1d case as a Functional Mock-up Unit for co-simulation with GT-SUITE / Simulink / OpenModelica / Dymola. aeolus1d.io.interop.fmu.write_fmu(case, out.fmu) packages the slave class + the case TOML into a single FMI archive via pythonfmu. CLI: aeolus1d export-fmu <case.toml> <out.fmu>. Inputs: per-BC stagnation / back pressures. Outputs: per-pipe end-state primitives + per-junction (p, T, m). 3 integration tests (skipped without pythonfmu).
  • Phase 45 — Catalyst / SCR / TWC closure stack. A 10-phase series that turns the Phase-15 thermal-only catalyst into a proper aftertreatment toolkit. 45a monolith spec resolver (CPSI / OFA / substrate → ρ_s·c_s, a_v, h, d_h); 45b variable-Nu Hawthorn h(Re, T) for flowing exhaust; 45c Pipe.surface_state registry for storage species; 45d NH3 surface coverage with adsorption / desorption / storage-SCR; 45e three-pathway SCR (Standard / Fast / NO₂-only) with the Olsson 2009 √(c_NO·c_NO₂) Fast kinetics; 45f NH3 oxidation side reaction (high-T selectivity knee); 45g TWC multi-step (CO + HC + NO with Voltz-Morgan-Liederman LH inhibition); 45h TWC Oxygen Storage Capacity (Ce/Zr) for λ-excursion buffering; 45i Thiele effectiveness factor for washcoat-diffusion-limited regime; 45j two-storage-site SCR (Cu²⁺ strong + Brønsted weak with two-peak TPD signature). Each closure ships standalone with closed-box scipy tests (saturation, equilibrium, conservation). 58 tests.
  • Phase 46 — Mechanism-driven catalyst (CHEMKIN-Surface + Cantera-YAML). End-to-end import of catalyst surface mechanisms into the solver. 46a CHEMKIN-Surface SITE / SDEN parser extension (mol/cm² → SI mol/m²); 46b Cantera YAML ideal-surface phase importer (handles both legacy cgs and modern SI units); 46c SurfaceMechanismKinetics — combined gas + surface stoichiometry matrix, per-cell Arrhenius rates, f(T, Y_gas, theta, rho) → (dY_gas, dtheta, s_E) callable; 46d [pipe.catalyst] kind = "mechanism" loader auto-detects format from extension and resolves per-species M_w from the species DB; 46e timeloop integration via the Strang-split chemistry sub-step (per-cell scipy LSODA on the joined (Y_gas, theta, E) state). The user can hand the solver a Cu-zeolite .inp or Cantera YAML and the rest is automatic. 35 tests.

Milestone 6 (closed) — full GUI control via sub-component dropdowns

The current Inspector exposes a flat form per component. The schema is heavily polymorphic via Literal["..."] "kind" fields and nested dataclasses (init.kind, area_profile, friction.kind, heat.kind, porous.kind, catalyst.kind, species.reaction, eos.kind, all the junction kinds, all the BC kinds, cylinder lift schedules, turbo maps), and most of those polymorphism axes are not reachable from the GUI today — power users hand-edit the TOML. This milestone closes that gap with a sub-component dropdown pattern, where each polymorphic field becomes a kind picker and the form below reshapes to the chosen kind.

  • Phase 47a — Form-factory v2: Literal-kind dispatch. The form factory recognises Literal["a", "b"] fields, renders them as a QComboBox, and reshapes the rest of the form when the kind changes. Generic; every later 47.x rides on it. Preserves shared fields across kind changes (e.g. T_solid_init survives a catalyst kind flip).
  • Phase 47b — Pipe init editor. All five [pipe.init].kind options (uniform / riemann / parabolic / regional / mode), including the regions: list[RegionSpec] table widget. Closes the original backlog item (regional init editor was hand-edit-only).
  • Phase 47c — Closure editors (friction / heat / porous / catalyst). Each becomes a kind dropdown + the right form. Catalyst gains the monolith form (cpsi / ofa / substrate combo) and the Phase-46d mechanism path picker.
  • Phase 47d — Area-profile + EOS editors. Two more polymorphic fields — parabolic / tanh-step area profile; ideal / Redlich-Kwong / Peng-Robinson / tabulated EOS.
  • Phase 47e — Species editor. Primary Y + the extras: list[SpeciesExtraSpec] table + reaction-callable selector (built-in demos / mechanism path).
  • Phase 47f — BC editor. All nine BC kinds, including transient_* with an inline (t, value) table widget shared with Phase 47h's cylinder-lift table.
  • Phase 47g — Junction editor. Plenum / Corberán / EHR picker; ports: list[PortSpec] table; per-port K_loss for Corberán; EHR neck geometry.
  • Phase 47h — Engine-side editors. Cylinder lift schedules; turbo map CSV importer + plot; valve orifice rigorous-form parameters; shaft-inertia ODE.
  • Phase 47i — Case header full editor. Reach scheme / limiter / RK / solver / chem_solver / chem_rtol / hllc_low_mach from the GUI (the Runner-dock strip exposes the common ones today; this is the rest).
  • Phase 47j — Round-trip + persistence audit. ✅ Every TOML field reachable from the GUI round-trips Case → save_case → load_case → Case cleanly. aeolus1d_gui/io/case_io.py now emits [pipe.friction] / [pipe.heat] / [pipe.porous] / [pipe.catalyst] (all three kinds, including surface_init inline-table) / [pipe.species] (extras list + reaction-kind parameters) / [pipe.eos], the missing transient-BC tables (p0_t/v, T0_t/v, p_back_t/v), the stagnation_inlet p_back field, and the case-header chem_solver / chem_rtol / chem_atol / hllc_low_mach. 28-test parametrised round-trip audit covers every polymorphic axis the GUI exposes (tests/integration/test_toml_roundtrip_47j.py).
  • Phase 47k — Component / sub-component palette polish. ✅ Palette now covers all nine BC kinds — including transient_stagnation_inlet / transient_pressure_outlet that were previously only reachable from hand-edited TOML. Form-factory kind dropdowns surface the dataclass docstring as a hover tooltip, so inspectors get an in-context summary of what each polymorphic kind models without leaving the form. Closes the Phase 47 series.
  • Phase 48a — schema_help registry + ? button infrastructure. ✅ New aeolus1d.io.schema_help module ships a ClassName.field_name → markdown registry seeded with the most-edited fields (CaseHeader, GasSpec, PipeSpec, UniformInit, closures, EOSSpec, BCSpec, JunctionSpec, PortSpec, SpeciesSpec). FormBinder renders a ? button next to every field that has a registry entry; clicking the button opens a modeless info window rendering the markdown body. Fields without an entry stay row-clean — no stub popups.
  • Phase 48b — Inspector field-help coverage sweep. ✅ Registry now covers 147 fields across 24 schema classes — every editable field on TurboSpec, TurboMapSpec, ShaftSpec, CylinderSpec, CylinderPortSpec, the four init kinds (UniformInit, RiemannInit, RegionSpec/RegionalInit, ModeInit), both area profiles (ParabolicProfile, TanhStepProfile), the full CatalystSpec (including all four kinds), the BC transient + valve fields, the PortSpec extras (neck, py_junction_bc, lodi), and the three polymorphic kind selectors (CatalystSpec.kind, JunctionSpec.kind, TurboMapSpec.kind). Tests pin both registry coverage and that every key points at a real dataclass field.
  • Phase 48c — Optimize-dock parameter-path help. ✅ New schema_help.help_for_path(case, path) walks a dotted- bracketed path against a loaded Case (reusing the optimizer's _parse_path / _resolve_index) and looks up the leaf field's help body. OptimizeDock now loads the full Case alongside the OptimizeSpec on file-open and populates each parameter-table path-cell's tooltip with the resolved help body — refreshes live on cell-edit, falls back to clear hints when the path is invalid or no case is loaded.
  • Phase 49a — SCR ads/des kinetics + chemistry-step surface integration. ✅ New aeolus1d.closures.scr_kinetics ships make_nh3_ads_des and make_nh3_scr_kinetics — surface- coverage closures with explicit per-rate scalar parameters reachable to the optimizer via pipes[0].species.A_ads etc. Pipe.py_reaction_surface is the new closure slot; _chemistry_step_surface_lsoda packs the per-cell state vector as (Y, *Y_extras, *theta_keys, E) and integrates via LSODA cell-by-cell. SpeciesSpec.reaction = "nh3_ads_des" / "nh3_scr_kinetics" are the new kinds; the loader auto- initialises Pipe.surface_state["theta_NH3"]. End-to-end optimize-recovery test pins that 50 % off A_ads, A_des estimates converge to within 10 % of truth in 50 evaluations via Nelder-Mead probe-match against a synthesised Y_NH3(t) trace. Tracer + probe extended with var = "Y" so any primary-species fit re-uses the same path.
  • Phase 49b — SCR rate-fitting example pack. ✅ Bundle of ready-to-run files at examples/optimize/scr/: two TOML optimization cases (scr_ads_des.toml for the ads/des-only fit recovering A_ads, A_des; scr_full.toml for the full SCR mechanism at T = 700 K recovering A_ads, A_des, A_scr with a Y_NO probe alongside Y_NH3). Synthetic measurement-CSV generator (generate_measurement.py) seeded at known truth params so users can verify the full workflow before swapping in real experimental data — CSV format documented in examples/optimize/scr/README.md. Three plot scripts: plot_scr_optimize_fit.py (ads/des) + plot_scr_full_optimize_fit.py (full SCR with Y_NH3 + Y_NO probes and 3 rate trajectories) + plot_scr_csv_preview.py (sanity-check a measurement CSV before fitting). Tracer + probe extended with var = "Y_extra" + extra_index so secondary-species traces (NO, etc.) feed the optimizer.

Implicit time integration stack (closed) — Phases I-1..I-9 + quasi-1D / var-step BDF2 / per-cell γ / real-gas EOS / py-closures

Single-pipe + network implicit time schemes that lift the explicit acoustic CFL ceiling. Motivation: the convective-only CFL analyzer (commit 8fb4d70a) flagged that low-Mach engine flows pay an ~20× acoustic penalty under the explicit RK schemes; the implicit family unlocks that headroom. Architecture + phase plan is in docs/implicit-solver-plan.md; it lives outside the Milestone-7 list because it's a solver- core addition rather than an engine-feature one. Empirically, on a stagnation-driven two-pipe + plenum case at moderate Mach, the BDF2 + adaptive-Δt run at CFL≈20 takes 31× fewer steps than the explicit RK2 reference at CFL=0.5 with mass-flow agreement to ~3%.

  • Phase I-1 — implicit-Euler scaffolding. Single pipe, first-order spatial (HLLC), uniform area + uniform γ, no closures. Block-tridiagonal Jacobian via central-difference perturbation of HLLC face fluxes; block-Thomas direct linear solve; Newton iteration with backtracking line search + positivity guard + Δt-halving off-ramp.
  • Phase I-2 — closures into the implicit residual. Native friction / heat / porous source terms folded into R(U); their Jacobian ∂S/∂U (FD per cell) augments the diagonal blocks. Sparsity stays block-tridiagonal because closures are local. Validated against analytic exponential-decay / Newton-cooling solutions at CFL > 1.
  • Phase I-3 — MUSCL with frozen reconstruction. 2nd-order spatial via primitive-slope edge reconstruction; slopes are computed once per Newton iteration and held constant (frozen) so the limiter doesn't widen the Jacobian band. Block-tridiagonal sparsity preserved, one extra Newton iteration on average vs. cell-center.
  • Phase I-4 — network coupling (block-split). Multiple pipes via plenum-ghost junctions. Pipe-side fully implicit (junction end is Dirichlet — chain rule = 0, ghost held fixed at start-of-step plenum state); plenum mass + energy update from boundary fluxes after pipe Newton solves. The follow-up commit added BC_NUMERICAL — a numerical chain rule via FD on apply_ghosts that captures arbitrary far BCs (stagnation_inlet, pressure_outlet, transient BCs). Without this, the BC chain was transmissive-approximated and Newton halved Δt back to the acoustic limit.
  • Phase I-5 — BDF2. 2nd-order time discretization on the same Newton + linear-solve plumbing — generic F(U) = U - U_hist - dt_eff · R(U) form unifies IE (U_hist = U_n, dt_eff = dt) and BDF2 (U_hist = β₁ U_n + β₂ U_{n-1}, dt_eff = h_n (1+ω)/(1+2ω), ω = h_n / h_{n-1}). The variable-step coefficients are Hairer-Wanner II §III.5: at ω = 1 they reduce to the classic fixed-step (4/3, -1/3, 2 h / 3); at ω → 0 they collapse smoothly to implicit Euler — that's what makes a Δt-halving graceful, no IE fall-back needed. First step is IE startup; from step 2 on, BDF2 runs through every dt change including post-halving and adaptive-Δt jumps.
  • Phase I-6 — adaptive Δt (PI controller). Newton-iteration- driven step-size control: easy step (≤1 iter) two times in a row → grow Δt; hard step (n_iters > target) → shrink via factor = (target/n_iters)^p; Newton failure → shrink hard. Final Δt clipped to [cfl_min, cfl_max] · dt_acoustic_unit. Available on all three implicit drivers.
  • Phase I-7 — TOML / CLI surface. rk = "implicit_euler" / rk = "bdf2" extends the existing schema field; new newton_tol, newton_max_iters, adaptive_dt, adaptive_cfl_min, adaptive_cfl_max, adaptive_target_iters knobs. The runner in aeolus1d.io.case dispatches implicit-rk single-pipe cases through the implicit drivers and implicit-rk network cases through advance_network_implicit. Save → reload round- trips the implicit fields; explicit cases stay unpolluted. Runner-dock RK combo + CFL Analyzer headroom banner make the implicit family discoverable from the GUI.
  • Phase I-9 — Radau IIA(3). Three-stage 5th-order fully- implicit Runge-Kutta. L-stable; stiffly accurate (so Y_3 = U^{n+1} directly). Per step solves a coupled (3 stages × N cells × 3 conservatives) system via simplified Newton (Jacobian frozen at predictor) and block-Thomas on 9×9 blocks. Same Kronecker-product trick that lets Radau inherit the existing block-tridiagonal sparsity from the IE/BDF2 path: J = I − Δt · A ⊗ J_spatial. Both the single-pipe driver and the network block-split path (advance_network_implicit) accept rk = "radau_iia". Compared to BDF2 at the same large CFL, Radau IIA's error is ~1000× smaller on smooth flow — the 5th-order vs 2nd- order temporal accuracy gap. Demo on a smooth pipe (t_end = 0.005 s, N = 80): explicit RK2 cfl=0.5 takes 299 steps; Radau IIA cfl=10 takes 15 steps + matches the explicit trace. Demo on a two-pipe + plenum network (stagnation → pipe → volume → pipe → outlet, t_end = 0.005 s): explicit RK2 cfl=0.5 takes 286 steps; Radau IIA cfl=10 takes 15 steps with mass-flow agreement to 0.5 %.
  • Quasi-1D in the implicit drivers. First-order + frozen- slope MUSCL flux divergence + geometric pressure source on the variable-area path. The driver gate that previously rejected not pipe.is_uniform_area is gone for IE / BDF2 / Radau, both single-pipe and network. Area-aware face fluxes scale by A_face[f] / A_cell[*]; the geometric source p · (A_R - A_L) / (A_c · dx) adds a closed-form contribution to the diagonal Jacobian blocks. For is_uniform_area=True the new code path reduces bit- identically to the previous one. Unlocks the implicit-CFL headroom for engine cases with CD nozzles, runner area variation, and port flow.
  • Per-cell γ in the implicit drivers. Multi-species reacting flow runs with pipe.use_per_cell_gamma=True and a fresh pipe.gamma_pad array refreshed once per outer step from the species composition. The implicit drivers now accept this — the HLLC face flux + frozen- slope MUSCL edge reconstruction route through hllc_flux_var_gamma (with γ held constant during Newton); the closure source uses cell-centre γ for pressure recovery; the geometric source's closed-form Jacobian contribution uses the cell-centre γ. For uniform γ_pad the new var-γ branches are bit-exact to the scalar fast path. Unlocks BDF2 / Radau headroom for the multi-species reacting cases that previously had to fall back to explicit RK.
  • Python closure overrides in the implicit drivers. Pipes with py_friction / py_heat / py_porous callables set now run on IE / BDF2 / Radau. Each override is folded into _cell_closure_source (called per cell with single- element arrays so the override's array-API signature is preserved) and the existing diagonal-block FD augmentation picks up the contribution automatically. Override takes precedence over the matching native kind (no double-count). Unlocks ML-driven closure development under implicit time.
  • Real-gas EOS in the implicit drivers. Cubic-EOS (Peng-Robinson, Redlich-Kwong) and tabulated CoolProp- style EOS now route through IE / BDF2 / Radau. Face flux uses hllc_flux_eos_jit; cell-centre pressure (closure source + geometric source) routes through primitives_from_U_eos_jit (which inverts e(T, ρ) for T via Newton + bracketing); MUSCL frozen-slope edge reconstruction uses EOS recovery for the (ρ, u, p) primitives. The geometric-source Jacobian drops to FD per cell (closed-form ∂p/∂U is ideal- gas-specific). Real-gas + per-cell γ is rejected — γ is irrelevant on the EOS path (pressure comes from (T, ρ)). Unlocks BDF2 / Radau headroom for CH₄ / NG / supercritical / cryogenic cases that previously had to use explicit RK because the implicit driver silently used ideal-gas (γ-1) ρ e.
  • Variable-step BDF2 coefficients. Hairer-Wanner II §III.5 variable-step BDF2 replaces the fixed-step (4/3, -1/3) history blend. Δt-halvings on Newton failure no longer trigger an implicit-Euler fall-back — the variable-step formula handles arbitrary step ratios and reduces smoothly to IE in the ω → 0 limit. Adaptive Δt + BDF2 now keeps formal 2nd-order accuracy across every dt change (the pre-existing fixed-step path lost order on every halving). Single-pipe + network drivers both updated.
  • Phase I-Y — neck pipes, LODI, and py_junction_bc in the network driver. Three previously-blocked port configurations now work under advance_network_implicit:
  • Neck-pipe junctions (neck_cells > 0) — neck pipes are collected from their owning junction and advanced as additional BC_DIRICHLET pipes alongside the user pipes. Their ghosts are filled once per step by apply_port_bcs (main ghost ← neck interior; neck far ghost ← plenum; neck near ghost ← main interior) and held frozen during Newton.
  • LODI ports — the Poinsot-Lele characteristic ghost is filled once per step from the start-of-step pipe interior via _fill_lodi_ghost; O(Δt) split, consistent with the existing plenum coupling.
  • py_junction_bc callbacks — Python port callbacks are called once per step (same O(Δt) split). 3 new integration tests cover all three paths.

  • Species transport + chemistry in the implicit path. All three implicit schemes (IE, BDF2, Radau IIA) now carry multi-species and chemistry via Lie-split operator passes appended after each Newton step:

  • Step 5 — catalyst Lie split. catalyst_rates_into overlays gas-solid heat exchange onto the energy conserved variable and advances T_solid with an explicit Euler step (operator-split from the pipe Newton).

  • Step 6 — DPF Lie split. cake_into adds cake back- pressure to momentum + energy; regen_into advances soot loading and cake temperature (T_cake) and overlays exothermic heat release onto the energy equation.
  • Step 7 — species transport Lie split. HLLC-upwind explicit sub-stepping of d(ρY)/dt. At implicit CFL > 1 the species advection CFL would violate stability; n_sp = ceil(Δt / Δt_acoustic) micro-steps at sp_dt = Δt / n_sp each cap the species CFL at ≤ 1. Both primary species (Y) and extra species (Y_extras) are advanced; DPF soot loading picks up the convective sink.
  • Step 8 — chemistry Lie split. LSODA / Bader-Deuflhard per-cell ODE integration of py_reaction / py_reaction_ multi / py_reaction_surface callbacks — identical to the Strang-split explicit path. Surface-coverage reactions (e.g. twc_brinkmeier) route through _surface_mechanism_step_on_pipe.

This completes the root-cause fix for avl_cruise_m cases using rk = "bdf2" + reaction = "twc_brinkmeier" where species concentrations were silently frozen (no transport, no chemistry) under the implicit path. All four passes are guarded by lazy imports — pure-gas cases pay zero overhead.

Chemistry stays operator-split — see the doc for the explicit boundary between this stack and aeolus1d.solver.chem_bs_jit.

  • Phase I-Z (planned) — avl_cruise_m implicit end-to-end smoke test. Run the full local_cases/avl_cruise_m case with rk = "bdf2" and confirm CO/NO species concentrations evolve, TWC conversion is non-zero, and T_solid advances from its initial value. Validates the root-cause fix for the silent species freeze.
  • Phase I-Strang (planned) — Strang split for species + chemistry in the implicit path. Upgrade the current Lie split (T(Δt) ∘ C(Δt) — first-order) to a symmetric Strang split (T(Δt/2) ∘ C(Δt) ∘ T(Δt/2)) for the species transport and chemistry operator passes. Recovers second-order accuracy in time, matching the BDF2 / Radau IIA temporal order for the reacting species sub-system. Catalyst and DPF remain Lie-split (their coupling is O(Δt) by the block-split component architecture anyway).

Bundled example: examples/40_implicit_adaptive_pipe.toml (BDF2 + adaptive Δt on a porous-plug pipe; reaches t_end in ~40 steps with the controller picking CFL up to ~20×).

Architecture refactor — closure-kind registry (in progress)

Phased migration of closure-kind dispatch from hand-maintained parallel lists (schema Literal[...] arms, KIND_FIELDS dicts, _attach_closures elif chains, kernel integer constants) to a single aeolus1d.closures.registry source of truth. The Numba @njit dispatch chains stay (Numba can't compile a registry lookup), but everything above the kernel — schema, TOML loader / saver, inspector form factory — reads the same registered ClosureKind records. Plan: docs/architecture-refactor-plan.md Phase 1.

  • Phase A1 — friction family migrated. register_closure declarations in src/aeolus1d/closures/friction.py, _attach_closures switched to code_for_name("friction", ...), and test_closure_registry.py pins the schema-vs-registry consistency (Literal arms, KIND_FIELDS, code round-trip, params-slot layout).
  • Phase A2 — heat family migrated. Same pattern as A1: constant_twall / none registered with the user-facing field tuple ("h", "T_wall", "D_h") (the trailing params[3] = R_gas slot stays an injected derived value, not a registered field). HeatSpec carries the matching Phase-A2 comment; consistency tests mirror the friction ones.
  • Phase A3 — porous family migrated. darcy_forchheimer / none registered with the six user-facing fields (alpha, beta, mu_ref, viscosity_model, sutherland_T_ref, sutherland_S). The viscosity_model Literal is encoded into a float-coded slot at params[3] by _attach_closures and the trailing params[6] = R_gas slot is injected from gas.R_gas; both encodings are hidden from the registry, which carries only the schema- surface tuple. PorousSpec docstring grew a Phase-A3 note pinning the schema-vs-registry contract.
  • Phase A4 — catalyst family migrated. Four schema kinds registered: none / two_phase_heat (Phase 15) / two_phase_heat_var_nu (Phase 45b Hawthorn correlation) / mechanism (Phase 46d, mechanism-driven surface kinetics). The mechanism kind is registry-only: a new CATALYST_MECHANISM = 3 integer code lives in the registry so the inspector form factory and TOML round-trip know about its fields, but the @njit catalyst dispatch has no branch for it — _attach_closures routes mechanism through _attach_mechanism_catalyst and overwrites the kernel-visible catalyst_kind back to none so the catalyst path stays a no-op. The Phase-45b kind's registered fields are the user-input monolith parameters (cpsi, ofa, substrate, k_gas, nu, mu_gas, Pr, T_solid_init); the resolved (rho_s_cs, a_v, d_h, …) block stays a kernel-side detail produced by resolve_monolith_var_nu_params.
  • Phase A5 — eos family migrated. Three schema-visible kinds registered: ideal / redlich_kwong / peng_robinson. The kernel-only EOS_TABULATED = 3 code (Phase 29b) is the inverse of catalyst's mechanism — selected via direct programmatic construction in closures/eos_tabulated.py, never via TOML — so it is intentionally absent from the registry. A dedicated test_eos_tabulated_is_kernel_only consistency test pins this contract. The schema-side fields tuples describe the user-input layouts (ideal: R_gas / cv; redlich_kwong: + Tc / Pc; peng_robinson: + ω); the kernel params block is constructed by make_ideal_params / make_rk_params / make_pr_params inside _attach_closures.
  • Phase A6 — reaction family migrated. The largest closure family (9 kinds across 5 implementation modules) is structurally different from the per-pipe closures: Pipe has no reaction_kind: int attribute, and _attach_reactions dispatches on the schema string directly to build Python callables stored on pipe.py_reaction / py_reaction_multi / py_reaction_surface. So the registry's value here is consistency-pinning + driving the inspector form factory, not collapsing an integer dispatch chain — the registry codes are synthetic (0..8). Registered kinds: none / arrhenius_scr_demo / arrhenius_scr_multi _demo / nh3_ads_des / nh3_scr_kinetics / hc_oxidation / doc_oxidation / twc / twc_brinkmeier. Registrations centralised in closures/reaction_demo.py for simplicity (build.py eager-imports it); a future migration may distribute them to each kind's implementation module. Replaces the schema's 38-line hand-maintained KIND_FIELDS dict with a closure_kind_fields("reaction") lookup pinned by the consistency test.
  • Phase A7 — bc family migrated. Same string-keyed pattern as reaction: _make_bc (in aeolus1d.io.build) returns Python callables stored as pipe.bc_left / pipe.bc_right — no pipe.bc_kind: int. Registry codes are synthetic (0..12). Thirteen BC kinds registered: wall / transmissive / stagnation_inlet / pressure_outlet / nonreflecting / transient_stagnation_inlet / transient_pressure_outlet / valve_inlet (Phase 25 Chalet-Chesse-Hetet) / valve_orifice (Phase 25b Winterbone-Pearson) / mass_flow_inlet / velocity_inlet (Phase 51c-fix) / transient_mass_flow_inlet / transient_velocity_inlet. Registrations live in aeolus1d/bc/__init__.py so they fire automatically whenever any aeolus1d.bc submodule is imported (Python runs the parent package's __init__ first regardless of which submodule was asked for). A test_bc_dispatch_is_string_keyed consistency test pins the synthetic-codes contract in the same shape as the reaction test.
  • Phase A8 — junction family migrated. Three kinds registered: plenum (Milestone-1 0-D Volume + plenum-ghost coupling), corberan (Phase 14b algebraic CP+h N-branch junction with no chamber state), ehr (Phase 14d Extended Helmholtz Resonator with per-port u_n state). Class-keyed dispatch like reaction / bc: the loader instantiates a different Python class per kind (Junction / CorberanJunction / EHRJunction) and appends them to the network's junction list — registry codes 0..2 are synthetic. Registrations live in aeolus1d/components/__init__.py; a dedicated test_junction_dispatch_is_class_keyed consistency test pins the three-distinct-classes contract.

Phase 1 of the docs/architecture-refactor-plan.md is now structurally complete — all 8 closure families (friction / heat / porous / catalyst / eos / reaction / bc / junction) read from a single aeolus1d.closures.registry source of truth. Future closure additions edit one schema Literal[...] arm + one register_closure call instead of the historical 3-7-place edit pattern. - Phase A9 — py_junction_bc family migrated. The Phase-17 ML-seam selector on PortSpec (three kinds: none / resistive_demo / learned_npz) follows the same string- keyed pattern as reaction / bc / junction. _resolve_py _junction_bc returns a Python callable per kind (None for none, an inline closure for resistive_demo, a deserialised LearnedJunctionBC for learned_npz). Registry codes 0..2 are synthetic. Registrations live in aeolus1d/components/__init__.py next to the junction family. - Phase A10 — turbo family migrated (turbo_map + turbo). turbo_map (TurboMapSpec.kind) selects the (PR, N_corr) → (ṁ, η) lookup the Compressor / Turbine reads — three kinds (synthetic_centrifugal / synthetic_turbine / learned_npz) with kind-specific KIND_FIELDS pinned by the registry. turbo (TurboSpec.kind) selects between Compressor / Turbine instances — class-keyed, no KIND_FIELDS dict (every field on TurboSpec is relevant to both kinds). Both register synthetic codes; both fire from aeolus1d/components/__init__.py. - Phase A11 — pipe-state tagged unions (init + area_profile). Each kind has its own dataclass (rather than a single dataclass with a multi-arm Literal), so the schema shape differs from the prior families: every dataclass carries a single-arm Literal[<name>] discriminator and the registered fields are the dataclass's fields minus kind. Two families register: init (4 kinds — uniform / riemann / regional / mode) dispatched by _apply_init via isinstance; area_profile (3 kinds — parabolic / tanh_step / tabulated) dispatched by _area_profile_array the same way. Registrations live in aeolus1d/io/__init__.py. The consistency tests cross-check each registered kind against its dataclass's actual dataclasses.fields() so drift between the registry and the dataclass field list fails loudly.

This closes the registry-coverage migration: every TOML- visible discriminated kind across the schema (closure families, BCs, junctions, components, pipe-state tagged unions) is now registered. The schema's hand-maintained KIND_FIELDS dicts are pinned against the registry by the consistency suite, and adding a new kind to any family means one schema edit + one register_closure call. - Phase A12 — Snapshot value type + tracer rewrite. Phase 2 of the architecture refactor. aeolus1d.io.snapshot now exports a frozen Snapshot dataclass and a take_snapshot(t, pipes, junctions) walker — the single source of truth for "what is the saveable state of a network at one moment in time". The writer takes one Snapshot for the final-state HDF5 group; the tracer holds a list[Snapshot] and stacks per-field slices on demand for its accessor methods (pipe_rho / pipe_Y / pipe_surface_state / junction_p ...). The tracer's previous nine parallel per-field accumulator dicts (samples_rho / samples_u / samples_p / samples_Y / samples_Y_extras / samples_surface / samples_vol_p / samples_vol_T / samples_vol_m) are gone, replaced by snapshots: list[Snapshot]; samples_t survives as a backward-compat property returning [s.t for s in self.snapshots]. preload_from _prior reconstructs prior snapshots from the PriorSeries and prepends them. The HDF5 layout is byte-identical (every existing trace plot / resume test passes unchanged); the win is in the abstraction — adding a new field of state is now one edit in snapshot.py and both the writer and the tracer pick it up automatically. - Phase A13 — Versioned-checkpoint upgrade framework. Phase 3 of the architecture refactor. The previously inline schema-version check in aeolus1d.io.resume is promoted to a real subsystem in aeolus1d.io.schema_upgrades: a CURRENT_SCHEMA_VERSION constant, a MIN_SUPPORTED_VERSION floor, an UpgradeStep frozen dataclass, and a resolve_upgrades(from, to) resolver that walks the UPGRADES chain. Today the chain is empty (the format has been stable at v1 since the resume- checkpoint refactor in commit f6557bfc), but the framework is in place — the next breaking layout change appends one UpgradeStep and a unit test, and the resume reader picks up the migration automatically. The resume reader now also accumulates version_notes on the ResumeState (no-attr fallback, applied upgrade steps); the CLI runner echoes them alongside the per-pipe / per- junction seed notes so the operator sees what the upgrade reader did. test_schema_upgrades.py pins the framework's edge cases (identity, gap, downgrade-rejection, future-version-refusal, well-formed-chain). - Phase A14a — _cmd_run decomposition. First half of Phase 4 of the architecture refactor (the second half — MainWindow decomposition — is GUI-side and tracked separately). The historical 295-line _cmd_run god- function is now a thin orchestrator (74 lines) over six composable helpers: _read_resume (resume HDF5 + version notes); _build_live_control (LiveControl construction); _build_eager_network (Clock + build_network + resume seed); _build_tracer (Tracer construction + prior- series preload); _StepHooks class composing the per-step on_step callback (progress emit, info-every-N, record_if_due, live-flush HDF5 rewrites); and _print_run_summary / _write_run_outputs for finalisation. The trace / no-trace / resume duplication is collapsed: any path that needs an eager network (trace OR resume) goes through _build_eager_network + dispatch_advance; the no-trace no-resume path keeps going through run_case for backward simplicity. test_cli_live_flush's upper-bound on flush count was loosened from 25 to 30 — the refactor's per-step dispatch has marginally different overhead which shifts how many trace samples land within each 5 ms wall-clock window; still well below the wall-clock-only ceiling (~40). - Phase A15 — Component protocol load-bearing in consumer signatures. Phase 5 of the architecture refactor. The aeolus1d.components.protocols.Component Protocol existed pre-fork as documentation-only; the timeloops and IO layer still typed their junctions parameters as list[Junction] even though the runtime list contains a mix of Junction / Cylinder / Compressor / Turbine / Shaft. This commit replaces those consumer-side annotations with list[Component] (or Iterable[Component]) across the network timeloops (advance_network / advance_network_strang / advance_network_implicit + their helpers), the IO layer (writer / cgns_writer / tracer / snapshot / resume / case / build return type), and the public run_case / build_network signatures. The Protocol's volume attribute also tightened from object to "Volume" (forward-typed via TYPE_CHECKING) so consumer code accessing concrete Volume attrs (.set_state / .m / .E_internal / ...) type-checks under mypy strict — the existing test_strict_modules_pass_mypy test catches drift. Producer-side (_build_* factories that construct concrete Junction / Compressor / etc. instances) keeps the concrete annotations. - Phase A16 — Test tiers (smoke / fast). Phase 7 of the architecture refactor. Two nested marker tiers give fast feedback for refactor loops: smoke (≤ 10 s wall, ~8 canonical happy-path tests — Sod, mass-flow inlet steady-state, two-pipes plenum init, CLI run-and-write- HDF5, closure-registry, Component protocol, checkpoint versioning, mypy strict) and fast (≤ 60 s wall, every test under tests/unit/ and aeolus1d_gui/tests/ plus the smoke tier — ~410 tests). The default unfiltered pytest continues to run the full ~6 minute suite for pre-commit / CI. Markers are auto-applied via two pytest_collection_modifyitems hooks (one in the repo- root conftest.py for path-based fast tagging, one in tests/conftest.py for the curated smoke ID list) — no per-test @pytest.mark.* decorators required. The "Test tiers" section in docs/architecture.md documents the curation criteria and budget targets. - Phase A17 — checkpoint_and_stop live-control field. Phase 8 of the architecture refactor's "Done when" criterion: cancel-with-checkpoint produces a resumable HDF5 within ≤ 1 s. The pre-existing file-based MVP IPC in aeolus1d.runtime.live_control (Phase-A8) already handled CFL / paused / adaptive_cfl_max live updates; this commit adds the missing checkpoint_and_stop field. LiveControl.should_stop() returns True after the next poll once the GUI / external script writes {"checkpoint_and_stop": true, "_revision": ...} to the control file. The three timeloop families (advance_network / explicit single-pipe advance / advance_network_implicit) all break out of their inner step-loops on the flag, returning the partial-run t_final; the CLI's _write_run_outputs then writes a resumable HDF5 with the partial state. wait_while _paused also short-circuits on the flag so a paused run doesn't deadlock against a "stop now" click. Four new unit tests pin the contract.

The full QLocalSocket / named-pipe upgrade described in the original Phase-8 plan is still queued — file-watch latency (~100 ms in practice) hasn't been a bottleneck yet, and the file-based protocol exercises the same semantics so a future socket variant slots in without rewiring the CLI / timeloop sides. - Phase A18 — concepts.md update + architecture.md cross-link. Phase 9 of the architecture refactor. The pre-fork docs/concepts.md was already a complete end-to-end runtime narrative (TOML → Case → network → timeloop → HDF5 → resume → GUI), but it was stale w.r.t. the recent A12-A17 phases. This commit refreshes the Snapshot / schema-upgrade / Component / live-control / test-tier sections so each names the load-bearing module + recent phase and points new contributors at the right edit sites. Adds a "Testing" section pointing at the smoke / fast tiers, and updates the "adding new things" cookbook for the CLI subcommand registry. Cross-link from docs/architecture.md to concepts.md makes the 30-minute onboarding read the canonical entry point. - Phase A19 — GitLab CI pipeline + expanded strict-mypy list. Phase 6 of the architecture refactor's CI half. The pre-fork pytest gate (test_strict_modules_pass _mypy) ran mypy strict against 7 modules at test time; this commit adds a real .gitlab-ci.yml so a regression fails the pipeline (not just the test suite) and grows the strict-typed list from 7 → 12 modules. Five new modules added: aeolus1d.bc (the BC family registrations), aeolus1d.components (junction + py_junction_bc + turbo registrations), aeolus1d.components .protocols (the Component Protocol that Phase A15 made load-bearing), aeolus1d.io (the package re-exports), and aeolus1d.io.schema_upgrades (Phase-A13's upgrade framework). All five passed mypy strict already; only aeolus1d.io needed an __all__ annotation on aeolus1d.io.case for explicit re-export of build_network / load_case.

Three CI stages: lint (ruff + mypy strict against the 12 modules — runs in ~30 s and fails first on a typo), test (smoke + fast + unit suite, ~60 s wall), and integration (the slow integration / validation suites + the GUI suite). The integration stage gates on CI_PIPELINE_SOURCE == "merge_request_event" / CI_COMMIT_REF_PROTECTED so feature-branch pushes don't pay for the ~6-minute integration run on every commit; the merge request itself fires it. The GUI suite further gates on changes: against aeolus1d_gui/** so a pure-solver MR doesn't pay for the Qt install. - Phase A42 — Per-sample cylinder V + compressor / turbine / shaft partition in the HDF5 reader. Two follow-ups to Phase-A41:

Per-sample V. The tracer's junction_state_dict already recorded V on every snapshot, but the writer only emitted it on the final-state group. Tracer.junction_V(idx) is added (mirroring junction_p / T / m) and the writer now creates a time_series/junctions/junction_NN/V dataset alongside the existing p / T / m series. The GUI's network view derives rho = m / V per sample for cylinder colour instead of the previous m / V_final approximation — cylinders' kinematic V(θ) swing now feeds the rho display exactly. Plain plenums emit a constant V series (tiny storage cost; simpler reader path); older v2 traces without the dataset fall back to the constant-V denominator.

Compressor / turbine / shaft partition. Phase-A41 already tagged every junction-list entry with kind and id, but the GUI reader only routed "cylinder" to its own dict — "compressor" / "turbine" / "shaft" fell through to run.junctions keyed by their indexed name. Reader now grows RunResult.compressors / .turbines / .shafts (and matching TimeSeries partitions) populated by a ComponentSnapshot dataclass (kind / id / p / T / rho / m / V); the writer's _component_kind_and_id was already emitting the right tags so no schema bump is needed. Compressor / turbine / shaft state is dummy-Volume constants today (their volume is a placeholder for timeloop compatibility), but the partition exists so future per-kind diagnostics (last_mdot on the turbo, omega on the shaft) can land here without a second schema bump. test_h5_reader.py gains an end-to-end test_writer_partitions_compressor_turbine_shaft against examples/38_engine_turbo_scr.toml pinning the new partition. - Phase A41 — Cylinder identification in HDF5 + GUI body colour from cylinder volume state. First real bump of the Phase-A13 versioned-checkpoint upgrade framework: schema v1 → v2. The writer (aeolus1d.io.writer) now tags each junctions/junction_NN final-state and trace group with kind ("junction" / "cylinder" / "compressor" / "turbine" / "shaft") and id (user-facing) string attrs so downstream readers can partition cylinders out of the flat junction list. Cylinder gains a runtime id: str field copied from CylinderSpec.id by the build path. Schema log gets a v2 entry; UPGRADES registers a UpgradeStep(1, 2, …) describing the additive layout change (v1 readers see the unchanged dataset shape; resume from a v1 checkpoint loses the cylinder ↔ junction distinction in the GUI display only — simulator state is unaffected). The GUI h5_reader grows a CylinderSnapshot dataclass and a RunResult.cylinders[cyl_id] partition (plus TimeSeries.cylinders); the network-view dock picks up each cylinder body + piston tick's colour from volume.p / T / rho directly when v2 data is available, and falls back to the port-pipe-average heuristic for legacy v1 files. End-to-end round-trip verified on examples/31_cylinder.toml: file lands at v2, cylinder state visible as run.cylinders["cyl1"], time-series pressure swing tracked across the slider sweep. - Phase A40 — Cruise M import warning sweep. Cuts the reference-fixture import warnings from 17 to 6, each remaining one describing a real gap (now in BACKLOG.md). Three groups:

Group A — silence info-confirmations and superseded-by- another-path warnings: * Drop "GeneralSpeciesTransport: applied" — confirmation, not a warning. * Drop "twc_brinkmeier: 11 reactions modelled — R1 ... R11" — closure documentation; visible by inspecting the catalyst's reaction list. * Silence ClassicSpeciesTransport when GeneralSpeciesTransport was successfully applied (Phase A21 path supersedes it). * Silence EnthalpyOutputAdaption (~100 J/kg accuracy delta with no actionable user remediation). * Silence empty gas_composition.PassiveSpecies / ParticulateMatterSpecie / Species containers (no data lost). * Silence empty PassiveSpeciesTransport / ParticulateMatterSpeciesTransport blocks via a new _is_empty_species_transport helper that walks the block's names container.

Group B — replace heuristic noise with hard data: * Drop the "Ω_OSC = storage × 1000 × ρ_washcoat = ..." confirmation; the math is a well-defined unit conversion, not a heuristic. * Suppress per-AUCI-block "AUCI rate constants for X — using sibling catalogue Y" lines when a .ucp is present — the per-catalyst-pipe note already names the catalogue.

Group C — surface concrete data instead of generic warnings: * TransferModelMultiplierConstDatabus warning now lists the actual multiplier values (e.g. [0.9435, 0.9793, 1, 1.55]) so the user can find them in the GUI and apply manually. New _extract_multiplier_const_value helper pulls each block's numeric const value via direct XML path navigation.

Group D — the remaining warnings each describe a real gap with no quick fix, queued in BACKLOG.md under "Interop — Cruise M import follow-ups (Phase A40+)": * Particulate-Matter species column + DPF component (Milestone 8 Phase 8f, 2-3 weeks). * RamFlow ram-air pressure correction (1 day, low priority). * Per-reaction rate-multiplier wiring (1-2 days, requires .proj reaction-to-multiplier graph mapping). * AUCI Darcy-Forchheimer α / β extraction (depends on user-supplied datasheet or a SiederTate-style transport .ucp catalogue).

_surface_untranslated_elements moved to after the catalyst-processing block so it can read the general_transport_applied / passive_species_applied flags set by the GeneralSpeciesTransport / passive-species paths and silence superseded warnings cleanly. test_untranslated_elements_audit_trail updated to assert the right reduced set; test_brinkmeier_reaction_list_in_report flipped to assert the message is no longer present (it was pinning the old confirmation text).

46 Cruise M + ucp + tune integration tests stay green.

  • Phase A36 — Cruise M .ucp catalogue importer. Adds end-to-end parsing of an AVL Cruise M .ucp archive (the source-code-shipped catalogue describing one mechanism's reactions + parameters). When a .ucp sits next to the .proj, the AUCI comparison dialog gains an extra column with exact, name-based matches to aeolus1d's Brinkmeier reactions — no more heuristic Eₐ-distance guessing.

New module aeolus1d.io.interop.ucp:

  • parse_ucp(path) -> UcpCatalogue walks the .ucp's sources/auci_cat.c and emits one UcpReaction per USR_Init_Reaction_RXXXX block with name + chemistry comment + heat of reaction (J/mol) + parameter list (label
    • unit + value).
  • match_brinkmeier_to_ucp(catalogue) resolves aeolus1d's 11 Brinkmeier reactions to .ucp ones via a normalised-name alias table ("CO Oxidation" ↔ "CO_oxidation", "Water-Gas Shift Reaction" ↔ "WGS", etc.). Reactions without a .ucp counterpart stay unmatched.
  • find_ucp_alongside(proj_path) — auto-detect a sibling .ucp in the same folder.

Wiring:

  • compare_auci_rate_constants auto-detects + parses any sibling .ucp; CatalystComparison gains an ucp_catalogue field carrying the parsed object.
  • auto_match_brinkmeier_via_ucp(catalogue) performs the name-based match: K1 → A_brnk[i], E1 → Ea_brnk[i] (×R for K → J/mol), HReac → dH_brnk[i] (sign-flipped to aeolus1d's exothermic-positive convention). All matches tagged "strong" confidence — no heuristics involved.
  • AuciCompareDialog adds a .ucp column when the catalogue is attached, with its own Select button. Source token "ucp" flows through collect_selections and apply_selections like the AUCI-block sources.

On the reference fixture (TWC_SC630_wct.ucp, 24 reactions, 85 parameters total), 8 aeolus1d Brinkmeier reactions get exact .ucp matches: H2 / CO / CH4 / C3H8 oxidations, NO+CO redox, WGS, CH4 + C3H8 steam reforming. The 3 NO+H2 / NO+C3H6 / NO+C3H8 reactions stay at aeolus1d literature defaults (no .ucp equivalent in this catalogue).

Nine new tests in test_ucp.py pin: parser walks auci_cat.c correctly (synthetic minimal fixture); rejects non-ZIP / missing-source inputs; alias resolution; find_ucp_alongside finds sibling .ucp; against the real fixture: 24 reactions parsed, 8 expected aeolus1d indices matched, compare_auci_rate_constants attaches the catalogue, auto_match_brinkmeier_via_ucp pulls real K1 / E1·R / -HReac values for CO Oxidation. 506 GUI tests + 29 tune+ucp tests stay green.

  • Phase A35 — GUI AUCI comparison dialog + auto-pop on Cruise M import. AuciCompareDialog lives in aeolus1d_gui/shell/auci_compare_dialog.py and renders Phase A33's Comparison as a per-row interactive table: one row per aeolus1d closure parameter, two Select push- buttons (aeolus1d default / project value) in an exclusive QButtonGroup, and a combobox listing every AUCI scalar across every attached block. The user clicks the project Select button on the rows they want to override; on Apply, the dialog emits selections_applied with the same dict shape apply_selections consumes, and the host writes those values onto the imported case.

Wiring: MainWindow._show_auci_compare is called after the interop-report dialog closes whenever the user imports a .proj via the AVL Cruise M action. Inert (returns silently) when the comparison has no catalysts or no AUCI scalars to choose from. apply_selections(case, dict) is the new bridge between GUI dict shapes and the existing per-element / scalar override machinery in avl_cruise_m_tune; the file-based apply_overrides delegates here too.

Default state is "every row uses aeolus1d literature defaults", so a user who hits Cancel (or Apply without flipping any rows) leaves the case exactly as read() produced it. Multi-block catalysts (the typical case — a brick with rate constants in one .ucp block and ageing data in another) surface every block's scalars in each row's combobox so all candidate values are reachable.

Six new tests in test_auci_compare_dialog.py pin: per-row layout (two Select buttons + combobox), default-state empty selections, project pick + click writes the right scalar value into the dict, every block's scalars surface in the combobox, the no-pick guard prevents locking a row into a "no value" state, and the dialog's signal payload feeds cleanly back through apply_selections.

  • Phase A33 — AUCI rate-constant comparison + override workflow. Replaces the read-only "AUCI rate constants … NOT IMPORTED" warning with an actionable file-based pipeline. aeolus1d.io.interop.avl_cruise_m_tune exposes:

  • compare_auci_rate_constants(proj) — returns a structured Comparison object with one CatalystComparison per catalyst pipe carrying both sides side-by-side: aeolus1d's closure parameters with their literature defaults / units / physics descriptions, and the project's AUCI scalars from the .ucp dump (value + unit).

  • format_comparison(cmp) — renders both tables as readable text for stdout / logs.
  • write_overrides_template(cmp, path) — emits an editable TOML file with one [catalyst.<id>] section per catalyst, every closure parameter listed as a commented-out line at the aeolus1d default, every AUCI scalar listed below as a candidate. The user uncomments + edits the lines they want to override (typically copying values from the AUCI candidate list by matching unit + position).
  • apply_overrides(case, path) — reads the user-edited template and writes the chosen values onto the matching SpeciesSpec fields. Supports both per-element form (A_brnk_0 = 1.5e+10 overriding one Brinkmeier reaction's pre-exponential) and full-array form (A_brnk = [..] replacing all eleven). Validates pipe ids, parameter names, array bounds, and finiteness; raises ValueError with a concrete fix message on bad input.

CLI subcommand: aeolus1d auci-compare PROJ.proj [--template OVERRIDES.toml] prints the comparison and optionally writes a template file. Use to audit the project values before transcribing, then feed the user-edited template to apply_overrides after import.

Currently supports the twc_brinkmeier reaction kind (44 parameters: 11 A + 11 Eₐ + 11 ΔH array entries plus 11 OSC / inhibition / initial-condition scalars). Adding more closures = adding entries to _CLOSURE_PARAMS in avl_cruise_m_tune.py. Ten new tests in test_avl_cruise_m_tune.py pin: empty-input handling, format rendering both sides with units, template TOML round- trips through tomli cleanly with every line commented, per-element + full-array overrides apply correctly, and validation rejects unknown pipes / bad indices / unknown parameters. Two integration tests exercise the local Cruise M fixture end-to-end (skipped when absent).

  • Phase 8f — Diesel Particulate Filter (Milestone 8). Six- step build-out of the DPF physics chain on top of the existing catalyst infrastructure. 8f.1 declares ParticulateMatterSpeciesTransport columns from the Cruise M importer as inert Y_extras slots so PM mass fractions ride along with the catalyst gas. 8f.2 adds the FilterSpec schema arm + the filter closure-kind registry entry + a per-pipe soot_loading solid-phase buffer + the standalone filter_trap_rates kernel. 8f.3 wires the trap kernel into the species-advection step: every _species_rate call returns a 3-tuple now, with the third element a per-cell d_soot/dt that the timeloop forward-Euler-integrates onto soot_loading at RK stage 1 only (per-stage rate evaluation but single-stage integration so the slow solid-phase ODE state isn't double-counted). 8f.4 overlays a Darcy cake-layer pressure-drop -cake_resistance · μ_ref · u · soot_loading on the momentum residual at every stage. 8f.5 adds the exothermic regeneration kernel: parallel active-O₂ + passive-NO₂ Arrhenius pathways consume soot_loading and release regen_dH × |rate| of heat into the gas-energy residual. 8f.6 auto-attaches a [pipe.filter] block to any Cruise M pipe whose topology element is a particulate- filter class (no-op when the .proj has no DPF, as in the bundled TWC-only fixture). Bundled example 41_dpf_regen.toml runs the full chain — trap + cake + active-O₂ regen — on a single pipe at 1100 K with a pre-loaded 1 g/m³ cake: PM1 trapping peaks the soot at ~1.2 g/m³ at t ≈ 10 ms, then the O₂ Arrhenius regen drains it to ~0.15 g/m³ by t = 80 ms (regen not fully exhausted). The original example used a low initial load (1 g/m³) to keep the gas-direct exotherm within ~100 K; Phase 8g (below) adds the solid-phase thermal buffer that makes higher loads tractable. 14 DPF + 3 importer-side unit/integration tests cover the chain.

  • Phase 8g — DPF cake solid-phase thermal mass. Routes regen heat release through a per-cell T_cake solid-temperature state (mirroring the T_solid catalyst path) rather than depositing it directly into the gas energy residual. A ThermalSpec sub-spec added to FilterSpec carries three parameters: Omega_cake (ρ_s·c_p_s, J/(m³·K)), h_eff (volumetric Newton cooling coefficient, W/(m³·K)), and T_cake_init (K); omitting the [pipe.filter.thermal] TOML block leaves the Phase 8f.5 direct-to-gas path unchanged. The FilterState dataclass gains T_cake, Omega_cake, h_eff fields. The dpf_regen_thermal_apply_cells Numba kernel splits heat into cake and adds the Newton coupling flux to R_E; regen_into now dispatches to this kernel when T_cake.size > 0 and always returns a (d_soot_dt, dT_cake_dt) tuple (dT_cake_dt = None on the direct path). network_timeloop.py snapshots T_cake and advances it with the same SSP-RK2/RK3 combiners used for T_solid. Bundled example 41_dpf_regen.toml updated: soot_loading_init = 8 g/m³, t_end = 60 ms; the thermal block uses Omega_cake = 5000 J/(m³·K) (thin reactive cake layer for visualization) so ⟨T_cake⟩ rises ~21 K over the run while the gas is heated via Newton coupling. 8 new unit tests pin the thermal kernel, the Newton cooling term, the tuple return contract, and the end-to-end T_cake rise.

  • Phase 8h — DPF regen mass balance: O₂ depletion, CO₂ production, gas-density increase. Operator-split post-RK step: for each soot-consuming cell the consumed mass (soot_n − soot_final) adjusts rho·Y_O₂ (−32/12 per kg soot), rho·Y_CO₂ (+44/12 per kg soot), and p.U[:,0] (+consumed). Deferring to post-RK prevents the SSP combiner from diluting corrections injected mid-stage. FilterState gains a co2_col field (wired by build.py); example 41_dpf_regen.toml adds a CO₂ extra column. 5 new unit tests cover wiring, end-to-end depletion/production, and the stoichiometric mass balance (≤5 % relative error over 3 ms).

  • Phase 7g — Per-cycle engine performance metrics. New aeolus1d.analysis.engine_metrics module (and aeolus1d.analysis top-level re-export). CycleAccumulator is passed as the on_step= callback to advance_network; it accumulates (V, p) samples each macrostep, detects 720°-cycle boundaries by monitoring last_theta_deg wraparound, and on each cycle close computes: IMEP (Pa) = ∫p dV / V_disp via np.trapz; ISFC (g/kWh) = m_fuel_per_cycle / W_cycle × 3.6 × 10⁹; indicated thermal efficiency η_i = W_cycle / (m_fuel × LHV). The p–V trace (V_array, p_array, theta_array) for each cycle is stored for indicator diagrams. cycle_summary() aggregates mean/std across cycles (skipping the first partial cycle by default). EngineMetrics dataclass carries all per-cycle outputs. Works with single- and multi-injector cylinders (Phase 7a). LHV is read from Cylinder.combustion.LHV automatically or passed explicitly. Demo: examples/plots/plot_47_engine_metrics.py — 80 mm × 90 mm SI engine at 1800 rpm, 60 mg/cycle iso-octane, 5 cycles, prints per-cycle IMEP / ISFC / η_i table and saves 4-panel figure.

  • Phase 7d-2 — EGR cooler ε-NTU heat-exchanger component. New HeatSpec.kind = "egr_cooler" implements a counter-flow effectiveness-NTU heat exchanger for exhaust-gas recirculation cooling. The closure computes bulk gas heat capacity rate C_hot = ṁ·cp from the mean mass flux, derives NTU = UA/C_min and the counter-flow effectiveness ε, then distributes the total heat duty Q = ε·C_min·(T_mean − T_coolant) to cells in proportion to local temperature excess. C_coolant = 0 selects the infinite-capacity limit (T_coolant fixed). The @njit kernel is bypassed; cooling is applied through the existing py_heat seam. Example examples/44_egr_cooler.toml demonstrates a light-duty diesel EGR cooler (UA = 400 W/K, C_coolant = 950 W/K, T_coolant = 360 K) cooling 700 K exhaust gas. 6 unit tests cover py_heat wiring, UA = 0 no-op (vs kind = "none"), gas cooling, finite vs infinite capacity comparison, and TOML round-trip.

  • Phase 7c — Multi-cylinder engine network: MultiCycleAccumulator. New MultiCycleAccumulator class (aeolus1d.analysis) wraps N CycleAccumulator instances and is passed as a single on_step= callback to advance_network for multi-cylinder engines. After the run, per-cylinder cycle lists are available via multi.accumulators[i].cycles and fleet-aggregated statistics (mean IMEP, mean η_i, etc.) via multi.fleet_summary(). New example 50_two_cyl_engine.toml demonstrates a 360° parallel twin (80×90 mm, ω = 1800 rpm) with shared intake and exhaust manifold plena and two anti-phase cylinders (θ = 0°, 360°). 11 new unit tests cover the API, phase-offset wiring, fleet aggregation, and an end-to-end run.

  • Phase 7b-3 — HNCO hydrolysis kinetics. New hnco_hydrolysis reaction kind (src/aeolus1d/closures/hnco_hydrolysis.py) models the kinetically rate-limited step in urea decomposition: HNCO + H₂O → NH₃ + CO₂ (first-order Arrhenius, Ea ≈ 80 kJ/mol, literature defaults: Kleemann et al. 2000). Primary species = HNCO; extras[0] = NH₃ product; extras[1] = CO₂ (optional). The UWS injector gains a new hnco_frac parameter (0–1): when > 0, that fraction of urea goes through the thermolysis intermediate (CO(NH₂)₂ → NH₃ + HNCO) and the HNCO is tracked as primary Y so the hydrolysis closure can consume it kinetically. hnco_frac = 0 (default) preserves the original instant- conversion model exactly. Example examples/48_hnco_hydrolysis.toml demonstrates a mixer+injector feeding a hydrolysis pipe with hnco_frac = 1.0. 13 unit tests cover the rate law, molar balance, temperature dependence, UWS yield accounting, and TOML round-trip.

  • Phase 7b-2 — SCRF + two-stage SCR + standalone DOC. Three aftertreatment configurations added in this phase. (i) SCRF (SCR-on-Filter): a single pipe combines DPF mechanics (filter.kind = "diesel_particulate_filter") + washcoat thermal mass (catalyst.kind = "two_phase_heat") + SCR kinetics (species.reaction = "nh3_scr_kinetics"); examples/43_scrf.toml shows a mixer+UWS injector feeding the combined brick. (ii) V2O5 slip catalyst (nh3_slip_catalyst): new reaction kind backed by make_nh3_slip_catalyst in src/aeolus1d/closures/scr_kinetics.py. Same Langmuir-Eley-Rideal mechanism as the Cu-zeolite closure but with V2O5-WO3/TiO2 defaults: weaker NH₃ binding (Ea_des 60 kJ/mol), faster NH₃ oxidation (Ea_ox 80 kJ/mol), smaller storage (W_cat 50 mol/m³). Placed downstream of the primary Cu-zeolite SCR to oxidise NH₃ slip above 400 °C. examples/45_two_stage_scr.toml demonstrates mixer → SCR1 Cu-zeolite → SCR2 V2O5 slip in series. (iii) Standalone DOC: examples/46_doc.toml shows a Pt/Pd cordierite DOC with reaction = "doc_oxidation" (HC primary + CO extras[0]); CO lights off first (~150 °C), HC at ~250 °C. 9 unit tests cover closure shape, θ-depletion physics, build wiring, TOML round-trip, end-to-end NO reduction (two-stage SCR), and HC/CO consumption (DOC).

  • Phase 7b-1 — Urea-water-solution (UWS) spray injector. New UreaSpec on PipeSpec adds an operator-split post-RK NH₃ source at a configurable axial location (x_frac, n_spray_cells). Net decomposition model: CO(NH₂)₂ + H₂O → 2 NH₃ + CO₂ (0.567 kg NH₃

  • 0.733 kg CO₂ per kg urea). UWS concentration, injection-rate schedule (mdot_t/v), and conversion efficiency (eta_conv) are all configurable; NH₃ and CO₂ column indices are auto-detected from extra-species names. Example examples/42_scr_urea_injector.toml shows a mixer+SCR two-pipe layout. 9 unit tests cover build wiring, column detection, NH₃ and CO₂ mass-balance closure (≤5 % error), zero- injection invariance, and TOML round-trip.

  • Phase 8l — biased_plenum junction: declarative flow + AdBlue maldistribution. New build-time junction kind replaces N per-pipe inlet [[bc]] blocks with a single [[junction]] carrying two Weltens uniformity-index targets: flow_ui (velocity) and adblue_ui (NH3 mass fraction). At case-build time (aeolus1d.io.build._expand_biased_plenum) the junction is expanded into N transient_mass_flow_inlet BCs; no runtime junction object is created. Distribution follows a symmetric linear profile q_i = q̄·(1 + α·xᵢ); α is solved analytically to match the UI target exactly. The correlation field controls coupling direction: "aligned" (slow channels carry less NH3, i.e. poor mixing) or "opposed" (slow channels carry more NH3, over-dosed). Extra species (e.g. NO) are passed uniformly via distributor_extras. The scr_uniformity example is updated to the new junction kind. 23 new unit tests cover distribution math, build-time expansion, validation, schema/registry consistency, TOML round-trip, and an end-to-end 2-pipe flow check.

  • Phase VGT — schema wiring for VGT / wastegate actuators. TurboSpec gains vane_open: float = 1.0 and bypass_open: float = 0.0, forwarded to Turbine at build time. Both fields were already live on the Turbine class (Phase 7c); this phase closes the TOML ↔ object gap so operators can drive VGT position and wastegate fraction from a TOML case without bespoke Python. Parser and builder updated; defaults preserve full backward-compatibility. 2 new schema tests added to test_turbo_vgt_wastegate.py.

  • Phase VGT-ctrl — VGT boost-pressure PI controller. New aeolus1d.control.BoostController drives Turbine.vane_open via a velocity-form P/PI controller to maintain a target static pressure in a designated manifold pipe. Velocity form avoids integral windup and preserves whatever vane_open the network starts at when error is zero. Usable as an on_step= callback to advance_network. Diagnostic last_p_meas and last_error fields updated each step. 5 integration tests cover default state, underboost / overboost response, output clamping, and integral accumulation.

  • Phase 8m — DPF active regeneration heater. FilterState gains heater_power_w (W, default 0 — off) and soot_trigger_kg_m3 (kg/m³, default 8.0) fields. When heater_power_w > 0 and mean soot loading ≥ trigger, the heater deposits power uniformly into T_cake as dT/dt += P / (N · V_cell · Ω_cake). Below threshold, or with heater disabled, the thermal path is unchanged. last_heater_active (bool) on FilterState is set each step for diagnostic / control use. Schema, parser, and builder updated. New example 50_dpf_active_regen.toml. 8 new unit tests cover heater-off no-op, threshold activation, analytical dT rate, schema defaults, and TOML round-trip.

  • Phase 8m-2 (planned) — DPF regen implicit-path integration test. Targeted test mirroring test_implicit_solver_catalyst_ T_solid_updates: build a DPF pipe with initial soot loading, run advance_network_implicit (BDF2), and assert that soot_loading decreases and T_cake rises under the regen Lie split (Step 6). Completes the implicit-path coverage for the DPF cake/regen operator splits added alongside the catalyst and species work.

  • Phase 7d-3 — SI knock: timing-retard feedback. The Livengood-Wu knock integral (Phase 7d) now drives a timing-retard accumulator: each knock event adds knock_retard_deg_per_event (°CA, default 2.0) to last_timing_retard_deg, capped at knock_retard_max_deg (default 10.0). A quiet cycle (no knock between IVC and the next IVC) recovers knock_advance_deg_per_cycle (default 0.5°) back toward zero. The accumulated retard is subtracted from theta_start_deg in both the two-zone update and the Wiebe burn-rate computation. CylinderSpec and build.py wired accordingly. 4 new integration tests pin retard accumulation, cap enforcement, quiet-cycle recovery, and schema defaults.

  • Phase 8n — DPF soot-oxidation heat-release stability clamp. The regen soot burn and its exotherm are integrated with forward Euler (p.U[:, 2] += dt·R), so at high soot load a coarse step can release more heat than the gas / cake can absorb and overshoot the temperature. The opt-in [pipe.filter.regen] hrr_clamp caps the total oxidation rate at the value whose one-step heat release equals the local sensible energy — ρ·cp·T/(dH·dt) on the direct-to-gas path and Ω_cake·T_cake/(dH·dt) on the cake-thermal path — scaling the slow + fast soot fractions together so mass and heat stay consistent. Both Numba kernels take a trailing dt; the wrapper forwards it only when regen_hrr_clamp is set (0 otherwise), and all three drivers (explicit, QS, implicit) pass their step dt. hrr_clamp = false (new RegenSpec / FilterState field, default) is bit-identical. Complements the global [case] dt_hrr_fraction adaptive-Δt controller (example 55_dpf_regen_hrr.toml) — the clamp is a per-cell rate backstop, the controller a step-size governor. 8 new unit tests in test_dpf_hrr_clamp.py cover the energy cap analytically (direct + cake), inert-when-dt=0, inert-when-well-resolved, mass↔energy consistency, and schema / build wiring.

  • Phase 8k — DPF passive NO₂ regen: Langmuir NO back-inhibition. The NO₂-pathway rate law gains a Langmuir competitive-adsorption term: r_NO₂ = A_no2 · C_NO₂ / (1 + K_NO_inh · C_NO) · soot · exp(-Ea/RT). As NO accumulates from the regen reaction it occupies active carbon sites, progressively inhibiting further NO₂ attack. K_NO_inh = 0 (new RegenSpec / FilterState field, default 0.0) recovers the Phase 8j rate exactly. Both Numba kernels (direct-to-gas and thermal) and the post-RK mass-balance block are updated. New example 49_dpf_no_inhibition.toml demonstrates K_NO_inh = 50 m³/mol with an O₂-pathway-disabled DPF at T = 850 K. 9 new unit tests cover the inhibition factor analytically, thermal vs. direct path agreement, schema defaults, TOML wiring, and an end-to-end run.

  • Phase 8j — DPF regen NO₂-pathway split in post-RK mass balance. The Phase 8h post-RK correction now apportions consumed soot between the O₂ path (C + O₂ → CO₂, −32/12 kg O₂ per kg soot) and the NO₂ path (C + 2NO₂ → CO₂ + 2NO, −92/12 kg NO₂, +60/12 kg NO per kg soot). The split fraction is evaluated from the current gas state using the same Arrhenius rates as the regen kernel; f_o2 = r_O₂/(r_O₂ + r_NO₂) (defaults to 1.0 where both rates vanish). FilterState gains a no_col field wired by build.py when a "NO" extra- species column is present; example 41_dpf_regen.toml adds a NO column. 5 new unit tests cover field default, wiring, end-to-end NO₂ depletion/NO production, and the stoichiometric mass balance via the pure NO₂ path (≤5 % relative error over 3 ms).

  • Phase 8i — DPF filter state in snapshot (soot loading + T_cake persistence). pipe_state_dict now emits filter/soot_loading and (when thermal mass is active) filter/T_cake so both the tracer and the HDF5 writer pick them up automatically. The writer creates a filter/ subgroup in each DPF pipe group. resume.py reads filter/soot_loading and filter/T_cake from the prior HDF5 and restores them via seed_pipes_from_prior. PriorSeries gains a pipe_filter dict so preload_from_prior reconstructs filter-state snapshots across resume boundaries. Tracer adds pipe_filter_state(pipe_id, key) and pipe_filter_keys(pipe_id) accessors (shape (n_samples, N)). 9 new integration tests cover snapshot content, tracer time-series fidelity, HDF5 group structure, and round-trip resume.

  • Phase A32 — Cruise M importer emits passthrough for pipe-pipe joints with matched areas. Wires Phase A31's zero-volume junction into the AVL Cruise M import path. Cruise M chains pipes via direct (Outlet → Inlet) connections with no explicit volume; the historic translation always inserted a 1 L plenum at every such joint, injecting unphysical mixing capacitance. The new logic compares each joint's pipe cross-sections (π/4 · d² per pipe) and emits kind = "passthrough" when they agree within 5 %, kind = "plenum" otherwise (with a warning suggesting the user model the contraction explicitly via a variable-area pipe or a Corberán junction with K_loss).

Schema-side wiring: NeutralJunction gains an optional kind: NeutralJunctionKind field (Literal["plenum", "passthrough"], default "plenum") so any neutral- schematic importer can opt in. schematic_to_case propagates the kind into JunctionSpec; schematic_from_case round-trips it. GT-Power / AVL Boost / WAVE keep emitting plenum (their formats don't carry a passthrough analogue), so behavior on those paths is unchanged.

The reference Cruise M aftertreatment fixture (~286 mm catalyst-pipe diameters on both sides of the joint) imports with one passthrough junction instead of one plenum, and the 35 Cruise M integration tests stay green. The junction-transfer test fixture's window grew from 5 ms to 50 ms — passthrough has faster-but-oscillatory startup dynamics than plenum's slow capacitance-driven equilibration, so the steady-state assertions need a longer settling window.

  • Phase A31 — direct pipe-to-pipe passthrough junction. New JunctionSpec.kind = "passthrough" (registry code 3) realises a zero-volume face-to-face coupling between two pipes: each pipe's port-end ghost is filled from the OTHER pipe's port-end interior cell, and HLLC at the boundary then resolves a single conservative flux. No chamber state, no capacitance, no Helmholtz mode — equivalent to one merged pipe of the same total length when areas + γ + R_gas match.

Motivation: AVL Cruise M (and most production 1-D engine tools) connect components port-to-port without an intermediate plenum. The historic aeolus1d translation inserted a plenum junction on every such joint, which injected an unphysical mixing volume that smeared transient response. The passthrough junction is the right model whenever a Cruise M / GT-Suite source connection has no explicit volume between the two endpoints (catalyst → exhaust pipe being the canonical case).

Constraints (validated in __post_init__): exactly 2 ports; one LEFT end + one RIGHT end (head-to-tail arrangement); port areas match within 5 %; no neck / py_junction_bc / lodi (those are alternative port-BC strategies for plenum junctions). Same-end pairings or large area discontinuities raise — use a regular plenum junction or a variable-area pipe.

Implementation lives in aeolus1d/components/junction_passthrough.py (single new module, ~150 lines) with a PassthroughJunction dataclass that reuses the existing Port type. apply_port_bcs does mutual U_pad ghost-cell exchange + writes species_inflow_* for every species column on the partner pipe via the shared species_coupling helpers, so multi-species transport carries through the joint identically to an interior face. accumulate_rates returns (0, 0); a dummy Volume is kept solely for timeloop compatibility (same trick as CorberanJunction). apply_fuel_step is a no-op.

Nine new tests in test_junction_passthrough.py pin: validation (3-port rejection, same-end rejection, area mismatch, neck/lodi/py-bc rejection), ghost-cell exchange exactness, zero accumulate-rates, equivalence to a single merged pipe at the joint face (HLLC bit-exact), species inflow propagation, and end-to-end TOML round-trip via [[junction]] kind = "passthrough".

  • Phase A30 — drop bundled vendor case fixtures. AVL Cruise M .proj archives and GT-Power XML cases carry vendor / engine-program IP and shouldn't ship in the open- source repo. This commit removes tests/fixtures/interop/MFR_transition_SC630.proj and tests/fixtures/interop/gt_power_sample.xml from git tracking, adds defensive .gitignore patterns for vendor case extensions (*.proj / *.gtm / *.gtsub / *.cmp / *.trb / *.glx / *.rlt) plus bracketed paths under tests/fixtures/interop/ and a /production_case/ folder convention, and re-points the GT-Power integration tests at an in-test synthetic XML generated from _build_schematic() (round-trippable + hand-injected <TurboCharger/> element so the lossy- warning detection still has something to fire on). The Cruise M integration tests keep their EXAMPLE path constant pointing at a generic aftertreatment_case.proj slot and skip cleanly when absent — drop your own .proj at that path to enable the run. Documentation, source-code docstrings, and test-file comments are scrubbed of brand-identifying naming.

  • Phase A29 — multi-injector cylinders (Milestone 8 Phase 8c). First feature-side commit of Milestone 8. CylinderSpec gains an injectors: list[InjectorSpec] field; each InjectorSpec carries one injector's (cad, rate, m_fuel_per_cycle, T_fuel, cv_fuel, fuel_species_name) so a cylinder can model split injection (pilot + main of the same fuel) and dual-fuel cases (e.g. H₂ pilot + NH₃ main, or diesel pilot + NG main). The runtime caches each injector's cycle integral _norm at construction; _fuel_injection_rates(θ) walks the list and sums per-injector rates into the cylinder's total (ṁ_fuel, Ė_fuel), with each injector caching last_mdot / last_Edot so the Phase 25n species-coupling write at apply_port_bcs can attribute the right per-injector mass-flow to the right fuel_species_name channel on the connected pipe. Wiebe combustion uses the SUMMED m_fuel_per_cycle across injectors as the heat-release source.

Single-injector back-compat: existing Cylinder(...) callers using the legacy inject_cad / inject_rate / m_fuel_per_cycle keyword args continue to work — __post_init__ synthesises a one-element injectors list from the legacy fields when injectors is empty but m_fuel_per_cycle > 0. The runtime always operates in terms of the multi-injector list after the migration.

Eight new tests in test_cylinder_multi_injector.py pin: cycle-integral caching, single-injector firing at peak, rates summing when both fire simultaneously, dual-fuel per-injector cv·T contributions, Wiebe heat release using summed m_fuel_per_cycle, legacy single-injector auto-migration, and end-to-end TOML round-trip with [[cylinder.injector]] blocks.

  • Phase A26 — Cruise M K + C_F porous-coefficient fallback. The Phase-A22 importer auto-attached the catalyst's pressure-drop closure when the AUCI catalogue carried α (1/m) + β (1/m²); files using the alternative permeability form (K m² + C_F dimensionless) were silently dropped. Phase A26 adds a fallback extractor consulted only when α/β are absent — when both K (with a unit token matching m^2 / permeability~m2) and C_F ([-] / dimensionless~no) are found, the catalyst pipe gets the Phase-A25 darcy_forchheimer_kc closure auto-attached + a heuristic-extraction warning asking the user to verify against their datasheet. Sanity bounds (K ∈ [1e-14, 1e-3] m²) reject false-positive matches against unrelated m²-valued AUCI blocks (orifice areas, port cross-sections). Three new tests pin: the unit-match path, the implausible-K rejection, and the "neither form found" fallback warning that names both α/β and K/C_F so the user sees both manual-configuration options.

  • Phase A25 — Darcy-Forchheimer K + C_F porous-media form. New darcy_forchheimer_kc porous closure kind, an alternative storage convention to the existing α / β form. Takes permeability K (m²) and dimensionless Forchheimer coefficient C_F directly — the same SI units a catalyst datasheet quotes — so a user with measured (K, C_F) can specify them in TOML without algebraic conversion at authoring time. The kernel evaluates |dP/dx| = (μ(T) / K) · u + (ρ · C_F / √K) · |u| · u with no α/β intermediate; physically equivalent to β = 1/K, α = 2 · C_F / √K on the existing kind. Both forms share the same params-array layout (slots 0/1 carry α/β or K/C_F depending on kind); the kernel dispatches on the integer kind code. Sutherland viscosity branch fires identically for both. The parser's legacy- field rejection now scopes itself to the α/β kind so K and C_F are accepted as natural fields when kind = "darcy_forchheimer_kc" (still rejected with the migration recipe under "darcy_forchheimer"). Six new unit tests pin the kernel formula, the equivalence to α/β, the Sutherland branch, K → 0 short-circuit, sign convention for reversed flow, and TOML round-trip.

  • Phase A24 — Cruise M untranslated-element audit trail. Final pass on the AVL Cruise M importer's audit completeness. Every Cruise M object type the importer doesn't translate now surfaces an informational warning in the InteropReport so the user sees the full audit trail of what's in the .proj file vs what landed in the aeolus1d Case. Covered types: ClassicSpeciesTransport (superseded by GeneralSpeciesTransport), PassiveSpeciesTransport / PassiveSpecies (no aeolus1d auto-import; users add via Y_extras + reaction='none' manually), ParticulateMatterSpeciesTransport / ParticulateMatterSpecie / ParticulateMatterSpecies (no aeolus1d PM tracking; queued as Phase 7b DPF work), RamFlow (vehicle aero, no analogue), AuciTransferHeatProperty / AuciTransferMassProperty (Cruise M correlation choice, aeolus1d has Hawthorn baked in), TransferModelMultiplier ConstDatabus (per-reaction rate multipliers), EnthalpyOutputAdaption (T ↔ h reference state). The test fixture's import now surfaces 18 warnings covering every non-translated object type so the user can audit the delta before running.

  • Phase A23 — surface AUCI rate-constant table. Cruise M's per-reaction (A, Eₐ/R) pairs live in cruise.props.auci.AmountWithCustomUnit blocks nested inside each AuciCatalytic. The values are in the XML; the labels (which scalar maps to which reaction) live in a binary .ucp file the importer can't read. Pre-A23 the importer warned "tune in the inspector" and left users to start from closure literature defaults entirely blind. Phase A23 emits a structured warning per catalyst listing every (value, unit) pair so the user can manually transcribe the project-tuned rate constants into the inspector's A_* / Ea_* / dH_* fields. Doesn't auto-populate aeolus1d closure parameters: the order of scalars within each AuciCatalytic depends on the .ucp mechanism's internal slot layout. Manual transcription is the right answer until either AVL ships the .ucp schema or someone reverse-engineers the mapping for a specific catalyst family.

  • Phase A22 — Cruise M boundary pressures + AUCI α/β. Two more gaps closed: BoundaryConditions::pressure schedules — pre-A22 the importer hardcoded p_back = 1.0e5 Pa for every pressure_outlet BC and used 1 atm for pipe init. Now walks every BoundaryConditions object: the last entry's pressure populates the outlet p_back; the first entry seeds every pipe and junction's init pressure (so over-pressured aftertreatment cases — turbo exhaust at ~1.2-3 bar — start at the correct operating point instead of 1 atm). Tabulated pressure schedules collapse to the first time-point with a warning; full transient_pressure_outlet support is queued. AUCI Darcy-Forchheimer α/β: _extract_porous_alpha_beta finds AmountWithCustomUnit blocks with unit_userstring of 1/m (Forchheimer) or 1/m^2 (Darcy) anywhere in the AUCI catalogue. When both are found, every catalyst pipe gets PorousSpec(kind='darcy _forchheimer', alpha=.., beta=.., mu_ref=1.8e-5) auto- attached + a confirmation warning. When absent, the existing "α/β not translated" warning fires per catalyst pipe so users see what's missing.

  • Phase A21 — Cruise M GeneralSpeciesTransport decode (per-species inlet schedule). First Milestone 7 commit. The pre-A21 importer used the static cruise.props.gas_composition.GeneralSpecies::gas_fraction per species as the inlet composition, planted as a single- point schedule on the transient mass-flow inlet BC. That threw away the file's actual cruise.props.boundary .GeneralSpeciesTransport table — a per-species, per-time mole-fraction schedule that's the real boundary composition. A representative aftertreatment .proj carries a 2-point table of post-combustion exhaust gas (~14.3% CO2, ~11.9% H2O, ~72.4% N2, ~0.76% O2 — a real lean-burn diesel exhaust composition) that the pre-A21 importer silently replaced with pure air (~76.85% N2 + 23.14% O2 + zero everything- else) — so the imported case ran the catalyst on the wrong feed entirely.

This commit reads GeneralSpeciesTransport (codomain table + names list + input_type molar/mass + residual_specie + setting wet/dry) and converts the mole-fraction table to mass fractions via Y_i = M_i·X_i / Σ_j M_j·X_j with the residual-species fixup X_resid = 1 − Σ(others) applied first. Each species in the catalogue receives its own Y(t) schedule on the BC's primary Y_t / Y_v (or Y_extras_specs). Wet basis maps directly; dry basis is flagged for the user to re-normalise. Falls back to the static catalogue gas_fraction when the file has no transport block (back-compat).

Two new test_avl_cruise_m_import.py tests pin the decoded exhaust composition + the mass-fraction sum-to-1 invariant. The pre-existing static-air test was rewritten to assert the now-correct exhaust values.

Other Cruise M file fields documented as still being missed by the importer (queued as follow-ups, all flagged in the InteropReport): ClassicSpeciesTransport (lumped fuel + combust-products schedules — superseded by GeneralSpeciesTransport when the latter is present); BoundaryConditions::pressure schedule; AUCI SelectionAmountIfun / ValueSets (per-reaction calibrated rate constants); AUCI Darcy- Forchheimer α/β; PassiveSpeciesTransport / ParticulateMatterSpeciesTransport; EnthalpyOutputAdaption / Settings.

  • Phase A20 — CLI subcommand registry full migration. Phase 10 of the architecture refactor's burn-down. The pre-fork aeolus1d.cli.registry framework + 3 demo subcommands (examples, inspect-gtm, analyze-cfl) shipped earlier; this commit migrates the remaining 10 inline _cmd_* blocks (run, bench, gui, modes, import, export, export-fmu, optimize, audio, parallel) to the same @register_subcommand pattern. Each kept its existing _cmd_* body verbatim; the only changes are the new _add_<name>_args(p) helper that holds the argparse declarations, plus the decorator. The build_parser() function collapses from 540 → 30 lines — a single loop over all_subcommands() builds the subparser tree. cli/main.py drops from 1397 → 1135 lines. Adding a new CLI subcommand is now one new module + one @register_subcommand call, no edits to main.py's parser builder.

Architecture refactor — Phase 8f post-mortem (Phases B1–B7, ongoing)

Five targeted cleanup phases extracted from the Phase-8f post-mortem audit (recorded 2026-05-01). All five landed in the same work session:

  • Phase B1 — layout-heuristic library. The two _layout_* helpers in aeolus1d.io.interop.avl_cruise_m._core that derive a sensible intake/exhaust port order from the raw GtmPipe graph were extracted into a standalone aeolus1d.io.interop.avl_cruise_m.layout module with 11 unit tests.
  • **Phase B2 — species column cache + filter needs_quasi1d
  • regen-default leak.** Column-index lookups are now cached once at build time (no repeated index() in the hot path). The needs_quasi1d property gained a filter_kind != 0 branch so porous-filter pipes route through the quasi-1D kernel correctly. The regen-kinetics default-leak (missing else branch in _attach_closures) was closed.
  • Phase B3 — avl_cruise_m.py → subpackage. The 1 400- line monolith was split into _core, units, types, layout, __init__ (re-exports), with matching import and test updates.
  • Phase B4 — per-stage closure hooks. cake_into and regen_into moved from inline code in network_timeloop.py to aeolus1d.closures.filter; catalyst_rates_into moved to aeolus1d.closures.catalyst. network_timeloop.py now calls the module-level helpers.
  • Phase B5 — Pipe state composition. The 19 flat closure- state fields on Pipe (14 filter_* / soot_loading fields, 3 catalyst_* / T_solid fields, 2 porous_* fields) replaced by composable sub-state dataclasses (FilterState, CatalystState, PorousState) hung off pipe.filter / pipe.catalyst / pipe.porous; None means the family is inactive. Numba kernels keep receiving unpacked arrays — only the Python-side carrier changed. pipe.filter is None replaces the previous 9-field zero dance; getattr fallbacks and the missing-else risk in _attach_closures are gone. Every closure helper, all implicit-solver kernels, io/snapshot.py, io/resume.py, io/build.py, bench/porous_plug.py, and 13 test files updated.
  • Phase B6 — Nested FilterSpec sub-specs. FilterSpec's 11 flat physics fields (4 trap + 2 cake + 5 regen) replaced by nested dataclasses: TrapSpec (pm_species_name, filtration_efficiency, soot_loading_init, filter_density), CakeSpec (resistance, mu_ref), RegenSpec (A_o2, Ea_o2, A_no2, Ea_no2, dH). TOML now uses [pipe.filter.trap] / [pipe.filter.cake] / [pipe.filter.regen] sub-tables; sub-specs grow independently (Phase 8g's ThermalSpec slotted in without touching the parent). Registry KIND_FIELDS updated to ["trap", "cake", "regen"]; schema_help keys migrated from FilterSpec.* to TrapSpec.* / CakeSpec.* / RegenSpec.*. avl_cruise_m/_core.py importer and all test fixtures updated.
  • Phase B7 — Registry-driven kind codec for friction / heat / porous. The [pipe.friction] / [pipe.heat] / [pipe.porous] blocks used to be threaded through three hand-maintained surfaces — a parser _parse_<family> in aeolus1d.io.parsers that hard-listed fields and defaults, a writer _<family>_lines in aeolus1d_gui.io.case_io that hard-listed the same fields a second time, and a schema KIND_FIELDS ClassVar mirroring the registry's fields_for(family, kind). New aeolus1d.io._kind_codec module exposes parse_kind_block + emit_kind_block, both registry-driven; the three parser / writer pairs collapse to one-line dispatchers. Defaults now live exclusively on the dataclass (no more T_coolant=300.0, mu_ref=1.8e-5 duplicated in three files). Side-effect: the porous emitter now round-trips the Phase-A25 darcy_forchheimer_kc kind the previous emitter dropped on the floor, and emits Sutherland fields under viscosity_model = "constant" too — a few extra default-valued lines in saved TOMLs but bit-identical on reload. Catalyst (B7.4) and filter (B7.5) sub-phases deferred: catalyst's Phase-45a monolith resolver derives rho_s_cs/a_v/h from cpsi post-parse and filter's nested TrapSpec/CakeSpec/RegenSpec from Phase B6 don't fit the flat-field codec without recursion.
  • Phase B7.4 — Catalyst migrated via codec, monolith path declarative. _parse_catalyst keeps three branches — mechanism (carries a surface_init dict the codec _coerce doesn't descend into), cpsi > 0 (Phase-45a monolith path that crosses kind boundaries), and the plain explicit form. The plain form now goes through parse_kind_block so its 4-line float(d.get(...)) block disappears; the cpsi-driven path uses parse_flat_spec to read the 12+ Phase-A39 / A42 / A44 metadata fields without re-listing each default, then dataclasses.replace overrides rho_s_cs / a_v / h with resolve_monolith_params outputs. The back-compat ofa = 0.76 parser default (the dataclass default is 0.0) is the only schema/parser drift this phase doesn't reconcile — fixing it would change Cruise-M-derived defaults, deferred. _catalyst_lines migrates the two_phase_heat branch to emit_kind_block; two_phase_heat_var_nu keeps its skip-default conditional emit pattern, and mechanism keeps the _toml_inline_dict_str_float path for surface_init.
  • Phase B7.5 — Filter migrated via flat-spec codec extension. aeolus1d.io._kind_codec gains parse_flat_spec and emit_flat_spec helpers for nested sub-specs that don't carry a kind discriminator (TrapSpec / CakeSpec / RegenSpec / ThermalSpec). _parse_filter collapses the four parallel float(d.get("Omega_cake", 1.5e6)) defaults into parse_flat_spec(<sub_spec>, d.get(...)) calls — every sub-spec default now lives only on its dataclass. _filter_lines uses emit_flat_spec(skip_defaults=True) for cake / regen / thermal so fresh TOMLs stay minimal; trap's filtration_efficiency headline parameter stays unconditionally emitted for round-trip stability with pre- B7.5 saved cases. Catalyst (B7.4) remains deferred — its monolith resolver and surface_init dict need a different shape than flat-spec recursion.
  • Phase B8.1 — IE + BDF2 share the Newton outer driver. The implicit Euler and BDF2 single-pipe drivers in aeolus1d/solver/implicit_euler.py repeated the entire outer loop — BC setup, Δt selection (CFL / fixed / adaptive), Jacobian chain-cache seeding, the inner Newton + halving retry, history bookkeeping — differing only in the residual's (U_hist, dt_eff) term. New private _advance_single_pipe_implicit(stencil, …) owns the shell; the IE-vs-BDF2 difference collapses to a 4-line StencilStrategy.fn callable that returns (U_hist, dt_eff). Both public entry points (advance_single_pipe_implicit_euler / advance_single_pipe_bdf2) become thin wrappers binding the right stencil. Pure-Python orchestration; the @njit _newton_step / Jacobian kernels stay untouched. Sets up Phase B8.2 — Radau IIA drops onto the same shared kernel by registering its own 2-stage stencil.
  • Phase B8.2 — Radau IIA joins the shared driver. advance_single_pipe_radau_iia had a near-identical 70-line outer loop to the IE/BDF2 drivers, differing only in the Newton step (Radau's _newton_step_radau solves the 3-stage coupled system, while IE/BDF2 share a single-state _newton_step). StencilStrategy gains a newton_fn attribute so the per-step Newton kernel is also pluggable. Radau's stencil reuses _ie_stencil_fn (its residual uses U_n directly with no history fold-in) and binds _newton_step_radau for the kernel. The Radau outer driver collapses to a thin wrapper around _advance_single_pipe_implicit; implicit_radau.py drops 54 lines, the IE/BDF2/Radau outer loop now lives in exactly one place. compose and pipe_min_dt imports are now unused in implicit_radau.py and removed. Pure-Python orchestration; @njit kernels untouched.
  • Phase B8.3 — Residual-builder helpers absorb the per-branch duplication. _build_pipe_residual_fn in aeolus1d/solver/network_timeloop.py walks an irregular scheme × EOS × geometry × γ-mode tree; before B8.3 each quasi-1D leaf re-listed the same 4-line (fk, hk, pk, pp) closure-arg bundle and the same 2-line py-overlay check. Three new module-level helpers absorb the pattern: _slice_pad_to_ng1(pipe, pad) (Ng>=2 → Ng=1 first-order layout slice, formerly a nested closure), _closure_args (pipe) (returns (fk, hk, pk, pp) honouring py-side overrides), and _wrap_with_py_overlay(pipe, fn) (decorator that hoists the has_py_closure branch out of the per-step hot path). Each branch now reads at a glance — kernel + arg-tuple, no scaffolding. network_timeloop.py drops 42 lines net (1716 → 1674); the @njit kernels themselves are unchanged. The dispatch shape stays a nested-if cascade: the cube has irregular gaps (no real-gas-+-per-cell-γ kernel, no MUSCL-+-per-cell-γ-+-quasi1d kernel) that a flat dispatch table would have to encode anyway.
  • Phase B11 — Interop registry + single version constant. aeolus1d.io.interop had a hardcoded if-elif dispatcher in _tool() and three parallel name-tuple constants (_TOOLS_SCHEMATIC, _TOOLS_DIRECT, _TOOLS_EXPORT). Adding a new importer / exporter required editing all four sites. Phase B11 promotes the dispatch table to a single _INTEROP_TOOLS dict keyed by kind, where each entry binds a module path, an import direction ("schematic" / "direct"), and an export-support flag. _tool() becomes a one-line lazy importlib.import_module call; the three derived tuples are list comprehensions over the dict. Adding a new tool now needs one row in the table plus a sibling module — no surgery anywhere else. Version: aeolus1d .__version__ reads from importlib.metadata.version ("aeolus1d") with a development fallback for editable installs without metadata; aeolus1d_gui.__version__ simply re-exports it. pyproject.toml is the only place a release bump touches.
  • Phase B9 — Unified network ghost-fill helper. The explicit timeloop and the implicit network driver carried near-identical ghost-fill helpers — fill_ghosts (a closure inside advance_network / advance_network_strang, applied every RK stage) and _fill_component_ghosts (called once per implicit step before Newton starts). Both walked the same list of components and called apply_port_bcs; the explicit helper additionally applied user-far-BCs (the implicit driver re-fires its far-BCs inside each Newton iteration via a per-pipe closure for the FD chain rule, so it owns those separately). New aeolus1d.bc.base.fill_network_ghosts (junctions, far_bcs=()) is the single source: pass the junction list and (optionally) the far-BC list and it does the right thing. Both call sites delegate; the redundant _apply_far_bcs helper is gone. The BCAdapter dataclass / chain-rule unification described in the original B9 plan stays deferred — the BC-kind strings in implicit_kernels._ghost_state_and_chain are real mathematical chain-rule dispatchers (transmissive = identity, closed = mirror, dirichlet = zero, numerical = FD-approximated) and unifying them across explicit / implicit drivers is API surgery that needs a dedicated session.
  • Phase B10 — Two-pass case validator (initial framework). Before B10, every spec-level mistake — a typo in a junction port's pipe reference, a pipe with neither diameter nor area set, a BC pointing at a non-existent pipe — surfaced as a ValueError raised deep inside the builder pipeline, one at a time. New aeolus1d.io.validate runs the parsed Case through a list of Validator(case) -> list[Issue] callables BEFORE construction starts and aggregates every issue into a single CaseValidationError. Five validators ship in this commit (pipe area exclusivity, unique pipe ids, junction port refs, BC pipe refs, turbo inlet/outlet refs); more migrate from build.py over time. The inline raise ValueError sites stay as defense-in-depth; the validator pass is purely additive — it catches problems earlier with clearer messages. Wired into load_case so every TOML loaded through the public API gets the check automatically. TOML source-position plumbing ((line, col) per issue) is deferred — stdlib tomllib doesn't expose positions, and a third-party parser swap is invasive enough to deserve its own phase. 10 new unit tests pin the contracts; all 46 bundled example TOMLs still pass the validator.
  • Phase B10.2 — Five more validators. Builds on the B10 framework: (1) _validate_turbo_shaft_refs — turbo's optional shaft reference must point at a defined shaft when set; (2) _validate_cylinder_port_refs — cylinder valve ports' pipe ids must exist (multi-cylinder TOMLs copy-paste between cylinder blocks and miss-rename); (3) _validate_unique_bc_ends — two [[bc]] blocks targeting the same (pipe, end) pair silently shadow; (4) _validate_junction_port_count — kind-scoped to biased_plenum (it's a flow distributor, needs ≥ 2 outlet ports); plenum / EHR allow single-port stagnation reservoirs; (5) _validate_unique_junction_port_ends — two ports in one junction targeting the same (pipe, end) would clobber the same boundary ghost. 6 more unit tests pin the contracts; all 46 bundled examples still pass.
  • Phase B7.4-fix — Schema/parser ofa default reconciled. Phase B7.4 noted a drift: CatalystSpec.ofa defaulted to 0.0 in the dataclass while the parser patched a sane 0.76 cordierite value post-hoc when cpsi > 0. Programmatic CatalystSpec(kind="two_phase_heat_var_nu", cpsi=400) constructions therefore crashed inside resolve_monolith_params ("ofa must be in (0,1)"). Phase B7.4-fix moves the 0.76 onto the dataclass as the single source of truth and drops the parser's post-codec dataclasses.replace patch — the parser now reads ofa through parse_flat_spec like any other field. Tests: 174 catalyst / codec / round-trip / case-loader / closure-registry / case-validation tests + 27 GUI form-factory tests green; Cruise-M imports unaffected (they emit ofa explicitly).
  • Phase B9.2 — BCAdapter dataclass. The implicit drivers historically take apply_ghosts (the no-arg ghost-fill closure) and bc_kind (the chain-rule discriminator) as two parallel kwargs that must be paired correctly — a Dirichlet bc_kind requires an apply_ghosts that pre-fills the ghosts before Newton starts; a transmissive bc_kind expects no-op fill semantics. New aeolus1d.bc.base.BCAdapter bundles the pair into one frozen dataclass with a BCAdapter.from_legacy(apply_ghosts, bc_kind) ctor and .apply_ghosts / .bc_kind legacy-name properties. The driver public APIs stay on the two-kwarg shape for back-compat; BCAdapter is the recommended way to construct, store, and pass these values around in new call sites. The plan's bigger goal — replacing the explicit / implicit / network triple BC abstraction with one — needs more invasive surgery (the chain-rule strings drive real Newton dispatch in implicit_kernels._ghost_state_and_chain and aren't just labels) and stays open as a follow-up.
  • Phase B10.3 — Three more validators, one promoted to upfront check. _validate_unique_shaft_ids flags duplicate shaft ids that would silently shadow each other through every Turbo's shaft lookup. _validate_bc_vs_junction_conflict flags a [[bc]] block sharing a (pipe, end) with a junction port — both writing the same ghost cell each step is silently inconsistent physics. _validate_pipe_end_coverage is the marquee migration: every pipe end must have something writing its ghost cells (BC, junction port, turbo inlet/outlet, cylinder valve port); previously a late-fail check inside build_network that surfaced after several builders had already half-wired against the broken case. The CaseValidationError exception is now a subclass of ValueError so legacy pytest.raises(ValueError, ...) contracts in test_case_loader keep working — they catch the error one phase earlier (at load_case) but with the same exception type. 7 new unit tests; all 46 example TOMLs still pass.
  • Phase B10.4 — Five more validators. Continues filling out the validator pass: _validate_reaction_pipe_refs (typos in [[reaction]] blocks' pipe field used to surface late as KeyError in _attach_reactions); _validate_reaction _species_init_length (the species and init lists in one reaction must have equal length — silent indexing mismatch otherwise); _validate_unique_reaction_pipes (two reactions on one pipe shadow each other through pipe.py_reaction_multi); _validate_unique_cylinder_ids and _validate_unique_turbo_ids (mirror the unique-pipe-id / unique-shaft-id pattern for the remaining named-component spec lists). 5 new unit tests; all 46 bundled example TOMLs still load through load_case cleanly. Validator pass tally now: 18 validators wired into aeolus1d.io.validate, 25 unit-test contracts pinning them.
  • Phase B10.5 — Controller validators. _validate_controller_refs flags unknown turbine_id / manifold_pipe_id references in boost_pressure controllers (Phase VGT-ctrl.2 — closed-loop boost-pressure PI), including the common "user pointed at the compressor on the same shaft by mistake" case (which now surfaces as "not a defined turbine" rather than as a generic KeyError mid-build). _validate_unique_controller_ids mirrors the unique-id pattern for the controller list. 4 new unit tests; all 46 bundled examples still load through load_case cleanly. Validator pass tally: 20 validators, 29 test contracts.
  • Phase B10.6 — biased_plenum input validator. _validate_biased_plenum_inputs catches two Phase-8l-distributor mistakes upfront: mdot=0 (the schema default — a missing-from-TOML user mistake produces a zero-flow distributor that runs to a stalled steady state) and Y_primary × ψ_max > 1 (over-aggressive adblue maldistribution under high Y_primary that produces impossible mass fractions). The ψ_max derivation reproduces the build-time _bp_species_bias algebra without invoking the kernel side, so the validator runs at parse time. 3 new unit tests pin the contract; all 46 bundled examples still pass. Validator pass tally: 21 validators, 32 test contracts.
  • Phase B10.7 — Init-block T-or-rho validator. _validate_init_blocks checks that every UniformInit block (on pipes and on stateful junction kinds — plenum / ehr / biased_plenum) sets either T or rho so the EOS can recover the third primitive from p. Stateless junction kinds (corberan / passthrough) ignore the init block entirely and are skipped — the original kind-aware scoping came from running the new validator against the bundled examples and discovering 7 of them legitimately leave init empty on stateless junctions. The build helpers used to raise this late in two separate places (_state_tuple for pipes, _junction_init_T for junctions); the validator catches both cases upfront. 3 new unit tests; all 46 bundled examples still pass. Validator pass tally: 22 validators, 35 test contracts.
  • Phase B10.8 — TOML source-position lookup. Validator issue paths (pipe[3], junction[0].ports[1]) are precise but make the user count brackets to find the offending block. New aeolus1d.io.source_pos.TomlSourceMap bridges the gap with a post-hoc text scan: it walks the raw TOML once, records the line number of every [[<key>]] array-of-tables header, and exposes a lookup(path) that resolves a validator path to the matching source line. Sub-paths (junction[1].ports[2]) work by scoping the sub-block search to the line range of the parent block. Issue gains optional source_line / source_file fields populated automatically by load_case; CaseValidationError.__str__ renders them inline so an error now reads [error] pipe[3] (cases/foo.toml:42): pipe 'tube' must specify ... instead of just [error] pipe[3]: .... The pipe-end-coverage validator's path was upgraded from pipe[id='tube'] to pipe[idx] so it matches the source-map's index-based lookup. 8 new test_source_pos.py tests pin the lookup contracts; all 46 bundled example TOMLs still load cleanly. Validator pass tally: 22 validators, 43 test contracts.
  • Phase B12 — _inv3x3 always-inlined. The 3×3 Cramer's- rule inverse is called once per cell inside the block-Thomas forward sweep (and again in the final back-substitute), and it allocated a fresh (3, 3) array each time. Adding inline="always" to the @njit decorator lets LLVM fold the allocation into the caller's stack frame and elides the function-call boundary entirely. Bench effect (median of three runs on this WSL2 box): _block_thomas_solve_3x3_jit ~58 μs / ~9 Mcell/s vs the prior ~97 μs / ~5 Mcell/s — about 1.7× on the linear-solve component of every implicit Newton step. 47 implicit / BDF2 / Radau / quasi1d / real-gas / per-cell-γ / MUSCL unit tests stay green; algebra is unchanged, only the call-site lowering moves.
  • Phase B13 — Analytic HLLC Jacobian for ideal-gas uniform-γ. The implicit driver's hot kernel — the per-face FD perturbation inside _assemble_jacobian_fo_ideal_uniform_jit — used to fire 12 hllc_flux calls per face for the central-difference ∂F/∂(UL, UR). New aeolus1d.solver.hllc_jac ships a closed-form HLLC Jacobian for the ideal-gas uniform-γ case the most-called assembly path takes. Branches on the wave pattern (supersonic-left / supersonic-right / subsonic-left-star / subsonic-right-star); inside the star-state branches it computes analytic derivatives of SL, SR, S*, and the Toro Ch. 10 star state, then assembles the 3×3 blocks via the Leibniz chain rule. tests/unit/test_hllc_jac.py pins agreement with central-FD to ~1e-5 relative across 8 named wave patterns (Sod, reverse-Sod, near-rest, supersonic-L, supersonic-R, hot-cold, subsonic-L, subsonic-R) plus 50 random subsonic states.

Bench effect (median of three warm runs on the dev WSL2 box): _assemble_jacobian_fo_ideal_uniform_jit ~766 μs vs ~1325 μs pre-B13 — about 1.7× on assembly and ~2.4× on the full Newton step (~618 μs vs ~1490 μs). Combined with B12's inlined _inv3x3, the implicit Newton step is ~2.8× faster than the pre-B12 baseline (618 μs vs 1755 μs).

Variable-γ, real-gas EOS, and frozen-MUSCL variants stay on the FD path for now — their HLLC algebra is identical but the primitive-recovery layer differs; lifting the analytic approach to those is incremental follow-up. The test_uniform_gamma_pad_matches_scalar_jacobian consistency test was loosened from 1e-7 to 1e-5 rtol because the two paths now have intentionally different precision floors (analytic = machine precision; FD ≈ O(1e-6)); 1e-5 still catches any algebra bug on the O(1)..O(100) Jacobian blocks.

  • Phase B14 — chem_solver = "auto" default. The case-header chem_solver field now defaults to "auto" instead of "lsoda", and the per-cell chemistry sub-stepper routes to the @njit Bader-Deuflhard kernel ("bs_jit") when the pipe matches its single-species Arrhenius contract, falling back to LSODA otherwise. The routing predicate (_bs_jit_applicable) and the solver == "auto" branch were already in place — this phase flips the default so the order-of-magnitude speed-up ships unopted-in. Behaviour is unchanged for any pipe that doesn't match the JIT contract (multi-species, surface mechanisms, custom py_reaction callable). examples/15_strang_chemistry.toml pins chem_solver = "lsoda" explicitly so it stays the LSODA reference demo; examples/16_bs_chemistry.toml keeps its explicit "bs_jit" pin and notes that "auto" would already pick that kernel. The GUI's TOML emitter (case_io.save_case) only writes chem_solver to the file when it differs from the new "auto" default.

  • Phase B15 — Analytic HLLC Jacobian for variable-γ. Lifts the B13 closed-form approach off the FD path for the var-γ ideal-gas branch — the case multi-species engine runs hit whenever the species composition gives a non-uniform γ across cells. New _hllc_flux_jac_analytic_var_gamma_impl in aeolus1d.solver.hllc_jac mirrors the B13 wave-pattern branching (supersonic-left / supersonic-right / subsonic-left-star / subsonic-right-star); the wave-speed and star-state algebra is γ-agnostic so only primitive recovery, sound speed, and the Euler flux Jacobian on each side carry γ. Reduces bit-exactly to the uniform-γ analytic when γ_L == γ_R. _assemble_jacobian_general_jit picks up the analytic on the no-EOS no-MUSCL branch for both var_gamma_flag = 0 (uniform-γ ideal-gas, variable-area case the B13 hot path doesn't catch) and var_gamma_flag = 1.

Tests in tests/unit/test_hllc_jac.py pin agreement with central-FD at ~1e-5 across 8 named γ-jump wave patterns (subsonic both signs, Sod-with-γ-jump both directions, slow flow with γ jump, hot-cold with γ jump, supersonic-L/R) plus 50 random subsonic states with γ ∈ [1.20, 1.45] on each side. A reduction test pins variable-γ with γ_L == γ_R == 1.4 matching the uniform-γ analytic to machine precision (rtol 1e-12) on 20 random subsonic states.

Bench effect (median over 200 calls, N=200, dev WSL2 box): _assemble_jacobian_general_jit for the var-γ ideal-gas uniform-area first-order configuration drops from ~292 μs (FD) to ~177 μs (analytic) — about 1.65× on assembly.

Real-gas EOS (Peng-Robinson / Redlich-Kwong / tabulated) stays on the FD path for now: its non-ideal pressure recovery makes ∂p/∂U an EOS-specific derivative, so lifting analytic to PR / RK is its own incremental phase. Frozen-MUSCL var-γ also stays on FD — the slope reconstruction layer wraps every primitive, so the analytic chain rule has to thread through MUSCL too. Both follow-ups are mechanical.

  • Phase B15.1 — Analytic HLLC Jacobian for frozen-MUSCL. Lifts the B13/B15 closed-form HLLC Jacobian onto the frozen-slope MUSCL face flux that the implicit drivers use inside Newton iterations (rk = "implicit_euler", "bdf2", "radau_iia"). New _hllc_flux_jac_analytic_frozen_muscl_impl and _hllc_flux_jac_analytic_frozen_muscl_var_gamma_impl in aeolus1d.solver.hllc_jac build edge primitives from cell-centre states + the frozen primitive slopes σ, convert to edge conservative U_edge, call the existing B13 / B15 analytic HLLC Jacobian at the edge states, and chain through ∂U_edge/∂U_cell (a 3×3 algebraic transform via primitive recovery + edge-conservative composition). The slopes are constants during the Newton step, so their derivative drops out and the sparsity pattern matches first-order HLLC. The shared chain helper _muscl_edge_chain_into is used on both the uniform-γ and var-γ paths.

Tests pin agreement with central-FD at ~1e-5 across 30 random-state batteries each for uniform-γ and var-γ (γ ∈ [1.20, 1.45] on each side, σ ~ N(0, 0.05) per component). Two reduction tests pin (a) frozen-MUSCL with σ ≡ 0 matching first-order analytic to machine precision (∂U_edge/∂U_cell collapses to identity), (b) var-γ with γ_L = γ_R matching uniform-γ to machine precision.

Bench effect (median over 200 calls, N=200, dev WSL2 box): _assemble_jacobian_general_jit for the frozen-MUSCL ideal-gas uniform-area configuration drops from ~415 μs (FD) to ~241 μs (analytic) — about 1.72× for uniform γ and ~1.63× for variable γ.

  • Phase B15.2 — Analytic HLLC Jacobian for real-gas EOS. Lifts most of the FD perturbation cost off the EOS path (Peng-Robinson / Redlich-Kwong / tabulated). New _hllc_flux_jac_analytic_eos_impl recovers (ρ, u, p, a) on each side via primitives_from_U_eos_jit and gets the primitive gradients ∂(u, p, a)/∂U by central-FD on the EOS dispatch — six EOS evaluations per side, far cheaper than the twelve full HLLC + EOS evaluations the FD-on-flux path issues. The HLLC wave-pattern algebra around the gradients stays analytic, mirroring the B15 var-γ structure.

The supersonic branches (where the flux equals one side's Euler flux) need ∂F_euler/∂U for non-ideal pressure; the new EOS-agnostic helper _euler_flux_jac_with_pgrad_into takes (u, p, ∂p/∂U) directly and reduces bit-exactly to the existing closed-form _euler_flux_jac_into when the ideal-gas relation holds (g_p = (γ-1)·(0.5·u², -u, 1)).

Pure-analytic PR / RK would need EOS-specific second derivatives of p(T, ρ) for the sound-speed chain rule (since a² = ∂p/∂ρ|_T + T·(∂p/∂T)²/(ρ²·c_v)). That's a larger increment — the FD-on-EOS-primitives approach captures most of the speedup with O(eps²) ≈ 1e-12 relative gradient error, which the implicit Newton step tolerates without losing convergence behaviour.

Tests pin EOS-IDEAL falling back to the uniform-γ analytic to ~1e-7 rel / 1e-5 abs (the FD on closed-form ideal pressure has ~1e-12 relative error, but with Jacobian entries reaching O(1e6) at moderate flow speeds, the absolute scale is O(1e-6)). Six PR / RK named states cover sub-supercritical methane and CO2 + a low-density methane case in the near-ideal regime; 20 random PR-methane subsonic states (ρ ∈ [5, 50] kg/m³, T ∈ [220, 450] K, u ∈ [-100, 100] m/s) all pass at rtol 1e-4 / atol 1e-5 vs FD-on-flux.

Bench effect (median over 200 calls, N=200, dev WSL2 box): _assemble_jacobian_general_jit for the first-order EOS-PR uniform-area configuration drops from ~595 μs (FD) to ~386 μs (analytic + FD-on-EOS-primitives) — about 1.54×.

  • Phase B15.3 — Analytic HLLC Jacobian for frozen-MUSCL real-gas EOS. Closes the last FD-on-flux path in the implicit Jacobian assembly by composing B15.1's edge-chain helper with B15.2's EOS primitive-FD. New _hllc_flux_jac_analytic_frozen_muscl_eos_impl in aeolus1d.solver.hllc_jac recovers cell-centre primitives via primitives_from_U_eos_jit (with FD-derived gradients), builds edge primitives + edge conservative state via the same ideal-gas- style E_edge = p_edge/(γ-1) + ½·ρ·u² fallback the explicit EOS+MUSCL kernel uses (see hllc_eos_njit.py), runs the B15.2 EOS-aware analytic HLLC Jacobian at the edge state, and chains through ∂U_edge/∂U_cell via a new EOS-aware variant _muscl_edge_chain_eos_into that takes the FD-derived (g_u, g_p) instead of the closed-form ideal-gas formulas.

Subtle point: at σ ≡ 0 the frozen-MUSCL EOS analytic does NOT match the first-order EOS analytic for non-ideal EOS — the edge conservative state still goes through the ideal-gas-energy fallback (E_edge = p_eos/(γ-1) + ½·ρ·u²), so U_edge ≠ U_cell even at σ = 0. The two paths legitimately linearise different residuals; the test battery uses central-FD on the matching frozen-MUSCL EOS flux as the oracle in both regimes.

Tests pin agreement with central-FD on the frozen-MUSCL EOS flux at rtol 1e-4 / atol 1e-5 across two batteries: 15 σ ≡ 0 PR-methane states (sanity-check the σ → 0 limit) and 15 random PR-methane states with non-trivial σ slopes (σ ~ N(0, 0.02) for ρ/u, scaled by 1e5 for the pressure component to match production magnitudes).

Bench effect (median over 200 calls, N=200, dev WSL2 box): _assemble_jacobian_general_jit for the frozen-MUSCL EOS-PR uniform-area configuration drops from ~2102 μs (FD) to ~889 μs (analytic) — about 2.37×, the largest B15-series speedup because the frozen-MUSCL EOS path runs the heaviest FD inner loop (12 full HLLC + EOS evaluations per face).

All FD-on-flux paths in the implicit Jacobian assembly are now retired. The dispatch matrix (uniform-γ × variable-γ ideal-gas × real-gas EOS × first-order × frozen-MUSCL) is fully analytic up to FD-on-EOS-primitives where the EOS pressure law isn't closed-form in U.

  • Phase B15.4 — Pure-analytic EOS primitive gradients. Replaces the FD-on-EOS-primitives loop in _eos_primitive_grads_into (six EOS dispatches per side) with a closed-form chain. The new helpers eos_dpdT_dpdrho_jit, eos_dedT_dedrho_jit, and eos_dadT_dadrho_jit (in aeolus1d.closures.eos_njit) supply the EOS first partials in (T, ρ); the chain computes ∂T/∂U via the implicit function theorem on e_eos(T, ρ) = e_int(U) and ∂p/∂U, ∂a/∂U by chain rule. Newton-solve count drops from 6 calls per side to 1; the ideal-gas branch is fully closed-form (no Newton). ∂a/∂U for non-ideal EOS still uses central-FD on eos_sound_speed_jit over (T, ρ) — closed-form would need second derivatives of the cubic EOS, which are messy enough that FD-on-(T, ρ) is the pragmatic choice (no Newton solve in the perturbations because T is already known). Surfacing this required tightening eos_T_from_rho_e_jit's Newton residual tolerance from 1e-6 to 1e-10 (and step-size from 1e-8 to 1e-12). The previous tolerance left ~2.4e-4 K of T noise on each call; the FD-on-flux test oracle had this noise in every perturbed evaluation. Under B15.2/B15.3, the analytic side shared the same noisy Newton calls so the noise cancelled — the test was effectively a self-consistency check, not an absolute agreement test. The pure-analytic B15.4 chain doesn't share that noise, so it became visible (3 tests failed at ~1e-4 vs rtol 1e-4 before tightening; ~6e-8 after). Newton converges linearly with the approximate cv-only Jacobian, so the tighter tolerance adds ~3-4 iterations per Newton call — but B15.4 cuts Newton calls from 6 to 1 per side, so total Newton work per face still drops ~3.3×. Bench effect: hllc_flux_jac_analytic_eos per-face drops from 3.24 μs (B15.3) to 2.38 μs (B15.4) — a further 1.36× on top of B15.3's 2.37×.

  • Phase B18 — analytic porous-source Jacobian for Darcy-Forchheimer. Closes the analytic-Jacobian story across all three native closure families. The combined Darcy + Forchheimer source decomposes into a linear (Darcy) and a quadratic-in-|u| (Forchheimer) term:

    s_mom = -B · u - C · ρ · |u| · u s_E = u · s_mom

where (B, C) come from the per-kind reparameterisation:

  DARCY_FORCHHEIMER     B = β · μ_ref     C = α / 2
  DARCY_FORCHHEIMER_KC  B = μ_ref / K     C = C_F / √K

Differentiating w.r.t. U = (ρ, ρu, ρE) gives a 3×3 block whose nonzero entries are 4 polynomial expressions in (ρ, u, B, C). New aeolus1d.closures.porous.porous_jacobian writes them into a caller-provided buffer; the analytic_supported(kind, viscosity_model) predicate gates the dispatch on constant-viscosity only — Sutherland-viscosity (μ depends on T which depends on (ρ, p) which depends on U) carries an extra chain-rule layer that today's FD path covers. _cell_closure_source gains an include_porous flag symmetric to B16's include_friction and B17's include_heat. Bench (uniform-state pipe, N=200, porous-only closures): _augment_jacobian_with_closures drops from 3612 µs (FD) to 801 µs (analytic) — a 4.51× speedup, in the same ballpark as B16 (4.46×) / B17 (4.82×). The Solver Internals dock now shows darcy_forchheimer ✓ for porous pipes with constant viscosity.

  • Phase B23 — closed-form ∂a/∂T, ∂a/∂ρ for tabulated EOS. Closes the last FD branch in aeolus1d.closures.eos_njit.eos_dadT_dadrho_jit. The tabulated EOS layout already stores precomputed ∂p/∂T and ∂p/∂ρ panels at every grid point; B23 differentiates the bilinear interpolant of those panels (gradients are piecewise-constant within a cell, computable in closed form) and chains through the standard sound-speed identity

    a² = ∂p/∂ρ|_T + T · (∂p/∂T|_v)² / (ρ²·cv)

to assemble ∂a/∂T and ∂a/∂ρ directly. No new table data needed; reuses the existing 4-panel layout.

New helper:

  • _bilinear_grad(params, table_offset, T, rho) — returns (f, ∂f/∂T, ∂f/∂ρ) from a single bilinear- interpolation pass. Same boundary-clamping rule as the plain _bilinear so out-of-range queries pin to the boundary cell's slope.
  • eos_dadT_dadrho_tab_jit(T, rho, params) — the closed-form helper assembling a² and its (T, ρ) partials.

Note. The bilinearly-interpolated cross-derivatives ∂(∂p/∂T)/∂ρ and ∂(∂p/∂ρ)/∂T are not exactly equal (the two panels were tabulated independently and bilinear interpolation doesn't preserve Schwarz symmetry), but this matches what the prior FD-on-a(T, ρ) reference was implicitly doing — the FD step samples the same interpolant.

Bench. Microbench at (T=450 K, ρ=50 kg/m³, CO2 RK table): closed-form ~257 ns/call vs FD ~1.1 µs/call → 4.28× speedup, comparable to the cubic-EOS B20 result.

Tests (2 new in test_eos_tabulated.py):

  • test_closed_form_dadT_dadrho_tab_matches_high_precision_fd — closed form must match FD to ~10% on quarter-cell steps (bilinear has piecewise-constant slopes within a cell; FD step picks up cell-boundary jumps as truncation noise). Also pins the eos_dadT_dadrho_jit dispatch routes through the closed-form helper for EOS_TABULATED.
  • test_closed_form_dadT_dadrho_tab_handles_negative_a2 — degenerate-state contract pin (a² ≤ 0 returns (0, 0) matching eos_sound_speed_tab_jit).

All 26 test_eos_njit + 47 test_eos_tabulated + 34 test_hllc_jac tests pass.

This closes the last item from the closure-Jacobian programme: every shipped closure × kind × viscosity × EOS combination now has a closed-form analytic Jacobian path.

  • Phase B22 — Wiebe combustion + knock-controller TOML support. Closes the gap noted in the Phase 7d-4 plot README: the WiebeCombustion heat-release model and the knock_* closed-loop retard controller were programmatic-only. With B22, both are first-class TOML citizens.

Schema additions:

  • New WiebeSpec dataclass in aeolus1d.io.schema — theta_start_deg, duration_deg, a, m, m_fuel, LHV, two_zone (Phase 25j), gamma_burnt / gamma_unburnt (Phase 25k two-zone γ).
  • New CylinderSpec.combustion: WiebeSpec | None field (default None keeps the cylinder air-breathing — bit-equivalent to omitting the block, so existing TOMLs round-trip identically).
  • Existing knock_A / knock_n / knock_B / knock_retard_deg_per_event / knock_retard_max_deg / knock_advance_deg_per_cycle fields on CylinderSpec were dataclass-only — the parser ignored them and the saver never wrote them. B22 wires both directions.

Parser. _parse_cylinder reads the [cylinder.combustion] block (when present) into a WiebeSpec and the flat knock_* fields with defaults matching the Cylinder dataclass.

Build path. _build_cylinders translates WiebeSpec into a live WiebeCombustion(...) instance, attaches it to the Cylinder via the combustion= constructor keyword, and applies all six knock-controller fields post-construction. End-to-end: a TOML case now produces a cylinder that can fire combustion and run the closed-loop retard controller without any Python wiring.

Saver. case_io.save_case emits [cylinder.combustion] + the knock fields ONLY when they differ from defaults, so pre-B22 TOMLs keep round-tripping bit-identically. The injection-schedule + knock fields are now emitted as flat [[cylinder]] fields before the [cylinder.init] sub-table opens — TOML doesn't care about leading whitespace, only the most recent [...] header, so putting flat fields after a sub-table silently absorbed them under the wrong scope.

Bundled case. New examples/52_knock_event.toml — boosted single-cylinder SI engine (CR=10, 1.4-bar boost, 50 mg/cycle) that drives the Livengood-Wu integral past threshold and exercises the closed-loop retard controller. Slotted into EXAMPLE_CATEGORIES under "Engines & cylinders" so the GUI runner can find it. The existing plot_52_knock_margin.py keeps its inline-engine setup but the plots README now points to 52_knock_event.toml as the equivalent TOML case.

Schema-help entries (15 new) on the new fields so the GUI Inspector tooltips match the rest of the schema.

Tests (5 new in aeolus1d_gui/tests/test_io_gui5.py):

  • test_wiebe_combustion_round_trip — end-to-end save→load with all 9 WiebeSpec fields (including two_zone = true and per-zone γ).
  • test_no_combustion_block_omitted_in_save — pre-B22 TOMLs (combustion = None) must NOT emit the block.
  • test_knock_fields_round_trip_when_overridden — knock fields must round-trip when any default is overridden.
  • test_knock_fields_omitted_in_save_when_default — pre-B22 TOMLs must NOT emit the knock fields.
  • test_wiebe_spec_builds_real_wiebe_combustion_instance — build path attaches the runtime WiebeCombustion.

All 11 round-trip tests + 25 cylinder integration tests pass.

  • Phase B21 — shared T_app chain across real-gas closures. Closes a redundancy in the closure-Jacobian augmenter: _real_gas_heat_jacobian (Phase B17-real-gas) and _real_gas_porous_sutherland_jacobian (Phase B18-ext-real- gas) each independently Newton-solved T_eos and recomputed the EOS partials (dedT, dedrho, dpdT, dpdrho, p_eos) on the same cell. When both fired on the same pipe, the augmenter paid for the EOS chain twice.

The fix factors the shared EOS-chain prelude into a tiny helper, _real_gas_T_app_chain, that writes a 4-element bundle (T_app, ∂T_app/∂U[0], ∂T_app/∂U[1], ∂T_app/∂U[2]) given U_cell. Two new _from_chain variants of the heat and porous helpers consume that bundle and skip the prelude:

  • _real_gas_heat_jacobian_from_chain
  • _real_gas_porous_sutherland_jacobian_from_chain

The original helpers stay unchanged (so single-closure pipes, which are the common case, are bit-identical to before). The augmenter only invokes the shared-chain path when pipe.has_real_gas AND analytic-heat AND analytic- Sutherland-porous all fire on the same pipe.

Tests (2 new in test_implicit_closures_phase2.py):

  • test_shared_real_gas_chain_matches_separate_helpers — feeding the shared chain into both _from_chain variants must give bit-identical Jacobian blocks to running each independent helper (same algebra, prelude factored out).
  • test_augment_jacobian_real_gas_heat_plus_porous_sutherland — end-to-end on a Peng-Robinson pipe with both heat AND porous-Sutherland active. The shared-chain dispatch routes correctly and the result still matches FD on the full closure source (rtol = 2e-3, same tolerance as other real-gas tests).

Bench. N=200 PR pipe with combined heat + Sutherland porous augmenter call: ~18.5 ms/call. The EOS-prelude reuse saves one Newton solve + 3 EOS-helper calls per cell, on the order of ~25% of the per-cell EOS work — modest in absolute terms but meaningful when scaled across a long real-gas catalyst-with-heat run.

All 33 closure phase-2 tests pass.

  • Phase B20 — closed-form ∂a/∂T, ∂a/∂ρ for cubic EOS. Closes the FD branch in aeolus1d.closures.eos_njit.eos_dadT_dadrho_jit for Peng-Robinson and Redlich-Kwong. Tabulated EOS stays on central-FD (no closed form available without a 2-D second-derivative table).

Math. Define M = ∂p/∂ρ|_T + T·(∂p/∂T|_v)²/(ρ²·cv) = a²; chain a = √M so ∂a/∂T = ∂M/∂T / (2a) and ∂a/∂ρ = ∂M/∂ρ / (2a). The cubic EOS gives all five required partials of p(T, v) in closed form:

∂p/∂v|_T   ∂p/∂T|_v   ∂²p/∂v²   ∂²p/∂T²   ∂²p/∂v∂T

With those, ∂M/∂T and ∂M/∂ρ are pure algebra (see the function docstring). For Peng-Robinson the temperature chain runs through the second derivative α''(T) = (κ/(2T)) · (κ/Tc + s/√(T·Tc)) of the Soave-style cubic-correction factor; for Redlich-Kwong it collapses to -3·a₀ / (4·T^(5/2)·v(v+b)).

Where it's used. eos_dadT_dadrho_jit is on the real-gas HLLC analytic-Jacobian hot path (hllc_jac._eos_primitive_gradients → called per face, per Newton iteration). Closing the FD here removes 4 eos_sound_speed_jit evaluations per call.

Bench. Microbench at (T=400 K, ρ=5 kg/m³, CH4-PR) gives closed-form ~1.03 µs/call vs FD ~4.25 µs/call → 4.11× speedup. The actual inline-from-@njit speedup will be larger since the Python-call overhead dominates the microbench; the closed-form path stays inside the JIT context end-to-end.

Tests (7 new in test_eos_njit.py):

  • 6 parametrised cases spanning RK / PR × CH4 / N2 × gas-like / liquid-like states; closed form must match high-precision FD (rtol = 5e-5, eps = 1e-6).
  • 1 pin test for the unchanged ideal-gas branch (∂a/∂ρ == 0 exactly; ∂a/∂T = γR/(2a)).

All 26 test_eos_njit tests + 34 test_hllc_jac tests (the actual downstream caller) pass.

  • Phase B19 — analytic-vs-FD closure-Jacobian bench. Consolidated bench in tests/benchmarks/test_bench_closures.py that times _augment_jacobian_with_closures analytic vs FD across every closure-Jacobian kind shipped in B16 → B18-ext-real- gas. Marked @pytest.mark.slow so it doesn't run in default CI; useful as a regression detector and as a one-stop shop for telemetry on the closure-Jacobian work.

Bench scenarios (7 new tests, N=200 cells, all pinning that the analytic branch is at least as fast as FD with a 25% noise tolerance):

  • DARCY_CONST friction (B16)
  • COLEBROOK friction (B16-ext)
  • HEAT_CONSTANT_TWALL ideal-gas (B17)
  • HEAT_CONSTANT_TWALL real-gas (B17-real-gas)
  • DARCY_FORCHHEIMER constant-μ (B18)
  • DARCY_FORCHHEIMER Sutherland-μ ideal-gas (B18-ext)
  • DARCY_FORCHHEIMER Sutherland-μ real-gas (B18-ext- real-gas)

How "FD" is forced. Each scenario builds a parallel pipe that stashes a pass-through Python callable in the matching pipe.py_* attribute. That trips the analytic- gate off (the override always wins) and routes the call through _cell_closure_source's FD central-difference branch — same code path the augmenter runs when no analytic kernel is available. Pass-through unwraps the single-element arrays to scalars, calls the @njit kernel, and re-wraps so the FD path's bookkeeping is unchanged.

Reference numbers (current laptop, single thread, N=200):

  • Friction DARCY_CONST: ~7.5× (~0.77 ms → ~5.8 ms)
  • Friction COLEBROOK: ~7.5× (~0.80 ms → ~6.0 ms)
  • Heat CONSTANT_TWALL ideal: ~7.0× (~0.74 ms → ~5.2 ms)
  • Heat CONSTANT_TWALL real-gas: ~3.2× (~1.68 ms → ~5.4 ms)
  • Porous const-μ: ~7.5× (~0.81 ms → ~6.1 ms)
  • Porous Sutherland ideal: ~7.0× (~0.91 ms → ~6.4 ms)
  • Porous Sutherland real-gas: ~2.7× (~2.28 ms → ~6.1 ms)

The real-gas branches are slower than the ideal-gas ones because of the EOS T-solve and the extra eos_dpdT_dpdrho_jit calls; those still beat the FD baseline, just with a smaller margin.

  • Phase B18-ext-real-gas — porous Jacobian under Sutherland μ(T) on real-gas pipes. Closes the gap noted in B18-ext: Sutherland-viscosity is now analytic on real-gas pipes too, via the EOS-derivative chain (Phase B15.4 helpers).

Subtle point. porous_source under Sutherland reads T = p/(ρR) — the ideal-gas relation applied to the EOS-recovered pressure, NOT the EOS-converged temperature. Same convention as heat_source uses, and the same chain rule that powered Phase B17-real-gas:

∂T_app/∂U  =  (1/(ρR)) · [∂p/∂U − (p/ρ) · ∂ρ/∂U]
∂p/∂U      =  dpdT · ∂T_eos/∂U  +  dpdρ · ∂ρ/∂U
∂T_eos/∂U  via implicit-function theorem on
            ``e_eos(T, ρ) = e_int(U)``    (Phase B15.4).

Then the Sutherland correction layered on top of the const-μ block is identical in shape to B18-ext — only the source of ∂T/∂U swaps from the ideal-gas formula to the real-gas ∂T_app/∂U.

New helper:

  • _real_gas_porous_sutherland_jacobian(U_cell, eos_kind, eos_params, porous_kind, porous_params, out) in aeolus1d.solver.implicit_kernels. Module-private, lives next to _real_gas_heat_jacobian so the EOS-derivative helpers (eos_T_from_rho_e_jit, eos_dedT_dedrho_jit, eos_dpdT_dpdrho_jit, eos_pressure_jit) are imported in one place. Reuses porous_jacobian_analytic_supported_sutherland_ideal for the (kind, viscosity_model) predicate — no new predicate function needed.

Wiring. _augment_jacobian_with_closures now branches the Sutherland-analytic path on pipe.has_real_gas: ideal → porous_jacobian_sutherland_ideal (Phase B18-ext); real- gas → _real_gas_porous_sutherland_jacobian (this phase). Constant-μ stays on the unchanged porous_jacobian.

Tests (2 new in test_implicit_closures_phase2.py):

  • test_real_gas_porous_sutherland_jacobian_matches_fd — Peng-Robinson + DARCY_FORCHHEIMER + Sutherland: the new helper matches FD on _cell_closure_source to rtol = 2e-3 (Newton noise from the EOS T-solve, same tolerance as Phase B17-real-gas).
  • test_augment_jacobian_with_real_gas_porous_sutherland_matches_fd — end-to-end on a porous-only Peng-Robinson pipe.

Dock changes:

  • test_porous_sutherland_real_gas_marked_analytic flips the previous test_porous_sutherland_real_gas_falls_back_to_fd expectation: real-gas + Sutherland now shows ✓ with a B18-ext-real-gas tooltip.
  • test_summary_counts_analytic_vs_fd_terms updated to use egr_cooler heat (which routes through py_heat and stays FD by construction) as the FD-term anchor — porous Sutherland is no longer a reliable "always FD" hook on either EOS branch.

All 31 closure phase-2 tests + 17 dock tests pass.

  • Phase B18-ext — porous Jacobian under Sutherland μ(T) on ideal-gas pipes. Closes the gap noted in B18: the Sutherland-viscosity branch of porous_source now has a closed-form ∂S/∂U on ideal-gas pipes, instead of falling through to FD.

Math. The combined Darcy-Forchheimer source

s_mom = -B(μ) · u  -  C · ρ · |u| · u
s_E   = u · s_mom

(with B = β·μ for DARCY_FORCHHEIMER and B = μ/K for DARCY_FORCHHEIMER_KC; C = α/2 or C_F/√K respectively) gains a μ-via-T chain on top of the Phase B18 const-μ block:

Δ(∂s_mom/∂U[k])  =  -∂B/∂U[k] · u
Δ(∂s_E/∂U[k])    =  -∂B/∂U[k] · u²
∂B/∂U[k]         =  (∂B/∂μ) · (dμ/dT) · ∂T/∂U[k]
dμ/dT            =  μ · (1.5/T  -  1/(T+S))      (Sutherland)
∂T/∂U[k]         =  (γ−1)/R · (∂(E_total − u²/2)/∂U[k]) / ρ

Real-gas + Sutherland still routes to FD — the EOS-derivative chain through ∂p_eos/∂U is left for a future B18-ext- real-gas, by analogy with how B17-real-gas extended B17.

New API:

  • porous_jacobian_sutherland_ideal(rho, u, p, gamma, kind, params, out) — the new njit kernel. Constant-μ keeps using the unchanged porous_jacobian.
  • porous_jacobian_analytic_supported_sutherland_ideal( kind, viscosity_model) — companion predicate; returns True only for the two Darcy-Forchheimer kinds with VISCOSITY_SUTHERLAND. Caller is responsible for the ideal-gas check (not pipe.has_real_gas).

Wiring. _augment_jacobian_with_closures now picks Sutherland-analytic when not pipe.has_real_gas, else falls through to const-μ-analytic-or-FD.

Tests (4 new in test_implicit_closures_phase2.py):

  • test_analytic_porous_sutherland_ideal_darcy_forchheimer — closed-form ∂S/∂U matches FD on the source for DARCY_FORCHHEIMER across u ∈ {5, 20, −10}.
  • test_analytic_porous_sutherland_ideal_kc — same for the DARCY_FORCHHEIMER_KC reparameterisation.
  • test_porous_sutherland_ideal_predicate — pins the new and old predicates' separation of concerns.
  • test_augment_jacobian_with_porous_sutherland_ideal_matches_fd — end-to-end on a porous-only ideal-gas pipe.

Dock changes:

  • test_porous_sutherland_ideal_marked_analytic — flips darcy_forchheimer + Sutherland + ideal from fd to ✓ with the B18-ext tooltip.
  • test_porous_sutherland_real_gas_falls_back_to_fd — pins that the real-gas variant still shows fd until a future B18-ext-real-gas.
  • test_summary_counts_analytic_vs_fd_terms — updated to use a real-gas + Sutherland combo to keep one FD term in the summary (ideal-gas + Sutherland is now analytic).

Speedup (N=200, DARCY_FORCHHEIMER + Sutherland on ideal gas): the closure-augmented Jacobian assembly drops from ~8.2 ms (py_porous-override FD path) to ~0.93 ms (analytic), comparable in scale to Phase B18's 4.51× win. All 29 closure phase-2 tests + 17 dock tests pass.

  • Phase B16-real-gas — friction Jacobian on real-gas pipes. Drops a (formerly implicit) restriction: the analytic friction Jacobian (Phases B16 / B16-ext) is now documented and tested as enabled on real-gas (Peng-Robinson, etc.) pipes too.

Why no code change was required. Friction reads only (ρ, u) — never p, never T. The closed-form ∂S_friction/∂U for DARCY_CONST and COLEBROOK therefore differentiates conservative state directly, without going through the EOS. Compare with the heat path, which did need a separate _real_gas_heat_jacobian branch in Phase B17-real-gas because s_E reads T = p/(ρR) with p taken from the EOS — friction has no such coupling.

_augment_jacobian_with_closures had no has_real_gas guard around use_analytic_friction to begin with, so this phase's contribution is:

  • Doc note on :py:func:aeolus1d.closures.friction.friction_jacobian_analytic_supported explicitly stating EOS-agnostic.
  • Tooltip on the Solver Internals dock's friction column header noting the same.
  • Two new tests in tests/unit/test_implicit_closures_phase2.py: test_analytic_friction_jacobian_works_on_real_gas_pipe (closed form is bit-identical at matching ρ, u between real-gas and ideal-gas pipes) and test_augment_jacobian_with_real_gas_friction_matches_fd (end-to-end Peng-Robinson + DARCY_CONST agrees with FD).
  • Dock test test_friction_marked_analytic_under_real_gas pins that DARCY_CONST and COLEBROOK both keep the ✓ marker on a real-gas pipe.

Pinning this contract guards against a future regression where someone might mistakenly add and not has_real_gas next to use_analytic_heat and copy-paste it onto the friction gate.

  • Phase B17-gui-live — Solver Internals runtime ↔ spec sync. Closes the Python-override gap noted in B17-gui's tooltip. The dock now ships two modes:

  • Spec mode (default) — reads the static TOML spec, fast, no side effects. The previous behaviour.

  • Live mode — opt-in via a "Refresh with built network" button. Calls aeolus1d.io.case.build_network(case) to materialise the live Pipe instances, scans them for py_friction / py_heat / py_porous Python-callable overrides, and forces the corresponding cells to FD even when the static kind would qualify for analytic.

A "mode: spec" / "mode: live (N pipe(s) with py overrides)" status label tells the user which path is active. Live-mode overrides clear automatically on every caseReplaced / caseStructureChanged / caseFieldChanged signal — a structural edit invalidates the previous build, so the dock reverts to spec mode and the user clicks the button again to re-check.

Defensive: build_network failures surface in the mode label (mode: build failed: <error>) rather than crashing the GUI — the dock is read-only diagnostic and must never take down the rest of the app.

Tests pin the contract: button exists + starts in spec mode, override map flips a friction cell from ✓ to fd, summary counts py-override term forced to FD, show_for_case resets stale overrides, build_network failures land in the mode label.

  • Phase B17-gui — Solver Internals dock. New SolverInternalsDock tabbed alongside the CFL Analyzer. Per-pipe table showing, for each closure family (friction / heat / porous), the configured kind plus a green ✓ marker when the implicit driver will use the analytic ∂S/∂U augmentation (Phase B16 / B16-ext / B17) or a grey fd marker when it falls back to FD.

Cells dispatch on the schema kinds + per-pipe eos.kind for the heat path: HEAT_CONSTANT_TWALL shows analytic on ideal-gas pipes, FD on real-gas pipes (a tooltip on the cell explains why). Top-line summary aggregates analytic vs FD term counts across the whole case.

Refreshes on caseReplaced / caseStructureChanged / caseFieldChanged so flipping a friction kind in the Inspector immediately re-routes the cell. Read-only — like the CFL analyzer dock, this is a diagnostic display, not an editor.

Documented limitation: the dock reads the static spec, so it can't see py_friction / py_heat / py_porous Python overrides which force the FD path at solver time. Tooltip on the friction / heat headers calls this out so users tuning a hot loop know to look at the runtime pipe.py_* if perf doesn't match the dock's prediction.

  • Phase B16-ext — analytic friction Jacobian for COLEBROOK. Extends Phase B16's closed-form ∂S_friction/∂U to the second native friction kind. The Colebrook-Haaland branch needs the chain rule routed through f(Re) where Re = |u|·D_h/ν, so the Jacobian is no longer the simple polynomial in (ρ, u) that DARCY_CONST gave. Implementation:

  • New _colebrook_haaland_f_and_df helper returns (f, df/dRe) analytically (Haaland argument's log derivative + the inv_sqrt_f → f cube-rule).

  • New _colebrook_f_and_df_du routes through the same laminar / transition / turbulent branching the source kernel uses, returning (f, df/d|u|) at the matching Re. The transition blend across [2300, 4000] picks up derivative contributions from BOTH branches plus the dα/dRe term — all closed-form.
  • friction_jacobian rewritten in a unified form ∂S/∂U from a generic g(u) = A(|u|)·|u|·u and dg/du (with dA/d|u| either 0 for DARCY_CONST or the Colebrook chain). Bit-equivalent to B16's polynomial form when dA/d|u| = 0.
  • friction_jacobian_analytic_supported now returns True for both DARCY_CONST and COLEBROOK; FRICTION_NONE stays the only kind that falls through to FD.

Tests pin the analytic ↔ FD agreement on all three Colebrook regimes (turbulent, laminar, transition) plus the |u| = 0 short-circuit (where sign(u) would otherwise blow up).

  • Phase B17-real-gas — analytic heat Jacobian on real-gas pipes. Closes the EOS gap in B17. The constant-T_wall heat source uses T = p / (ρ · R_gas) — an ideal-gas relation applied to the EOS-recovered pressure, NOT the EOS-converged T. So the differentiated quantity is the "apparent T":

    T_app = p_eos(T_eos(U), ρ(U)) / (ρ · R_gas)

Chain rule:

  ∂T_app/∂U = (1 / (ρ · R)) · [ ∂p/∂U − (p / ρ) · ∂ρ/∂U ]
  ∂p/∂U     = dpdT · ∂T_eos/∂U + dpdρ · ∂ρ/∂U
  ∂T_eos/∂U via the implicit function theorem on
             ``e_eos(T, ρ) = e_int(U)`` (Phase B15.4)

Then ∂s_E/∂U = -K · ∂T_app/∂U. The chain reuses Phase B15.4's eos_T_from_rho_e_jit, eos_dedT_dedrho_jit, eos_dpdT_dpdrho_jit, eos_pressure_jit helpers — no new EOS algebra.

New _real_gas_heat_jacobian in aeolus1d.solver.implicit_kernels writes the 3×3 (only row 2 is nonzero). The augmentation path drops the previous not pipe.has_real_gas gate and dispatches: has_real_gas → real-gas helper, otherwise the ideal-gas heat_jacobian from B17.

Solver Internals dock now shows constant_twall ✓ for real-gas pipes too, with a tooltip referencing the EOS-derivative chain (Phase B17-real-gas) so users see which branch is in use.

Tests pin the analytic Jacobian against FD on a Peng- Robinson + HEAT_CONSTANT_TWALL pipe at supercritical conditions (rtol = 2e-3 — looser than B17 because the FD oracle picks up Newton-iteration noise from the EOS T-solve inside primitives_from_U_eos_jit).

  • Phase B17 — analytic heat-source Jacobian for HEAT_CONSTANT_TWALL. Mirrors B16 on the heat side. Closed-form ∂S_heat/∂U for the constant-T_wall convective source under ideal-gas pressure recovery:

    s_E = K · (T_wall − T) with K = h · 4 / D_h T = (γ−1) / R_gas · (E_total − u² / 2)

Differentiating w.r.t. U = (ρ, ρu, ρE) gives three nonzero entries on row 2 (energy):

  ∂s_E/∂ρ      = +K · (γ−1)/R_gas · (E_total − u²) / ρ
  ∂s_E/∂(ρu)   = +K · (γ−1)/R_gas · u / ρ
  ∂s_E/∂(ρE)   = -K · (γ−1)/R_gas · 1 / ρ

Rows 0 / 1 are zero — constant-T_wall heat has no mass / momentum source. New aeolus1d.closures.heat.heat_jacobian writes those entries into a caller-provided buffer; the analytic_supported predicate gates the dispatch and the augmentation path additionally checks not pipe.has_real_gas (the closed-form holds only for ideal-gas pressure recovery — real-gas pipes need the EOS-derivative chain from B15.4 which the heat kernel doesn't yet consume). _cell_closure_source gains an include_heat flag symmetric to B16's include_friction; both default True for backward compat. Bench (uniform-state pipe, N=200, heat-only closures): _augment_jacobian_with_closures drops from 3213 µs (FD) to 667 µs (analytic) — a 4.82× speedup, comparable to B16's 4.46× on friction. New stacked test confirms B16 + B17 compose correctly when both are active on the same pipe.

  • Phase B16 — analytic friction-source Jacobian for DARCY_CONST. Replaces the FD-on-closure-source pass inside _augment_jacobian_with_closures with a closed-form ∂S_friction/∂U for the most common friction kind (constant Darcy factor):

    s_mom = -A·ρ·|u|·u with A = f / (2·D_h) s_E = u · s_mom = -A·ρ·u²·|u|

Differentiating w.r.t. U = (ρ, ρu, ρE) gives a 3×3 block whose nonzero entries are 4 polynomial expressions in (ρ, u). The new aeolus1d.closures.friction.friction_jacobian helper writes these into a caller-provided buffer; an analytic_supported predicate gates the dispatch so COLEBROOK / Python-override paths fall through to the existing FD treatment unchanged. _cell_closure_source gains an include_friction flag so the FD pass for OTHER closures (heat, porous, py overrides) skips friction when the analytic path adds it, avoiding double-counting. Bench (uniform-state pipe, N=200, friction-only closures): _augment_jacobian_with_closures drops from 3352 µs (FD) to 751 µs (analytic) — a 4.46× speedup. Tests pin the analytic ↔ FD agreement to FD-truncation tolerance and confirm the off-diagonal blocks stay untouched.

  • Phase VGT-3D-train — bundled training script. New scripts/train_learned_vgt_map.py closes the data-to- trained-model loop. Reads a manufacturer CSV / NPZ (Phase VGT-3D-loader schema), samples it via LearnedTurboMap. from_table, trains a 3-D MLP, and saves the result to NPZ for use via the existing [turbo.map] kind = "learned_npz" path. CLI knobs: --input (omit for the bundled synthetic VGT factory), --output, --n-pr / --n-n-corr / --n-vane (training-grid size), --hidden / --epochs / --lr / --seed (MLP). Prints a spot-check residual table at probe points so the user can sanity-check the surrogate before dropping it into a TOML. End-to-end-tested against the synthetic factory, a hand-written 3-D CSV, and a 2-D NPZ source (preserves backward compat with the 2-D learned-map path).

  • Phase VGT-3D-learned — 3-D MLP surrogate for VGT / IGV maps. LearnedTurboMap extended from a fixed 2-input network to a polymorphic 2-or-3-input network. The dispatch is data-driven — X.shape[1] at fit time determines the network topology; the trained instance carries has_vane_axis (mirrors TurboMap) and the same lookup(pr, n_corr, vane=...) semantics. The from_table factory auto-detects whether the source TurboMap has a vane axis and samples a 3-D (n_pr × n_n_corr × n_vane) grid accordingly. Backward compat: legacy 2-D NPZ files (including the bundled turbo_map_centrifugal.npz) load via the existing load_npz path; the new vane_min / vane_max fields are optional and default to 0 / 0 when absent. The 2-D inference path stays bit-exact. Lets users train a small NN surrogate on real measured VGT data, drop the NPZ into a TOML via kind = "learned_npz", and have the controller actuate a learned off-design map instead of the hand-coded synthetic one.

  • Phase VGT-3D-help — Inspector help-text for the new TurboMap kinds. Schema-help entries populated for the new TurboMapSpec kinds added by Phase VGT-3D / IGV-3D / loader / loader-gui: synthetic_vgt_turbine, synthetic_igv_centrifugal, from_npz, from_csv, plus polymorphic-npz_path and the csv_path / vane_levels fields. The Inspector ? button on each field now surfaces a focused blurb explaining what the kind does, when to pick it (e.g. from_npz for measured tabulated data vs learned_npz for a trained NN), and what the file format expects (CSV header row, NPZ array names, optional vane / vane_grid triggering the 3-D path). Closes the in-app onboarding gap left by Phase VGT-3D-loader-gui's silent field rendering.

  • Phase VGT-3D-loader-gui — TOML kinds + Inspector path picker. Two new TurboMapSpec kinds, from_csv and from_npz, route through import_turbo_map_csv / import_turbo_map_npz so users can drop a manufacturer CSV / NPZ map into a TOML without writing Python:

    [turbo.map] kind = "from_npz" npz_path = "/path/to/measured_map.npz"

Distinct from the existing learned_npz kind, which keeps routing through the trained-NN LearnedTurboMap loader; the _make_turbo_map dispatch picks the right path. Empty paths raise on build with a kind-specific error message. GUI Inspector: the form_factory now detects str-typed fields whose name ends in _path and renders them as a QLineEdit + Browse… button, with the QFileDialog filter set per field (*.csv for csv_path, *.npz for the rest). Plain string fields (id, etc.) keep their bare-line-edit rendering. Round-trip: aeolus1d_gui/io/case_io.py::save_case emits csv_path / npz_path for the new kinds and the vane_levels field for synthetic_vgt_turbine / synthetic_igv_centrifugal.

  • Phase VGT-3D-loader-demo — bundled round-trip example. examples/plots/plot_56_vgt_3d_loader.py (no TOML — pure Python). Builds a synthetic 3-D TurboMap via synthetic_vgt_turbine, serialises it to a tempdir as both CSV and NPZ, reloads each, and overlays source vs loaded mdot_corr(vane) / η_is(vane) at a fixed operating point. Diff panels (log-y) confirm the loaders are bit-exact at the grid nodes and machine-precision-equivalent off them. Footprint comparison: a 6 × 6 × 5 map is ~9 kB as CSV vs ~5 kB as NPZ.

  • Phase VGT-3D-loader — CSV / NPZ ingest for 3-D maps. Extends the existing 2-D loaders in aeolus1d.components.turbo_csv to detect the new third axis automatically. CSV: a vane column triggers the 3-D path; rows must form a complete (N_corr × PR × vane) cross-product or the loader raises with the same regular grid shape guard as the 2-D path. NPZ: an optional vane_grid array triggers the 3-D path; tables must be 3-D (N_n, N_pr, N_vane). save_turbo_map_npz emits vane_grid only when has_vane_axis is True so legacy 2-D NPZ files don't grow a spurious key, keeping pre-Phase-VGT-3D readers compatible. Lets users supply manufacturer-measured 3-D maps directly without authoring Python.

  • Phase IGV-3D — Compressor with optional vane axis. Symmetric extension to Phase VGT-3D on the compressor side. The Compressor dataclass gains a vane_open: float = 1.0 field (mirrors Turbine.vane_open); when the configured map has a vane axis (map.has_vane_axis), Compressor.apply_port_bcs passes the vane through the trilinear lookup so closing IGVs throttles ṁ_corr and shifts η_is per the map's off-design correction. 2-D maps stay vane-blind on the Compressor side (unlike Turbine which has a legacy linear-scaling fallback) — Compressor never had a 2-D linear-vane path so backward compat means vane_open is silently ignored when the map is 2-D. New TurboMap.synthetic_igv_centrifugal(vane_levels=...) factory: mdot ∝ vane^0.6 (slightly more direct than the VGT's vane^0.7 — IGV throttling is geometrically tighter)

  • Gaussian on (V − 0.8) for η so peak efficiency is near fully-open with a mild penalty as IGVs close. Closure-kind registry under family turbo_map registered with code 4. TOML knob: [turbo.map] kind = "synthetic_igv_centrifugal"; the build path now also threads tspec.vane_open into the Compressor constructor (was Turbine-only). Schema-help / unit-tests follow the VGT-3D template.

  • Phase VGT-3D-demo — bundled example for the 3-D VGT map. examples/55_vgt_3d_turbo.toml plus examples/plots/plot_55_vgt_3d.py exercises the new 3-D map end-to-end. Top half of the plot is the static envelope: mdot_corr(vane) and η_is(vane) at three (PR_inv, N_corr) operating points, showing the vane^0.7 mass-flow throttle + the η peak at vane ≈ 0.7. Bottom half loads the TOML, advances the network with an open-loop vane ramp 1.0 → 0.2 wired into on_step, and scatters the measured (mdot, η) trajectory against the static envelope — the run traces η climbing through the 0.7 peak and falling on the closed-vane off-design side. TOML is GUI-registered under Turbomachinery alongside 51_vgt_boost_ctrl and 54_adaptive_picard_turbo. (The boost controller is in the TOML but left at the rig's vane_max in this open-loop demo; closed-loop controller dynamics are exercised by examples 51 / 54.)

  • Phase VGT-3D — TurboMap with a vane axis. TurboMap gains an optional third vane_grid axis. When populated, mdot_corr / eta_is are 3-D (N_corr_grid, pr_grid, vane_grid) tables and lookup(pr, n_corr, vane=...) does trilinear interpolation on (PR, N_corr, vane_open). The Turbine code path branches on map.has_vane_axis: when True, the legacy linear ṁ *= vane_open rescaling drops out — the map carries the full off-design correction so a partly-closed vane lifts swallowing capacity and drops η_is together (rather than cleanly re-scaling mass flow as the 2-D

  • linear-scaling path does, which couldn't represent the efficiency penalty at all). 2-D maps stay exact for backward compat. New TurboMap.synthetic_vgt_turbine(vane_levels=...) factory builds a plausible map: mdot ∝ vane^0.7 (closing vanes throttles flow with diminishing returns) and a 1-D Gaussian on (V − 0.7) for η so peak efficiency lives at moderate vane (the rest is off-design). 5-level default axis is [0.1, 0.325, 0.55, 0.775, 1.0]. TOML knob: [turbo.map] kind = "synthetic_vgt_turbine" with optional vane_levels (default 5). Closure-kind registry under family turbo_map registered with code 3 alongside the existing synthetic_centrifugal (0), synthetic_turbine (1), learned_npz (2). The BoostController's actuation now produces a true off- design correction instead of a degenerate scale factor.

  • Phase I-Y-cli3 — total Newton-iter badge alongside the Picard badge. The CLI's --info-every-n line now emits newton = T after the existing picard = N/M token, where T is LiveControl.last_picard_total_newton_iters (sum across pipes × Picard iterations on the most recent step). The runner dock parses the new token into a Σ N: T badge in the header bar — symmetric to the Picard badge, hidden when T is 0 (no warm-start tracking populated). Reset on each new Run click alongside the Picard badge. Watching the two badges together tells users whether a saturated Picard run is also doing materially more Newton work — that's the signal to raise adaptive_picard_max.

  • Phase I-Y-cli2 — GUI runner-dock Picard badge. The CLI's --info-every-n line now emits picard = N/M (or just picard = N when no adaptive controller) where N is last_picard_iters and M is last_picard_max_iter_managed. RunnerDock parses the token off incoming stdout chunks and surfaces a Picard: N/M badge in the header bar next to the status label — blue when there's headroom (N < M), red when the adaptive controller has saturated at the cap (N == M). The badge resets on each new Run click and findall-handles multiline stdout chunks (always picks the latest match). No-op for explicit drivers and pre-Phase-I-Y LiveControl.

  • Phase I-Y-adaptive-demo — bundled example for the adaptive Picard controller. examples/54_adaptive_picard_turbo.toml

  • examples/plots/plot_54_adaptive_picard.py exercises the controller on a turbine + compressor + boost-PI network under implicit-Euler at CFL = 2. The initial transient saturates the budget (used == 4 == max) → controller grows to 6. As ghost-coupling residual drops, the controller ratchets back through 5 → 4 → 3 → 2 (the plenum-only natural floor where used + 1 > current) and eventually 2 → 1 once the network reaches steady state with zero Newton work per step. 4-panel plot: managed budget, picard_iters_used, ghost-change residual (log scale), and total Newton iters per step. The TOML is GUI-registered under Turbomachinery alongside 51_vgt_boost_ctrl so it shows up in the example dropdown.

  • Phase I-Y-adaptive — adaptive picard_max_iter controller. Mirrors the existing adaptive_dt pattern. When adaptive_picard = true is set in [case], the implicit-network driver adjusts the effective picard_max_iter per step: saturated = (used >= max AND ghost_change > picard_tol) → grow by grow_step; used + shrink_step <= max for adaptive_picard_n_easy consecutive steps → shrink by shrink_step. Bounded by [adaptive_picard_min, adaptive_picard_max]. "Easy" is defined by budget headroom, not by an absolute iteration count — the convergence-check iteration always fires for networks with far-BC ghosts that depend on the pipe interior, so the typical converged step has picard_iters_used == 2. The headroom-based shrink lets the controller settle cleanly even in that regime: a plenum-only network started at picard_max_iter = 4 walks down to 2 (the floor where used + 1 > current) and holds. Surfaces via LiveControl.last_picard_max_iter_managed alongside the existing last_picard_iters / last_picard_total_newton_iters / last_picard_ghost_change fields, so users can watch the controller settle. TOML knobs: adaptive_picard_min (default 1), adaptive_picard_max (default 8), adaptive_picard_n_easy (default 5). GUI Inspector entries populated via schema_help.

  • Phase I-Y-cache — hoist BC chain cache out of the Picard loop. _seed_jac_chain_cache is documented as "FD-perturbed against the start-of-step interior state and reused across Newton iterations" — that reference state is invariant across Picard re-iterations, so the seed call (one per pipe per Picard iter under Phase I-Y-warm) collapses to one per pipe per step. Modest perf: advance_network_implicit on a 2-pipe Compressor case at CFL 2, 2 ms run drops from 645 ms → 609 ms (~6 %) at picard_max_iter = 10 and from 280 ms → 277 ms (~1 %) at picard_max_iter = 3. Single-iter runs see no change. Cache also stays correct under Newton-failure recovery and dt-halving (the chain depends only on the seed-time interior state and BC kind, not on dt).

  • Phase I-Y-warm — warm-started Newton across Picard iterations. The Phase I-Y outer loop used to reset pipe.U to U_n (the start-of-step snapshot) on every Picard re-iteration, so each re-iter's Newton walked the full distance from U_n to the new converged state. Phase I-Y-warm removes that reset: the second / third / Nth Picard iter starts Newton from the previous iter's converged solution, which is by construction close to the new answer (the only thing that changed between iterations is the refreshed ghost coupling). The implicit equation solved is unchanged — U_hist (the residual's history term) still rolls back to U_n — so the fixed point is bit-identical; only the path Newton takes to reach it gets shorter. Newton-failure recovery inside the dt-halving retry loop still resets pipe.U[:] = U_n, so a divergent re-iter falls back to cold-start automatically without any extra plumbing. LiveControl gains a third diagnostic field, last_picard_total_newton_iters, summed across pipes × Picard iters per step. On a typical Compressor case at CFL 2, 10 Picard iters at picard_tol = 1e-10 use 22 total Newton iters (avg 2.2/iter) vs the cold-start estimate of ~40 (avg ~4/iter) — a 1.8× reduction in inner-loop work when Picard saturates a tight-tol budget. Plenum-only networks (where Picard exits at iter 2 anyway) see no perf change.

  • Phase I-Y-cli — CLI / GUI surfacing of Picard diagnostics. The CLI's --info-every-n line now appends picard = N when running under the implicit-network driver (where LiveControl.last_picard_iters gets written each step); explicit-driver runs leave it off so the historic line shape is preserved. _StepHooks reads last_picard_iters from the LiveControl instance the CLI already builds for --live-control, so no new wiring at the user level — enabling --info-every-n (or the GUI's "Info every N steps" checkbox) is enough. GUI: the RunnerDock checkbox tooltip notes the new token and points users at picard_max_iter in [case] when the budget saturates. Console output flows through the existing QProcessRunner stdout pipe — no GUI parser changes needed for the basic display.

  • Phase I-Y-toml — TOML schema for the Picard knobs. The Phase I-Y outer-iteration controls (picard_max_iter / picard_tol) now land on CaseHeader, propagate into the NewtonConfig the implicit network driver reads, and round-trip through the GUI's case_io.save_case (emitted only when non-default so explicit-only cases stay clean). Defaults match the dataclass-level Phase I-Y defaults (picard_max_iter = 3, picard_tol = 1e-6); set picard_max_iter = 1 to recover Phase I-X block-Gauss-Seidel-without-iteration. Schema-help and field- units entries are populated so the GUI Inspector explains the knobs in-app.

  • Phase I-Y-diag — Picard diagnostics on the implicit network step report. Surfaces the most recent step's Picard iteration count + final ghost-state change so outer callers can see whether a network actually needs picard_max_iter > 1 (the Phase I-Y default of 3) without instrumenting the driver. Two paths: LiveControl gains two new diagnostic fields, last_picard_iters and last_picard_ghost_change, which the implicit network driver writes after each completed step. Outer code (CLI run-progress printer, GUI runner-dock console, custom drivers) can poll these between steps. Newton-failure messages now include a Picard (Phase I-Y): on iteration N of M (tol = ...) line in the solver-state block, plus a new remedy suggestion [F] **Raise picard_max_iter.** when the failure happened on the last Picard iteration of a multi-iter budget. Lets users tell at a glance whether the failure was a Newton-blowup inside the first Picard iter (so the inner step is the problem) versus a Picard-non-convergence (so the outer coupling is the problem).

  • Phase 7d-4-demo / 8m-ctrl-demo — bundled examples + plot scripts. Two new visualisations exercise the knock-margin overlay and the DPF heater controller end-to-end: examples/plots/plot_52_knock_margin.py builds a boosted single-cylinder SI engine (CR=10, 1.4 bar boost, 1800 rpm) that drives the Livengood-Wu integral past 1.0 in early cycles; the closed-loop retard controller walks θ_spark from 8° BTDC to TDC over the first ~5 cycles, then saturates. Plots show knock_integral_max, knock_fired, timing_retard_deg, and IMEP per cycle — proof-of-concept for the Phase 7d-4 overlay. examples/53_dpf_heater_ctrl.toml + examples/plots/ plot_53_dpf_heater_ctrl.py exercises the new DPFHeaterController (velocity-form PI) on a single DPF pipe: cake starts at 750 K, target 870 K, P_max 3 kW; the pure-I controller saturates power, drives T_cake across T_target around 530 ms with a ~60 K overshoot, then backs off as the error sign flips. The TOML is GUI-registered under Aftertreatment alongside 50_dpf_active_regen so it shows up in the example dropdown.

  • Phase 8m-ctrl — DPF active-regen heater controller. Velocity-form PI controller in aeolus1d.control.DPFHeaterController, mirroring the BoostController pattern. Drives mean filter.T_cake toward T_target by modulating pipe.filter.heater_power_w; output clamped to [P_min, P_max]. Pass an instance as the on_step= callback to advance_network (or chain it via the case's [[controller]] block). The closure-side soot-eligibility gate (heater engages only when mean soot ≥ soot_trigger_kg_m3) is unchanged — the controller modulates the power level; soot trigger is a safety / readiness gate. Set soot_trigger_kg_m3 = 0.0 for an always-on configuration (useful for soak / pre-heat sequences). TOML: [[controller]] blocks accept kind = "dpf_heater" with required fields dpf_pipe_id, T_target, Kp and optional Ki / P_min / P_max. Validation rejects controllers pointing at a pipe without a [pipe.filter] block of kind = "diesel_particulate_filter". No-ops cleanly when the cake-thermal path is inactive (Phase 8m's direct-to-gas heat path doesn't track T_cake).

  • Phase 7d-4 — knock-margin overlay on the cycle accumulator. Three new fields on EngineMetrics expose the cylinder's per-cycle knock state through CycleAccumulator: knock_integral_max (peak Livengood-Wu integral reached during the cycle), knock_fired (latched True if cylinder.last_knock_fired was True at any point — the per-step flag itself resets at IVC), and timing_retard_deg (the spark retard at cycle close, i.e. what the controller will apply to the next cycle). cycle_summary reports the aggregates knock_fired_rate (frac of cycles that knocked), max_knock_integral_max (worst-case peak across the window — comfortably under 1.0 means the sweep had knock margin), and mean_timing_retard_deg. Fleet-level MultiCycleAccumulator.fleet_summary aggregates these as means except max_knock_integral_max, which uses max-across-cylinders so the worst-cylinder peak isn't lost to averaging. Defensive getattr on the cylinder fields keeps mock / stub cylinders without knock plumbing working (they report all knock fields as zero / False).

  • Phase I-Y — Picard outer iteration over pipe ↔ component ghost coupling. Phase I-X dispatched every component through the apply_port_bcs / accumulate_rates protocol but held the ghosts FIXED across the pipe Newton solve (block Gauss- Seidel without iteration). For Compressor / Turbine ports the ghost depends on the pipe interior — perturbations to the interior during Newton don't propagate back to the ghost, giving an O(Δt) splitting error on the turbo coupling. Same splitting applies to LODI ports (characteristic ghost from the pipe outgoing invariant) and py_junction_bc callbacks. Phase I-Y wraps the pipe Newton loop and the post-Newton apply_port_bcs refresh in a Picard iteration: for picard_iter in range(picard_max_iter): restart pipes from start-of-step → solve every pipe Newton → refresh ghosts → check max-relative-change on ghosts vs picard_tol → break or continue. At convergence the ghosts and pipe interiors satisfy the same implicit relation, removing the O(Δt) error. Two new NewtonConfig knobs: picard_max_iter (default 3) and picard_tol (default 1e-6). picard_max_iter == 1 reproduces Phase I-X exactly. The default picard_max_iter == 3 is enough for typical engine configs — plenum-only networks exit at iter 2 (ghosts are invariant under iteration because they read the frozen volume state); turbo configs typically exit at iter 2-3. Cylinder ports are unaffected because their ghost reads only the Cylinder volume state (frozen during the step), not the pipe interior.

  • Turbo-map visualizer. Post-Phase-9 follow-on. New module aeolus1d.analysis.turbo_map_viz renders manufacturer-style compressor / turbine maps for any TurboMap: iso-speed lines (each row of the table → one curve), iso-efficiency contours via ax.tricontour, the surge / choke envelope from the existing turbo_map_envelope helper, and an optional operating-point overlay. Compressor maps orient (ṁ_corr on x, PR on y); turbine maps swap to (PR on x, ṁ_corr on y) per industry convention. CLI: aeolus1d turbo-map <path>.cmp|.trb|.npz [--out img.png] registered in cli/main.py (auto-classifies by extension, toggles for --no-envelope / --no-eta, configurable DPI). New files: src/aeolus1d/analysis/turbo_map_viz.py, tests/unit/test_turbo_map_viz.py (7 tests covering minimal fixture rendering, axis orientation per kind, overlay toggles, operating-point scatter, real-F2 map render smoke-tests, unknown-extension error).

  • Phase 9h Tier-1 — 3-zone cylinder wall heat transfer. New aeolus1d physics. HeatTransferModel gains optional T_head / T_liner / T_piston fields (default 0.0 = fall back to bulk T_wall); when any are non-zero the cylinder's wall-heat block dispatches to a 3-zone Q̇ split:

    Q_head   = h · A_bore · (T_head   − T_gas)
    Q_piston = h · A_bore · (T_piston − T_gas)
    Q_liner  = h · π·B·L(θ) · (T_liner − T_gas)
    
    Gas-side h stays the bulk Woschni / Hohenberg correlation (Tier-2 zone-specific h is deferred). When all 3 zone Ts collapse onto T_wall the formula degenerates to the pre-9h single-zone Q̇ bit-exact — 15 pre-existing cylinder + cylinder- heat-transfer tests pass unchanged. Schema integration: CylinderSpec gets optional T_wall, T_head, T_liner, T_piston (all default 0.0 = legacy adiabatic). _build_cylinders constructs a WoschniHeatTransfer with the zone Ts when any are set. NeutralCylinder (and schematic_to_case) pass the zone Ts through. GT-Suite importer wiring: each emitted F2 cylinder receives diesel-grade zone defaults (T_head=600 K, T_liner=450 K, T_piston=550 K) sourced from Heywood §12.4 / GT-Suite production norms. Users edit per operating point. New file: tests/unit/test_cylinder_three_zone_wall.py (6 tests covering the helpers, bit-exact single-zone fallback, analytic 3-zone sum, dispatch asymmetry, and per-zone T_wall fallback).

  • Phase 9g — GT-Suite controller sub-assembly translation. F2-class cases reference 6-8 .gtsub controllers (wastegatecontrol / throttleCTRL / EXHBCTRL / PRAIL / SOI / deltactrl4 / HTCcalc / FMEPcalc) via <object t_id="314"> placeholders whose <va a_id="3079"> text names the sub-assembly file. Phase 9g maps the one with a direct aeolus1d analog (wastegate → WastegateController) and silent-skips the rest with a per-controller InteropReport lossy entry naming the dropped file + pointing the user at Phase 7e FMI co-simulation.

The wastegate translation picks the schematic's first turbine + the compressor's outlet pipe as the (turbine_id, manifold_pipe_id) target, with default gains Kp=1e-5 / Ki=1e-6 and a 2-bar boost setpoint placeholder. When the schematic carries no turbine + compressor pair the wastegate drops cleanly with an explanatory lossy entry.

On the F2 driver case 1 controller wires (the wastegate) and 7 silent-skip with named drops (deltactrl4, HTCcalc, FMEPcalc, SOI, PRAIL, throttleCTRL, EXHBCTRL). HTCcalc + FMEPcalc are physics models that already inform the Phase 9d thermal + the deferred cylinder-side phases; the rest belong on the control-system side. New file: tests/unit/test_gt_suite_native_phase_9g.py (3 tests covering the wastegate emission, the per-file lossy entries, and the no-turbo bail).

  • Phase 9d — GT-Suite thermal lookup-table dispatch. GT-Suite's thermal lookup tables (template id 16175, sample names avgport / inport / steelwall / cachcond / laddror / avgsam / bellow) carry HTC + coolant-T refs; Phase 9d splits the lookup-table family into lookup_table_thermal (16175) + lookup_table_material (16176), pulls the HTC value from each thermal table, and heuristically matches pipes to tables by longest-substring name match (so cachcond wins over cac on a cachcondPipe-1). Matched pipes get a constant_twall HeatSpec attached. NeutralHeatExchanger.kind was extended with constant_twall; schematic_to_case dispatches on kind to emit either an ε-NTU CAC/EGR closure or a direct h+T_wall spec on the target pipe. On the F2 driver case 10 of the 59 pipes now carry a real thermal closure — the 6 intake-port pipes (h=3000 W/m²K via inport) and 4 bellow pipes (h=35 W/m²K via bellow). The wall T defaults to 350 K because the GT [twater] / [TAMB_eng] parameter-block refs aren't yet resolved; a warning surfaces in the InteropReport. New file: tests/unit/test_gt_suite_native_phase_9d.py (6 tests covering classifier split, HTC extraction, low-h vs high-h match, no-match default, warning emission, longest substring greediness).

  • Phase 9f — GT-Suite .cmp / .trb ASCII map readers. New module aeolus1d.io.interop.gt_maps parses GT-Wide ASCII compressor (.cmp v7.x) and turbine (.trb v2016) maps into populated TurboMap objects. The reader handles the GT KEY> value header block, the DATA> marker, and the named numeric blocks (ARPM, APPR, PMAX, PMIN, MASS, EFFI), each terminated by <. Per-iso-speed actual PR ranges (where PR_actual(i,j) = pmin[i] + (pmax[i] - pmin[i]) * appr[j]) get re-sampled onto a regular (N_corr, PR) grid via piecewise-linear interpolation with edge clipping for transient overshoot. Turbine maps with PRMIN> no correctly default PMIN to 1.0 per iso-speed line. The GT-Suite native importer (read()) is wired to call the new readers: compressor parts pick up their map filename from <va a_id="17">, turbine parts from <va a_id="259">. Both use the GT <![CDATA[<NAME.ext>]]> bracket convention — Phase 9f also fixes _va_text to fall back from the text= attribute to the element's CDATA body so the references parse. Each parsed map is serialised next to the source as NAME.aeolus1d.npz and the resulting NeutralTurbo gets map_kind="from_npz" with the NPZ path; aeolus1d's existing from_npz TurboMap loader picks it up at build time. On the F2 driver case both turbos (CHP and THP-1) now carry their real BWT91L82 / F1.5-C9 maps (60×60 grids, peak efficiency ~0.65 / ~0.80) instead of synthetic defaults. New file: tests/unit/test_gt_maps_phase_9f.py (7 tests, including 3 skip-on-missing tests against the real F2 maps).

  • Phase 9i — GT-Suite end-to-end integration test + demo script. Caps the Phase 9j/9a/9b/9c/9e arc with a single end-to-end test pinning the user-visible goal: a GT-Suite native .gtm loads, builds a network, and steps the timeloop.

  • tests/integration/test_gt_suite_native_full_engine.py — skip-on-missing fixture test that asserts the F2 case imports to a Case with the expected counts (59 pipes, 109 junctions, 2 ambient BCs, 6 cylinders with 1 intake + 2 exhaust ports each, 2 turbos + 1 shared shaft) AND advances at least 5 timeloop steps without raising. Past step ~5 the default cylinder fuelling (0 kg/cycle) leaves the engine air-breathing and the turbo on a synthetic map; continued stability needs Phases 9b-2 / 9f / 9h's fidelity work.
  • examples/plots/import_gt_suite_native.py — runnable demo that points at any GT-Suite .gtm and prints the template inspection table + the post-import Case structure (cylinder geometry derivation, turbo wiring, interop-report warnings). The F2-BWC9-100 driver case isn't bundled — carries vendor IP — but the script works against it. Status: the F2 case loads and runs. Cylinder bore / stroke / con-rod / clearance reflect the cylgeom values exactly (130 mm bore, 22.9 CR, ~97 cm³ clearance); each cylinder has the intake + 2 exhaust pipes wired; ambient BCs sit at both breathing-path termini; the time loop advances. Subsequent phases (9d lookup-tables, 9f .cmp / .trb maps, 9g controllers, 9h cylinder wall) tighten fidelity.

  • Phase 9e — GT-Suite cylinder port wiring through valveconns. Originally scoped as "orifice / valveconn → BC translation" for an EngExhBrake top-level part, but the F2 case turned out to carry the exhaust-brake purely as .gtsub controller logic — the bigger user-visible gap was that the 6 imported cylinders had no ports and were isolated chambers. Phase 9e closes that gap. GT-Suite represents each intake / exhaust valve as a valveconn part (t_id=992 intake / 997 exhaust on F2) that links to both a cylinder and an adjacent manifold flowsplit. Sometimes the chain has further hops (intake side: cyl → valveconn → flowsplit → orifice → pipe). The new Pass 2f walks each cylinder's links → valveconn neighbours, then BFSes from the valveconn (depth-limited to 4 hops, visited- set to break cycles) through junction-like parts until it hits a pipe. The first pipe found becomes a NeutralCylinderPort on the cylinder. CAD windows default to typical 4-stroke diesel values (intake 350°–580°, exhaust 130°–380°); valve geometry + lift schedules land in Phase 9e-2 from the GT valve table. Junction ports that double-up on a cylinder-owned pipe-end are stripped; the wall-fallback synthesis pass also counts cylinder ports as "referenced" so no spurious wall BC overlaps a cylinder port. On the F2 driver case each cylinder picks up 3 ports — 1 intake (the shared InPort flowsplit feeds both IV valveconns at the same pipe-end) + 2 exhaust (separate ExP* pipes per EV valveconn, matching the twin-exhaust diesel geometry). Schematic counts go from junctions=127, bcs=2 (post-9c) to junctions=109, bcs=2 with 18 cylinder ports across the 6 cylinders. build_network succeeds; the case is now a fully wired fired-engine model. New file: tests/unit/test_gt_suite_native_phase_9e.py (5 tests). Touched: tests/unit/test_gt_suite_native.py for the post-9e production-fixture counts (109 junctions

  • per-cylinder intake / exhaust port assertions).

  • Phase 9c — GT-Suite EndEnvironment auto-attach. Two template-id corrections + a new orifice-chain wiring pass give the GT-Suite native importer real ambient BCs:

  • t_id=4 (sample names AMBIN / ambut / ambient-exh) is reclassified from param to endenvironment — this is the actual GT-Suite "Ambient" template that the F2 case uses at the intake-side air-filter inlet and exhaust-side outlet.
  • t_id=326 (sample names CDEXHFWDF2 / CDEXHREVF2) is demoted from endenvironment to a new endenvironment_ref class — these are internal state-reference objects pointed at from exhaust valveconn <va> blocks, never instantiated as part-level BCs.
  • New Pass 2e walks each ambient part's link to the adjacent orifice-passthrough junction, picks the first wired pipe-end on that junction, and stamps a stagnation_inlet / pressure_outlet BC there. Kind heuristic: part name containing out / exh or ending in ut → exhaust outlet; otherwise intake inlet. The orifice's port to the BC'd pipe-end is then stripped (otherwise the case validator rejects the BC + junction-port collision); orphan orifices that lose all ports get pruned.
  • The wall-fallback synthesis at Pass 3 also gained a guard so pipe-ends with a real BC don't get a duplicate wall BC.

On the F2 driver case the schematic now reads pipes=59 / junctions=127 / bcs=2 (down from 159 junctions because Phase 9c pruned 32 orphan orifices) with the 2 BCs being a real stagnation_inlet on befLF (intake) and a pressure_outlet on the exhaust terminus. Ambient pressures / temperatures default to 101.325 kPa / 298.15 K — the GT-Suite [PAMB] / [TAMB] parameter-block resolution lands in a later phase. A warning surfaces in the InteropReport. New file: tests/unit/test_gt_suite_native_phase_9c.py (6 tests). Touched: tests/unit/test_gt_suite_native.py for the post-9c production-fixture counts (159 → 127 junctions; 0 → 2 BCs) + the t_id=326 reclassification.

  • Phase 9b — GT-Suite engine breathing-path templates (MVP). gt_suite_native.read() now classifies and emits the engine-side templates that the F2-BWC9-100 driver case ships: EngCylinder (t_id=12), CylGeometry sub-template (t_id=112), Compressor (t_id=2 'CHP'), Turbine (t_id=24 'TurbineLeft'), generic Shaft (t_id=163), turbocharger shaft TCshaft (t_id=164), turbine volute (t_id=17), and the previously misclassified exhaust valveconn (t_id=997 -> valveconn). A new _resolve_cylinder_geometry extracts bore / stroke / con-rod length / compression ratio from the cylgeom sub- object; clearance volume is derived as V_disp / (CR - 1) and stamped onto every cylinder part. Compressor / turbine parts wire to their pipes through the GT-Suite orifice-passthrough chain (pipe ↔ orifice ↔ turbo); a Pass-2b traversal walks turbo-neighbour links, hopping through the intermediate junction's first port to find the upstream / downstream pipe id. Orphan turbos that can't fully resolve both pipes drop with a lossy-report entry rather than producing an invalid Case. Single-shaft cases auto-wire every turbo to the only available NeutralShaft. On the F2 fixture this lifts the schematic from pipes=59 / junctions=145 / bcs=12 to pipes=59 / junctions=159 / bcs=0 plus cylinders=6 / turbos=2 / shafts=1. The 6 cylinders parse with bore=130 mm, stroke=160 mm, con-rod=255 mm, V_clear≈97 cm³ (CR=22.9) — exact match to the case decode. The compressor wires CHPUT→CHPIN, the turbine wires THPvolute1→THPvolute2, both share HPSHAFT (inertia=2.988e-4 kg·m²). New file: tests/unit/test_gt_suite_native_phase_9b.py (6 tests). Touched: tests/unit/test_gt_suite_native.py for the post-9b production-fixture counts. Deferred to later phases: cylinder port wiring through valveconns (9b-2), <va>-block crank-speed parsing (9b-2), injector schedule + Wiebe params (9b-2 / 9g for SOI controller), turbo-map .cmp / .trb loading (9f), and proper EndEnvironment BCs (9c).

  • Phase 9a — GT-Suite .gtsub sub-assembly recursion. The GT-Suite native importer's inspect() and read() now walk the top-level .gtm's <object t_id="314"> sub-assembly placeholders, follow the <va a_id="3079" text="&lt;NAME.gtsub&gt;"/> filename pointer to the sibling file, and recursively in-line every referenced .gtsub before the part walk. IDs from each sub-tree are namespaced by file stem (o60003 → o{stem}_60003, preserving the leading prefix so the existing strip-alias logic still resolves sub-internal cross-references) and parts / links merge into the parent's <assembly_table><main><part_table> / <port_link_table>. Cycle detection via a visited-paths set; missing-sibling and cycle-detection conditions surface as warnings on the InteropReport. On the F2-BWC9-100 driver case this lifts inspect.n_objects from 323 to 452 (+129 from 8 sub-assemblies) and surfaces control_signal / controller / param objects that were previously invisible — the breathing-path part counts stay the same because the F2 sub-assemblies are all ECU-style control logic, not flow primitives. Phase 9b lights up cylinder / turbo / shaft templates next. New file: tests/unit/test_gt_suite_native_phase_9a.py (9 tests covering merge, namespacing, missing-file warning, cycle detection, no-op on plain .gtm, and the _namespace_id helper). Touched tests/unit/test_gt_suite_native.py to refresh the headline count assertions for the post-recursion model.

  • Phase 9j — Neutral schematic carries engine-side components. The aeolus1d.io.interop.schematic neutral layer gained NeutralCylinder (+ NeutralCylinderPort), NeutralShaft, NeutralTurbo, NeutralHeatExchanger, and NeutralController dataclasses, plus matching Schematic.cylinders / .shafts / .turbos / .heat_exchangers / .controllers lists. The schematic_to_case arm now materialises each onto the corresponding Case field (Cylinder → CylinderSpec + optional WiebeSpec; turbo pair → TurboSpec with TurboMapSpec dispatch; heat-exchanger → stamped onto the named pipe's HeatSpec; controller → ControllerSpec). The reverse direction (schematic_from_case) round-trips them too — turbos / cylinders / shafts / controllers are no longer reported as "lossy components" on export. CLI registration for the gt_suite_native import direction (already wired at cli/main.py:904 and the _INTEROP_TOOLS registry) becomes end-to-end usable as soon as the GT importer (Phase 9a / 9b) starts populating the new lists. Unblocks Phase 9b (engine breathing-path templates). New file: tests/unit/test_schematic_phase_9j.py (8 tests). Touched: src/aeolus1d/io/interop/schematic.py, src/aeolus1d/io/interop/__init__.py.

  • Phase 27c — Packed multi-pipe parallel residual kernel. All fast pipes (first-order ideal-gas, uniform area, no closures, no species) are packed into contiguous arrays U_pad_all[P, max_N+2, 3] and evaluated in a single @njit(parallel=True) kernel with prange over the pipe index. The classification and pre-allocation happen once before the time loop; each stage packs the ghost-filled U_pads (~4 μs for 13 pipes), then fires a single Numba dispatch instead of N separate dispatches. Kernel-level speedup 4.2× on 13 pipes × 30 cells (11.6 μs packed vs 48.6 μs sequential); full wall-time speedup is hardware-dependent (bare-metal Linux benefits most; WSL2 thread wake-up cost ~160 μs per prange call masks the gain for small networks). Slow pipes (EOS / quasi-1D / MUSCL / py-closures) still run through the sequential dispatch. All 519 unit tests pass bit-exact. New files: src/aeolus1d/solver/packed_pipes.py, tests/unit/test_packed_pipes.py.

Milestone 7 (planned) — feature parity with commercial 1D engine tools

Distilled from a 2026-04-27 literature scan of GT-POWER / AVL BOOST / Ricardo WAVE feature sets, recent SCR / aftertreatment research, and the alternative-fuels surge (H2 / NH3 / dual-fuel). Nine sub-phases ranked by impact relative to effort, with concrete scope notes in docs/roadmap.md. The Phase 25k → 25q chain (cylinder injection → blowback → manifold transport → plenum tracking → turbo pass-through → catalyst HC oxidation) and Phase 26 case-level [[species]] metadata table laid the architectural groundwork; Milestone 7 is mostly content + presets on top of the existing closure seams.

  • Phase 7a — Alternative-fuel kinetics. H2 / NH3 / NH3-H2 blend Wiebe parameter presets, multi-injector cylinders for dual-fuel cases, NH3 NOx-production kinetics for engines using ammonia as the primary fuel. Effort: ~1 week (presets + multi-injector); ~2 weeks (NH3 NOx kinetics).
  • Phase 7b — Aftertreatment depth. (Substantially complete: Phase 7b-1 UWS injector, Phase 7b-2 SCRF + V2O5 slip catalyst + standalone DOC example all shipped.) Remaining: HNCO hydrolysis kinetics, active DPF regeneration (burner BC), cold-start warm-up model. Effort: ~1-2 weeks per remaining component.
  • Phase 7c — Turbocharger depth. Variable Geometry Turbine (vane position as 3rd map axis), wastegate (parallel-pipe bypass with controlled-area orifice), twin / sequential topology. Compressor surge / choke detection. Effort: ~1 week (VGT); ~1 week (wastegate); rest is example work.
  • Phase 7d — SI knock + EGR cooler. Livengood-Wu integral for knock-onset prediction sitting on Phase 25j two-zone state; effectiveness-NTU EGR cooler component. Effort: ~2 days (knock); ~1 week (EGR cooler). Multi-mode combustion (HCCI / RCCI) is its own multi-week project.
  • Phase 7e — Co-simulation interfaces. FMI 2.0 / 3.0 export (was the deferred Phase 40 — promoted) so an aeolus1d case becomes an FMU consumable by GT-SUITE / Simulink / Modelica tooling; Simulink coupling; structural import from GT-POWER / AVL BOOST / WAVE via a neutral schematic. Effort: ~2 weeks (FMI); per-format import is its own multi-week project.
  • Phase 7f — Hardware-in-the-loop / real-time execution. Fixed-step integrator + UDP / shared-memory I/O for ECU test rigs. Phase 27c (packed multi-pipe Numba kernel) is now shipped — per-step deadline budget can be re-evaluated. Effort: ~6-10 weeks total. Only commit if there's explicit demand for HIL.
  • Phase 7g — Engine-side post-processing. (Shipped: IMEP / ISFC / η_i + p–V trace via CycleAccumulator; see milestone 6 entry above.) Remaining: knock-margin overlay (pairs with Phase 7d Livengood-Wu), specific-emissions integrators, drive-cycle integrators (WLTC / RDE / FTP-75). Effort: ~1 week per remaining item.
  • Phase 7h — Physics-informed ML closures. Wrapper preset on the existing py_reaction_multi seam that loads a trained PINN (PyTorch / JAX / ONNX) as a rate function; worked example + training script for a small surrogate of a published mechanism. Bigger value paired with the JAX-mirror speculative item — backprop through the solver enables joint PINN-plus-physics calibration. Effort: ~2 weeks.
  • Phase 7i — 1D / 3D CFD coupling. Hybrid loops where one component (cylinder, manifold collector, EGR mixer) runs in 3D CFD while the rest stays in 1D. Coupling via ZMQ or direct memory; 1D side provides BCs, 3D side returns aggregated quantities. Effort: ~6-10 weeks. Specialist — only commit when a specific case demands it.

Milestone 8 (planned) — Cruise M import completeness

A coverage-driven companion to Milestone 7, organised around AVL Cruise M element types rather than feature requests. The acceptance criterion for each phase is "this Cruise M element type imports cleanly into a runnable aeolus1d case." Splits into three buckets:

  • Importer-side gaps (Phases 8a + 8b) — aeolus1d already has the physics; the importer just doesn't recognise the topology elements. ~5 days total.
  • Aeolus feature gaps (Phases 8c-8i) — physics aeolus1d doesn't have today: multi-injector cylinders, ε-NTU heat exchanger (intercooler / EGR cooler), urea injector + spray decomposition, DPF + soot loading, VGT + wastegate, combined SCRF, particulate-matter species column. ~12-16 weeks total.
  • Out-of-scope items via FMI (Phase 8z) — vehicle dynamics, transmission, hybrid drive, cooling system, and ECU controller logic don't belong in a 1D gas-dynamics solver. Phase 7e's FMI export hands these off to a coupled-simulation tool.

Full plan in docs/milestone-8-cruise-import-completeness.md . Cross- references each Milestone 7 phase that addresses the same physics from a feature-driven angle so the two milestones can be scheduled together.

Milestone 9 (planned) — GT-Suite import completeness

A coverage-driven companion to Milestone 8, this time targeting Gamma Technologies' GT-POWER / GT-SUITE — the biggest commercial competitor to AVL Cruise M. Calibrated against a representative production heavy-duty diesel GT-Suite case (~940 parts in the resolved model; ~260 in the top-level .gtm; ~70% live in .gtsub sub-templates the importer doesn't yet follow). Vendor case files carry IP and are kept outside the repo (gitignored); the integration tests skip when the local fixture is absent.

  • GT-Suite-specific importer phases (9a-9h) — .gtsub sub-assembly resolution, engine-breathing template parsing (cylinder, valve, injector, compressor, turbine, shaft), EndEnvironment auto-attach, lookup-table → closure dispatch, orifice / valveconn translation, binary .cmp / .trb turbo-map readers, controller silent-skip with warnings, cylinder wall thermal conduction (Tier-1 3-zone or Tier-2 1-D radial). ~6-9 weeks total.
  • Shared physics with Milestone 8 — multi-injector (8c), heat-exchanger ε-NTU (8d), DPF (8f), VGT (8g), particulate-matter species (8i). Each scheduled once and counted against both milestones.
  • Out-of-scope items via FMI — same Phase 8z / 7e hand-off as Milestone 8.

Full plan in docs/milestone-9-gt-suite-import-completeness.md .