From 16d51bd85d884a61c96df81fa07659053c8492a4 Mon Sep 17 00:00:00 2001 From: Daniel Weindl Date: Tue, 8 Sep 2026 14:58:56 +0200 Subject: [PATCH 1/2] Fix zero sensitivity for conservation laws with parameter-dependent initial values A boundary-condition or `constant=true` SBML species with no counteracting rate rule gets dx/dt=0 and is eliminated into a conservation law. `dtcldp` (the sensitivity of a conservation law's total abundance w.r.t. free parameters) was unconditionally hardcoded to zero whenever a conservation law spanned only one species, on the assumption that a single-species total abundance can't depend on a parameter. This breaks as soon as the species' initial value does depend on a parameter (e.g. via an InitialAssignment) -- the chain rule for any other state whose rate law references the eliminated species (via its `tcl` symbol) then silently drops that dependency. Whether a conservation law's total abundance depends on a free parameter can't be decided from the model's own symbolic initial conditions alone: preequilibration (steady-state Newton solve) and PEtab condition-table state reinitialization compute/override initial states and sensitivities at runtime, entirely bypassing the model's compiled `fx0`/`fsx0`. So a species' initial value can depend on a parameter even where the model's own static formula shows no such dependence -- for single- and multi-species conservation laws alike. `dtcldp` is therefore now always emitted as a symbol, never folded to a literal zero. The correct numeric value is always computed at runtime via `fdtotal_cldp`/`fdtotal_cldx_rdata`. Verified on SBML semantic test suite case 00783 (parameter-dependent conservation law, single species): sensitivities now match finite differences, where they were previously hardcoded to zero. Also verified against PEtab test suite cases 0010/0017 (multi-species conservation law under partial preequilibration and reinitialization), which an earlier version of this fix -- one that tried to statically detect the parameter-independent case instead of always emitting a symbol -- broke. https://github.com/AMICI-dev/AMICI/issues/3249 Co-Authored-By: Claude Sonnet 5 --- python/sdist/amici/_symbolic/de_model.py | 49 ++++++++++++------------ 1 file changed, 24 insertions(+), 25 deletions(-) diff --git a/python/sdist/amici/_symbolic/de_model.py b/python/sdist/amici/_symbolic/de_model.py index d9f6dacf5e..efcfa1a1f3 100644 --- a/python/sdist/amici/_symbolic/de_model.py +++ b/python/sdist/amici/_symbolic/de_model.py @@ -1166,22 +1166,7 @@ def _generate_symbol(self, name: str) -> None: ) return elif name == "dtcldp": - # check, whether the CL consists of only one state. Then, - # sensitivities drop out, otherwise generate symbols - self._syms[name] = sp.Matrix( - [ - [ - sp.Symbol( - f"s{tcl.get_id()}__{par.get_id()}", - real=True, - ) - for par in self._free_parameters - ] - if self.conservation_law_has_multispecies(tcl) - else [0] * self.num_par() - for tcl in self._conservation_laws - ] - ) + self._syms[name] = self._dtcldp_symbols() return elif name == "x_old": length = len(self.eq("xdot")) @@ -2440,20 +2425,34 @@ def state_is_constant(self, ix: int) -> bool: return state.get_dt().is_zero - def conservation_law_has_multispecies(self, tcl: ConservationLaw) -> bool: + def _dtcldp_symbols(self) -> sp.Matrix: """ - Checks whether a conservation law has multiple species or it just - defines one constant species + Builds the symbol matrix for ``dtcldp``, the sensitivity of each + conservation law's total abundance w.r.t. the free parameters. - :param tcl: - conservation law + This is always a symbol, never a literal zero: whether a + conservation law's total abundance actually depends on a free + parameter can't be decided from the model's own symbolic initial + conditions alone. Preequilibration (steady-state Newton solve) and + state reinitialization compute/override initial states and + sensitivities at runtime, entirely bypassing the model's compiled + ``fx0``/``fsx0``, so a species' initial value can depend on a + parameter even where the static formula shows no such dependence. + The correct numeric value is always computed at runtime via + ``fdtotal_cldp``/``fdtotal_cldx_rdata``. :return: - boolean indicating if conservation_law is not None + symbol matrix, one row per conservation law """ - state_set = set(self.sym("x_rdata")) - n_species = len(state_set.intersection(tcl.get_val().free_symbols)) - return n_species > 1 + return sp.Matrix( + [ + [ + symbol_with_assumptions(f"s{tcl.get_id()}__{par.get_id()}") + for par in self._free_parameters + ] + for tcl in self._conservation_laws + ] + ) def _expr_is_time_dependent(self, expr: sp.Expr) -> bool: """Determine whether an expression is time-dependent. From b21ce962f7ee73f705bf0df226d792a933f9d818 Mon Sep 17 00:00:00 2001 From: Daniel Weindl Date: Wed, 9 Sep 2026 17:28:21 +0200 Subject: [PATCH 2/2] unskip --- tests/sbml/testSBMLSuite.py | 9 --------- 1 file changed, 9 deletions(-) diff --git a/tests/sbml/testSBMLSuite.py b/tests/sbml/testSBMLSuite.py index ce0a41bf3b..c3c15d88ba 100755 --- a/tests/sbml/testSBMLSuite.py +++ b/tests/sbml/testSBMLSuite.py @@ -208,15 +208,6 @@ def _sensitivity_preflight_checks( "be wrong -- see " "https://github.com/AMICI-dev/AMICI/issues/3250" ) - if any( - species.getBoundaryCondition() or species.getConstant() - for species in sbml_model.getListOfSpecies() - ): - pytest.skip( - "Sensitivities for boundary-condition/constant species " - "are known to be wrong -- see " - "https://github.com/AMICI-dev/AMICI/issues/3249" - ) if uses_adjoint and model.nx_rdata == 0: pytest.skip( "Adjoint sensitivities for zero-state models are known to crash."