diff --git a/CHANGELOG.md b/CHANGELOG.md
index cba09d146c..ce826758df 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -4,6 +4,27 @@ See also our [versioning policy](https://amici.readthedocs.io/en/latest/versioni
## v0.X Series
+### v0.30.1 (2025-02-18)
+
+Bugfix-only release.
+
+* Removed `eqx.debug.nan`, fixes #2629
+ by @FFroehlich in https://github.com/AMICI-dev/AMICI/pull/2630
+* Fixes an SBML import issue that led to incorrect results for models with
+ species-dependent initial assignments (fixes #2642)
+ by @FFroehlich in https://github.com/AMICI-dev/AMICI/pull/2643
+* Fixed `CVodeGetSensDky` error message
+ by @dweindl in https://github.com/AMICI-dev/AMICI/pull/2644
+* Disabled `CvodeF` checkpointing to prevent certain rare crashes when forward
+ integration takes exactly `maxsteps` integration steps, plus some additional
+ yet unclear condition.
+ by @dweindl in https://github.com/AMICI-dev/AMICI/pull/2645
+* Fixed rare crashes due to uncaught exceptions in `~FinalStateStorer`
+ by @dweindl in https://github.com/AMICI-dev/AMICI/pull/2647
+
+**Full Changelog**: https://github.com/AMICI-dev/AMICI/compare/v0.30.0...v0.30.1
+
+
### v0.30.0 (2024-12-10)
*Please note that the amici JAX model generation introduced in v0.29.0 is
diff --git a/include/amici/forwardproblem.h b/include/amici/forwardproblem.h
index 40a4147f4e..cad434a3e1 100644
--- a/include/amici/forwardproblem.h
+++ b/include/amici/forwardproblem.h
@@ -441,19 +441,37 @@ class FinalStateStorer : public ContextManager {
/**
* @brief destructor, stores simulation state
*/
- ~FinalStateStorer() {
+ ~FinalStateStorer() noexcept(false) {
if (fwd_) {
- fwd_->final_state_ = fwd_->getSimulationState();
- // if there is an associated output timepoint, also store it in
- // timepoint_states if it's not present there.
- // this may happen if there is an error just at
- // (or indistinguishably before) an output timepoint
- auto final_time = fwd_->getFinalTime();
- auto const timepoints = fwd_->model->getTimepoints();
- if (!fwd_->timepoint_states_.count(final_time)
- && std::find(timepoints.cbegin(), timepoints.cend(), final_time)
- != timepoints.cend()) {
- fwd_->timepoint_states_[final_time] = fwd_->final_state_;
+ try {
+ // This may throw in `CVodeSolver::getSens`
+ // due to https://github.com/LLNL/sundials/issues/82.
+ // Therefore, this dtor must be `noexcept(false)` to avoid
+ // programm termination.
+ fwd_->final_state_ = fwd_->getSimulationState();
+ // if there is an associated output timepoint, also store it in
+ // timepoint_states if it's not present there.
+ // this may happen if there is an error just at
+ // (or indistinguishably before) an output timepoint
+ auto final_time = fwd_->getFinalTime();
+ auto const timepoints = fwd_->model->getTimepoints();
+ if (!fwd_->timepoint_states_.count(final_time)
+ && std::find(timepoints.cbegin(), timepoints.cend(), final_time)
+ != timepoints.cend()) {
+ fwd_->timepoint_states_[final_time] = fwd_->final_state_;
+ }
+ } catch (std::exception const&) {
+ // We must not throw in case we are already in the stack
+ // unwinding phase due to some other active exception, otherwise
+ // this will also lead to termination.
+ //
+ // In case there is another active exception,
+ // `fwd_->{final_state_,timepoint_states_}` won't be set,
+ // and we assume that they are either not accessed anymore, or
+ // that there is appropriate error handling in place.
+ if(!std::uncaught_exceptions()) {
+ throw;
+ }
}
}
}
diff --git a/python/sdist/amici/jax/petab.py b/python/sdist/amici/jax/petab.py
index b5834223fb..b2b42f5c2a 100644
--- a/python/sdist/amici/jax/petab.py
+++ b/python/sdist/amici/jax/petab.py
@@ -500,7 +500,7 @@ def run_simulation(
simulation_condition[0], p
)
return self.model.simulate_condition(
- p=eqx.debug.backward_nan(p),
+ p=p,
ts_init=jax.lax.stop_gradient(jnp.array(ts_preeq)),
ts_dyn=jax.lax.stop_gradient(jnp.array(ts_dyn)),
ts_posteq=jax.lax.stop_gradient(jnp.array(ts_posteq)),
diff --git a/python/sdist/amici/sbml_import.py b/python/sdist/amici/sbml_import.py
index 557ad02d0f..29a9608c4b 100644
--- a/python/sdist/amici/sbml_import.py
+++ b/python/sdist/amici/sbml_import.py
@@ -2292,9 +2292,29 @@ def _make_initial(
sym_math, rateof_to_dummy = _rateof_to_dummy(sym_math)
- for species_id, species in self.symbols[SymbolId.SPECIES].items():
- if "init" in species:
- sym_math = smart_subs(sym_math, species_id, species["init"])
+ # we can't rely on anything else being properly initialized at this point, so we need to
+ # compute all initial values from scratch, recursively
+ for var in sym_math.free_symbols:
+ element_id = str(var)
+ # already recursive since _get_element_initial_assignment calls _make_initial
+ if (
+ ia := self._get_element_initial_assignment(element_id)
+ ) is not None:
+ sym_math = sym_math.subs(var, ia)
+ elif (species := self.sbml.getSpecies(element_id)) is not None:
+ # recursive!
+ init = self._make_initial(get_species_initial(species))
+ sym_math = sym_math.subs(var, init)
+ elif var in self.symbols[SymbolId.SPECIES]:
+ sym_math = sym_math.subs(
+ var, self.symbols[SymbolId.SPECIES][var]["init"]
+ )
+ elif (
+ element := self.sbml.getElementBySId(element_id)
+ ) and self.is_rate_rule_target(element):
+ # no need to recurse here, as value is numeric
+ init = sp.Float(element.getValue())
+ sym_math = sym_math.subs(var, init)
sym_math = smart_subs(sym_math, sbml_time_symbol, sp.Float(0))
diff --git a/python/tests/sbml_models/regression_2642.xml b/python/tests/sbml_models/regression_2642.xml
new file mode 100644
index 0000000000..63dd9602d4
--- /dev/null
+++ b/python/tests/sbml_models/regression_2642.xml
@@ -0,0 +1,130 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/python/tests/test_sbml_import.py b/python/tests/test_sbml_import.py
index 87b18da92d..9047b76b4a 100644
--- a/python/tests/test_sbml_import.py
+++ b/python/tests/test_sbml_import.py
@@ -852,3 +852,26 @@ def test_import_same_model_name():
assert model_module_1c.get_model().getParameters()[0] == 1.0
assert model_module_1c.get_model().module is model_module_1c
+
+
+@skip_on_valgrind
+def test_regression_2642():
+ sbml_file = Path(__file__).parent / "sbml_models" / "regression_2642.xml"
+ sbml_importer = amici.SbmlImporter(sbml_file)
+ model_name = "regression_2642"
+ with TemporaryDirectory(prefix="regression_2642") as outdir:
+ sbml_importer.sbml2amici(
+ model_name=model_name,
+ output_dir=outdir,
+ )
+ module = amici.import_model_module(
+ module_name=model_name, module_path=outdir
+ )
+ model = module.getModel()
+ solver = model.getSolver()
+ model.setTimepoints(np.linspace(0, 1, 3))
+ r = amici.runAmiciSimulation(model, solver)
+ assert (
+ len(np.unique(r.w[:, model.getExpressionIds().index("binding")]))
+ == 1
+ )
diff --git a/python/tests/valgrind-python.supp b/python/tests/valgrind-python.supp
index eea8347d27..6b66869a3f 100644
--- a/python/tests/valgrind-python.supp
+++ b/python/tests/valgrind-python.supp
@@ -1006,3 +1006,33 @@
fun:_ZN4rpds*
...
}
+
+{
+ Python
+ Memcheck:Leak
+ match-leak-kinds: definite
+ fun:realloc
+ fun:_PyUnicodeWriter_Finish
+ obj:/usr/bin/python3.*
+ ...
+}
+
+{
+ Something matplotlib. Cannot reproduce.
+ Memcheck:Leak
+ match-leak-kinds: definite
+ fun:_Znwm
+ fun:_ZL14PyFT2Font_initN8pybind116objectElSt8optionalISt6vectorIP9PyFT2FontSaIS4_EEEi.lto_priv.0
+ fun:_ZZN8pybind1112cpp_function10initializeIZNOS_6detail8initimpl7factoryIPFP9PyFT2FontNS_6objectElSt8optionalISt6vectorIS6_SaIS6_EEEiEPFNS2_9void_typeEvESD_SG_E7executeINS_6class_IS5_JEEEJNS_3argENS_5arg_vENS_7kw_onlyESN_SN_PKcEEEvRT_DpRKT0_EUlRNS2_16value_and_holderES7_lSC_iE_vJSY_S7_lSC_iEJNS_4nameENS_9is_methodENS_7siblingENS2_24is_new_style_constructorESM_SN_SO_SN_SN_SQ_EEEvOSR_PFT0_DpT1_EDpRKT2_ENUlRNS2_13function_callEE1_4_FUNES1F_
+ fun:_ZN8pybind1112cpp_function10dispatcherEP7_objectS2_S2_
+ fun:cfunction_call
+ fun:_PyObject_MakeTpCall
+ fun:_PyObject_VectorcallTstate
+ fun:_PyObject_VectorcallTstate
+ fun:method_vectorcall
+ fun:slot_tp_init
+ fun:type_call
+ fun:pybind11_meta_call
+ fun:_PyObject_MakeTpCall
+ fun:_PyEval_EvalFrameDefault
+}
diff --git a/scripts/run-valgrind-py.sh b/scripts/run-valgrind-py.sh
index 7ca73c427c..de140de686 100755
--- a/scripts/run-valgrind-py.sh
+++ b/scripts/run-valgrind-py.sh
@@ -23,7 +23,7 @@ if [ $# -eq 0 ]
# for does not match any known type: falling back to type probe function.
else
# Run whatever was passed as arguments
- command=($@)
+ command=("$@")
fi
diff --git a/src/solver_cvodes.cpp b/src/solver_cvodes.cpp
index efe50eb8a9..8d34fedfb4 100644
--- a/src/solver_cvodes.cpp
+++ b/src/solver_cvodes.cpp
@@ -667,7 +667,7 @@ void CVodeSolver::getSensDky(realtype const t, int const k) const {
int status
= CVodeGetSensDky(solver_memory_.get(), t, k, sx_.getNVectorArray());
if (status != CV_SUCCESS)
- throw CvodeException(status, "CVodeGetSens");
+ throw CvodeException(status, "CVodeGetSensDky");
}
void CVodeSolver::getDkyB(realtype const t, int const k, int const which)
@@ -716,7 +716,7 @@ void CVodeSolver::adjInit() const {
status = CVodeAdjReInit(solver_memory_.get());
} else {
status = CVodeAdjInit(
- solver_memory_.get(), static_cast(maxsteps_),
+ solver_memory_.get(), static_cast(maxsteps_ + 1),
static_cast(interp_type_)
);
setAdjInitDone();
diff --git a/version.txt b/version.txt
index c25c8e5b74..1a44cad74d 100644
--- a/version.txt
+++ b/version.txt
@@ -1 +1 @@
-0.30.0
+0.30.1