Steady-state d(func)/dp still finite-differences: route it through the compiled bngsim_codegen_output_sens

[wshlavacek]

Follow-up from #73, which fixed the correctness of this block but left it on finite differences.

compute_ss_output_sensitivity (src/steady_state.cpp:1284) builds the expression output sensitivity as

d(func_m)/dp = Σ_i (∂func_m/∂x_i)·dY_ss_i/dp  +  ∂func_m/∂p

with both terms finite-differenced: n_species perturbations for the state chain and one per parameter for the explicit term (src/steady_state.cpp:1357). The compiled bngsim_codegen_output_sens (#198) already computes exactly this quantity analytically, in the same .so as the RHS, and is what the CVODES run() path uses. The steady-state path should prefer it and keep FD as the fallback — the same shape ∂f/∂p already has in compute_ss_sensitivity since #63.

This is not a "resolve one more symbol" change. Measured below.

The blocker: the symbol is not in the artifact the solver gets

bngsim_codegen_output_sens is emitted only when the model carries _want_output_sens, which Simulator.__init__ sets from its constructor sensitivity_params (python/bngsim/_simulator.py). steady_state() takes its own sensitivity_params as a method argument, so the ordinary call never sets the flag:

construction	_want_output_sens	output_sens symbol in artifact
Simulator(m, method="ode") then .steady_state(sensitivity_params=[...])	False	absent
Simulator(m, method="ode", sensitivity_params=["kon"]) then .steady_state(...)	True	present

Measured on tests/data/ss_expr_sens_derived.net with BNGSIM_CODEGEN_THRESHOLD=1, checking the emitted C source / nm -gU on the .so, with ss.rhs_backend == "codegen-so" in both rows. So the first row — the documented usage in the steady_state() docstring and every test in test_steady_state_output_sensitivities.py — is precisely the case with no symbol to call.

Four pieces

1. Python re-prep. The #205 dance compute_all_sensitivities already implements (_simulator.py, the n_functions > 0 and not _want_output_sens branch): set the flag, drop the plain artifact, _auto_codegen_for_sensitivity, restore the previous artifact if regeneration yields nothing. steady_state() calls _auto_codegen_for_sensitivity today but without the flag dance, so it can inherit a plain-RHS .so that shadows the sensitivity one.
2. Symbol resolve in SteadyStateRhs — try_symbol<CodegenOutputSensFn>, alongside the existing sens_fn_ / jac_fn_, plus a has_analytical_output_sens() predicate.
3. Transpose. The ABI wants const double *const *state_sens, one contiguous n_species column per parameter. result.sensitivity is species-major (sensitivity[i*np + p]), so it needs materializing into np column buffers. plist[c] is the differentiated parameter index (>= n_params marks an IC column, not applicable here).
4. NaN sentinels. The codegen writes NaN for functions it marks unsupported and the Result rejects them at selection time. The steady-state path has no equivalent rejection, so it must either adopt one or fall back per-function to FD — the FD block stays either way, which is why this is a preferred path, not a replacement.

Oracle

Already in the tree. TestDerivedParameterChain::test_matches_long_cvode_run (added in #73) pins the steady-state answer against the last point of a converged CVODES forward-sensitivity run, which goes through bngsim_codegen_output_sens — i.e. against exactly the evaluator this issue wires in. ss_expr_sens_derived.net additionally has a full closed form (A ⇌ B with _rateLaw1 = chi*kon, flux() = _rateLaw1*A_tot), so agreement can be checked without an integrator. ss_output_sens.net and ss_birthdeath.net cover the reduced- and full-space solves.

Worth adding: a case where the codegen declines a function, to pin the fallback rather than discovering it as a NaN in a gradient.

Expected payoff

Accuracy first, speed second. The FD explicit term is a sqrt(eps) one-sided quotient — the same class of noise #63 removed from ∂f/∂p and #71 removed from the sympy path. The state-chain term additionally costs n_species full observable+function re-evaluations, which on a large model with many functions is the dominant cost of the output-sensitivity block; the compiled evaluator does it in one call per output row.

Not yet measured on the corpus — worth doing before committing to a size, following the #66 precedent where the budgeted work turned out to be mostly already implemented.

Rationale already in-tree

src/steady_state.cpp:1267 carries a comment block with this analysis, so the reasoning does not need re-deriving from scratch when someone picks this up.