Implicit-solver reference — Newton paths + knobs¶

A practical reference for the implicit-solver knobs in NewtonConfig / [case]. This doc is the what to set companion to numerics.md, which explains why each knob is the shape it is.

If you're tuning convergence on a stiff case, chord-Jacobian knobs is the section. If you're benchmarking, the [cost-shape comparison](

cost-shape-comparison) gives the empirical baseline.¶

Which Newton-solve path do I want?¶

There is one implicit Newton path: the direct block-tridiagonal Thomas solve with chord-Jacobian reuse. It is the right cost shape for the ≤ a few-thousand-cells × 3-unknowns systems aeolus1d solves — the matrix-free JFNK variants that once shipped here were measured slower than direct at every N tested and were removed in v0.9 (see the changelog). So there is no path to pick: tune the chord-Jacobian aging if a stiff case needs it, and otherwise use the defaults.

The reasoning, with bench data, is in numerics.md.

NewtonConfig knob reference¶

Most flags are settable both on NewtonConfig (Python API) and as [case] keys (TOML), but the names do not always match, and a few knobs are Python-API-only. The mapping lives in aeolus1d.io.case._implicit_configs_from_case:

tol ↔ TOML newton_tol and max_iters ↔ TOML newton_max_iters (renamed on the case header).
picard_max_iter, picard_tol, chord_max_age, chord_rate_threshold, parallel_pipes and parallel_pipes_auto_threshold share the same name in both.
backtrack_max, backtrack_factor, eps_rel_jac and dt_halvings_max have no [case] key — they always take the NewtonConfig default and can only be changed from the Python API.

Newton iteration¶

Flag	TOML key	Default	Purpose
`tol`	`newton_tol`	`1.0e-8`	`‖F(U)‖ / ‖U‖` convergence test for the outer Newton. Loosen for transient stability; tighten for accuracy-critical runs.
`max_iters`	`newton_max_iters`	`8`	Cap on inner Newton iterations per step. Hitting the cap triggers a Δt-halving retry (up to `dt_halvings_max` times).
`backtrack_max`	(API only)	`5`	Backtracking line-search retries per Newton iteration.
`backtrack_factor`	(API only)	`0.5`	Backtracking step shrink factor (`α ← α · factor`).
`eps_rel_jac`	(API only)	`1.0e-6`	FD perturbation size on the flux Jacobian (legacy; analytic paths bypass this).
`dt_halvings_max`	(API only)	`5`	Δt halvings before the step gives up and raises.

Picard outer loop (network only)¶

Flag	Default	Purpose
`picard_max_iter`	`3`	Cap on Picard outer iterations per step. `= 1` is block-Gauss-Seidel without iteration.
`picard_tol`	`1.0e-6`	Ghost-state convergence threshold.

Chord-Jacobian (Phase Perf-3)¶

Flag	Default	Purpose
`chord_max_age`	`2`	Hard cap on Newton iterations between Jacobian rebuilds. `= 0` reproduces the legacy full-Newton path bit-exactly.
`chord_rate_threshold`	`0.5`	If `F_norm_new / F_norm_old > threshold`, force a Jacobian refresh on the next iteration. Loosen (e.g. `0.7`) for stiffer cases where the chord ages faster than expected.

Parallel pipes (Phase 27c)¶

Flag	Default	Purpose
`parallel_pipes`	`False`	Run per-pipe Newton solves on threads. Default off because the Newton-path `@njit` kernels are not compiled `nogil` — the threads serialize on the GIL, so this path adds `ThreadPoolExecutor` coordination overhead without parallel execution (audit 2026-06-10 measured the serial path ~2× faster at typical pipe counts). Kept as an explicit opt-in for its failure-isolated first-pass semantics (on any pipe's Newton failure the state is restored and the loop falls back to the sequential Δt-halving path) and as the seam for a future `nogil` kernel build.
`parallel_pipes_auto_threshold`	`0`	Auto-enable the parallel pass when `parallel_pipes` is False but the network has at least this many pipes. `0` (default) disables the auto path. The previous default of `12` silently routed every ≥12-pipe implicit network onto the GIL-serialized threaded path above — pure overhead — so it was zeroed (audit 2026-06-10 batch B item 6).

Adaptive Δt (implicit schemes only)¶

Inert when rk is an explicit scheme; when False, Δt = cfl · dt_acoustic_unit per step.

Flag	Default	Purpose
`adaptive_dt`	`False`	Newton-iteration-driven Δt control. When on, the implicit driver scales Δt step-by-step inside `[adaptive_cfl_min, adaptive_cfl_max] · dt_acoustic_unit` toward `adaptive_target_iters` Newton iterations.
`adaptive_cfl_min`	`0.5`	Lower CFL bound for the adaptive Δt window.
`adaptive_cfl_max`	`50.0`	Upper CFL bound — implicit schemes unlock CFL ≫ 1.
`adaptive_target_iters`	`4`	Target Newton iteration count the controller steers toward (grows Δt below it, shrinks above).
`adaptive_error_control`	`False`	True local-error (step-doubling) Δt control: with `adaptive_dt` also on, Δt is driven by a Richardson one-step-vs-two-half-steps estimate held near `adaptive_atol + adaptive_rtol·\\|U\\|` instead of the Newton-difficulty heuristic. Only valid for one-step stencils (`implicit_euler` / `radau_iia`); raises on BDF2.
`adaptive_atol`	`1.0e-6`	Absolute tolerance for the step-doubling error estimate.
`adaptive_rtol`	`1.0e-4`	Relative tolerance for the step-doubling error estimate.

Adaptive Picard (network only)¶

Flag	Default	Purpose
`adaptive_picard`	`False`	Adaptive `picard_max_iter` controller. When on, the implicit-network driver grows the effective cap when the Picard budget saturates without converging and shrinks it after `adaptive_picard_n_easy` steps that exit at iter 1, bounded by `[adaptive_picard_min, adaptive_picard_max]`.
`adaptive_picard_min`	`1`	Lower bound on the adaptive Picard cap.
`adaptive_picard_max`	`8`	Upper bound on the adaptive Picard cap.
`adaptive_picard_n_easy`	`5`	Consecutive iter-1 exits before the cap shrinks.

Cost-shape comparison¶

Bench numbers from scripts/bench_long_pipe.py on a single 80 / 500-cell pipe with implicit BDF2 + adaptive Δt + porous-Forchheimer plug.

N	direct + chord
80	0.010 s
500	0.227 s

Per Newton-iter cost shape: the direct + chord path costs ~1 block-Thomas solve per Newton iter (with the Jacobian cached on most iters via chord aging — see Chord-Jacobian).

TOML examples¶

Default direct + chord (most users)¶

[case]
rk = "bdf2"
adaptive_dt = true
# everything else uses the defaults

Disable chord (debugging, or extra-tight convergence)¶

[case]
rk = "bdf2"
chord_max_age = 0              # full Newton on every iter

QS-driver knobs¶

The quasi-steady driver (solver = "quasi_steady") does not run the Newton path described above — it integrates only the slow solid states in time and resolves the gas-phase fields per macro step via a SIMPLE-style pressure/velocity march plus an outer junction Picard. The Newton / chord flags in the tables above are therefore inactive under QS; the knobs in this section take over. For the AT-targeted long-cycle wins (AVL Cruise-M, drive-cycle runs) these are the dials you want.

A quick defaults snapshot (post-2026-05-31 main):

Knob	Default	What flipped recently
`qs_picard_relax`	`0.7`	Secant-Newton damping; lower (`0.3`–`0.1`) on sharp inlet transients. MR !130 added the knob; MR !131 root-fixed the underlying stale-secant 2-cycle with a sign-bracket-bisection guard, so the knob is now a belt-and-braces lever rather than a required workaround.
`qs_picard_tol`	`1.0e-6`	Picard convergence threshold on the per-junction relative ṁ imbalance. Loosen to `1e-4` for bare-duct AT pipes whose α-noise floor lives higher.
`qs_picard_max_iter`	`100`	Cap on Picard outer iterations per macro step. Raise on larger networks where the FD-Newton's off-diagonal coupling needs more sweeps.
`qs_steady_energy`	`True`	Gas-phase energy as a steady advection-reaction balance (only slow solid states integrate in time). Removes the regen-ignition Δt collapse on long drive cycles. `False` falls back to the legacy explicit enthalpy advection.
`qs_chem_fused_rhs`	`True`	Flipped True 2026-05-31 (MR !134). @njit fused LSODA RHS for Brinkmeier-TWC surface chemistry. 2.08× speedup on AVL Cruise-M `MFR_transition_SC630` (97.7 s → 47.0 s on a 12 s sim; 8.15× → 3.92× realtime), bit-equivalent. Anchored by `tests/integration/test_qs_chem_fused_rhs.py`. Non-Brinkmeier closures (generic builder, NH3, `py_reaction`) auto-fall-back to the Python rhs — no per-case opt-out.
`chem_every_n`	`1`	Lie super-cycle for the surface-chemistry sub-step (the single canonical knob shared by every driver; the QS-specific `qs_chem_every_n` twin was collapsed into it in v0.9). `> 1` fires chemistry once every N macro steps over the accumulated Σ Δt; transport / steady-energy still run every step. Cuts chemistry call rate ~N× when the kinetic timescale ≫ macro Δt (warm catalyst, steady AT).
`qs_chem_settle_tol`	`0.0`	HRR-activity gate that suspends the super-cycle during transients. Skipping chemistry on non-Nth steps is safe only when the catalyst is warm/near-steady; during light-off the kinetics are faster than Δt and an unconditional super-cycle is catastrophic (measured: `chem_every_n=4` took 2.25M steps at light-off). Set `1e-2` (1 %/step) for a safe gate; only consulted when `chem_every_n > 1`.
`qs_chem_coverage_settle_tol`	`0.0`	Coverage-aware companion to `qs_chem_settle_tol`. On a hot catalyst the gas energy can be settled (HRR ≈ 0) while the oxygen-storage θ is still relaxing; the energy gate then super-cycles through it and degrades θ / NOx fidelity. Forces a chemistry step while the lagged max `\|Δθ\|` exceeds the tol. Fires on EITHER gate.
`chem_interval`	`0.0`	Wall-clock-time (simulated seconds) gate on the chemistry super-cycle, shared by every driver. Companion to `chem_every_n` — the step-counter gate fails when macro Δt grows (idle / cold-start ramp / adaptive HRR controller releasing) so `N · Δt` spans more than the chemistry timescale. The interval gate decouples that. Typical: `1e-3`–`1e-1` s on transient AT runs. Fires on EITHER.
`qs_cfl_adv`	`0.5`	Advective CFL ceiling on the QS macro step: `Δt ≤ qs_cfl_adv · min(dx / \|u\|)`. Bigger means fewer macro steps and proportionally fewer chemistry sub-step calls — a ~4× wall lever on chemistry-dominated cases. The species step is mass-conserving first-order upwind, stable up to CFL = 1; `1.0`–`2.0` typically works for reactive-AT cases. Validate species conservation when pushing past 1.0.

A few non-obvious interaction notes:

qs_picard_relax is not the primary mitigation any more — MR !131 fixed the stale-secant 2-cycle that motivated the knob with a sign-bracket-bisection Newton guard. Lower it only as a secondary lever on cases that still report a stuck relative imbalance after the bisection kicks in. The permanent diagnostic for that failure mode lives at scripts/diagnose_qs_picard_2cycle.py + scripts/analyze_qs_picard_trace.py.
chem_every_n without qs_chem_settle_tol / qs_chem_coverage_settle_tol / chem_interval is the legacy unconditional super-cycle — kept for backward compatibility but unsafe through light-off and θ-relaxation windows. Pair it with at least one settle/interval gate on any new QS case.
qs_chem_fused_rhs = True only matters on pipes whose reaction closure publishes py_reaction_surface._jit_params (Brinkmeier-TWC today). Other closures fall back to the Python rhs transparently; the flag is a no-op for them. Flip to False only if a regression surfaces against tests/integration/test_qs_chem_fused_rhs.py.

QS junction-pressure solve method¶

The decoupled per-junction scalar-secant Picard (the default) is exact for serial chains but diverges on parallel branches because it ignores the off-diagonal ∂net_j/∂p_k coupling. Two opt-in solvers handle looped / parallel-branch networks (e.g. the Cruise-T cas1 dual SCR-bank import); both are off by default so validated cases stay bit-exact.

Knob	Default	Purpose
`qs_junction_coupled_newton`	`False`	Solve the whole junction-pressure vector at once with a coupled finite-difference (GGA-style nodal) Newton. Handles parallel SCR banks and split/merge collectors that the decoupled Picard can't.
`qs_junction_loop_flow`	`False`	LOOP-FLOW / null-space solve (Hardy Cross; Abraham & Stoianov 2016): solves flows in the cycle space of the incidence matrix (well-conditioned SPD loop-Newton), then recovers junction pressures. Robust on looped networks with extreme resistance contrast where the nodal coupled-Newton stalls. Takes precedence over `qs_junction_coupled_newton`. EXPERIMENTAL.
`qs_flow_reg_velocity`	`0.0`	Low-flow head-loss regularization threshold `u_reg` [m/s]. A purely quadratic (Forchheimer) pipe has `dΔp/dṁ → 0` as ṁ → 0, so the junction Jacobian goes singular near zero flow; setting `u_reg > 0` adds a small linear drag below this velocity for a finite non-zero slope (negligible at physical flow). `0.0` is bit-identical.
`qs_flow_reg_alpha_floor`	`1.0`	Forchheimer-equivalent α floor [1/m] applied during regularization so genuinely zero-drag pipes also get a non-singular zero-flow slope.

QS gas-energy HRR coupling¶

Knob	Default	Purpose
`qs_surface_q_rxn_into_energy`	`True`	Include the lagged surface-chemistry reaction HRR (`pipe.q_rxn`) as a source in the steady-energy solve. Set False only as a diagnostic to A/B the HRR coupling.
`qs_surface_q_rxn_into_substrate`	`False`	Route the lagged `q_rxn` into the catalyst substrate `T_solid` (bounded by the substrate heat capacity, so the gas-side capacitance limiter isn't needed for this term) instead of folding it straight into the gas-energy source. Only applies to two-phase-catalyst pipes; a `q_rxn` pipe without one falls back to the gas overlay with a one-time warning.
`qs_dt_init`	`1.0e-4`	Bootstrap macro step [s] before the HRR Δt controller engages.
`qs_dt_max`	`0.0`	Cap on the macro step [s]; `0.0` = uncapped.

The full, exhaustive list of QS-driver flags (with per-knob doc comments) lives in src/aeolus1d/io/schema.py on CaseHeader — search for the qs_ prefix. See the v0.9 changelog for the recent default flips and the algorithmic story behind MR !131 (sign-bracket Picard guard) and MR !132 / !134 (fused Brinkmeier RHS + default-on).