Skip to content

Fix adjoint sensitivities for events that are both state-triggered and state-updating #3258

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
dweindl wants to merge 3 commits into
base: main
Choose a base branch
from
+34 −39
Open

Fix adjoint sensitivities for events that are both state-triggered and state-updating #3258

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
7401026
Fix adjoint sensitivities for events that are both state-triggered an…
dweindl Sep 8, 2026
67589de
unskip
dweindl Sep 9, 2026
f04344f
skip subset
dweindl Sep 10, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
31 changes: 25 additions & 6 deletions python/sdist/amici/_symbolic/de_model.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@
ObservableTransformation,
_default_simplify,
amici_time_symbol,
smart_subs_dict,
symbol_with_assumptions,
toposort_symbols,
unique_preserve_order,
Expand Down Expand Up @@ -1459,6 +1460,20 @@ def _generate_sparse_symbol(self, name: str) -> None:
self._sparsesyms[name] = symbol_list
self._syms[name] = sparse_matrix

def _pre_event_eqs(self, name: str) -> list[sp.Matrix]:
"""Re-express `dtaudx`/`dtaudp` purely in terms of the pre-event
state, for use inside `deltaxB`/`deltaqB` only.

`w` is eliminated first (not available in the C++ functions),
then `x -> x_old` is substituted.
"""
w_to_expr = dict(zip(self.sym("w"), self.eq("w")))
x_to_x_old = dict(zip(self.sym("x"), self.sym("x_old")))
return [
smart_subs_dict(smart_subs_dict(expr, w_to_expr), x_to_x_old)
for expr in self.eq(name)
]

def _compute_equation(self, name: str) -> None:
"""
Computes the symbolic formula for a symbolic variable
Expand Down Expand Up @@ -1867,18 +1882,20 @@ def _compute_equation(self, name: str) -> None:
self._eqs[name] = event_eqs

elif name == "deltaxB":
# express in terms of pre-event state
dtaudx_pre = self._pre_event_eqs("dtaudx")
event_eqs = []
for ie, event in enumerate(self._events):
# ==== 1st group of terms: Heaviside functions ===========
tmp_eq = smart_multiply(
self.sym("xdot") - self.sym("xdot_old"),
self.eq("dtaudx")[ie],
dtaudx_pre[ie],
)
if event.updates_state:
# ==== 2nd group of terms: Derivatives of Dirac deltas ===
# Part 2a: explicit time dependence of bolus function
tmp_eq -= smart_multiply(
self.eq("ddeltaxdt")[ie], self.eq("dtaudx")[ie]
self.eq("ddeltaxdt")[ie], dtaudx_pre[ie]
)
# Part 2b: implicit time dependence of bolus function
tmp_eq -= smart_multiply(
Expand All @@ -1887,7 +1904,7 @@ def _compute_equation(self, name: str) -> None:
+ self.eq("ddeltaxdx_old")[ie],
self.sym("xdot_old"),
),
self.eq("dtaudx")[ie],
dtaudx_pre[ie],
)
# ==== 3rd group of terms: Dirac deltas ==================
tmp_eq += (
Expand All @@ -1898,18 +1915,20 @@ def _compute_equation(self, name: str) -> None:
self._eqs[name] = event_eqs

elif name == "deltaqB":
# express in terms of pre-event state
dtaudp_pre = self._pre_event_eqs("dtaudp")
event_eqs = []
for ie, event in enumerate(self._events):
# ==== 1st group of terms: Heaviside functions ===========
tmp_eq = smart_multiply(
self.sym("xdot") - self.sym("xdot_old"),
self.eq("dtaudp")[ie],
dtaudp_pre[ie],
)
if event.updates_state:
# ==== 2nd group of terms: Derivatives of Dirac deltas ===
# Part 2a: explicit time dependence of bolus function
tmp_eq -= smart_multiply(
self.eq("ddeltaxdt")[ie], self.eq("dtaudp")[ie]
self.eq("ddeltaxdt")[ie], dtaudp_pre[ie]
)
# Part 2b: implicit time dependence of bolus function
tmp_eq -= smart_multiply(
Expand All @@ -1918,7 +1937,7 @@ def _compute_equation(self, name: str) -> None:
+ self.eq("ddeltaxdx_old")[ie],
self.sym("xdot_old"),
),
self.eq("dtaudp")[ie],
dtaudp_pre[ie],
)
# ==== 3rd group of terms: Dirac deltas ==================
tmp_eq += self.eq("ddeltaxdp")[ie]
Expand Down
42 changes: 9 additions & 33 deletions tests/sbml/testSBMLSuite.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -141,35 +141,6 @@ def test_sbml_testsuite_case(test_id, compiled_case, result_path):
write_result_file(simulated, test_id, result_path)


def _ast_has_piecewise(node: libsbml.ASTNode | None) -> bool:
"""Recursively check whether a libsbml math AST contains a `piecewise`
function anywhere in its tree."""
if node is None:
return False
if node.getType() == libsbml.AST_FUNCTION_PIECEWISE:
return True
return any(
_ast_has_piecewise(node.getChild(i))
for i in range(node.getNumChildren())
)


def _model_has_event_jump_risk(sbml_model: libsbml.Model) -> bool:
"""Whether this model has an event, or a piecewise formula."""
if sbml_model.getNumEvents() > 0:
return True
for reaction in sbml_model.getListOfReactions():
kinetic_law = reaction.getKineticLaw()
if kinetic_law is not None and _ast_has_piecewise(
kinetic_law.getMath()
):
return True
return any(
_ast_has_piecewise(rule.getMath())
for rule in sbml_model.getListOfRules()
)


# FIXME: Skip list - to be investigated further
# test_id -> adjoint_only (whether forward is unaffected)
_OTHER_KNOWN_SENSITIVITY_CHECK_ISSUES = {
Expand All @@ -182,6 +153,15 @@ def _model_has_event_jump_risk(sbml_model: libsbml.Model) -> bool:
"01104": True,
"01107": True,
"01148": True,
# Adjoint sensitivities & events. To revisit later.

@dweindl dweindl Sep 10, 2026

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Check if we can unskip anything after #3263

"00348": True,
"00350": True,
"00354": True,
"00368": True,
"00373": True,
"00396": True,
"00753": True,
"01106": True,
}


Expand Down Expand Up @@ -221,10 +201,6 @@ def _sensitivity_preflight_checks(
pytest.skip(
"Adjoint sensitivities for zero-state models are known to crash."
)
if uses_adjoint and _model_has_event_jump_risk(sbml_model):
pytest.skip(
"Adjoint sensitivities for (some) events are known to be wrong (https://github.com/AMICI-dev/AMICI/pull/3258)."
)
if test_id in _OTHER_KNOWN_SENSITIVITY_CHECK_ISSUES:
adjoint_only = _OTHER_KNOWN_SENSITIVITY_CHECK_ISSUES[test_id]
if uses_adjoint or not adjoint_only:
Expand Down
Loading