From a55d71083e48e95c837ce38dce0e3cf5adcc2f5d Mon Sep 17 00:00:00 2001 From: Daniel Weindl Date: Wed, 1 Oct 2025 10:31:32 +0200 Subject: [PATCH] `ReturnDataView.{preeq,posteq}_status` to `list[SteadyStateStatus]`; squeeze numsteps * `ReturnDataView.posteq_numsteps` and `ReturnDataView.posteq_numsteps` now return a one-dimensional array of shape `(num_timepoints,)` instead of a two-dimensional array of shape `(1, num_timepoints)`. * `ReturnDataView.posteq_status` and `ReturnDataView.preeq_status` now return `list[SteadyStateStatus]` instead of an `ndarray[int]` of shape `(1, 3)`. Closes #1944. --- CHANGELOG.md | 7 +++- doc/examples/example_errors.ipynb | 15 +++----- .../ExampleEquilibrationLogic.ipynb | 34 +++++++++---------- python/sdist/amici/numpy.py | 20 +++++++---- .../test_compare_conservation_laws_sbml.py | 15 ++++++-- python/tests/test_preequilibration.py | 30 +++++++++------- 6 files changed, 69 insertions(+), 52 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 38807fd6f8..1938928aad 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,7 +9,12 @@ See also our [versioning policy](https://amici.readthedocs.io/en/latest/versioni BREAKING CHANGES * The MATLAB interface has been removed. - +* `ReturnDataView.posteq_numsteps` and `ReturnDataView.posteq_numsteps` now + return a one-dimensional array of shape `(num_timepoints,)` instead of a + two-dimensional array of shape `(1, num_timepoints)`. +* `ReturnDataView.posteq_status` and `ReturnDataView.preeq_status` now + return `list[SteadyStateStatus]` instead of an `ndarray[int]` of shape + `(1, 3)`. ## v0.X Series diff --git a/doc/examples/example_errors.ipynb b/doc/examples/example_errors.ipynb index 87da94c912..dbc00e9be6 100644 --- a/doc/examples/example_errors.ipynb +++ b/doc/examples/example_errors.ipynb @@ -32,7 +32,6 @@ " SensitivityMethod,\n", " SensitivityOrder,\n", " runAmiciSimulation,\n", - " SteadyStateStatus,\n", ")\n", "from amici.petab.petab_import import import_petab_problem\n", "from amici.petab.simulations import simulate_petab, RDATAS, EDATAS\n", @@ -921,7 +920,7 @@ "outputs": [], "source": [ "rdata = res[RDATAS][0]\n", - "list(map(SteadyStateStatus, rdata.preeq_status.flatten()))" + "rdata.preeq_status" ] }, { @@ -968,9 +967,7 @@ ")\n", "\n", "rdata = res[RDATAS][0]\n", - "print(\n", - " f\"preeq_status={list(map(SteadyStateStatus, rdata.preeq_status.flatten()))}\"\n", - ")\n", + "print(f\"preeq_status={rdata.preeq_status}\")\n", "print(f\"{rdata.preeq_numsteps=}\")\n", "\n", "# hard to reproduce on GHA\n", @@ -1015,9 +1012,7 @@ ")\n", "\n", "rdata = res[RDATAS][0]\n", - "print(\n", - " f\"preeq_status={list(map(SteadyStateStatus, rdata.preeq_status.flatten()))}\"\n", - ")\n", + "print(f\"preeq_status={rdata.preeq_status}\")\n", "print(f\"{rdata.preeq_numsteps=}\")\n", "assert all(rdata.status == amici.AMICI_SUCCESS for rdata in res[RDATAS])" ] @@ -1056,9 +1051,7 @@ ")\n", "\n", "rdata = res[RDATAS][0]\n", - "print(\n", - " f\"preeq_status={list(map(SteadyStateStatus, rdata.preeq_status.flatten()))}\"\n", - ")\n", + "print(f\"preeq_status={rdata.preeq_status}\")\n", "print(f\"{rdata.preeq_numsteps=}\")\n", "# hard to reproduce on GHA\n", "if os.getenv(\"GITHUB_ACTIONS\") is None:\n", diff --git a/doc/examples/example_steady_states/ExampleEquilibrationLogic.ipynb b/doc/examples/example_steady_states/ExampleEquilibrationLogic.ipynb index f50207fc51..8d104ccf8b 100644 --- a/doc/examples/example_steady_states/ExampleEquilibrationLogic.ipynb +++ b/doc/examples/example_steady_states/ExampleEquilibrationLogic.ipynb @@ -332,23 +332,21 @@ "source": [ "The `x` value contains state values corresponding to the time points that were set for the model simulation. In this case, it is the state values at $t=\\infty$, i.e., the steady state. One can also see that `xdot` values are equal to zero or very small, which means that the steady-state condition ($\\dot{x} \\approx 0$) is fulfilled. The actual weighted root-mean-square of the right-hand side when steady state was reached is reported in `posteq_wrms`.\n", "\n", - "The fields `posteq_status` and `posteq_numsteps` in `rdata` tell us how post-equilibration worked. In both arrays each entry is, respectively, the status or the number of steps taken for\n", + "The fields `posteq_status` and `posteq_numsteps` in `rdata` tell us how post-equilibration worked. In both lists each entry is, respectively, the status or the number of steps taken for\n", "\n", " 1. the Newton's method run;\n", " 2. the numerical integration run;\n", " 3. the second Newton's method run, starting from the simulation results.\n", "\n", - "The status is encoded as an integer flag with the following meanings:\n", + "The status is of type `SteadyStateStatus` with the following meanings:\n", "\n", - " * ` 1`: Successful run\n", - " * ` 0`: Did not run\n", - " * `-1`: Error: No further specification is given, the error message should give more information.\n", - " * `-2`: Error: The method did not converge to a steady state within the maximum number of steps (Newton's method or simulation).\n", - " * `-3`: Error: The Jacobian of the right hand side is singular (only Newton's method)\n", - " * `-4`: Error: The damping factor in Newton's method was reduced until it met the lower bound without success (Newton's method only)\n", - " * `-5`: Error: The model was simulated past the timepoint `t=1e100` without finding a steady state. Therefore, it is likely that the model has no steady state for the given parameter vector.\n", - "\n", - " Don't check the integer values directly, but use the `SteadyStateStatus` enum instead, which is more readable and stable:\n" + " * `SteadyStateStatus.success`: Successful run\n", + " * `SteadyStateStatus.not_run`: Did not run\n", + " * `SteadyStateStatus.failed`: Error: No further specification is given, the error message should give more information.\n", + " * `SteadyStateStatus.failed_convergence`: Error: The method did not converge to a steady state within the maximum number of steps (Newton's method or simulation).\n", + " * `SteadyStateStatus.failed_factorization`: Error: The Jacobian of the right hand side is singular (only Newton's method)\n", + " * `SteadyStateStatus.failed_damping`: Error: The damping factor in Newton's method was reduced until it met the lower bound without success (Newton's method only)\n", + " * `SteadyStateStatus.failed_too_long_simulation`: Error: The model was simulated past the timepoint `t=1e100` without finding a steady state. Therefore, it is likely that the model has no steady state for the given parameter vector.\n" ] }, { @@ -440,7 +438,7 @@ "solver_reduced.setMaxSteps(100)\n", "rdata_reduced = runAmiciSimulation(model_reduced, solver_reduced)\n", "assert rdata_reduced.status == AMICI_SUCCESS\n", - "assert rdata_reduced.posteq_status[0][0] == SteadyStateStatus.success\n", + "assert rdata_reduced.posteq_status[0] == SteadyStateStatus.success\n", "print(\n", " \"Simulation status:\",\n", " simulation_status_to_str(rdata_reduced[\"status\"]),\n", @@ -673,7 +671,7 @@ "solver_reduced.setMaxSteps(1000)\n", "rdata_reduced = runAmiciSimulation(model_reduced, solver_reduced)\n", "assert rdata_reduced.status == AMICI_SUCCESS\n", - "assert list(rdata_reduced.posteq_status[0]) == [\n", + "assert rdata_reduced.posteq_status == [\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", " SteadyStateStatus.not_run,\n", @@ -903,7 +901,7 @@ "solver.setSensitivityOrder(SensitivityOrder.first)\n", "rdata = runAmiciSimulation(model, solver, edata)\n", "assert rdata.status == AMICI_SUCCESS\n", - "assert list(rdata.preeq_status[0]) == [\n", + "assert rdata.preeq_status == [\n", " SteadyStateStatus.failed_factorization,\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", @@ -970,7 +968,7 @@ "rdata_reduced = runAmiciSimulation(model_reduced, solver_reduced, edata)\n", "\n", "assert rdata_reduced.status == AMICI_SUCCESS\n", - "assert list(rdata_reduced.preeq_status[0]) == [\n", + "assert rdata_reduced.preeq_status == [\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", " SteadyStateStatus.not_run,\n", @@ -995,7 +993,7 @@ "rdata_reduced = runAmiciSimulation(model_reduced, solver_reduced, edata)\n", "\n", "assert rdata_reduced.status == AMICI_SUCCESS\n", - "assert list(rdata_reduced.preeq_status[0]) == [\n", + "assert rdata_reduced.preeq_status == [\n", " SteadyStateStatus.not_run,\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", @@ -1022,7 +1020,7 @@ "solver_reduced.setSensitivityOrder(SensitivityOrder.first)\n", "rdata_reduced = runAmiciSimulation(model_reduced, solver_reduced, edata)\n", "assert rdata_reduced.status == AMICI_SUCCESS\n", - "assert list(rdata_reduced.preeq_status[0]) == [\n", + "assert rdata_reduced.preeq_status == [\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", " SteadyStateStatus.not_run,\n", @@ -1054,7 +1052,7 @@ "solver.setSensitivityOrder(SensitivityOrder.first)\n", "rdata = runAmiciSimulation(model, solver, edata)\n", "assert rdata.status == AMICI_SUCCESS\n", - "assert list(rdata.preeq_status[0]) == [\n", + "assert rdata.preeq_status == [\n", " SteadyStateStatus.failed_factorization,\n", " SteadyStateStatus.success,\n", " SteadyStateStatus.not_run,\n", diff --git a/python/sdist/amici/numpy.py b/python/sdist/amici/numpy.py index a90b32c769..f7b85d8b96 100644 --- a/python/sdist/amici/numpy.py +++ b/python/sdist/amici/numpy.py @@ -14,7 +14,14 @@ import numpy as np import sympy as sp from sympy.abc import _clash -from . import ExpData, ExpDataPtr, Model, ReturnData, ReturnDataPtr +from . import ( + ExpData, + ExpDataPtr, + Model, + ReturnData, + ReturnDataPtr, + SteadyStateStatus, +) StrOrExpr = str | sp.Expr @@ -290,11 +297,9 @@ def __init__(self, rdata: ReturnDataPtr | ReturnData): "w": [rdata.nt, rdata.nw], "xdot": [rdata.nx_solver], "preeq_numlinsteps": [rdata.newton_maxsteps, 2], - "preeq_numsteps": [1, 3], - "preeq_status": [1, 3], + "preeq_numsteps": [3], "posteq_numlinsteps": [rdata.newton_maxsteps, 2], - "posteq_numsteps": [1, 3], - "posteq_status": [1, 3], + "posteq_numsteps": [3], "numsteps": [rdata.nt], "numrhsevals": [rdata.nt], "numerrtestfails": [rdata.nt], @@ -309,7 +314,7 @@ def __init__(self, rdata: ReturnDataPtr | ReturnData): def __getitem__( self, item: str - ) -> np.ndarray | ReturnDataPtr | ReturnData | float: + ) -> np.ndarray | ReturnDataPtr | ReturnData | float | int | list: """ Access fields by name. @@ -322,6 +327,9 @@ def __getitem__( if item == "status": return int(super().__getitem__(item)) + if item in ("preeq_status", "posteq_status"): + return list(map(SteadyStateStatus, super().__getitem__(item))) + if item == "t": item = "ts" diff --git a/python/tests/test_compare_conservation_laws_sbml.py b/python/tests/test_compare_conservation_laws_sbml.py index 2956213f16..44701f70d5 100644 --- a/python/tests/test_compare_conservation_laws_sbml.py +++ b/python/tests/test_compare_conservation_laws_sbml.py @@ -4,7 +4,8 @@ import amici import numpy as np import pytest -from numpy.testing import assert_allclose, assert_array_equal +from amici import SteadyStateStatus +from numpy.testing import assert_allclose from amici.testing import skip_on_valgrind @@ -229,9 +230,17 @@ def test_compare_conservation_laws_sbml(models, edata_fixture): ) assert rdata["status"] == amici.AMICI_SUCCESS # check that steady state computation succeeded only by sim in full model - assert_array_equal(rdata["preeq_status"], np.array([[-3, 1, 0]])) + assert rdata["preeq_status"] == [ + SteadyStateStatus.failed_factorization, + SteadyStateStatus.success, + SteadyStateStatus.not_run, + ] # check that steady state computation succeeded by Newton in reduced model - assert_array_equal(rdata_cl["preeq_status"], np.array([[1, 0, 0]])) + assert rdata_cl["preeq_status"] == [ + SteadyStateStatus.success, + SteadyStateStatus.not_run, + SteadyStateStatus.not_run, + ] # compare state sensitivities with edata and preequilibration for field in ["x", "x_ss", "sx", "llh", "sllh"]: diff --git a/python/tests/test_preequilibration.py b/python/tests/test_preequilibration.py index 49422fbe0a..660acb2424 100644 --- a/python/tests/test_preequilibration.py +++ b/python/tests/test_preequilibration.py @@ -13,7 +13,7 @@ skip_on_valgrind, ) from amici.gradient_check import check_derivatives -from amici import SensitivityMethod, SensitivityOrder +from amici import SensitivityMethod, SensitivityOrder, SteadyStateStatus pytestmark = pytest.mark.filterwarnings( # https://github.com/AMICI-dev/AMICI/issues/18 @@ -517,28 +517,32 @@ def test_steadystate_computation_mode(preeq_fixture): # assert successful simulation assert rdatas[mode]["status"] == amici.AMICI_SUCCESS + assert rdatas[amici.SteadyStateComputationMode.integrationOnly][ + "preeq_status" + ] == [ + SteadyStateStatus.not_run, + SteadyStateStatus.success, + SteadyStateStatus.not_run, + ] - assert np.all( - rdatas[amici.SteadyStateComputationMode.integrationOnly][ - "preeq_status" - ][0] - == [0, 1, 0] - ) assert ( rdatas[amici.SteadyStateComputationMode.integrationOnly][ "preeq_numsteps" - ][0][0] + ][0] == 0 ) - assert np.all( - rdatas[amici.SteadyStateComputationMode.newtonOnly]["preeq_status"][0] - == [1, 0, 0] - ) + assert rdatas[amici.SteadyStateComputationMode.newtonOnly][ + "preeq_status" + ] == [ + SteadyStateStatus.success, + SteadyStateStatus.not_run, + SteadyStateStatus.not_run, + ] assert ( rdatas[amici.SteadyStateComputationMode.newtonOnly]["preeq_numsteps"][ 0 - ][0] + ] > 0 )