From f6a7e798e0e5300fa0242bf5171067402a56b753 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Sun, 5 Jul 2026 11:41:14 +0100 Subject: [PATCH 01/19] #42 Support no-ODE models via placeholder state variable A model with no reactions, ODEs or rate rules has no state variables of its own, but Chaste's ODE system integrates a vector of state variables. Synthesise a single placeholder state variable with dy/dt = 0 so the generated ODE system always has something for the solver to integrate; the model's real outputs (constants, time-dependent assignment rules, event-driven changes) are recomputed each step as before. Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/_config.py | 4 ++ chaste_sbml/_model_builder.py | 26 +++++++-- chaste_sbml/tests/test_model_builder.py | 77 ++++++++++++++++++++++++- 3 files changed, 102 insertions(+), 5 deletions(-) diff --git a/chaste_sbml/_config.py b/chaste_sbml/_config.py index e9d70ebd..c3c16ec0 100644 --- a/chaste_sbml/_config.py +++ b/chaste_sbml/_config.py @@ -16,6 +16,10 @@ DERIVATIVE_PREFIX = "d_" DERIVATIVE_SUFFIX = "_dt" +# Base id for the synthetic zero-derivative state variable added to models with no ODEs, so the +# generated ODE system always has at least one variable for the solver to integrate. +PLACEHOLDER_STATE_ID = "placeholder" + class EquationType(Enum): """Enumeration of equation types for code generation.""" diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index 49f809b3..aca17321 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -33,6 +33,7 @@ DERIVATIVE_SUFFIX, INITIAL_ASSIGNMENT_PREFIX, NON_DIM_UNITS, + PLACEHOLDER_STATE_ID, PREFIX_SEP, EquationType, EventType, @@ -327,6 +328,19 @@ def _add_state_variable( return state_var + def _add_placeholder_state_variable(self) -> None: + """Add a synthetic state variable with a zero derivative for a model with no ODEs. + + Chaste's ODE system integrates a vector of state variables, so a model with no continuous + dynamics still needs at least one. The placeholder is constant (``dy/dt = 0``), starts at + zero and is not an SBML output, so it never appears in the expected results; it exists only + to give the solver a trivial state to advance in time while the model's real outputs are + recomputed each step. + """ + placeholder_id = self._names.reserve(PREFIX_SEP.join([CHASTE_PREFIX, PLACEHOLDER_STATE_ID])) + state_var = self._add_state_variable(placeholder_id, "Placeholder state variable (model has no ODEs)", 0.0) + self._add_equation(var=state_var.derivative_id, math=parseL3Formula("0.0"), eq_type=EquationType.DERIVATIVE) + def _add_stoichiometry_variable(self, species_reference: "SpeciesReference") -> None: """Add a stoichiometry variable to the template variables. @@ -449,9 +463,9 @@ def _update_ode(species_ref: "SpeciesReference", rxn: "SbmlReaction", add: bool) math = rule.math odes[var] = math - if not odes: - raise NotImplementedError("Models without ODEs are not supported.") - + # An empty odes dict is allowed: a model with no reactions or rate rules has no + # continuous dynamics, and build() synthesises a placeholder state variable so the + # generated ODE system still has something for the solver to integrate. self._odes = odes def _total_time_derivative(self, ast_node: "ASTNode") -> str: @@ -1207,7 +1221,11 @@ def build(self) -> None: self._format_parameters() if not self._state_variables: - raise NotImplementedError("Models without state variables are not supported.") + # The model has no continuous dynamics (no reactions, ODEs or rate rules). Add a + # placeholder state variable with a zero derivative so the generated ODE system still + # has something for the solver to integrate; its outputs (constants, assignment rules + # of time, event-driven changes) are recomputed each step like any other model. + self._add_placeholder_state_variable() self._format_reactions() self._format_function_definitions() diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index 8822ef87..524c06b2 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -4,7 +4,14 @@ import pytest from libsbml import parseL3Formula -from chaste_sbml._config import ROOT_DIR, EquationType, VarType +from chaste_sbml._config import ( + CHASTE_PREFIX, + PLACEHOLDER_STATE_ID, + PREFIX_SEP, + ROOT_DIR, + EquationType, + VarType, +) from chaste_sbml._model_builder import ModelBuilder from chaste_sbml._names import NameConflictError, NameManager from chaste_sbml._records import ( @@ -171,6 +178,74 @@ def test_check_name_conflicts_flags_local_parameter_shadowing_time(): builder._check_name_conflicts() +PLACEHOLDER_ID = PREFIX_SEP.join([CHASTE_PREFIX, PLACEHOLDER_STATE_ID]) + + +def _no_ode_model(tmp_path): + """Write a no-ODE SBML model (constant species + time-dependent assignment rule) to disk. + + The model has no reactions, ODEs or rate rules, so it has no state variables of its own. + """ + doc = libsbml.SBMLDocument(3, 2) + model = doc.createModel() + model.setId("NoOde") + compartment = model.createCompartment() + compartment.setId("cell") + compartment.setConstant(True) + compartment.setSize(1.0) + compartment.setSpatialDimensions(3) + # A parameter that is an explicit function of time via an assignment rule. + parameter = model.createParameter() + parameter.setId("P") + parameter.setConstant(False) + parameter.setValue(0.0) + rule = model.createAssignmentRule() + rule.setVariable("P") + rule.setMath(parseL3Formula("2 * time")) + # A constant species. + species = model.createSpecies() + species.setId("S") + species.setCompartment("cell") + species.setConstant(True) + species.setBoundaryCondition(True) + species.setHasOnlySubstanceUnits(False) + species.setInitialConcentration(3.0) + + path = tmp_path / "NoOde.xml" + libsbml.writeSBMLToFile(doc, str(path)) + return path + + +def test_build_synthesises_placeholder_state_variable_for_no_ode_model(tmp_path): + """A model with no ODEs gets a single zero-derivative placeholder state variable.""" + builder = _build(_no_ode_model(tmp_path)) + + assert [s.id for s in builder._state_variables] == [PLACEHOLDER_ID] + placeholder = builder._state_variables[0] + assert placeholder.initial_value == 0.0 + + deriv_equations = [e for e in builder._equations if e.type == EquationType.DERIVATIVE] + assert len(deriv_equations) == 1 + assert deriv_equations[0].var == placeholder.derivative_id + assert deriv_equations[0].rhs == "0.0" + + +def test_build_no_ode_model_preserves_real_outputs(tmp_path): + """The placeholder does not displace the model's real outputs (assignment rule, species).""" + builder = _build(_no_ode_model(tmp_path)) + + # The time-dependent assignment-rule parameter is still emitted as a derived quantity, and + # the constant boundary species is still emitted (as a parameter). Neither is the placeholder. + assert "P" in {d.id for d in builder._derived_quantities} + assert "S" in {p.id for p in builder._parameters} + + +def test_build_does_not_add_placeholder_when_model_has_odes(): + """A model with genuine ODEs is left untouched: no placeholder state variable is added.""" + builder = _build(REFERENCE / "Goldbeter1991" / "Goldbeter1991.xml") + assert PLACEHOLDER_ID not in {s.id for s in builder._state_variables} + + def test_build_extracts_state_variables_and_derivatives(): """Each SBML species with an ODE becomes a state variable with a d__dt derivative.""" builder = _build(REFERENCE / "Goldbeter1991" / "Goldbeter1991.xml") From ef44f9968abcddecaadc25222a327515e4f2a840 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Sun, 5 Jul 2026 16:41:34 +0100 Subject: [PATCH 02/19] #42 Report renamed SBML ids under their original id + triage no-ODE cases When a real SBML id is a C++ keyword or reserved Chaste name (e.g. a parameter `time` or `true`), it is escaped for use as a C++ identifier (`time_`). Previously the escaped identifier was also used as the name reported to Chaste, so the variable could not be looked up by its real SBML id. Record the original id in NameManager and register it as the variable's name, keeping the escaped identifier only in the emitted code. Also re-triage the 253 no-ODE test-suite cases now that no-ODE models are supported: 161 pass, delays/random-event stay unsupported, and remaining failures (event ordering, n-ary relational, codegen gaps, missing MathML) are marked fail_known with notes. Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 506 +++++++++--------- chaste_sbml/_model_builder.py | 3 + chaste_sbml/_names.py | 17 + chaste_sbml/_records.py | 6 + chaste_sbml/templates/ode/ode.cpp | 6 +- chaste_sbml/tests/test_names.py | 21 + 6 files changed, 303 insertions(+), 256 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 1447d775..dc0abe4f 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -7,11 +7,11 @@ test,status,notes 00006,pass,ok 00007,pass,ok 00008,pass,ok -00009,unsupported,no odes +00009,pass,ok (no-ode placeholder) 00010,pass,ok 00011,pass,ok 00012,pass,ok -00013,unsupported,no odes +00013,pass,ok (no-ode placeholder) 00014,pass,ok 00015,pass,ok 00016,pass,ok @@ -27,8 +27,8 @@ test,status,notes 00026,pass,ok 00027,pass,ok 00028,pass,ok -00029,unsupported,no odes -00030,unsupported,no odes +00029,pass,ok (no-ode placeholder) +00030,pass,ok (no-ode placeholder) 00031,pass,ok 00032,pass,ok 00033,pass,ok @@ -116,7 +116,7 @@ test,status,notes 00115,pass,ok 00116,pass,ok 00117,pass,ok -00118,unsupported,no odes +00118,pass,ok (no-ode placeholder) 00119,pass,ok 00120,pass,ok 00121,pass,ok @@ -172,7 +172,7 @@ test,status,notes 00171,pass,ok 00172,pass,ok 00173,pass,ok -00174,unsupported,no odes +00174,pass,ok (no-ode placeholder) 00175,pass,ok 00176,pass,ok 00177,pass,ok @@ -211,7 +211,7 @@ test,status,notes 00210,pass,ok 00211,pass,ok 00212,pass,ok -00213,unsupported,no odes +00213,pass,ok (no-ode placeholder) 00214,pass,ok 00215,pass,ok 00216,pass,ok @@ -226,7 +226,7 @@ test,status,notes 00225,pass,ok 00226,pass,ok 00227,pass,ok -00228,unsupported,no odes +00228,pass,ok (no-ode placeholder) 00229,pass,ok 00230,pass,ok 00231,pass,ok @@ -577,7 +577,7 @@ test,status,notes 00576,unsupported,algebraic rule 00577,pass,ok 00578,pass,ok -00579,unsupported,no odes +00579,pass,ok (no-ode placeholder) 00580,pass,ok 00581,pass,ok 00582,pass,ok @@ -918,54 +918,54 @@ test,status,notes 00917,pass,ok 00918,pass,ok 00919,pass,ok -00920,unsupported,no odes -00921,unsupported,no odes -00922,unsupported,no odes -00923,unsupported,no odes -00924,unsupported,no odes -00925,unsupported,no odes +00920,pass,ok (no-ode placeholder) +00921,pass,ok (no-ode placeholder) +00922,pass,ok (no-ode placeholder) +00923,pass,ok (no-ode placeholder) +00924,pass,ok (no-ode placeholder) +00925,pass,ok (no-ode placeholder) 00926,pass,ok 00927,pass,ok 00928,pass,ok 00929,pass,ok -00930,unsupported,no odes -00931,unsupported,no odes +00930,fail_known,event priority/ordering +00931,fail_known,event priority/ordering 00932,unsupported,delay 00933,unsupported,delay -00934,unsupported,no odes -00935,unsupported,no odes -00936,unsupported,no odes -00937,unsupported,no odes +00934,fail_known,event priority/ordering +00935,fail_known,event priority/ordering +00936,unsupported,delay +00937,unsupported,delay function 00938,unsupported,delay 00939,unsupported,delay -00940,unsupported,no odes -00941,unsupported,no odes -00942,unsupported,no odes -00943,unsupported,no odes +00940,unsupported,delay function +00941,unsupported,delay function +00942,unsupported,delay function +00943,unsupported,delay function 00944,pass,ok 00945,pass,ok 00946,pass,ok 00947,pass,ok 00948,pass,ok -00949,unsupported,no odes -00950,unsupported,no odes -00951,unsupported,no odes -00952,unsupported,no odes -00953,unsupported,no odes -00954,unsupported,no odes -00955,unsupported,no odes -00956,unsupported,no odes -00957,unsupported,no odes -00958,unsupported,no odes -00959,unsupported,no odes -00960,unsupported,no odes -00961,unsupported,no odes -00962,unsupported,no odes -00963,unsupported,no odes -00964,unsupported,no odes -00965,unsupported,no odes -00966,unsupported,no odes -00967,unsupported,no odes +00949,pass,ok (no-ode placeholder) +00950,fail_known,codegen: inf constant +00951,fail_known,codegen: inf constant +00952,unsupported,random event execution +00953,fail_known,event priority/ordering +00954,pass,ok (no-ode placeholder) +00955,pass,ok (no-ode placeholder) +00956,pass,ok (no-ode placeholder) +00957,pass,ok (no-ode placeholder) +00958,pass,ok (no-ode placeholder) +00959,pass,ok (no-ode placeholder) +00960,pass,ok (no-ode placeholder) +00961,pass,ok (no-ode placeholder) +00962,unsupported,random event execution +00963,fail_known,event priority/ordering +00964,unsupported,random event execution +00965,unsupported,random event execution +00966,unsupported,random event execution +00967,fail_known,event priority/ordering 00968,pass,ok 00969,pass,ok 00970,pass,ok @@ -976,9 +976,9 @@ test,status,notes 00975,pass,ok 00976,pass,ok 00977,pass,ok -00978,unsupported,no odes -00979,unsupported,no odes -00980,unsupported,no odes +00978,fail_known,event priority/ordering +00979,pass,ok (no-ode placeholder) +00980,unsupported,delay 00981,unsupported,delay 00982,unsupported,delay 00983,unsupported,algebraic rule @@ -993,9 +993,9 @@ test,status,notes 00992,pass,ok 00993,unsupported,algebraic rule 00994,pass,ok -00995,unsupported,no odes +00995,pass,ok (no-ode placeholder) 00996,pass,ok -00997,unsupported,no odes +00997,pass,ok (no-ode placeholder) 00998,pass,ok 00999,pass,ok 01000,unsupported,delay @@ -1110,33 +1110,33 @@ test,status,notes 01109,pass,ok 01110,pass,ok 01111,pass,ok -01112,unsupported,no odes -01113,unsupported,no odes -01114,unsupported,no odes -01115,unsupported,no odes -01116,unsupported,no odes +01112,pass,ok (no-ode placeholder) +01113,pass,ok (no-ode placeholder) +01114,pass,ok (no-ode placeholder) +01115,pass,ok (no-ode placeholder) +01116,fail_known,codegen: piecewise arity 01117,pass,ok 01118,pass,ok -01119,unsupported,no odes +01119,unsupported,delay 01120,unsupported,delay 01121,pass,ok 01122,pass,ok 01123,pass,ok -01124,unsupported,no odes -01125,unsupported,no odes -01126,unsupported,no odes +01124,pass,ok (no-ode placeholder) +01125,pass,ok (no-ode placeholder) +01126,pass,ok (no-ode placeholder) 01127,pass,ok -01128,unsupported,no odes -01129,unsupported,no odes -01130,unsupported,no odes -01131,unsupported,no odes -01132,unsupported,no odes -01133,unsupported,no odes -01134,unsupported,no odes -01135,unsupported,no odes -01136,unsupported,no odes -01137,unsupported,no odes -01138,unsupported,no odes +01128,pass,ok (no-ode placeholder) +01129,pass,ok (no-ode placeholder) +01130,pass,ok (no-ode placeholder) +01131,pass,ok (no-ode placeholder) +01132,pass,ok (no-ode placeholder) +01133,pass,ok (no-ode placeholder) +01134,pass,ok (no-ode placeholder) +01135,pass,ok (no-ode placeholder) +01136,pass,ok (no-ode placeholder) +01137,pass,ok (no-ode placeholder) +01138,pass,ok (no-ode placeholder) 01139,pass,ok 01140,pass,ok 01141,pass,ok @@ -1144,23 +1144,23 @@ test,status,notes 01143,pass,ok 01144,pass,ok 01145,pass,ok -01146,unsupported,no odes -01147,unsupported,no odes +01146,pass,ok (no-ode placeholder) +01147,pass,ok (no-ode placeholder) 01148,pass,ok -01149,unsupported,no odes +01149,pass,ok (no-ode placeholder) 01150,pass,ok 01151,pass,ok -01152,unsupported,no odes -01153,unsupported,no odes -01154,unsupported,no odes -01155,unsupported,no odes -01156,unsupported,no odes -01157,unsupported,no odes -01158,unsupported,no odes +01152,pass,ok (no-ode placeholder) +01153,pass,ok (no-ode placeholder) +01154,pass,ok (no-ode placeholder) +01155,pass,ok (no-ode placeholder) +01156,unsupported,delay +01157,pass,ok (no-ode placeholder) +01158,pass,ok (no-ode placeholder) 01159,pass,ok 01160,pass,ok 01161,pass,ok -01162,unsupported,no odes +01162,pass,ok (no-ode placeholder) 01163,pass,ok 01164,pass,ok 01165,pass,ok @@ -1169,13 +1169,13 @@ test,status,notes 01168,pass,ok 01169,pass,ok 01170,pass,ok -01171,unsupported,no odes +01171,pass,ok (no-ode placeholder) 01172,pass,ok -01173,unsupported,no odes +01173,unsupported,delay function 01174,unsupported,algebraic rule 01175,pass,ok 01176,unsupported,delay -01177,unsupported,no odes +01177,unsupported,delay 01178,pass,ok 01179,pass,ok 01180,pass,ok @@ -1207,16 +1207,16 @@ test,status,notes 01206,pass,ok 01207,pass,ok 01208,pass,ok -01209,unsupported,no odes -01210,unsupported,no odes -01211,unsupported,no odes -01212,unsupported,no odes -01213,unsupported,no odes -01214,unsupported,no odes +01209,pass,ok (no-ode placeholder) +01210,pass,ok (no-ode placeholder) +01211,pass,ok (no-ode placeholder) +01212,fail_known,nary relational operator +01213,unsupported,delay +01214,pass,ok (no-ode placeholder) 01215,pass,ok -01216,unsupported,no odes -01217,unsupported,no odes -01218,unsupported,no odes +01216,fail_known,nary relational operator +01217,pass,ok (no-ode placeholder) +01218,pass,ok (no-ode placeholder) 01219,pass,ok 01220,pass,ok 01221,pass,ok @@ -1232,22 +1232,22 @@ test,status,notes 01231,pass,ok 01232,pass,ok 01233,pass,ok -01234,unsupported,no odes -01235,unsupported,no odes -01236,unsupported,no odes -01237,unsupported,no odes -01238,unsupported,no odes -01239,unsupported,no odes -01240,unsupported,no odes -01241,unsupported,no odes -01242,unsupported,no odes -01243,unsupported,no odes +01234,pass,ok (no-ode placeholder) +01235,pass,ok (no-ode placeholder) +01236,pass,ok (no-ode placeholder) +01237,pass,ok (no-ode placeholder) +01238,fail_known,missing mathml (event trigger) +01239,fail_known,missing mathml (event trigger) +01240,pass,ok (no-ode placeholder) +01241,fail_known,missing mathml (event delay) +01242,pass,ok (no-ode placeholder) +01243,pass,ok (no-ode placeholder) 01244,unsupported,algebraic rule -01245,unsupported,no odes +01245,pass,ok (no-ode placeholder) 01246,pass,ok -01247,unsupported,no odes -01248,unsupported,no odes -01249,unsupported,no odes +01247,pass,ok (no-ode placeholder) +01248,pass,ok (no-ode placeholder) +01249,pass,ok (no-ode placeholder) 01250,pass,ok 01251,pass,ok 01252,pass,ok @@ -1269,27 +1269,27 @@ test,status,notes 01268,unsupported,delay 01269,unsupported,delay 01270,unsupported,delay -01271,unsupported,no odes -01272,unsupported,no odes -01273,unsupported,no odes -01274,unsupported,no odes -01275,unsupported,no odes -01276,unsupported,no odes -01277,unsupported,no odes -01278,unsupported,no odes -01279,unsupported,no odes -01280,unsupported,no odes -01281,unsupported,no odes -01282,unsupported,no odes -01283,unsupported,no odes -01284,unsupported,no odes -01285,unsupported,no odes -01286,unsupported,no odes -01287,unsupported,no odes +01271,fail_known,missing mathml (function definition) +01272,pass,ok (no-ode placeholder) +01273,pass,ok (no-ode placeholder) +01274,fail_known,codegen: implies operator +01275,fail_known,codegen: implies operator +01276,fail_known,codegen: implies operator +01277,pass,ok (no-ode placeholder) +01278,pass,ok (no-ode placeholder) +01279,fail_known,codegen: single-arg min +01280,fail_known,codegen: implies operator +01281,fail_known,codegen: implies operator +01282,pass,ok (no-ode placeholder) +01283,pass,ok (no-ode placeholder) +01284,pass,ok (no-ode placeholder) +01285,pass,ok (no-ode placeholder) +01286,pass,ok (no-ode placeholder) +01287,unsupported,delay 01288,pass,ok -01289,unsupported,no odes +01289,pass,ok (no-ode placeholder) 01290,pass,ok -01291,unsupported,no odes +01291,pass,ok (no-ode placeholder) 01292,unsupported,algebraic rule 01293,pass,ok 01294,pass,ok @@ -1298,50 +1298,50 @@ test,status,notes 01297,pass,ok 01298,pass,ok 01299,unsupported,delay -01300,unsupported,no odes -01301,unsupported,no odes +01300,pass,ok (no-ode placeholder) +01301,pass,ok (no-ode placeholder) 01302,pass,ok -01303,unsupported,no odes -01304,unsupported,no odes -01305,unsupported,no odes -01306,unsupported,no odes +01303,pass,ok (no-ode placeholder) +01304,pass,ok (no-ode placeholder) +01305,unsupported,delay +01306,pass,ok (no-ode placeholder) 01307,pass,ok 01308,pass,ok 01309,pass,ok -01310,unsupported,no odes -01311,unsupported,no odes -01312,unsupported,no odes -01313,unsupported,no odes -01314,unsupported,no odes -01315,unsupported,no odes -01316,unsupported,no odes -01317,unsupported,no odes -01318,unsupported,no odes -01319,unsupported,no odes +01310,pass,ok (no-ode placeholder) +01311,pass,ok (no-ode placeholder) +01312,fail_known,codegen: expression called as function +01313,pass,ok (no-ode placeholder) +01314,pass,ok (no-ode placeholder) +01315,pass,ok (no-ode placeholder) +01316,pass,ok (no-ode placeholder) +01317,pass,ok (no-ode placeholder) +01318,unsupported,delay function +01319,unsupported,delay function 01320,unsupported,delay 01321,pass,ok 01322,pass,ok -01323,unsupported,no odes +01323,pass,ok (no-ode placeholder) 01324,unsupported,delay 01325,unsupported,delay 01326,unsupported,delay 01327,unsupported,delay -01328,unsupported,no odes -01329,unsupported,no odes -01330,unsupported,no odes -01331,unsupported,no odes -01332,unsupported,no odes -01333,unsupported,no odes -01334,unsupported,no odes -01335,unsupported,no odes -01336,unsupported,no odes -01337,unsupported,no odes +01328,unsupported,delay +01329,unsupported,delay +01330,pass,ok (no-ode placeholder) +01331,fail_known,event priority/ordering +01332,pass,ok (no-ode placeholder) +01333,fail_known,event priority/ordering +01334,pass,ok (no-ode placeholder) +01335,unsupported,delay +01336,pass,ok (no-ode placeholder) +01337,fail_known,event priority/ordering 01338,pass,ok 01339,pass,ok 01340,pass,ok -01341,unsupported,no odes -01342,unsupported,no odes -01343,unsupported,no odes +01341,pass,ok (no-ode placeholder) +01342,pass,ok (no-ode placeholder) +01343,pass,ok (no-ode placeholder) 01344,pass,ok 01345,pass,ok 01346,pass,ok @@ -1351,30 +1351,30 @@ test,status,notes 01350,unsupported,algebraic rule 01351,pass,ok 01352,pass,ok -01353,unsupported,no odes -01354,unsupported,no odes -01355,unsupported,no odes -01356,unsupported,no odes -01357,unsupported,no odes -01358,unsupported,no odes +01353,pass,ok (no-ode placeholder) +01354,pass,ok (no-ode placeholder) +01355,pass,ok (no-ode placeholder) +01356,pass,ok (no-ode placeholder) +01357,pass,ok (no-ode placeholder) +01358,unsupported,delay 01359,unsupported,algebraic rule 01360,pass,ok 01361,pass,ok -01362,unsupported,no odes -01363,unsupported,no odes -01364,unsupported,no odes -01365,unsupported,no odes -01366,unsupported,no odes -01367,unsupported,no odes +01362,pass,ok (no-ode placeholder) +01363,pass,ok (no-ode placeholder) +01364,pass,ok (no-ode placeholder) +01365,pass,ok (no-ode placeholder) +01366,pass,ok (no-ode placeholder) +01367,unsupported,delay 01368,unsupported,algebraic rule 01369,pass,ok 01370,pass,ok -01371,unsupported,no odes -01372,unsupported,no odes -01373,unsupported,no odes -01374,unsupported,no odes -01375,unsupported,no odes -01376,unsupported,no odes +01371,pass,ok (no-ode placeholder) +01372,pass,ok (no-ode placeholder) +01373,pass,ok (no-ode placeholder) +01374,pass,ok (no-ode placeholder) +01375,pass,ok (no-ode placeholder) +01376,unsupported,delay 01377,unsupported,algebraic rule 01378,pass,ok 01379,pass,ok @@ -1388,11 +1388,11 @@ test,status,notes 01387,pass,ok 01388,pass,ok 01389,pass,ok -01390,unsupported,no odes -01391,unsupported,no odes -01392,unsupported,no odes -01393,unsupported,no odes -01394,unsupported,no odes +01390,pass,ok (no-ode placeholder) +01391,pass,ok (no-ode placeholder) +01392,pass,ok (no-ode placeholder) +01393,pass,ok (no-ode placeholder) +01394,pass,ok (no-ode placeholder) 01395,pass,ok 01396,unsupported,fast reaction 01397,unsupported,fast reaction @@ -1414,7 +1414,7 @@ test,status,notes 01413,unsupported,delay 01414,unsupported,delay 01415,unsupported,delay -01416,unsupported,no odes +01416,unsupported,delay function 01417,unsupported,delay 01418,unsupported,delay 01419,unsupported,delay @@ -1459,13 +1459,13 @@ test,status,notes 01458,pass,ok 01459,pass,ok 01460,pass,ok -01461,unsupported,no odes +01461,pass,ok (no-ode placeholder) 01462,pass,ok 01463,pass,ok 01464,pass,ok 01465,pass,ok -01466,unsupported,no odes -01467,unsupported,no odes +01466,pass,ok (no-ode placeholder) +01467,pass,ok (no-ode placeholder) 01468,pass,ok 01469,pass,ok 01470,pass,ok @@ -1483,19 +1483,19 @@ test,status,notes 01482,unsupported,algebraic rule 01483,unsupported,algebraic rule 01484,unsupported,algebraic rule -01485,unsupported,no odes -01486,unsupported,no odes -01487,unsupported,no odes -01488,unsupported,no odes -01489,unsupported,no odes -01490,unsupported,no odes -01491,unsupported,no odes +01485,pass,ok (no-ode placeholder) +01486,fail_known,codegen: unconverted ^ operator +01487,pass,ok (no-ode placeholder) +01488,pass,ok (no-ode placeholder) +01489,pass,ok (no-ode placeholder) +01490,fail_known,codegen: unconverted ^ operator +01491,fail_known,codegen: unconverted ^ operator 01492,pass,ok 01493,pass,ok -01494,unsupported,no odes -01495,unsupported,no odes -01496,unsupported,no odes -01497,unsupported,no odes +01494,fail_known,nary relational operator +01495,pass,ok (no-ode placeholder) +01496,pass,ok (no-ode placeholder) +01497,fail_known,codegen: implies operator 01498,pass,ok 01499,unsupported,algebraic rule 01500,unsupported,algebraic rule @@ -1520,7 +1520,7 @@ test,status,notes 01519,unsupported,delay 01520,unsupported,delay 01521,unsupported,delay -01522,unsupported,no odes +01522,unsupported,delay function 01523,unsupported,delay 01524,unsupported,delay 01525,unsupported,delay @@ -1528,10 +1528,10 @@ test,status,notes 01527,pass,ok 01528,unsupported,delay 01529,unsupported,delay -01530,unsupported,no odes +01530,pass,ok (no-ode placeholder) 01531,partly_supported,clustered events 01532,unsupported,delay -01533,unsupported,no odes +01533,pass,ok (no-ode placeholder) 01534,unsupported,delay 01535,unsupported,delay 01536,unsupported,delay @@ -1586,23 +1586,23 @@ test,status,notes 01585,unsupported,delay 01586,unsupported,delay 01587,unsupported,delay -01588,unsupported,no odes +01588,unsupported,random event execution 01589,unsupported,algebraic rule -01590,unsupported,no odes +01590,unsupported,delay 01591,unsupported,random event execution 01592,unsupported,delay 01593,unsupported,delay -01594,unsupported,no odes -01595,unsupported,no odes -01596,unsupported,no odes -01597,unsupported,no odes -01598,unsupported,no odes +01594,unsupported,delay +01595,unsupported,delay +01596,pass,ok (no-ode placeholder) +01597,unsupported,delay +01598,unsupported,delay 01599,unsupported,random event execution -01600,unsupported,no odes -01601,unsupported,no odes -01602,unsupported,no odes -01603,unsupported,no odes -01604,unsupported,no odes +01600,unsupported,delay +01601,unsupported,delay +01602,unsupported,delay +01603,unsupported,delay +01604,unsupported,delay 01605,unsupported,random event execution 01606,unsupported,flux balance 01607,unsupported,flux balance @@ -1633,13 +1633,13 @@ test,status,notes 01632,pass,ok 01633,pass,ok 01634,pass,ok -01635,unsupported,no odes +01635,pass,ok (no-ode placeholder) 01636,pass,ok 01637,pass,ok -01638,unsupported,no odes +01638,pass,ok (no-ode placeholder) 01639,pass,ok 01640,pass,ok -01641,unsupported,no odes +01641,pass,ok (no-ode placeholder) 01642,pass,ok 01643,pass,ok 01644,pass,ok @@ -1656,13 +1656,13 @@ test,status,notes 01655,pass,ok 01656,pass,ok 01657,pass,ok -01658,unsupported,no odes -01659,unsupported,no odes -01660,unsupported,no odes -01661,unsupported,no odes -01662,unsupported,no odes -01663,unsupported,no odes -01664,unsupported,no odes +01658,pass,ok (no-ode placeholder) +01659,unsupported,delay +01660,unsupported,delay +01661,unsupported,delay +01662,pass,ok (no-ode placeholder) +01663,pass,ok (no-ode placeholder) +01664,pass,ok (no-ode placeholder) 01665,pass,ok 01666,pass,ok 01667,pass,ok @@ -1691,14 +1691,14 @@ test,status,notes 01690,unsupported,delay 01691,unsupported,delay 01692,unsupported,delay -01693,unsupported,no odes -01694,unsupported,no odes -01695,unsupported,no odes -01696,unsupported,no odes -01697,unsupported,no odes -01698,unsupported,no odes -01699,unsupported,no odes -01700,unsupported,no odes +01693,pass,ok (no-ode placeholder) +01694,pass,ok (no-ode placeholder) +01695,pass,ok (no-ode placeholder) +01696,pass,ok (no-ode placeholder) +01697,pass,ok (no-ode placeholder) +01698,pass,ok (no-ode placeholder) +01699,pass,ok (no-ode placeholder) +01700,pass,ok (no-ode placeholder) 01701,unsupported,delay 01702,unsupported,delay 01703,unsupported,delay @@ -1752,12 +1752,12 @@ test,status,notes 01751,pass,ok 01752,pass,ok 01753,pass,ok -01754,unsupported,no odes -01755,unsupported,no odes -01756,unsupported,no odes -01757,unsupported,no odes -01758,unsupported,no odes -01759,unsupported,no odes +01754,unsupported,delay +01755,unsupported,delay +01756,unsupported,delay +01757,unsupported,delay +01758,unsupported,delay +01759,unsupported,delay 01760,pass,ok 01761,pass,ok 01762,pass,ok @@ -1777,10 +1777,10 @@ test,status,notes 01776,pass,ok 01777,pass,ok 01778,pass,ok -01779,unsupported,no odes -01780,unsupported,no odes -01781,unsupported,no odes -01782,unsupported,no odes +01779,pass,ok (no-ode placeholder) +01780,pass,ok (no-ode placeholder) +01781,fail_known,nary relational operator +01782,pass,ok (no-ode placeholder) 01783,pass,ok 01784,pass,ok 01785,unsupported,algebraic rule @@ -1794,29 +1794,29 @@ test,status,notes 01793,unsupported,algebraic rule 01794,unsupported,algebraic rule 01795,pass,ok -01796,unsupported,no odes +01796,pass,ok (no-ode placeholder) 01797,pass,ok 01798,unsupported,delay 01799,pass,ok 01800,pass,ok 01801,pass,ok 01802,pass,ok -01803,unsupported,no odes -01804,unsupported,no odes -01805,unsupported,no odes -01806,unsupported,no odes +01803,pass,ok (no-ode placeholder) +01804,pass,ok (no-ode placeholder) +01805,pass,ok (no-ode placeholder) +01806,pass,ok (no-ode placeholder) 01807,pass,ok 01808,pass,ok 01809,pass,ok -01810,unsupported,no odes -01811,unsupported,no odes -01812,unsupported,no odes -01813,unsupported,no odes -01814,unsupported,no odes -01815,unsupported,no odes -01816,unsupported,no odes -01817,unsupported,no odes -01818,unsupported,no odes -01819,unsupported,no odes -01820,unsupported,no odes -01821,unsupported,no odes +01810,pass,ok (no-ode placeholder) +01811,fail_known,codegen: inf in expected results +01812,fail_known,codegen: string literal in math +01813,fail_known,codegen: string literal in math +01814,pass,ok (no-ode placeholder) +01815,pass,ok (no-ode placeholder) +01816,pass,ok (no-ode placeholder) +01817,pass,ok (no-ode placeholder) +01818,pass,ok (no-ode placeholder) +01819,fail_known,pi constant vs parameter collision +01820,pass,ok (no-ode placeholder) +01821,pass,ok (no-ode placeholder) diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index aca17321..b54d2341 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -179,6 +179,7 @@ def _add_derived_quantity( is_conversion=is_conversion, is_reaction=is_reaction, units=units, + name=self._names.sbml_name(id_), ) ) if not is_reaction: @@ -280,6 +281,7 @@ def _add_parameter( label=label, initial_value=initial_value, units=units, + name=self._names.sbml_name(id_), ) ) self._variable_types[id_] = VarType.PARAMETER @@ -321,6 +323,7 @@ def _add_state_variable( label=label, initial_value=initial_value, units=units, + name=self._names.sbml_name(id_), ) self._state_variables.append(state_var) diff --git a/chaste_sbml/_names.py b/chaste_sbml/_names.py index 31635fd8..bbc35dce 100644 --- a/chaste_sbml/_names.py +++ b/chaste_sbml/_names.py @@ -257,6 +257,10 @@ def __init__(self, sbml_model) -> None: """:param sbml_model: The libsbml Model whose ids are being turned into C++.""" self._sbml_model = sbml_model self._taken: set = set() + # C++ identifier (the renamed id) -> original SBML id, for ids renamed to dodge a C++ + # keyword or reserved name. Lets the generated code report a variable under its real SBML + # id even though its C++ identifier had to be escaped. + self._sbml_names: dict[str, str] = {} def resolve_real_id_conflicts(self) -> None: """Rename real SBML ids that are C++ keywords or reserved Chaste names. @@ -296,6 +300,19 @@ def resolve_real_id_conflicts(self) -> None: for referrer in model.getListOfAllElements(): referrer.renameSIdRefs(old_id, new_id) taken.add(new_id) + self._sbml_names[new_id] = old_id + + def sbml_name(self, cpp_id: str) -> str: + """Return the original SBML id for a (possibly renamed) C++ identifier. + + A variable is reported to Chaste under this name, so it is looked up by its real SBML id + even when the emitted C++ identifier had to be escaped (e.g. a parameter ``time`` emitted + as ``time_``). Ids that were never renamed are returned unchanged. + + :param cpp_id: The C++ identifier as it appears in the generated code. + :return: The original SBML id if ``cpp_id`` was a renamed id, else ``cpp_id``. + """ + return self._sbml_names.get(cpp_id, cpp_id) def reset(self) -> None: """Recompute the taken-name set from the (resolved) model, ready for a build. diff --git a/chaste_sbml/_records.py b/chaste_sbml/_records.py index 1058850b..d98b4050 100644 --- a/chaste_sbml/_records.py +++ b/chaste_sbml/_records.py @@ -45,6 +45,8 @@ class StateVariable: label: str initial_value: Optional[float] units: str + # SBML id reported to Chaste; equals ``id`` unless the id was escaped for C++ (see NameManager). + name: str = "" @dataclass @@ -57,6 +59,8 @@ class Parameter: label: str initial_value: Optional[float] units: str + # SBML id reported to Chaste; equals ``id`` unless the id was escaped for C++ (see NameManager). + name: str = "" @dataclass @@ -70,6 +74,8 @@ class DerivedQuantity: is_conversion: bool is_reaction: bool units: str + # SBML id reported to Chaste; equals ``id`` unless the id was escaped for C++ (see NameManager). + name: str = "" @dataclass diff --git a/chaste_sbml/templates/ode/ode.cpp b/chaste_sbml/templates/ode/ode.cpp index e22ee65d..93aadb4d 100644 --- a/chaste_sbml/templates/ode/ode.cpp +++ b/chaste_sbml/templates/ode/ode.cpp @@ -327,7 +327,7 @@ void CellwiseOdeSystemInformation<{{ ode_class_name }}>::Initialise() { // STATE VARIABLES {% for var in state_variables %} - this->mVariableNames.push_back("{{ var['id'] }}"); + this->mVariableNames.push_back("{{ var['name'] }}"); this->mVariableUnits.push_back("{{ var['units'] }}"); {% if var["initial_value"] is not none %} this->mInitialConditions.push_back({{ var['initial_value'] }}); @@ -339,14 +339,14 @@ void CellwiseOdeSystemInformation<{{ ode_class_name }}>::Initialise() // DERIVED QUANTITIES {% for dq in derived_quantities %} - this->mDerivedQuantityNames.push_back("{{ dq['id'] }}"); + this->mDerivedQuantityNames.push_back("{{ dq['name'] }}"); this->mDerivedQuantityUnits.push_back("{{ dq['units'] }}"); {% endfor %} // PARAMETERS {% for var in parameters %} - this->mParameterNames.push_back("{{ var['id'] }}"); + this->mParameterNames.push_back("{{ var['name'] }}"); this->mParameterUnits.push_back("{{ var['units'] }}"); {% endfor %} diff --git a/chaste_sbml/tests/test_names.py b/chaste_sbml/tests/test_names.py index a188e504..cb8be04e 100644 --- a/chaste_sbml/tests/test_names.py +++ b/chaste_sbml/tests/test_names.py @@ -166,6 +166,27 @@ def test_name_manager_resolve_leaves_safe_ids_untouched(): assert sbml_model.getElementBySId("cell") is not None +def test_name_manager_sbml_name_recovers_original_id(): + """After a keyword id is renamed, sbml_name maps the C++ id back to the original SBML id. + + A variable must be reported to Chaste under its real SBML id even though the emitted C++ + identifier had to be escaped, so it can be looked up by that id. + """ + doc = libsbml.SBMLDocument(3, 2) + sbml_model = doc.createModel() + parameter = sbml_model.createParameter() + parameter.setId("true") # 'true' is a C++ keyword, so it is renamed to 'true_' + parameter.setConstant(True) + parameter.setValue(1.0) + + manager = NameManager(sbml_model) + manager.resolve_real_id_conflicts() + + assert sbml_model.getParameter("true_") is not None + assert manager.sbml_name("true_") == "true" # renamed id maps back to the real SBML id + assert manager.sbml_name("k") == "k" # an un-renamed id is returned unchanged + + def test_generate_header_guard(): """Test header guard generation.""" test_cases = [ From c52f53e51e827013058d00b0530d949c96d3bdab Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Sun, 5 Jul 2026 17:05:10 +0100 Subject: [PATCH 03/19] #42 Emit n-ary relationals via the sm:: SbmlMath helpers SBML MathML relationals are n-ary: lt(a, b, c) means a < b < c (a2 operand) relational AST nodes to function-call form before stringifying, so they dispatch to the variadic sm::lt/gt/leq/geq/eq helpers with the correct semantics. Binary relationals stay as infix operators, so existing generated code is unchanged. Marks test-suite cases 01212, 01216, 01494, 01781 as passing. Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 8 +-- chaste_sbml/_expressions.py | 49 ++++++++++++++++++- chaste_sbml/tests/test_expressions.py | 32 ++++++++++++ 3 files changed, 84 insertions(+), 5 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index dc0abe4f..9726039c 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1210,11 +1210,11 @@ test,status,notes 01209,pass,ok (no-ode placeholder) 01210,pass,ok (no-ode placeholder) 01211,pass,ok (no-ode placeholder) -01212,fail_known,nary relational operator +01212,pass,ok (nary relational via sm:: helpers) 01213,unsupported,delay 01214,pass,ok (no-ode placeholder) 01215,pass,ok -01216,fail_known,nary relational operator +01216,pass,ok (nary relational via sm:: helpers) 01217,pass,ok (no-ode placeholder) 01218,pass,ok (no-ode placeholder) 01219,pass,ok @@ -1492,7 +1492,7 @@ test,status,notes 01491,fail_known,codegen: unconverted ^ operator 01492,pass,ok 01493,pass,ok -01494,fail_known,nary relational operator +01494,pass,ok (nary relational via sm:: helpers) 01495,pass,ok (no-ode placeholder) 01496,pass,ok (no-ode placeholder) 01497,fail_known,codegen: implies operator @@ -1779,7 +1779,7 @@ test,status,notes 01778,pass,ok 01779,pass,ok (no-ode placeholder) 01780,pass,ok (no-ode placeholder) -01781,fail_known,nary relational operator +01781,pass,ok (nary relational via sm:: helpers) 01782,pass,ok (no-ode placeholder) 01783,pass,ok 01784,pass,ok diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index e3a304fe..7078cfb5 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -9,7 +9,18 @@ import re from typing import TYPE_CHECKING, Optional -from libsbml import AST_FUNCTION_DELAY, AST_NAME, AST_NAME_AVOGADRO, formulaToL3String +from libsbml import ( + AST_FUNCTION, + AST_FUNCTION_DELAY, + AST_NAME, + AST_NAME_AVOGADRO, + AST_RELATIONAL_EQ, + AST_RELATIONAL_GEQ, + AST_RELATIONAL_GT, + AST_RELATIONAL_LEQ, + AST_RELATIONAL_LT, + formulaToL3String, +) from ._config import CHASTE_PREFIX, PREFIX_SEP, VarType @@ -87,6 +98,39 @@ def strip_ast_units(node: "ASTNode") -> None: strip_ast_units(node.getChild(i)) +# SBML MathML relational operators, and the SbmlMath (``sm::``) helper each maps to. +_NARY_RELATIONAL_NAMES = { + AST_RELATIONAL_LT: "lt", + AST_RELATIONAL_GT: "gt", + AST_RELATIONAL_LEQ: "leq", + AST_RELATIONAL_GEQ: "geq", + AST_RELATIONAL_EQ: "eq", +} + + +def rewrite_nary_relational(node: Optional["ASTNode"]) -> None: + """Rewrite n-ary (>2 operand) relational AST nodes into function-call form, in place. + + SBML MathML relationals are n-ary -- ``lt(a, b, c)`` means ``a < b < c`` (i.e. ``a < b`` and + ``b < c``), and ``eq(a, b, c)`` means all equal -- but ``formulaToL3String`` renders them as a + chained infix expression (``a < b < c``), which C++ mis-evaluates left-to-right as + ``(a < b) < c``. Converting each such node to a function call (rendered ``lt(a, b, c)``, later + mapped to the variadic ``sm::lt`` helper that implements the correct semantics) fixes this. + Binary relationals are left as infix operators, so existing generated code is unchanged. + ``neq`` is binary-only in MathML, so it is never n-ary and is not rewritten here. + + :param node: The root AST node to rewrite. + """ + if node is None: + return + name = _NARY_RELATIONAL_NAMES.get(node.getType()) + if name is not None and node.getNumChildren() > 2: + node.setType(AST_FUNCTION) + node.setName(name) + for i in range(node.getNumChildren()): + rewrite_nary_relational(node.getChild(i)) + + def search_ast_type(root: Optional["ASTNode"], node_type: int) -> bool: """Recursively search the AST for a node of a certain type. @@ -260,6 +304,9 @@ def formula_to_string( # constants below) so it stays distinct from a same-named parameter. avogadro_placeholder = f"{CHASTE_PREFIX}{PREFIX_SEP}avogadro" replace_avogadro_csymbol(math, avogadro_placeholder) + # Route n-ary relationals (a < b < c) to the sm:: helpers before stringifying, so they are + # not emitted as C++-misevaluating chained infix. Binary relationals stay as infix operators. + rewrite_nary_relational(math) formula = formulaToL3String(math) # Convert all integer literals to doubles to fix integer division. diff --git a/chaste_sbml/tests/test_expressions.py b/chaste_sbml/tests/test_expressions.py index 24ad0fbd..92c7e923 100644 --- a/chaste_sbml/tests/test_expressions.py +++ b/chaste_sbml/tests/test_expressions.py @@ -19,6 +19,7 @@ convert_infix_operator_to_function_syntax, formula_to_string, replace_avogadro_csymbol, + rewrite_nary_relational, search_ast_type, strip_ast_units, substitute_ast_names, @@ -75,6 +76,37 @@ def test_formula_to_string_functions_constants_and_power(): assert out == "k * std::pow(S, 2.0) + M_PI" +@pytest.mark.parametrize( + ("formula", "expected"), + [ + # n-ary relationals (>2 operands) become sm:: calls with the correct semantics... + ("lt(1, 2, 3)", "sm::lt(1.0, 2.0, 3.0)"), + ("gt(2, 1, 2)", "sm::gt(2.0, 1.0, 2.0)"), + ("leq(1, 2, 3, 4)", "sm::leq(1.0, 2.0, 3.0, 4.0)"), + ("geq(2, 1, 2)", "sm::geq(2.0, 1.0, 2.0)"), + ("eq(1, 1, 2)", "sm::eq(1.0, 1.0, 2.0)"), + # ...while binary relationals stay as infix operators (unchanged output). + ("a < b", "a < b"), + ("a >= b", "a >= b"), + ], +) +def test_formula_to_string_nary_relational_uses_sbmlmath(formula, expected): + """n-ary relationals route to the variadic sm:: helpers; binary ones stay infix.""" + variable_types = {"a": VarType.PARAMETER, "b": VarType.PARAMETER} + assert formula_to_string(parseL3Formula(formula), variable_types, state_variables=[]) == expected + + +def test_rewrite_nary_relational_leaves_binary_untouched(): + """rewrite_nary_relational converts only >2-operand relationals, leaving binary ones as-is.""" + binary = parseL3Formula("a < b") + rewrite_nary_relational(binary) + assert formulaToL3String(binary) == "a < b" + + nary = parseL3Formula("lt(a, b, c)") + rewrite_nary_relational(nary) + assert formulaToL3String(nary) == "lt(a, b, c)" + + def test_formula_to_string_rateof_state_variable_uses_derivative(): """rateOf(state var) resolves to that variable's derivative id.""" math = parseL3Formula("rateOf(X)") From af7440b645508c17b6d3d6dbd8f8e3a38812d597 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Sun, 5 Jul 2026 17:53:12 +0100 Subject: [PATCH 04/19] #42 Convert powers to pow() on the AST; tidy formula_to_string Replace the string-based infix-power rewrite (convert_infix_operator_to_ function_syntax / convert_infix_power_to_function_syntax) with an AST pass that converts power nodes to pow() before stringifying. It handles both power node types -- AST_POWER (L3 `^`/`pow`) and AST_FUNCTION_POWER (from MathML ``, used by every reference model) -- and, unlike the string rewrite, converts nested powers correctly (`a^b^c` no longer leaves an inner `^`, which was invalid C++). Also tidy formula_to_string: - hoist the constant/function lookup tables to module scope (were rebuilt on every call); - drop the misleading single-element `unsupported_functions` loop; - normalise a deepCopy so the caller's AST (and the libsbml model) is no longer mutated as a side effect. Regenerate the affected reference models: the only semantic change is pow((x), y) -> pow(x, y); the large line counts are clang-format re-aligning trailing comments on the now-shorter lines. Co-Authored-By: Claude Opus 4.8 --- .../Chen2000/Chen2000SbmlOdeSystem.cpp | 282 ++++++------- .../Chen2004/Chen2004SbmlOdeSystem.cpp | 2 +- .../Gardner1998/Gardner1998SbmlOdeSystem.cpp | 56 +-- .../Goldbeter1991SbmlOdeSystem.cpp | 42 +- .../TysonNovak2001SbmlOdeSystem.cpp | 220 +++++------ chaste_sbml/_expressions.py | 374 +++++++----------- chaste_sbml/tests/test_expressions.py | 66 ++-- 7 files changed, 473 insertions(+), 569 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/src/reference/Chen2000/Chen2000SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/Chen2000/Chen2000SbmlOdeSystem.cpp index ee1100fc..9681a281 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Chen2000/Chen2000SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Chen2000/Chen2000SbmlOdeSystem.cpp @@ -64,119 +64,119 @@ void Chen2000SbmlOdeSystem::EvaluateYDerivatives(double time, const std::vector< void Chen2000SbmlOdeSystem::Initialise(double time) { - COMpartment = 1.0; // - Cln2 = 0.0078; // - ks_n2 = 0.0; // - ks_n2_ = 0.05; // - kd_n2 = 0.1; // - Clb2_T = 0.2342; // - Hct1_T = 1.0; // - ks_b2 = 0.002; // - ks_b2_ = 0.05; // - kd_b2 = 0.01; // - kd_b2_ = 2.0; // - kd_b2__ = 0.05; // - Clb5_T = 0.0614; // - ks_b5 = 0.006; // - ks_b5_ = 0.02; // - kd_b5 = 0.1; // - kd_b5_ = 0.25; // - Bck2_0 = 0.0027; // - Jn3 = 6.0; // - Dn3 = 1.0; // - Cln3_max = 0.02; // - Sic1_T = 0.1231; // - ks_c1 = 0.02; // - ks_c1_ = 0.1; // - Clb2_Sic1 = 0.079; // - kas_b2 = 50.0; // - kdi_b2 = 0.05; // - Clb5_Sic1 = 0.0207; // - kas_b5 = 50.0; // - kdi_b5 = 0.05; // - kd2_c1 = 0.3; // - epsilonc1_n3 = 20.0; // - epsilonc1_k2 = 2.0; // - epsilonc1_b5 = 1.0; // - epsilonc1_b2 = 0.067; // - Cdc20_T = 0.8332; // - ks_20 = 0.005; // - ks_20_ = 0.06; // - Cdc20 = 0.6848; // - ka_20 = 1.0; // - ki_20 = 0.1; // - ki_20_ = 10.0; // - Hct1 = 0.9946; // - ka_t1 = 0.04; // - ka_t1_ = 2.0; // - ki_t1 = 0.0; // - ki_t1_ = 0.64; // - Ji_t1 = 0.05; // - Ja_t1 = 0.05; // - epsiloni_t1_n2 = 1.0; // - epsiloni_t1_b5 = 0.5; // - epsiloni_t1_b2 = 1.0; // - mass = 0.6608; // - mu = 0.005776; // - ORI = 0.0; // - ks_ori = 2.0; // - kd_ori = 0.06; // - epsilonori_b2 = 0.4; // - BUD = 0.0; // - ks_bud = 0.3; // - kd_bud = 0.06; // - epsilonbud_b5 = 1.0; // - SPN = 0.0; // - ks_spn = 0.08; // - kd_spn = 0.06; // - J_spn = 0.2; // - ka_sbf = 1.0; // - ki_sbf = 0.5; // - ki_sbf_ = 6.0; // - Ji_sbf = 0.01; // - Ja_sbf = 0.01; // - epsilonsbf_n3 = 75.0; // - epsilonsbf_b5 = 0.5; // - ka_mcm = 1.0; // - ki_mcm = 0.15; // - Ji_mcm = 1.0; // - Ja_mcm = 1.0; // - ka_swi = 1.0; // - ki_swi = 0.3; // - ki_swi_ = 0.2; // - Ji_swi = 0.1; // - Ja_swi = 0.1; // - kd1_c1 = 0.01; // - kd_20 = 0.08; // - Jd2_c1 = 0.05; // - Vd_b2 = kd_b2 * (Hct1_T - Hct1) + kd_b2_ * Hct1 + kd_b2__ * Cdc20; // - Clb2 = Clb2_T - Clb2_Sic1; // - Clb5 = Clb5_T - Clb5_Sic1; // - Sic1 = Sic1_T - (Clb2_Sic1 + Clb5_Sic1); // - Vd_b5 = kd_b5 + kd_b5_ * Cdc20; // - Bck2 = Bck2_0 * mass; // - Cln3 = Cln3_max * Dn3 * mass / (Jn3 + Dn3 * mass); // - Vd2_c1 = kd2_c1 * (epsilonc1_n3 * Cln3 + epsilonc1_k2 * Bck2 + Cln2 + epsilonc1_b5 * Clb5 + epsilonc1_b2 * Clb2); // - Vi_20 = sm::piecewise(ki_20_, ORI >= 1.0, ki_20, SPN >= 1.0, 0.1); // - Vi_t1 = ki_t1 + ki_t1_ * (Cln3 + epsiloni_t1_n2 * Cln2 + epsiloni_t1_b5 * Clb5 + epsiloni_t1_b2 * Clb2); // - Va_sbf = ka_sbf * (Cln2 + epsilonsbf_n3 * (Cln3 + Bck2) + epsilonsbf_b5 * Clb5); // - SBF = 2.0 * Va_sbf * Ji_sbf / ((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf) + std::sqrt(std::pow((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf), 2.0) - 4.0 * Va_sbf * Ji_sbf * (ki_sbf + ki_sbf_ * Clb2 - Va_sbf))); // - MBF = SBF; // - Mcm1 = 2.0 * ka_mcm * Clb2 * Ji_mcm / ((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2) + std::sqrt(std::pow((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2), 2.0) - 4.0 * (ki_mcm - ka_mcm * Clb2) * ka_mcm * Clb2 * Ji_mcm)); // - Swi5 = 2.0 * ka_swi * Cdc20 * Ji_swi / ((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20) + std::sqrt(std::pow((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20), 2.0) - 4.0 * (ki_swi + ki_swi_ * Clb2 - ka_swi * Cdc20) * ka_swi * Cdc20 * Ji_swi)); // - d_Cln2_dt = mass * (ks_n2 + ks_n2_ * SBF) - kd_n2 * Cln2; // - d_Clb2_T_dt = mass * (ks_b2 + ks_b2_ * Mcm1) - Vd_b2 * Clb2_T; // - d_Clb5_T_dt = mass * (ks_b5 + ks_b5_ * MBF) - Vd_b5 * Clb5_T; // - d_Sic1_T_dt = ks_c1 + ks_c1_ * Swi5 - Sic1_T * (kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Clb2_Sic1_dt = kas_b2 * Clb2 * Sic1 - Clb2_Sic1 * (kdi_b2 + Vd_b2 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Clb5_Sic1_dt = kas_b5 * Clb5 * Sic1 - Clb5_Sic1 * (kdi_b5 + Vd_b5 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Cdc20_T_dt = ks_20 + ks_20_ * Clb2 - kd_20 * Cdc20_T; // - d_Cdc20_dt = ka_20 * (Cdc20_T - Cdc20) - Cdc20 * (Vi_20 + kd_20); // - d_Hct1_dt = (ka_t1 + ka_t1_ * Cdc20) * (Hct1_T - Hct1) / (Ja_t1 + Hct1_T - Hct1) - Vi_t1 * Hct1 / (Ji_t1 + Hct1); // - d_mass_dt = mu * mass; // - d_ORI_dt = ks_ori * (Clb5 + epsilonori_b2 * Clb2) - kd_ori * ORI; // - d_BUD_dt = ks_bud * (Cln2 + Cln3 + epsilonbud_b5 * Clb5) - kd_bud * BUD; // - d_SPN_dt = ks_spn * Clb2 / (J_spn + Clb2) - kd_spn * SPN; // + COMpartment = 1.0; // + Cln2 = 0.0078; // + ks_n2 = 0.0; // + ks_n2_ = 0.05; // + kd_n2 = 0.1; // + Clb2_T = 0.2342; // + Hct1_T = 1.0; // + ks_b2 = 0.002; // + ks_b2_ = 0.05; // + kd_b2 = 0.01; // + kd_b2_ = 2.0; // + kd_b2__ = 0.05; // + Clb5_T = 0.0614; // + ks_b5 = 0.006; // + ks_b5_ = 0.02; // + kd_b5 = 0.1; // + kd_b5_ = 0.25; // + Bck2_0 = 0.0027; // + Jn3 = 6.0; // + Dn3 = 1.0; // + Cln3_max = 0.02; // + Sic1_T = 0.1231; // + ks_c1 = 0.02; // + ks_c1_ = 0.1; // + Clb2_Sic1 = 0.079; // + kas_b2 = 50.0; // + kdi_b2 = 0.05; // + Clb5_Sic1 = 0.0207; // + kas_b5 = 50.0; // + kdi_b5 = 0.05; // + kd2_c1 = 0.3; // + epsilonc1_n3 = 20.0; // + epsilonc1_k2 = 2.0; // + epsilonc1_b5 = 1.0; // + epsilonc1_b2 = 0.067; // + Cdc20_T = 0.8332; // + ks_20 = 0.005; // + ks_20_ = 0.06; // + Cdc20 = 0.6848; // + ka_20 = 1.0; // + ki_20 = 0.1; // + ki_20_ = 10.0; // + Hct1 = 0.9946; // + ka_t1 = 0.04; // + ka_t1_ = 2.0; // + ki_t1 = 0.0; // + ki_t1_ = 0.64; // + Ji_t1 = 0.05; // + Ja_t1 = 0.05; // + epsiloni_t1_n2 = 1.0; // + epsiloni_t1_b5 = 0.5; // + epsiloni_t1_b2 = 1.0; // + mass = 0.6608; // + mu = 0.005776; // + ORI = 0.0; // + ks_ori = 2.0; // + kd_ori = 0.06; // + epsilonori_b2 = 0.4; // + BUD = 0.0; // + ks_bud = 0.3; // + kd_bud = 0.06; // + epsilonbud_b5 = 1.0; // + SPN = 0.0; // + ks_spn = 0.08; // + kd_spn = 0.06; // + J_spn = 0.2; // + ka_sbf = 1.0; // + ki_sbf = 0.5; // + ki_sbf_ = 6.0; // + Ji_sbf = 0.01; // + Ja_sbf = 0.01; // + epsilonsbf_n3 = 75.0; // + epsilonsbf_b5 = 0.5; // + ka_mcm = 1.0; // + ki_mcm = 0.15; // + Ji_mcm = 1.0; // + Ja_mcm = 1.0; // + ka_swi = 1.0; // + ki_swi = 0.3; // + ki_swi_ = 0.2; // + Ji_swi = 0.1; // + Ja_swi = 0.1; // + kd1_c1 = 0.01; // + kd_20 = 0.08; // + Jd2_c1 = 0.05; // + Vd_b2 = kd_b2 * (Hct1_T - Hct1) + kd_b2_ * Hct1 + kd_b2__ * Cdc20; // + Clb2 = Clb2_T - Clb2_Sic1; // + Clb5 = Clb5_T - Clb5_Sic1; // + Sic1 = Sic1_T - (Clb2_Sic1 + Clb5_Sic1); // + Vd_b5 = kd_b5 + kd_b5_ * Cdc20; // + Bck2 = Bck2_0 * mass; // + Cln3 = Cln3_max * Dn3 * mass / (Jn3 + Dn3 * mass); // + Vd2_c1 = kd2_c1 * (epsilonc1_n3 * Cln3 + epsilonc1_k2 * Bck2 + Cln2 + epsilonc1_b5 * Clb5 + epsilonc1_b2 * Clb2); // + Vi_20 = sm::piecewise(ki_20_, ORI >= 1.0, ki_20, SPN >= 1.0, 0.1); // + Vi_t1 = ki_t1 + ki_t1_ * (Cln3 + epsiloni_t1_n2 * Cln2 + epsiloni_t1_b5 * Clb5 + epsiloni_t1_b2 * Clb2); // + Va_sbf = ka_sbf * (Cln2 + epsilonsbf_n3 * (Cln3 + Bck2) + epsilonsbf_b5 * Clb5); // + SBF = 2.0 * Va_sbf * Ji_sbf / ((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf) + std::sqrt(std::pow(ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf, 2.0) - 4.0 * Va_sbf * Ji_sbf * (ki_sbf + ki_sbf_ * Clb2 - Va_sbf))); // + MBF = SBF; // + Mcm1 = 2.0 * ka_mcm * Clb2 * Ji_mcm / ((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2) + std::sqrt(std::pow(ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2, 2.0) - 4.0 * (ki_mcm - ka_mcm * Clb2) * ka_mcm * Clb2 * Ji_mcm)); // + Swi5 = 2.0 * ka_swi * Cdc20 * Ji_swi / ((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20) + std::sqrt(std::pow(ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20, 2.0) - 4.0 * (ki_swi + ki_swi_ * Clb2 - ka_swi * Cdc20) * ka_swi * Cdc20 * Ji_swi)); // + d_Cln2_dt = mass * (ks_n2 + ks_n2_ * SBF) - kd_n2 * Cln2; // + d_Clb2_T_dt = mass * (ks_b2 + ks_b2_ * Mcm1) - Vd_b2 * Clb2_T; // + d_Clb5_T_dt = mass * (ks_b5 + ks_b5_ * MBF) - Vd_b5 * Clb5_T; // + d_Sic1_T_dt = ks_c1 + ks_c1_ * Swi5 - Sic1_T * (kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Clb2_Sic1_dt = kas_b2 * Clb2 * Sic1 - Clb2_Sic1 * (kdi_b2 + Vd_b2 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Clb5_Sic1_dt = kas_b5 * Clb5 * Sic1 - Clb5_Sic1 * (kdi_b5 + Vd_b5 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Cdc20_T_dt = ks_20 + ks_20_ * Clb2 - kd_20 * Cdc20_T; // + d_Cdc20_dt = ka_20 * (Cdc20_T - Cdc20) - Cdc20 * (Vi_20 + kd_20); // + d_Hct1_dt = (ka_t1 + ka_t1_ * Cdc20) * (Hct1_T - Hct1) / (Ja_t1 + Hct1_T - Hct1) - Vi_t1 * Hct1 / (Ji_t1 + Hct1); // + d_mass_dt = mu * mass; // + d_ORI_dt = ks_ori * (Clb5 + epsilonori_b2 * Clb2) - kd_ori * ORI; // + d_BUD_dt = ks_bud * (Cln2 + Cln3 + epsilonbud_b5 * Clb5) - kd_bud * BUD; // + d_SPN_dt = ks_spn * Clb2 / (J_spn + Clb2) - kd_spn * SPN; // mStateVariables.push_back(Cln2); mStateVariables.push_back(Clb2_T); @@ -392,34 +392,34 @@ std::vector Chen2000SbmlOdeSystem::RunModelEquations(double time, const kd_20 = GetParameter(69); Jd2_c1 = GetParameter(70); - Vd_b2 = kd_b2 * (Hct1_T - Hct1) + kd_b2_ * Hct1 + kd_b2__ * Cdc20; // - Clb2 = Clb2_T - Clb2_Sic1; // - Clb5 = Clb5_T - Clb5_Sic1; // - Sic1 = Sic1_T - (Clb2_Sic1 + Clb5_Sic1); // - Vd_b5 = kd_b5 + kd_b5_ * Cdc20; // - Bck2 = Bck2_0 * mass; // - Cln3 = Cln3_max * Dn3 * mass / (Jn3 + Dn3 * mass); // - Vd2_c1 = kd2_c1 * (epsilonc1_n3 * Cln3 + epsilonc1_k2 * Bck2 + Cln2 + epsilonc1_b5 * Clb5 + epsilonc1_b2 * Clb2); // - Vi_20 = sm::piecewise(ki_20_, ORI >= 1.0, ki_20, SPN >= 1.0, 0.1); // - Vi_t1 = ki_t1 + ki_t1_ * (Cln3 + epsiloni_t1_n2 * Cln2 + epsiloni_t1_b5 * Clb5 + epsiloni_t1_b2 * Clb2); // - Va_sbf = ka_sbf * (Cln2 + epsilonsbf_n3 * (Cln3 + Bck2) + epsilonsbf_b5 * Clb5); // - SBF = 2.0 * Va_sbf * Ji_sbf / ((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf) + std::sqrt(std::pow((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf), 2.0) - 4.0 * Va_sbf * Ji_sbf * (ki_sbf + ki_sbf_ * Clb2 - Va_sbf))); // - MBF = SBF; // - Mcm1 = 2.0 * ka_mcm * Clb2 * Ji_mcm / ((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2) + std::sqrt(std::pow((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2), 2.0) - 4.0 * (ki_mcm - ka_mcm * Clb2) * ka_mcm * Clb2 * Ji_mcm)); // - Swi5 = 2.0 * ka_swi * Cdc20 * Ji_swi / ((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20) + std::sqrt(std::pow((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20), 2.0) - 4.0 * (ki_swi + ki_swi_ * Clb2 - ka_swi * Cdc20) * ka_swi * Cdc20 * Ji_swi)); // - d_Cln2_dt = mass * (ks_n2 + ks_n2_ * SBF) - kd_n2 * Cln2; // - d_Clb2_T_dt = mass * (ks_b2 + ks_b2_ * Mcm1) - Vd_b2 * Clb2_T; // - d_Clb5_T_dt = mass * (ks_b5 + ks_b5_ * MBF) - Vd_b5 * Clb5_T; // - d_Sic1_T_dt = ks_c1 + ks_c1_ * Swi5 - Sic1_T * (kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Clb2_Sic1_dt = kas_b2 * Clb2 * Sic1 - Clb2_Sic1 * (kdi_b2 + Vd_b2 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Clb5_Sic1_dt = kas_b5 * Clb5 * Sic1 - Clb5_Sic1 * (kdi_b5 + Vd_b5 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // - d_Cdc20_T_dt = ks_20 + ks_20_ * Clb2 - kd_20 * Cdc20_T; // - d_Cdc20_dt = ka_20 * (Cdc20_T - Cdc20) - Cdc20 * (Vi_20 + kd_20); // - d_Hct1_dt = (ka_t1 + ka_t1_ * Cdc20) * (Hct1_T - Hct1) / (Ja_t1 + Hct1_T - Hct1) - Vi_t1 * Hct1 / (Ji_t1 + Hct1); // - d_mass_dt = mu * mass; // - d_ORI_dt = ks_ori * (Clb5 + epsilonori_b2 * Clb2) - kd_ori * ORI; // - d_BUD_dt = ks_bud * (Cln2 + Cln3 + epsilonbud_b5 * Clb5) - kd_bud * BUD; // - d_SPN_dt = ks_spn * Clb2 / (J_spn + Clb2) - kd_spn * SPN; // + Vd_b2 = kd_b2 * (Hct1_T - Hct1) + kd_b2_ * Hct1 + kd_b2__ * Cdc20; // + Clb2 = Clb2_T - Clb2_Sic1; // + Clb5 = Clb5_T - Clb5_Sic1; // + Sic1 = Sic1_T - (Clb2_Sic1 + Clb5_Sic1); // + Vd_b5 = kd_b5 + kd_b5_ * Cdc20; // + Bck2 = Bck2_0 * mass; // + Cln3 = Cln3_max * Dn3 * mass / (Jn3 + Dn3 * mass); // + Vd2_c1 = kd2_c1 * (epsilonc1_n3 * Cln3 + epsilonc1_k2 * Bck2 + Cln2 + epsilonc1_b5 * Clb5 + epsilonc1_b2 * Clb2); // + Vi_20 = sm::piecewise(ki_20_, ORI >= 1.0, ki_20, SPN >= 1.0, 0.1); // + Vi_t1 = ki_t1 + ki_t1_ * (Cln3 + epsiloni_t1_n2 * Cln2 + epsiloni_t1_b5 * Clb5 + epsiloni_t1_b2 * Clb2); // + Va_sbf = ka_sbf * (Cln2 + epsilonsbf_n3 * (Cln3 + Bck2) + epsilonsbf_b5 * Clb5); // + SBF = 2.0 * Va_sbf * Ji_sbf / ((ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf) + std::sqrt(std::pow(ki_sbf + ki_sbf_ * Clb2 + Va_sbf * Ji_sbf + (ki_sbf + ki_sbf_ * Clb2) * Ja_sbf - Va_sbf, 2.0) - 4.0 * Va_sbf * Ji_sbf * (ki_sbf + ki_sbf_ * Clb2 - Va_sbf))); // + MBF = SBF; // + Mcm1 = 2.0 * ka_mcm * Clb2 * Ji_mcm / ((ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2) + std::sqrt(std::pow(ki_mcm + ka_mcm * Clb2 * Ji_mcm + ki_mcm * Ja_mcm - ka_mcm * Clb2, 2.0) - 4.0 * (ki_mcm - ka_mcm * Clb2) * ka_mcm * Clb2 * Ji_mcm)); // + Swi5 = 2.0 * ka_swi * Cdc20 * Ji_swi / ((ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20) + std::sqrt(std::pow(ki_swi + ki_swi_ * Clb2 + ka_swi * Cdc20 * Ji_swi + (ki_swi + ki_swi_ * Clb2) * Ja_swi - ka_swi * Cdc20, 2.0) - 4.0 * (ki_swi + ki_swi_ * Clb2 - ka_swi * Cdc20) * ka_swi * Cdc20 * Ji_swi)); // + d_Cln2_dt = mass * (ks_n2 + ks_n2_ * SBF) - kd_n2 * Cln2; // + d_Clb2_T_dt = mass * (ks_b2 + ks_b2_ * Mcm1) - Vd_b2 * Clb2_T; // + d_Clb5_T_dt = mass * (ks_b5 + ks_b5_ * MBF) - Vd_b5 * Clb5_T; // + d_Sic1_T_dt = ks_c1 + ks_c1_ * Swi5 - Sic1_T * (kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Clb2_Sic1_dt = kas_b2 * Clb2 * Sic1 - Clb2_Sic1 * (kdi_b2 + Vd_b2 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Clb5_Sic1_dt = kas_b5 * Clb5 * Sic1 - Clb5_Sic1 * (kdi_b5 + Vd_b5 + kd1_c1 + Vd2_c1 / (Jd2_c1 + Sic1_T)); // + d_Cdc20_T_dt = ks_20 + ks_20_ * Clb2 - kd_20 * Cdc20_T; // + d_Cdc20_dt = ka_20 * (Cdc20_T - Cdc20) - Cdc20 * (Vi_20 + kd_20); // + d_Hct1_dt = (ka_t1 + ka_t1_ * Cdc20) * (Hct1_T - Hct1) / (Ja_t1 + Hct1_T - Hct1) - Vi_t1 * Hct1 / (Ji_t1 + Hct1); // + d_mass_dt = mu * mass; // + d_ORI_dt = ks_ori * (Clb5 + epsilonori_b2 * Clb2) - kd_ori * ORI; // + d_BUD_dt = ks_bud * (Cln2 + Cln3 + epsilonbud_b5 * Clb5) - kd_bud * BUD; // + d_SPN_dt = ks_spn * Clb2 / (J_spn + Clb2) - kd_spn * SPN; // std::vector derivatives(13); derivatives[0] = d_Cln2_dt; diff --git a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp index 9c1b0035..cdc4aae3 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp @@ -1682,7 +1682,7 @@ inline double Chen2004SbmlOdeSystem::BB_218(double A1, double A2, double A3, dou inline double Chen2004SbmlOdeSystem::GK_219(double A1, double A2, double A3, double A4) { - return 2.0 * A4 * A1 / ((A2 - A1) + A3 * A2 + A4 * A1 + sm::root(2.0, std::pow(((A2 - A1) + A3 * A2 + A4 * A1), 2.0) - 4.0 * (A2 - A1) * A4 * A1)); + return 2.0 * A4 * A1 / ((A2 - A1) + A3 * A2 + A4 * A1 + sm::root(2.0, std::pow((A2 - A1) + A3 * A2 + A4 * A1, 2.0) - 4.0 * (A2 - A1) * A4 * A1)); } inline double Chen2004SbmlOdeSystem::MichaelisMenten_220(double M1, double J1, double k1, double S1) diff --git a/chaste_sbml/SbmlRefModels/src/reference/Gardner1998/Gardner1998SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/Gardner1998/Gardner1998SbmlOdeSystem.cpp index 2bee481e..78bf29b0 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Gardner1998/Gardner1998SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Gardner1998/Gardner1998SbmlOdeSystem.cpp @@ -74,22 +74,22 @@ void Gardner1998SbmlOdeSystem::EvaluateYDerivatives(double time, const std::vect void Gardner1998SbmlOdeSystem::Initialise(double time) { - Cell = 1.0; // - C = 0.0; // - X = 0.0; // - M = 0.0; // - Y = 1.0; // - Z = 1.0; // - K6 = 0.3; // - V1p = 0.75; // - V3p = 0.3; // - C = C / Cell; // - X = X / Cell; // - M = M / Cell; // - Y = Y / Cell; // - Z = Z / Cell; // - V1 = C * V1p * std::pow((C + K6), -1.0); // - V3 = M * V3p; // + Cell = 1.0; // + C = 0.0; // + X = 0.0; // + M = 0.0; // + Y = 1.0; // + Z = 1.0; // + K6 = 0.3; // + V1p = 0.75; // + V3p = 0.3; // + C = C / Cell; // + X = X / Cell; // + M = M / Cell; // + Y = Y / Cell; // + Z = Z / Cell; // + V1 = C * V1p * std::pow(C + K6, -1.0); // + V3 = M * V3p; // // reaction1: { [[maybe_unused]] double vi = 0.1; @@ -105,7 +105,7 @@ void Gardner1998SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction2 = C * k1 * X * std::pow((C + K5), -1.0); + this->reaction2 = C * k1 * X * std::pow(C + K5, -1.0); } // reaction3: { @@ -121,7 +121,7 @@ void Gardner1998SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction4 = (1.0 + -1.0 * M) * V1 * std::pow((K1 + -1.0 * M + 1.0), -1.0); + this->reaction4 = (1.0 + -1.0 * M) * V1 * std::pow(K1 + -1.0 * M + 1.0, -1.0); } // reaction5: { @@ -130,7 +130,7 @@ void Gardner1998SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction5 = M * V2 * std::pow((K2 + M), -1.0); + this->reaction5 = M * V2 * std::pow(K2 + M, -1.0); } // reaction6: { @@ -138,7 +138,7 @@ void Gardner1998SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction6 = V3 * (1.0 + -1.0 * X) * std::pow((K3 + -1.0 * X + 1.0), -1.0); + this->reaction6 = V3 * (1.0 + -1.0 * X) * std::pow(K3 + -1.0 * X + 1.0, -1.0); } // reaction7: { @@ -147,7 +147,7 @@ void Gardner1998SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction7 = V4 * X * std::pow((K4 + X), -1.0); + this->reaction7 = V4 * X * std::pow(K4 + X, -1.0); } // reaction8: { @@ -259,8 +259,8 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con V1p = GetParameter(1); V3p = GetParameter(2); - V1 = C * V1p * std::pow((C + K6), -1.0); // - V3 = M * V3p; // + V1 = C * V1p * std::pow(C + K6, -1.0); // + V3 = M * V3p; // // reaction1: { [[maybe_unused]] double vi = 0.1; @@ -276,7 +276,7 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction2 = C * k1 * X * std::pow((C + K5), -1.0); + this->reaction2 = C * k1 * X * std::pow(C + K5, -1.0); } // reaction3: { @@ -292,7 +292,7 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction4 = (1.0 + -1.0 * M) * V1 * std::pow((K1 + -1.0 * M + 1.0), -1.0); + this->reaction4 = (1.0 + -1.0 * M) * V1 * std::pow(K1 + -1.0 * M + 1.0, -1.0); } // reaction5: { @@ -301,7 +301,7 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction5 = M * V2 * std::pow((K2 + M), -1.0); + this->reaction5 = M * V2 * std::pow(K2 + M, -1.0); } // reaction6: { @@ -309,7 +309,7 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction6 = V3 * (1.0 + -1.0 * X) * std::pow((K3 + -1.0 * X + 1.0), -1.0); + this->reaction6 = V3 * (1.0 + -1.0 * X) * std::pow(K3 + -1.0 * X + 1.0, -1.0); } // reaction7: { @@ -318,7 +318,7 @@ std::vector Gardner1998SbmlOdeSystem::RunModelEquations(double time, con // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction7 = V4 * X * std::pow((K4 + X), -1.0); + this->reaction7 = V4 * X * std::pow(K4 + X, -1.0); } // reaction8: { diff --git a/chaste_sbml/SbmlRefModels/src/reference/Goldbeter1991/Goldbeter1991SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/Goldbeter1991/Goldbeter1991SbmlOdeSystem.cpp index 4cc43c96..65deeb9a 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Goldbeter1991/Goldbeter1991SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Goldbeter1991/Goldbeter1991SbmlOdeSystem.cpp @@ -64,15 +64,15 @@ void Goldbeter1991SbmlOdeSystem::EvaluateYDerivatives(double time, const std::ve void Goldbeter1991SbmlOdeSystem::Initialise(double time) { - cell = 1.0; // - C = 0.01; // - M = 0.01; // - X = 0.01; // - VM1 = 3.0; // - VM3 = 1.0; // - Kc = 0.5; // - V1 = C * VM1 * std::pow((C + Kc), -1.0); // - V3 = M * VM3; // + cell = 1.0; // + C = 0.01; // + M = 0.01; // + X = 0.01; // + VM1 = 3.0; // + VM3 = 1.0; // + Kc = 0.5; // + V1 = C * VM1 * std::pow(C + Kc, -1.0); // + V3 = M * VM3; // // reaction1: { [[maybe_unused]] double vi = 0.025; @@ -96,7 +96,7 @@ void Goldbeter1991SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction3 = C * cell * vd * X * std::pow((C + Kd), -1.0); + this->reaction3 = C * cell * vd * X * std::pow(C + Kd, -1.0); } // reaction4: { @@ -104,7 +104,7 @@ void Goldbeter1991SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction4 = cell * (1.0 + -1.0 * M) * V1 * std::pow((K1 + -1.0 * M + 1.0), -1.0); + this->reaction4 = cell * (1.0 + -1.0 * M) * V1 * std::pow(K1 + -1.0 * M + 1.0, -1.0); } // reaction5: { @@ -113,7 +113,7 @@ void Goldbeter1991SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction5 = cell * M * V2 * std::pow((K2 + M), -1.0); + this->reaction5 = cell * M * V2 * std::pow(K2 + M, -1.0); } // reaction6: { @@ -121,7 +121,7 @@ void Goldbeter1991SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction6 = cell * V3 * (1.0 + -1.0 * X) * std::pow((K3 + -1.0 * X + 1.0), -1.0); + this->reaction6 = cell * V3 * (1.0 + -1.0 * X) * std::pow(K3 + -1.0 * X + 1.0, -1.0); } // reaction7: { @@ -130,7 +130,7 @@ void Goldbeter1991SbmlOdeSystem::Initialise(double time) // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction7 = cell * V4 * X * std::pow((K4 + X), -1.0); + this->reaction7 = cell * V4 * X * std::pow(K4 + X, -1.0); } d_C_dt = (reaction1 - reaction2 - reaction3) / cell; // d_M_dt = (reaction4 - reaction5) / cell; // @@ -184,8 +184,8 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c VM3 = GetParameter(1); Kc = GetParameter(2); - V1 = C * VM1 * std::pow((C + Kc), -1.0); // - V3 = M * VM3; // + V1 = C * VM1 * std::pow(C + Kc, -1.0); // + V3 = M * VM3; // // reaction1: { [[maybe_unused]] double vi = 0.025; @@ -209,7 +209,7 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction3 = C * cell * vd * X * std::pow((C + Kd), -1.0); + this->reaction3 = C * cell * vd * X * std::pow(C + Kd, -1.0); } // reaction4: { @@ -217,7 +217,7 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction4 = cell * (1.0 + -1.0 * M) * V1 * std::pow((K1 + -1.0 * M + 1.0), -1.0); + this->reaction4 = cell * (1.0 + -1.0 * M) * V1 * std::pow(K1 + -1.0 * M + 1.0, -1.0); } // reaction5: { @@ -226,7 +226,7 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction5 = cell * M * V2 * std::pow((K2 + M), -1.0); + this->reaction5 = cell * M * V2 * std::pow(K2 + M, -1.0); } // reaction6: { @@ -234,7 +234,7 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction6 = cell * V3 * (1.0 + -1.0 * X) * std::pow((K3 + -1.0 * X + 1.0), -1.0); + this->reaction6 = cell * V3 * (1.0 + -1.0 * X) * std::pow(K3 + -1.0 * X + 1.0, -1.0); } // reaction7: { @@ -243,7 +243,7 @@ std::vector Goldbeter1991SbmlOdeSystem::RunModelEquations(double time, c // Qualify with this-> so the assignment targets the reaction member even when a local // parameter shadows its name (an SBML local parameter may shadow the reaction ID, e.g. a // reaction J1 with a local parameter also named J1). - this->reaction7 = cell * V4 * X * std::pow((K4 + X), -1.0); + this->reaction7 = cell * V4 * X * std::pow(K4 + X, -1.0); } d_C_dt = (reaction1 - reaction2 - reaction3) / cell; // d_M_dt = (reaction4 - reaction5) / cell; // diff --git a/chaste_sbml/SbmlRefModels/src/reference/TysonNovak2001/TysonNovak2001SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/TysonNovak2001/TysonNovak2001SbmlOdeSystem.cpp index ba1391b3..cb94a31f 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/TysonNovak2001/TysonNovak2001SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/TysonNovak2001/TysonNovak2001SbmlOdeSystem.cpp @@ -109,83 +109,83 @@ void TysonNovak2001SbmlOdeSystem::EvaluateYDerivatives(double time, const std::v void TysonNovak2001SbmlOdeSystem::Initialise(double time) { - cell = 1.0; // - CycBt = 0.001; // - Cdc20a = 0.001; // - Cdh1 = 0.001; // - m = 0.5; // - Cdc20t = 0.001; // - IEP = 0.001; // - CKIt = 0.001; // - SK = 0.001; // - k1 = 0.04; // - k2p = 0.04; // - k2pp = 1.0; // - k2ppp = 1.0; // - k3p = 1.0; // - k3pp = 10.0; // - J3 = 0.04; // - k4 = 35.0; // - k5p = 0.005; // - k5pp = 0.2; // - J5 = 0.3; // - k6 = 0.1; // - n = 4.0; // - k7 = 1.0; // - J7 = 0.001; // - k8 = 0.5; // - J8 = 0.001; // - k9 = 0.1; // - k10 = 0.02; // - mu = 0.005; // - k11 = 1.0; // - k12p = 0.2; // - k12pp = 50.0; // - mmax = 10.0; // - k12ppp = 100.0; // - Keq = 1000.0; // - k13 = 1.0; // - k14 = 1.0; // - k15p = 1.5; // - k15pp = 0.05; // - k16p = 1.0; // - k16pp = 3.0; // - J15 = 0.01; // - J16 = 0.01; // - k4p = 2.0; // - J4 = 0.04; // - CycB = CycBt - 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow((std::pow((CycBt + CKIt + 1.0 / Keq), 2.0) - 4.0 * CycBt * CKIt), (1.0 / 2.0))); // - Trimer = 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow((std::pow((CycBt + CKIt + 1.0 / Keq), 2.0) - 4.0 * CycBt * CKIt), (1.0 / 2.0))); // - Mad = 1.0; // - TF = GK(k15p * m + k15pp * SK, k16p + k16pp * m * CycB, J15, J16); // - CycBt_synthesis = k1; // - CycBdegradation = k2p * CycBt; // - CycBdegradationviaCdh1 = k2pp * Cdh1 * CycBt; // - CycBtdegradationviaCdc20a = k2ppp * Cdc20a * CycBt; // - Cdh1synthesis = (k3p + k3pp * Cdc20a) * (1.0 - Cdh1) / (J3 + 1.0 - Cdh1); // - Cdh1degradation = (k4p * SK * Cdh1 + k4 * m * CycB * Cdh1) / (J4 + Cdh1); // - Cdc20tsynthesis = k5p + k5pp * std::pow((CycB * m / J5), n) / (1.0 + std::pow((CycB * m / J5), n)); // - Cdc20t_deg = k6 * Cdc20t; // - Cdc20activation = k7 * IEP * (Cdc20t - Cdc20a) / (J7 + Cdc20t - Cdc20a); // - Cdc20ainhibition = k8 * Mad * Cdc20a / (J8 + Cdc20a); // - Cdc20adegradation = k6 * Cdc20a; // - IEPsynthesis = k9 * m * CycB * (1.0 - IEP); // - IEPdegradation = k10 * IEP; // - growth = mu * m * (1.0 - m / mmax); // - CKItsynthesis = k11; // - CKIdegradation = k12p * CKIt; // - CKItphosphorilationviaSK = k12pp * SK * CKIt; // - eq_7 = k12ppp * m * CycB * CKIt; // - SKsynthesis = k13 * TF; // - SKdegradation = k14 * SK; // - d_CycBt_dt = CycBt_synthesis - CycBdegradation - CycBdegradationviaCdh1 - CycBtdegradationviaCdc20a; // - d_Cdc20a_dt = Cdc20activation - Cdc20ainhibition - Cdc20adegradation; // - d_Cdh1_dt = Cdh1synthesis - Cdh1degradation; // - d_m_dt = growth; // - d_Cdc20t_dt = Cdc20tsynthesis - Cdc20t_deg; // - d_IEP_dt = IEPsynthesis - IEPdegradation; // - d_CKIt_dt = CKItsynthesis - CKIdegradation - CKItphosphorilationviaSK - eq_7; // - d_SK_dt = SKsynthesis - SKdegradation; // + cell = 1.0; // + CycBt = 0.001; // + Cdc20a = 0.001; // + Cdh1 = 0.001; // + m = 0.5; // + Cdc20t = 0.001; // + IEP = 0.001; // + CKIt = 0.001; // + SK = 0.001; // + k1 = 0.04; // + k2p = 0.04; // + k2pp = 1.0; // + k2ppp = 1.0; // + k3p = 1.0; // + k3pp = 10.0; // + J3 = 0.04; // + k4 = 35.0; // + k5p = 0.005; // + k5pp = 0.2; // + J5 = 0.3; // + k6 = 0.1; // + n = 4.0; // + k7 = 1.0; // + J7 = 0.001; // + k8 = 0.5; // + J8 = 0.001; // + k9 = 0.1; // + k10 = 0.02; // + mu = 0.005; // + k11 = 1.0; // + k12p = 0.2; // + k12pp = 50.0; // + mmax = 10.0; // + k12ppp = 100.0; // + Keq = 1000.0; // + k13 = 1.0; // + k14 = 1.0; // + k15p = 1.5; // + k15pp = 0.05; // + k16p = 1.0; // + k16pp = 3.0; // + J15 = 0.01; // + J16 = 0.01; // + k4p = 2.0; // + J4 = 0.04; // + CycB = CycBt - 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow(std::pow(CycBt + CKIt + 1.0 / Keq, 2.0) - 4.0 * CycBt * CKIt, 1.0 / 2.0)); // + Trimer = 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow(std::pow(CycBt + CKIt + 1.0 / Keq, 2.0) - 4.0 * CycBt * CKIt, 1.0 / 2.0)); // + Mad = 1.0; // + TF = GK(k15p * m + k15pp * SK, k16p + k16pp * m * CycB, J15, J16); // + CycBt_synthesis = k1; // + CycBdegradation = k2p * CycBt; // + CycBdegradationviaCdh1 = k2pp * Cdh1 * CycBt; // + CycBtdegradationviaCdc20a = k2ppp * Cdc20a * CycBt; // + Cdh1synthesis = (k3p + k3pp * Cdc20a) * (1.0 - Cdh1) / (J3 + 1.0 - Cdh1); // + Cdh1degradation = (k4p * SK * Cdh1 + k4 * m * CycB * Cdh1) / (J4 + Cdh1); // + Cdc20tsynthesis = k5p + k5pp * std::pow(CycB * m / J5, n) / (1.0 + std::pow(CycB * m / J5, n)); // + Cdc20t_deg = k6 * Cdc20t; // + Cdc20activation = k7 * IEP * (Cdc20t - Cdc20a) / (J7 + Cdc20t - Cdc20a); // + Cdc20ainhibition = k8 * Mad * Cdc20a / (J8 + Cdc20a); // + Cdc20adegradation = k6 * Cdc20a; // + IEPsynthesis = k9 * m * CycB * (1.0 - IEP); // + IEPdegradation = k10 * IEP; // + growth = mu * m * (1.0 - m / mmax); // + CKItsynthesis = k11; // + CKIdegradation = k12p * CKIt; // + CKItphosphorilationviaSK = k12pp * SK * CKIt; // + eq_7 = k12ppp * m * CycB * CKIt; // + SKsynthesis = k13 * TF; // + SKdegradation = k14 * SK; // + d_CycBt_dt = CycBt_synthesis - CycBdegradation - CycBdegradationviaCdh1 - CycBtdegradationviaCdc20a; // + d_Cdc20a_dt = Cdc20activation - Cdc20ainhibition - Cdc20adegradation; // + d_Cdh1_dt = Cdh1synthesis - Cdh1degradation; // + d_m_dt = growth; // + d_Cdc20t_dt = Cdc20tsynthesis - Cdc20t_deg; // + d_IEP_dt = IEPsynthesis - IEPdegradation; // + d_CKIt_dt = CKItsynthesis - CKIdegradation - CKItphosphorilationviaSK - eq_7; // + d_SK_dt = SKsynthesis - SKdegradation; // mStateVariables.push_back(CycBt); mStateVariables.push_back(Cdc20a); @@ -406,38 +406,38 @@ std::vector TysonNovak2001SbmlOdeSystem::RunModelEquations(double time, k4p = GetParameter(34); J4 = GetParameter(35); - CycB = CycBt - 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow((std::pow((CycBt + CKIt + 1.0 / Keq), 2.0) - 4.0 * CycBt * CKIt), (1.0 / 2.0))); // - Trimer = 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow((std::pow((CycBt + CKIt + 1.0 / Keq), 2.0) - 4.0 * CycBt * CKIt), (1.0 / 2.0))); // - Mad = 1.0; // - TF = GK(k15p * m + k15pp * SK, k16p + k16pp * m * CycB, J15, J16); // - CycBt_synthesis = k1; // - CycBdegradation = k2p * CycBt; // - CycBdegradationviaCdh1 = k2pp * Cdh1 * CycBt; // - CycBtdegradationviaCdc20a = k2ppp * Cdc20a * CycBt; // - Cdh1synthesis = (k3p + k3pp * Cdc20a) * (1.0 - Cdh1) / (J3 + 1.0 - Cdh1); // - Cdh1degradation = (k4p * SK * Cdh1 + k4 * m * CycB * Cdh1) / (J4 + Cdh1); // - Cdc20tsynthesis = k5p + k5pp * std::pow((CycB * m / J5), n) / (1.0 + std::pow((CycB * m / J5), n)); // - Cdc20t_deg = k6 * Cdc20t; // - Cdc20activation = k7 * IEP * (Cdc20t - Cdc20a) / (J7 + Cdc20t - Cdc20a); // - Cdc20ainhibition = k8 * Mad * Cdc20a / (J8 + Cdc20a); // - Cdc20adegradation = k6 * Cdc20a; // - IEPsynthesis = k9 * m * CycB * (1.0 - IEP); // - IEPdegradation = k10 * IEP; // - growth = mu * m * (1.0 - m / mmax); // - CKItsynthesis = k11; // - CKIdegradation = k12p * CKIt; // - CKItphosphorilationviaSK = k12pp * SK * CKIt; // - eq_7 = k12ppp * m * CycB * CKIt; // - SKsynthesis = k13 * TF; // - SKdegradation = k14 * SK; // - d_CycBt_dt = CycBt_synthesis - CycBdegradation - CycBdegradationviaCdh1 - CycBtdegradationviaCdc20a; // - d_Cdc20a_dt = Cdc20activation - Cdc20ainhibition - Cdc20adegradation; // - d_Cdh1_dt = Cdh1synthesis - Cdh1degradation; // - d_m_dt = growth; // - d_Cdc20t_dt = Cdc20tsynthesis - Cdc20t_deg; // - d_IEP_dt = IEPsynthesis - IEPdegradation; // - d_CKIt_dt = CKItsynthesis - CKIdegradation - CKItphosphorilationviaSK - eq_7; // - d_SK_dt = SKsynthesis - SKdegradation; // + CycB = CycBt - 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow(std::pow(CycBt + CKIt + 1.0 / Keq, 2.0) - 4.0 * CycBt * CKIt, 1.0 / 2.0)); // + Trimer = 2.0 * CycBt * CKIt / (CycBt + CKIt + 1.0 / Keq + std::pow(std::pow(CycBt + CKIt + 1.0 / Keq, 2.0) - 4.0 * CycBt * CKIt, 1.0 / 2.0)); // + Mad = 1.0; // + TF = GK(k15p * m + k15pp * SK, k16p + k16pp * m * CycB, J15, J16); // + CycBt_synthesis = k1; // + CycBdegradation = k2p * CycBt; // + CycBdegradationviaCdh1 = k2pp * Cdh1 * CycBt; // + CycBtdegradationviaCdc20a = k2ppp * Cdc20a * CycBt; // + Cdh1synthesis = (k3p + k3pp * Cdc20a) * (1.0 - Cdh1) / (J3 + 1.0 - Cdh1); // + Cdh1degradation = (k4p * SK * Cdh1 + k4 * m * CycB * Cdh1) / (J4 + Cdh1); // + Cdc20tsynthesis = k5p + k5pp * std::pow(CycB * m / J5, n) / (1.0 + std::pow(CycB * m / J5, n)); // + Cdc20t_deg = k6 * Cdc20t; // + Cdc20activation = k7 * IEP * (Cdc20t - Cdc20a) / (J7 + Cdc20t - Cdc20a); // + Cdc20ainhibition = k8 * Mad * Cdc20a / (J8 + Cdc20a); // + Cdc20adegradation = k6 * Cdc20a; // + IEPsynthesis = k9 * m * CycB * (1.0 - IEP); // + IEPdegradation = k10 * IEP; // + growth = mu * m * (1.0 - m / mmax); // + CKItsynthesis = k11; // + CKIdegradation = k12p * CKIt; // + CKItphosphorilationviaSK = k12pp * SK * CKIt; // + eq_7 = k12ppp * m * CycB * CKIt; // + SKsynthesis = k13 * TF; // + SKdegradation = k14 * SK; // + d_CycBt_dt = CycBt_synthesis - CycBdegradation - CycBdegradationviaCdh1 - CycBtdegradationviaCdc20a; // + d_Cdc20a_dt = Cdc20activation - Cdc20ainhibition - Cdc20adegradation; // + d_Cdh1_dt = Cdh1synthesis - Cdh1degradation; // + d_m_dt = growth; // + d_Cdc20t_dt = Cdc20tsynthesis - Cdc20t_deg; // + d_IEP_dt = IEPsynthesis - IEPdegradation; // + d_CKIt_dt = CKItsynthesis - CKIdegradation - CKItphosphorilationviaSK - eq_7; // + d_SK_dt = SKsynthesis - SKdegradation; // std::vector derivatives(8); derivatives[0] = d_CycBt_dt; @@ -454,7 +454,7 @@ std::vector TysonNovak2001SbmlOdeSystem::RunModelEquations(double time, // MODEL FUNCTIONS inline double TysonNovak2001SbmlOdeSystem::GK(double A1, double A2, double A3, double A4) { - return 2.0 * A4 * A1 / ((A2 - A1) + A3 * A2 + A4 * A1 + sm::root(2.0, std::pow(((A2 - A1) + A3 * A2 + A4 * A1), 2.0) - 4.0 * (A2 - A1) * A4 * A1)); + return 2.0 * A4 * A1 / ((A2 - A1) + A3 * A2 + A4 * A1 + sm::root(2.0, std::pow((A2 - A1) + A3 * A2 + A4 * A1, 2.0) - 4.0 * (A2 - A1) * A4 * A1)); } template <> diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index 7078cfb5..2154f12d 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -12,8 +12,10 @@ from libsbml import ( AST_FUNCTION, AST_FUNCTION_DELAY, + AST_FUNCTION_POWER, AST_NAME, AST_NAME_AVOGADRO, + AST_POWER, AST_RELATIONAL_EQ, AST_RELATIONAL_GEQ, AST_RELATIONAL_GT, @@ -153,126 +155,126 @@ def search_ast_type(root: Optional["ASTNode"], node_type: int) -> bool: return False -def convert_infix_operator_to_function_syntax(formula: str, operator: str, function_name: str) -> str: - """Convert infix operator expressions to function syntax. +def rewrite_power(node: Optional["ASTNode"]) -> None: + """Rewrite power-operator (``a ^ b``) AST nodes into ``pow()`` function-call form, in place. - Example: with operator='^' and function_name='pow', rewrites - ``a ^ b`` to ``pow(a, b)``. + ``formulaToL3String`` renders a power as the infix ``a ^ b``, which is not valid C++ (``^`` is + bitwise xor there). Converting each power node to a ``pow`` call before stringifying (rendered + ``pow(a, b)``, later mapped to ``std::pow`` by the function-name token map) handles nesting + correctly -- a string-level rewrite of ``a ^ b ^ c`` would convert only the outer operator and + leave an inner ``^`` behind. Both power node types are handled: ``AST_POWER`` (from ``^``/``pow`` + in L3 infix, and from MathML ```` that libsbml normalises) and ``AST_FUNCTION_POWER`` + (from other MathML ```` forms) -- both render as ``a ^ b``. - This parser handles parenthesized operands (including nested parentheses) - and simple symbolic/numeric tokens. - - :param formula: Formula text in SBML infix style. - :param operator: Infix operator token to rewrite. - :param function_name: Function name used for the replacement. - :return: Formula text with infix operators rewritten as function calls. + :param node: The root AST node to rewrite. """ - if not operator: - raise ValueError("operator must be a non-empty string") - if not function_name: - raise ValueError("function_name must be a non-empty string") - - token_chars = set("ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789_:.eE") - - def read_left_operand(s: str, operator_index: int) -> tuple[int, int] | None: - i = operator_index - 1 - while i >= 0 and s[i].isspace(): - i -= 1 - if i < 0: - return None - - if s[i] == ")": - depth = 1 - j = i - 1 - while j >= 0: - if s[j] == ")": - depth += 1 - elif s[j] == "(": - depth -= 1 - if depth == 0: - break - j -= 1 - if depth != 0: - return None - # j points to the matching '('; also include any function name before it - k = j - 1 - while k >= 0 and s[k] in token_chars: - k -= 1 - return (k + 1, i + 1) - - j = i - while j >= 0 and s[j] in token_chars: - j -= 1 - - start = j + 1 - end = i + 1 - if start >= end: - return None - return (start, end) - - def read_right_operand(s: str, operator_index: int) -> tuple[int, int] | None: - i = operator_index + len(operator) - n = len(s) - while i < n and s[i].isspace(): - i += 1 - if i >= n: - return None - - if s[i] == "(": - depth = 1 - j = i + 1 - while j < n: - if s[j] == "(": - depth += 1 - elif s[j] == ")": - depth -= 1 - if depth == 0: - return (i, j + 1) - j += 1 - return None - - j = i - if s[j] in "+-": - j += 1 - - while j < n and s[j] in token_chars: - j += 1 - - if j <= i: - return None - return (i, j) - - max_rewrites = 500 - rewrites = 0 - search_start = 0 - while rewrites < max_rewrites: - operator_index = formula.find(operator, search_start) - if operator_index < 0: - break - - left = read_left_operand(formula, operator_index) - right = read_right_operand(formula, operator_index) - if left is None or right is None: - search_start = operator_index + len(operator) - continue - - l_start, l_end = left - r_start, r_end = right - - lhs = formula[l_start:l_end].strip() - rhs = formula[r_start:r_end].strip() - replacement = f"{function_name}({lhs}, {rhs})" - - formula = formula[:l_start] + replacement + formula[r_end:] - search_start = l_start + len(replacement) - rewrites += 1 - - return formula - - -def convert_infix_power_to_function_syntax(formula: str) -> str: - """Convert infix power expressions (a ^ b) to function syntax pow(a, b).""" - return convert_infix_operator_to_function_syntax(formula=formula, operator="^", function_name="pow") + if node is None: + return + if node.getType() in (AST_POWER, AST_FUNCTION_POWER): + node.setType(AST_FUNCTION) + node.setName("pow") + for i in range(node.getNumChildren()): + rewrite_power(node.getChild(i)) + + +# Placeholder id the avogadro csymbol is renamed to (see replace_avogadro_csymbol) so it stays +# distinct from a same-named parameter; mapped to sm::AVOGADRO in _CONSTANTS below. +AVOGADRO_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}avogadro" + +# SBML symbolic constants replaced with their C++ equivalents. (true/false are deliberately absent.) +_CONSTANTS = { + "avogadro": "sm::AVOGADRO", + AVOGADRO_PLACEHOLDER: "sm::AVOGADRO", + "exponentiale": "M_E", + "inf": "std::numeric_limits::infinity()", + "infinity": "std::numeric_limits::infinity()", + "nan": "NAN", + "notanumber": "NAN", + "pi": "M_PI", + "time": "time", + "t": "time", + "s": "time", +} + +# SBML functions whose name matches the C++ (std::) equivalent. +_UNCHANGED_FUNCTIONS = { + "acos", + "acosh", + "asin", + "asinh", + "atan", + "atanh", + "ceil", + "cos", + "cosh", + "exp", + "floor", + "pow", + "sin", + "sinh", + "sqrt", + "tan", + "tanh", +} + +# SBML functions with a different name in C++ (std::). +_RENAMED_FUNCTIONS = { + "abs": "fabs", + "arccos": "acos", + "arccosh": "acosh", + "arcsin": "asin", + "arcsinh": "asinh", + "arctan": "atan", + "arctanh": "atanh", + "ceiling": "ceil", + "ln": "log", + "power": "pow", + "rem": "fmod", +} + +# SBML functions with custom implementations in SbmlMath (sm::). +_CUSTOM_FUNCTIONS = { + "and": "and_", + "acot": "acot", + "acoth": "acoth", + "acsc": "acsc", + "acsch": "acsch", + "asec": "asec", + "asech": "asech", + "arccot": "acot", + "arccoth": "acoth", + "arccsc": "acsc", + "arccsch": "acsch", + "cot": "cot", + "coth": "coth", + "csc": "csc", + "csch": "csch", + "eq": "eq", + "factorial": "factorial", + "geq": "geq", + "gt": "gt", + "leq": "leq", + "log": "log", + "lt": "lt", + "max": "max", + "min": "min", + "neq": "neq", + "not": "not_", + "or": "or_", + "piecewise": "piecewise", + "plus": "plus", + "quotient": "quotient", + "root": "root", + "sec": "sec", + "sech": "sech", + "sqr": "sqr", + "times": "times", + "xor": "xor_", +} + +# TODO: From SBML Level 3 upwards, log defaults to base 10. +# SBML versions lower than 3 default to base e. +# See https://sbml.org/software/libsbml/5.18.0/docs/formatted/python-api/namespacelibsbml.html#a8e96a5a70569ae32655c6302638f6dc3 # noqa: B950 def formula_to_string( @@ -283,7 +285,7 @@ def formula_to_string( ) -> str: """Convert an AST math formula to an equivalent C++ string. - :param math: The AST math formula. + :param math: The AST math formula. It is not mutated (a copy is normalised internally). :param variable_types: Mapping of model id to :class:`VarType`, used to tell model variables apart from SBML constants/functions of the same spelling. :param state_variables: The model's state variables, used to resolve ``rateOf`` to a state @@ -293,20 +295,22 @@ def formula_to_string( ``rateOf`` applied to one is zero. :return: The equivalent C++ string. """ - unsupported_functions = ["delay"] - for func in unsupported_functions: - if search_ast_type(math, AST_FUNCTION_DELAY): - raise NotImplementedError(f"SBML function not supported: '{func}'.") + if search_ast_type(math, AST_FUNCTION_DELAY): + raise NotImplementedError("SBML function not supported: 'delay'.") + + # Normalise a copy so the caller's AST -- and the libsbml model it belongs to -- is not mutated. + math = math.deepCopy() strip_ast_units(math) # The avogadro csymbol and a parameter both named 'avogadro' are distinct AST nodes that - # render identically. Rename the csymbol to a placeholder (mapped to sm::AVOGADRO in the - # constants below) so it stays distinct from a same-named parameter. - avogadro_placeholder = f"{CHASTE_PREFIX}{PREFIX_SEP}avogadro" - replace_avogadro_csymbol(math, avogadro_placeholder) - # Route n-ary relationals (a < b < c) to the sm:: helpers before stringifying, so they are - # not emitted as C++-misevaluating chained infix. Binary relationals stay as infix operators. + # render identically. Rename the csymbol to a placeholder (mapped to sm::AVOGADRO in _CONSTANTS) + # so it stays distinct from a same-named parameter. + replace_avogadro_csymbol(math, AVOGADRO_PLACEHOLDER) + # Route operators that formulaToL3String would render as C++-invalid infix to function-call + # form before stringifying: n-ary relationals (a < b < c) to the sm:: helpers, and the power + # operator (a ^ b) to pow(). Binary relationals stay as infix operators. rewrite_nary_relational(math) + rewrite_power(math) formula = formulaToL3String(math) # Convert all integer literals to doubles to fix integer division. @@ -314,106 +318,6 @@ def formula_to_string( # nodes to AST_REAL. This should only need to apply to numbers used in a division. formula = re.sub(r"(?::infinity()", - "infinity": "std::numeric_limits::infinity()", - "nan": "NAN", - "notanumber": "NAN", - "pi": "M_PI", - "time": "time", - "t": "time", - "s": "time", - } - # skip: "true", "false" - - # SBML functions with same name as C++ equivalents - unchanged_functions = { - "acos", - "acosh", - "asin", - "asinh", - "atan", - "atanh", - "ceil", - "cos", - "cosh", - "exp", - "floor", - "pow", - "sin", - "sinh", - "sqrt", - "tan", - "tanh", - } - - # SBML functions with different names in C++ - renamed_functions = { - "abs": "fabs", - "arccos": "acos", - "arccosh": "acosh", - "arcsin": "asin", - "arcsinh": "asinh", - "arctan": "atan", - "arctanh": "atanh", - "ceiling": "ceil", - "ln": "log", - "power": "pow", - "rem": "fmod", - } - - # SBML functions with custom implementations - custom_functions = { - "and": "and_", - "acot": "acot", - "acoth": "acoth", - "acsc": "acsc", - "acsch": "acsch", - "asec": "asec", - "asech": "asech", - "arccot": "acot", - "arccoth": "acoth", - "arccsc": "acsc", - "arccsch": "acsch", - "arcsec": "asec", - "arcsech": "asech", - "cot": "cot", - "coth": "coth", - "csc": "csc", - "csch": "csch", - "eq": "eq", - "factorial": "factorial", - "geq": "geq", - "gt": "gt", - "leq": "leq", - "log": "log", - "lt": "lt", - "max": "max", - "min": "min", - "neq": "neq", - "not": "not_", - "or": "or_", - "piecewise": "piecewise", - "plus": "plus", - "quotient": "quotient", - "root": "root", - "sec": "sec", - "sech": "sech", - "sqr": "sqr", - "times": "times", - "xor": "xor_", - } - - # TODO: From SBML Level 3 upwards, log defaults to base 10. - # SBML versions lower than 3 default to base e. - # See https://sbml.org/software/libsbml/5.18.0/docs/formatted/python-api/namespacelibsbml.html#a8e96a5a70569ae32655c6302638f6dc3 # noqa: B950 - tokens = re.findall(r"\w+|\W+", formula) local_param_ids = {param.id for param in (local_parameters or [])} @@ -426,17 +330,17 @@ def formula_to_string( # model variable or a local parameter (e.g. a species named "s" or "t" must not be # replaced with the SBML time symbol, and a local parameter named "avogadro" must # stay the local, distinct from the avogadro csymbol handled above). - if token in constants and token not in variable_types and token not in local_param_ids: - cpp_token = f"{constants[token]}" + if token in _CONSTANTS and token not in variable_types and token not in local_param_ids: + cpp_token = f"{_CONSTANTS[token]}" - elif token in unchanged_functions: + elif token in _UNCHANGED_FUNCTIONS: cpp_token = f"std::{token}" - elif token in renamed_functions: - cpp_token = f"std::{renamed_functions[token]}" + elif token in _RENAMED_FUNCTIONS: + cpp_token = f"std::{_RENAMED_FUNCTIONS[token]}" - elif token in custom_functions: - cpp_token = f"sm::{custom_functions[token]}" + elif token in _CUSTOM_FUNCTIONS: + cpp_token = f"sm::{_CUSTOM_FUNCTIONS[token]}" cpp_tokens.append(cpp_token) cpp_formula = "".join(cpp_tokens) diff --git a/chaste_sbml/tests/test_expressions.py b/chaste_sbml/tests/test_expressions.py index 92c7e923..1fe22ef5 100644 --- a/chaste_sbml/tests/test_expressions.py +++ b/chaste_sbml/tests/test_expressions.py @@ -1,5 +1,6 @@ """Tests for SBML expression/AST -> C++ formula translation.""" +import libsbml import pytest from libsbml import ( AST_FUNCTION_COS, @@ -16,10 +17,10 @@ from chaste_sbml._config import VarType from chaste_sbml._expressions import ( collect_ast_names, - convert_infix_operator_to_function_syntax, formula_to_string, replace_avogadro_csymbol, rewrite_nary_relational, + rewrite_power, search_ast_type, strip_ast_units, substitute_ast_names, @@ -27,44 +28,43 @@ from chaste_sbml._records import StateVariable -def test_convert_infix_operator_to_function_syntax_power(): - """Converts infix power expressions to pow calls.""" - converted = convert_infix_operator_to_function_syntax( - formula="(a + b) ^ (1.0 / 2.0)", - operator="^", - function_name="pow", - ) - - assert converted == "pow((a + b), (1.0 / 2.0))" - - -def test_convert_infix_operator_to_function_syntax_custom_pair(): - """Supports arbitrary operator/function conversion pairs.""" - converted = convert_infix_operator_to_function_syntax( - formula="alpha @@ beta", operator="@@", function_name="combine" - ) +@pytest.mark.parametrize( + ("formula", "expected"), + [ + ("a ^ b", "std::pow(a, b)"), + ("2 ^ 3", "std::pow(2.0, 3.0)"), + ("(a + b) ^ 2", "std::pow(a + b, 2.0)"), + # Nested power fully converts (the old string rewrite left an inner ^ behind). + ("a ^ b ^ c", "std::pow(a, std::pow(b, c))"), + ], +) +def test_formula_to_string_power_uses_std_pow(formula, expected): + """Power operators (including nested) become std::pow calls.""" + variable_types = {name: VarType.PARAMETER for name in ("a", "b", "c")} + assert formula_to_string(parseL3Formula(formula), variable_types, state_variables=[]) == expected - assert converted == "combine(alpha, beta)" +def test_rewrite_power_converts_only_power_nodes(): + """rewrite_power converts AST_POWER nodes to pow() and leaves other nodes alone.""" + node = parseL3Formula("a ^ b ^ c") + rewrite_power(node) + assert formulaToL3String(node) == "pow(a, pow(b, c))" -def test_convert_infix_operator_to_function_syntax_nested_and_repeated(): - """Converts nested and repeated operator expressions in one formula.""" - converted = convert_infix_operator_to_function_syntax(formula="(x ^ y) ^ z", operator="^", function_name="pow") + unaffected = parseL3Formula("a * b + c") + rewrite_power(unaffected) + assert formulaToL3String(unaffected) == "a * b + c" - assert converted == "pow((pow(x, y)), z)" +def test_rewrite_power_converts_mathml_function_power(): + """A MathML (AST_FUNCTION_POWER, distinct from L3 AST_POWER) is also converted. -@pytest.mark.parametrize( - ("operator", "function_name", "error_message"), - [ - ("", "pow", "operator must be a non-empty string"), - ("^", "", "function_name must be a non-empty string"), - ], -) -def test_convert_infix_operator_to_function_syntax_invalid_args(operator: str, function_name: str, error_message: str): - """Validates required arguments for operator/function conversion.""" - with pytest.raises(ValueError, match=error_message): - convert_infix_operator_to_function_syntax(formula="a ^ b", operator=operator, function_name=function_name) + Guards against only handling the L3-parsed power node type: SBML MathML powers parse to a + different AST type that ``formulaToL3String`` still renders as ``^``. + """ + mathml = '' "a2" + node = libsbml.readMathMLFromString(mathml) + rewrite_power(node) + assert formulaToL3String(node) == "pow(a, 2)" def test_formula_to_string_functions_constants_and_power(): From 9173e6b8e99d4c77034a2f43af2a8e50ddd5316a Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Sun, 5 Jul 2026 22:14:13 +0100 Subject: [PATCH 05/19] #42 Handle SBML symbols and rateOf on the AST in formula_to_string Move the fragile string/regex handling of special symbols to the AST, before stringifying: - Symbols avogadro, pi, exponentiale and time are renamed to placeholders (replace_constant_symbols), so each stays distinct from a same-named parameter -- fixing a `` constant being confused with a parameter named `pi` (case 01819) and removing the hazard of a variable named `t`/`s` being treated as the time symbol. - rateOf is resolved on the AST (resolve_rate_of), robust to the argument's spelling and to nested parentheses that broke the old regex. - _CONSTANTS now maps the capitalised INF/NaN that formulaToL3String actually emits (the old map only had lowercase), and the int->double regex no longer mangles positive exponents (1e+5). Cases 00950/00951/01811 remain failing but are re-noted "non-finite expected results (INF/NaN)": their reference output is literally INF/NaN, which the finite-tolerance test comparison cannot validate -- a test-harness limitation, not a codegen bug. 01819 now passes. Behaviour is unchanged for models without these symbols (the reference generation tests are byte-identical), so no reference models are regenerated. Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 8 +- chaste_sbml/_expressions.py | 194 +++++++++++------- chaste_sbml/tests/test_expressions.py | 61 +++++- 3 files changed, 176 insertions(+), 87 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 9726039c..d986f2cd 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -948,8 +948,8 @@ test,status,notes 00947,pass,ok 00948,pass,ok 00949,pass,ok (no-ode placeholder) -00950,fail_known,codegen: inf constant -00951,fail_known,codegen: inf constant +00950,fail_known,non-finite expected results (INF/NaN) +00951,fail_known,non-finite expected results (INF/NaN) 00952,unsupported,random event execution 00953,fail_known,event priority/ordering 00954,pass,ok (no-ode placeholder) @@ -1809,7 +1809,7 @@ test,status,notes 01808,pass,ok 01809,pass,ok 01810,pass,ok (no-ode placeholder) -01811,fail_known,codegen: inf in expected results +01811,fail_known,non-finite expected results (INF/NaN) 01812,fail_known,codegen: string literal in math 01813,fail_known,codegen: string literal in math 01814,pass,ok (no-ode placeholder) @@ -1817,6 +1817,6 @@ test,status,notes 01816,pass,ok (no-ode placeholder) 01817,pass,ok (no-ode placeholder) 01818,pass,ok (no-ode placeholder) -01819,fail_known,pi constant vs parameter collision +01819,pass,ok (pi constant vs parameter distinguished) 01820,pass,ok (no-ode placeholder) 01821,pass,ok (no-ode placeholder) diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index 2154f12d..72a830c5 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -10,12 +10,17 @@ from typing import TYPE_CHECKING, Optional from libsbml import ( + AST_CONSTANT_E, + AST_CONSTANT_PI, AST_FUNCTION, AST_FUNCTION_DELAY, AST_FUNCTION_POWER, + AST_FUNCTION_RATE_OF, AST_NAME, AST_NAME_AVOGADRO, + AST_NAME_TIME, AST_POWER, + AST_REAL, AST_RELATIONAL_EQ, AST_RELATIONAL_GEQ, AST_RELATIONAL_GT, @@ -68,21 +73,75 @@ def substitute_ast_names(node: "ASTNode", replacements: dict) -> "ASTNode": return result -def replace_avogadro_csymbol(node: "ASTNode", placeholder: str) -> None: - """Recursively rename avogadro csymbol nodes to a placeholder identifier. +# Placeholder ids the SBML math symbols are renamed to. A parameter may share a symbol's spelling +# (e.g. a parameter named 'pi' or 'time'), and ``formulaToL3String`` renders the symbol and the +# parameter identically; renaming the symbol node to a unique placeholder keeps them distinct, so +# the symbol maps to its C++ value (below) while the parameter keeps its own name. +AVOGADRO_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}avogadro" +PI_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}pi" +E_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}exponentiale" +TIME_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}time" + +_SYMBOL_PLACEHOLDERS = { + AST_NAME_AVOGADRO: AVOGADRO_PLACEHOLDER, + AST_CONSTANT_PI: PI_PLACEHOLDER, + AST_CONSTANT_E: E_PLACEHOLDER, + AST_NAME_TIME: TIME_PLACEHOLDER, +} + - The avogadro csymbol and a parameter both named 'avogadro' are distinct AST nodes that - ``formulaToL3String`` renders identically. Renaming the csymbol to a placeholder keeps it - distinct so it can be mapped to ``sm::AVOGADRO`` while the parameter keeps its own name. +def replace_constant_symbols(node: "ASTNode") -> None: + """Recursively rename SBML math-symbol nodes (avogadro, pi, exponentiale, time) to placeholders. + + Each of these symbols renders identically to a same-named parameter, so it is renamed to a + unique placeholder that maps to its C++ value, keeping it distinct from a parameter of the same + spelling. Handling this at the AST level (rather than string-matching the rendered name) means a + real variable named ``t``, ``pi`` etc. is never mistaken for the symbol. :param node: The root AST node. - :param placeholder: The identifier to rename avogadro csymbol nodes to. """ - if node.getType() == AST_NAME_AVOGADRO: + placeholder = _SYMBOL_PLACEHOLDERS.get(node.getType()) + if placeholder is not None: node.setType(AST_NAME) node.setName(placeholder) for i in range(node.getNumChildren()): - replace_avogadro_csymbol(node.getChild(i), placeholder) + replace_constant_symbols(node.getChild(i)) + + +def resolve_rate_of( + node: "ASTNode", + variable_types: dict, + state_variables: list["StateVariable"], + local_param_ids: set, +) -> None: + """Recursively resolve ``rateOf(x)`` nodes to a state variable's derivative id, or ``0``. + + ``rateOf`` of a state variable is that variable's time derivative (its ``d__dt`` id); of + anything constant -- a parameter, or a kinetic-law local parameter that shadows a global -- it + is zero. Resolving on the AST (rather than by regex on the rendered string) is robust to the + argument's spelling and to surrounding parentheses. + + :param node: The root AST node. + :param variable_types: Mapping of model id to :class:`VarType`. + :param state_variables: The model's state variables (to look up a derivative id). + :param local_param_ids: Ids of local parameters in scope (constant, so their rate is zero). + """ + if node.getType() == AST_FUNCTION_RATE_OF and node.getNumChildren() >= 1: + var = node.getChild(0).getName() or "" + derivative = None + if var not in local_param_ids and variable_types.get(var, VarType.UNKNOWN) == VarType.STATE_VARIABLE: + derivative = next((sv.derivative_id for sv in state_variables if sv.id == var), None) + while node.getNumChildren() > 0: + node.removeChild(0) + if derivative is not None: + node.setType(AST_NAME) + node.setName(derivative) + else: + node.setType(AST_REAL) + node.setValue(0.0) + return + for i in range(node.getNumChildren()): + resolve_rate_of(node.getChild(i), variable_types, state_variables, local_param_ids) def strip_ast_units(node: "ASTNode") -> None: @@ -90,10 +149,11 @@ def strip_ast_units(node: "ASTNode") -> None: SBML allows numeric literals to carry units annotations (e.g. ````). ``formulaToL3String`` includes these in its output (e.g. ``0.00015 mole``), which is not - valid C++. Stripping them here is safe because the units carry no mathematical information. + valid C++. Stripping them here is safe because the units are not used. :param node: The root AST node. """ + # TODO: Handle units properly if node.isSetUnits(): node.unsetUnits() for i in range(node.getNumChildren()): @@ -113,13 +173,15 @@ def strip_ast_units(node: "ASTNode") -> None: def rewrite_nary_relational(node: Optional["ASTNode"]) -> None: """Rewrite n-ary (>2 operand) relational AST nodes into function-call form, in place. - SBML MathML relationals are n-ary -- ``lt(a, b, c)`` means ``a < b < c`` (i.e. ``a < b`` and - ``b < c``), and ``eq(a, b, c)`` means all equal -- but ``formulaToL3String`` renders them as a - chained infix expression (``a < b < c``), which C++ mis-evaluates left-to-right as - ``(a < b) < c``. Converting each such node to a function call (rendered ``lt(a, b, c)``, later - mapped to the variadic ``sm::lt`` helper that implements the correct semantics) fixes this. - Binary relationals are left as infix operators, so existing generated code is unchanged. - ``neq`` is binary-only in MathML, so it is never n-ary and is not rewritten here. + SBML MathML relationals are n-ary: ``lt(a, b, c)`` means ``a < b < c`` + (i.e. ``a < b`` and ``b < c``), and ``eq(a, b, c)`` means all equal. However, + ``formulaToL3String`` renders them as a chained infix expression (``a < b < c``), + which C++ mis-evaluates left-to-right as ``(a < b) < c``. Converting each + such node to a function call (rendered ``lt(a, b, c)``, later mapped to the + variadic ``sm::lt`` helper that implements the correct semantics) fixes this. + Binary relationals are left as infix operators, so existing generated code + is unchanged. ``neq`` is binary-only in MathML, so it is never n-ary and is + not rewritten here. :param node: The root AST node to rewrite. """ @@ -158,13 +220,13 @@ def search_ast_type(root: Optional["ASTNode"], node_type: int) -> bool: def rewrite_power(node: Optional["ASTNode"]) -> None: """Rewrite power-operator (``a ^ b``) AST nodes into ``pow()`` function-call form, in place. - ``formulaToL3String`` renders a power as the infix ``a ^ b``, which is not valid C++ (``^`` is - bitwise xor there). Converting each power node to a ``pow`` call before stringifying (rendered - ``pow(a, b)``, later mapped to ``std::pow`` by the function-name token map) handles nesting - correctly -- a string-level rewrite of ``a ^ b ^ c`` would convert only the outer operator and - leave an inner ``^`` behind. Both power node types are handled: ``AST_POWER`` (from ``^``/``pow`` - in L3 infix, and from MathML ```` that libsbml normalises) and ``AST_FUNCTION_POWER`` - (from other MathML ```` forms) -- both render as ``a ^ b``. + ``formulaToL3String`` renders a power as the infix ``a ^ b``, which is not + valid C++ (``^`` is bitwise xor). Converting each power node to a ``pow`` + call before stringifying (rendered ``pow(a, b)``, later mapped to ``std::pow`` + by the function-name token map) handles nesting correctly. Both power node + types are handled: ``AST_POWER`` (from ``^``/``pow`` in L3 infix, and from + MathML ```` that libsbml normalises) and ``AST_FUNCTION_POWER`` + (from other MathML ```` forms). :param node: The root AST node to rewrite. """ @@ -177,23 +239,21 @@ def rewrite_power(node: Optional["ASTNode"]) -> None: rewrite_power(node.getChild(i)) -# Placeholder id the avogadro csymbol is renamed to (see replace_avogadro_csymbol) so it stays -# distinct from a same-named parameter; mapped to sm::AVOGADRO in _CONSTANTS below. -AVOGADRO_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}avogadro" - -# SBML symbolic constants replaced with their C++ equivalents. (true/false are deliberately absent.) +# SBML symbolic constants replaced with their C++ equivalents, keyed by the +# placeholder each symbol is renamed to (see replace_constant_symbols) plus the +# infinity/NaN literals formulaToL3String emits. (avogadro/pi/exponentiale/time +# are handled via placeholders; true/false deliberately absent.) _CONSTANTS = { - "avogadro": "sm::AVOGADRO", AVOGADRO_PLACEHOLDER: "sm::AVOGADRO", - "exponentiale": "M_E", + PI_PLACEHOLDER: "M_PI", + E_PLACEHOLDER: "M_E", + TIME_PLACEHOLDER: "time", "inf": "std::numeric_limits::infinity()", + "INF": "std::numeric_limits::infinity()", "infinity": "std::numeric_limits::infinity()", "nan": "NAN", + "NaN": "NAN", "notanumber": "NAN", - "pi": "M_PI", - "time": "time", - "t": "time", - "s": "time", } # SBML functions whose name matches the C++ (std::) equivalent. @@ -285,51 +345,50 @@ def formula_to_string( ) -> str: """Convert an AST math formula to an equivalent C++ string. - :param math: The AST math formula. It is not mutated (a copy is normalised internally). - :param variable_types: Mapping of model id to :class:`VarType`, used to tell model variables - apart from SBML constants/functions of the same spelling. - :param state_variables: The model's state variables, used to resolve ``rateOf`` to a state - variable's derivative id. - :param local_parameters: Local parameters in scope (e.g. a reaction's kinetic-law - parameters). These shadow global symbols of the same name and are constant, so - ``rateOf`` applied to one is zero. + :param math: The AST math formula. It is not mutated. + :param variable_types: Mapping of model id to :class:`VarType`, used to tell + model variables apart from SBML constants/functions of the same spelling. + :param state_variables: The model's state variables, used to resolve ``rateOf`` + to a state variable's derivative id. + :param local_parameters: Local parameters in scope (e.g. a reaction's + kinetic-law parameters). These shadow global symbols of the same name + and are constant, so ``rateOf`` applied to one is zero. :return: The equivalent C++ string. """ if search_ast_type(math, AST_FUNCTION_DELAY): raise NotImplementedError("SBML function not supported: 'delay'.") - # Normalise a copy so the caller's AST -- and the libsbml model it belongs to -- is not mutated. + # Make a copy so the AST is not mutated. math = math.deepCopy() + local_param_ids = {param.id for param in (local_parameters or [])} strip_ast_units(math) - # The avogadro csymbol and a parameter both named 'avogadro' are distinct AST nodes that - # render identically. Rename the csymbol to a placeholder (mapped to sm::AVOGADRO in _CONSTANTS) - # so it stays distinct from a same-named parameter. - replace_avogadro_csymbol(math, AVOGADRO_PLACEHOLDER) - # Route operators that formulaToL3String would render as C++-invalid infix to function-call - # form before stringifying: n-ary relationals (a < b < c) to the sm:: helpers, and the power - # operator (a ^ b) to pow(). Binary relationals stay as infix operators. + # Modify the AST, renaming the avogadro/pi/exponentiale/time symbol nodes to + # placeholders (so each stays distinct from a same-named parameter), and + # resolve rateOf to a derivative id (or 0). + replace_constant_symbols(math) + resolve_rate_of(math, variable_types, state_variables, local_param_ids) + # Route operators that formulaToL3String would render as C++-invalid infix + # to function-call form: n-ary relationals (a < b < c) to the SBMLMath + # helpers (sm::lt(a, b, c)), and the power operator (a ^ b) to std::pow(). + # Binary relationals stay as infix operators (a < b). rewrite_nary_relational(math) rewrite_power(math) formula = formulaToL3String(math) - # Convert all integer literals to doubles to fix integer division. - # TODO: Perhaps instead of regex, traverse AST and convert some AST_INTEGER - # nodes to AST_REAL. This should only need to apply to numbers used in a division. - formula = re.sub(r"(?", "M_PI"), + ("", "M_E"), + ("", "std::numeric_limits::infinity()"), + ("", "NAN"), + ], +) +def test_formula_to_string_math_constants(mathml_symbol, expected): + """SBML math constants (incl. the capitalised INF/NaN libsbml emits) map to their C++ form.""" + math = libsbml.readMathMLFromString(f'{mathml_symbol}') + assert formula_to_string(math, {}, state_variables=[]) == expected + + +def test_formula_to_string_pi_constant_and_parameter_are_distinct(): + """The constant maps to M_PI while a same-named parameter keeps its own name.""" + math = libsbml.readMathMLFromString( + 'pi' + ) + assert formula_to_string(math, {"pi": VarType.PARAMETER}, state_variables=[]) == "pi + M_PI" + + +def test_formula_to_string_time_symbol_and_named_variable_are_distinct(): + """The time csymbol maps to `time`; a plain variable named `t` is left as-is (not time).""" + assert formula_to_string(parseL3Formula("time"), {}, state_variables=[]) == "time" + assert formula_to_string(parseL3Formula("t"), {"t": VarType.PARAMETER}, state_variables=[]) == "t" + + +def test_formula_to_string_rateof_local_parameter_is_zero(): + """rateOf of a kinetic-law local parameter (constant) is zero, not a global's derivative.""" + from chaste_sbml._records import LocalParameter + + out = formula_to_string( + parseL3Formula("rateOf(k)"), + variable_types={"k": VarType.STATE_VARIABLE}, + state_variables=[ + StateVariable(index=0, id="k", derivative_id="d_k_dt", label="", initial_value=None, units="") + ], + local_parameters=[LocalParameter(id="k", label="", value="1")], + ) + assert out == "0.0" + + def test_strip_ast_units_removes_units_recursively(): """Units annotations are stripped from a node and its children.""" node = parseL3Formula("2 * 3") From 67102bb5c19af1bdfe03b16cb22b8ff294ca62dd Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 05:11:00 +0100 Subject: [PATCH 06/19] #42 Precompile shared headers to speed up SbmlRefModels rebuilds Add target_precompile_headers for the heavy header every generated case model shares (AbstractSbmlOdeSystem.hpp) on the project library, and the shared CVODE/test-helper headers on each case test target. Cuts a single-case incremental rebuild's compile time noticeably. Serialization export wrappers and FakePetscSetup are left out (order-sensitive / per-TU). Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/CMakeLists.txt | 8 ++++++++ chaste_sbml/SbmlRefModels/test/CMakeLists.txt | 14 ++++++++++++++ 2 files changed, 22 insertions(+) diff --git a/chaste_sbml/SbmlRefModels/CMakeLists.txt b/chaste_sbml/SbmlRefModels/CMakeLists.txt index 9ae069b5..88a617a2 100644 --- a/chaste_sbml/SbmlRefModels/CMakeLists.txt +++ b/chaste_sbml/SbmlRefModels/CMakeLists.txt @@ -1,2 +1,10 @@ find_package(Chaste COMPONENTS cell_based) chaste_do_project(SbmlRefModels) + +# Precompile the heavy header every generated case model includes, to speed up rebuilds. The +# serialization export wrappers are deliberately excluded (they are order-sensitive and belong in +# each .cpp, not a shared PCH). +if (TARGET chaste_project_SbmlRefModels) + target_precompile_headers(chaste_project_SbmlRefModels PRIVATE + "${CMAKE_CURRENT_SOURCE_DIR}/src/AbstractSbmlOdeSystem.hpp") +endif() diff --git a/chaste_sbml/SbmlRefModels/test/CMakeLists.txt b/chaste_sbml/SbmlRefModels/test/CMakeLists.txt index cea25f65..abee9aea 100644 --- a/chaste_sbml/SbmlRefModels/test/CMakeLists.txt +++ b/chaste_sbml/SbmlRefModels/test/CMakeLists.txt @@ -1 +1,15 @@ chaste_do_test_project(SbmlRefModels) + +# Precompile the heavy Chaste headers shared by every generated case test (the CVODE solver and +# test helpers), to speed up rebuilds. FakePetscSetup and the per-test model header are left out. +get_property(_sbml_test_targets DIRECTORY PROPERTY BUILDSYSTEM_TARGETS) +foreach(_sbml_test_target ${_sbml_test_targets}) + get_target_property(_sbml_test_type ${_sbml_test_target} TYPE) + if (_sbml_test_type STREQUAL "EXECUTABLE") + target_precompile_headers(${_sbml_test_target} PRIVATE + + + + ) + endif() +endforeach() From abe0f40fa84cec048482329ca6fe9caf33b6f7d6 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 05:55:55 +0100 Subject: [PATCH 07/19] #42 Handle events/functions with missing MathML instead of crashing Degenerate SBML that omits a MathML element used to crash the generator (AttributeError on a null math node). Handle each per SBML semantics, so the model still simulates: - an event with no trigger, or a trigger with no MathML, can never fire, so the event is skipped; - an event delay with no MathML is treated as no delay (zero); - a function definition with no body is skipped (a valid model never calls it). Marks test-suite cases 01238, 01239, 01241, 01271 as passing. Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 8 +-- chaste_sbml/_model_builder.py | 34 +++++++---- chaste_sbml/tests/test_model_builder.py | 56 +++++++++++++++++++ 3 files changed, 83 insertions(+), 15 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index d986f2cd..3a7c84f1 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1236,10 +1236,10 @@ test,status,notes 01235,pass,ok (no-ode placeholder) 01236,pass,ok (no-ode placeholder) 01237,pass,ok (no-ode placeholder) -01238,fail_known,missing mathml (event trigger) -01239,fail_known,missing mathml (event trigger) +01238,pass,ok (missing mathml treated as no-op) +01239,pass,ok (missing mathml treated as no-op) 01240,pass,ok (no-ode placeholder) -01241,fail_known,missing mathml (event delay) +01241,pass,ok (missing mathml treated as no-op) 01242,pass,ok (no-ode placeholder) 01243,pass,ok (no-ode placeholder) 01244,unsupported,algebraic rule @@ -1269,7 +1269,7 @@ test,status,notes 01268,unsupported,delay 01269,unsupported,delay 01270,unsupported,delay -01271,fail_known,missing mathml (function definition) +01271,pass,ok (missing mathml treated as no-op) 01272,pass,ok (no-ode placeholder) 01273,pass,ok (no-ode placeholder) 01274,fail_known,codegen: implies operator diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index b54d2341..f1429ab9 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -548,12 +548,16 @@ def _format_equations(self) -> None: def _format_events(self) -> None: """Add events to template variables.""" for event in self._sbml_events: + trigger = event.getTrigger() + if trigger is None or trigger.getMath() is None: + # An event with no trigger, or a trigger with no MathML, can never fire. + continue self._reject_unsupported_delay(event) label = event.getName().strip() event_type = self._guess_event_type(label) # Compute the trigger formula first: it mutates the trigger AST in place (stripping # units, renaming the avogadro csymbol), and the distance below deep-copies that AST. - trigger_formula = self._formula_to_string(event.getTrigger().getMath()) + trigger_formula = self._formula_to_string(trigger.getMath()) trigger_distance = self._event_trigger_distance(event) assignments = self._event_assignments(event) self._compensate_compartment_resizes(event, assignments) @@ -573,16 +577,21 @@ def _format_events(self) -> None: def _reject_unsupported_delay(self, event) -> None: """Raise if the event has a non-zero delay (delays are not supported).""" - if event.isSetDelay(): - math = self._formula_to_string(event.getDelay().getMath()) - try: - delay = float(math) - except ValueError: - delay = 9999 - - if delay != 0.0: - # Delay of zero is equivalent to no delay - raise NotImplementedError("Events with delays are not supported.") + if not event.isSetDelay(): + return + delay_math = event.getDelay().getMath() + if delay_math is None: + # A delay element with no MathML is treated as no delay. + return + math = self._formula_to_string(delay_math) + try: + delay = float(math) + except ValueError: + delay = 9999 + + if delay != 0.0: + # Delay of zero is equivalent to no delay + raise NotImplementedError("Events with delays are not supported.") @staticmethod def _guess_event_type(label: str) -> EventType: @@ -753,6 +762,9 @@ def _event_priority(self, event) -> Optional[str]: def _format_function_definitions(self) -> None: """Add function definitions to template variables.""" for fd in self._sbml_function_definitions: + if fd.getBody() is None: + # A function definition with no MathML body cannot be emitted. + continue fd_id = fd.getId() label = fd.getName().strip() arg_list = get_function_definition_arguments(fd) diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index 524c06b2..a4b6f182 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -283,6 +283,62 @@ def test_build_extracts_cell_division_event(): assert event.assignments == [EventAssignment(index=3, lhs="m", rhs="m / 2.0", type=VarType.STATE_VARIABLE)] +def _event_builder(): + """A bare builder with an empty model, wired just enough for _format_events.""" + builder = _builder_without_init() + doc = libsbml.SBMLDocument(3, 2) + builder._sbml_model = doc.createModel() + builder._keep_doc = doc # keep the owning document (and thus the model) alive + builder._sbml_events = builder._sbml_model.getListOfEvents() + builder._sbml_function_definitions = builder._sbml_model.getListOfFunctionDefinitions() + builder._variable_types = {} + builder._state_variables = [] + builder._events = [] + builder._functions = [] + return builder + + +def test_format_events_skips_event_with_missing_trigger_math(): + """An event whose trigger has no MathML can never fire, so no event is emitted.""" + builder = _event_builder() + event = builder._sbml_model.createEvent() + event.createTrigger() # trigger element present but with no math + event.createEventAssignment().setVariable("p") + + builder._format_events() + + assert builder._events == [] + + +def test_format_events_skips_event_with_no_trigger(): + """An event with no trigger element at all can never fire, so no event is emitted.""" + builder = _event_builder() + builder._sbml_model.createEvent() # no trigger + + builder._format_events() + + assert builder._events == [] + + +def test_reject_unsupported_delay_allows_missing_delay_math(): + """A delay element with no MathML is treated as no delay, so it is not rejected.""" + builder = _event_builder() + event = builder._sbml_model.createEvent() + event.createDelay() # delay element present but with no math + + builder._reject_unsupported_delay(event) # should not raise + + +def test_format_function_definitions_skips_missing_body(): + """A function definition with no MathML body is skipped rather than crashing.""" + builder = _event_builder() + builder._sbml_model.createFunctionDefinition().setId("foo") # no body + + builder._format_function_definitions() + + assert builder._functions == [] + + def _distance_builder() -> ModelBuilder: """A bare builder wired just enough for _event_trigger_distance and trigger formatting.""" builder = _builder_without_init() From 97a3ba220dbe895ad237a689d3efe5575f893de5 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 07:04:18 +0100 Subject: [PATCH 08/19] #42 Classify event-modified constant species as parameters A species with no reactions, rules or ODE but modified by an event was classified as a derived quantity. Derived quantities are recomputed from the current member each step, so CalculateDerivedQuantitiesAndParameters reported the post-event value at every earlier time point. Model such a species as a (variable) parameter instead -- time-resolved per step, event assignments via the deferred parameter mechanism -- mirroring the boundary-species handling. Fixes cases 00930, 00931, and corrects reference model Chen2004 (BUB2/LTE1/MAD2 move to parameters; its OdeSystem test expectations are updated). The remaining event-ordering cases are marked unsupported (priority ordering, useValuesFromTriggerTime, persistence, competing simultaneous events). Co-Authored-By: Claude Opus 4.8 --- .../Chen2004/Chen2004SbmlOdeSystem.cpp | 377 ++++++++++-------- .../Chen2004/Chen2004SbmlOdeSystem.hpp | 6 +- .../test/data/sbml_test_suite_status.csv | 22 +- .../reference/TestChen2004SbmlOdeSystem.hpp | 8 +- chaste_sbml/_model_builder.py | 11 +- chaste_sbml/tests/test_model_builder.py | 38 ++ 6 files changed, 276 insertions(+), 186 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp index cdc4aae3..6355da1a 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.cpp @@ -11,7 +11,7 @@ namespace sm = sbmlmath; Chen2004SbmlOdeSystem::Chen2004SbmlOdeSystem() - : AbstractSbmlOdeSystem(36, 143, 4) + : AbstractSbmlOdeSystem(36, 146, 4) { mpSystemInfo.reset(new CellwiseOdeSystemInformation); @@ -30,8 +30,8 @@ Chen2004SbmlOdeSystem::Chen2004SbmlOdeSystem() mEventSatisfied = { true, true, true, true }; // From SBML trigger initialValue mEventTriggered.resize(4, false); - mEventAdjustedParameters.resize(143, false); - mEventAdjustedParameterValues.resize(143, 0.0); + mEventAdjustedParameters.resize(146, false); + mEventAdjustedParameterValues.resize(146, 0.0); mEventAdjustedStateVars.resize(36, false); mEventAdjustedStateValues.resize(36, 0.0); @@ -104,7 +104,6 @@ std::vector Chen2004SbmlOdeSystem::ComputeDerivedQuantities(double time, dqs.push_back(cell); dqs.push_back(BCK2); - dqs.push_back(BUB2); dqs.push_back(CDC14T); dqs.push_back(CDC15i); dqs.push_back(CDC6T); @@ -113,8 +112,6 @@ std::vector Chen2004SbmlOdeSystem::ComputeDerivedQuantities(double time, dqs.push_back(CLB5T); dqs.push_back(CLN3); dqs.push_back(IE); - dqs.push_back(LTE1); - dqs.push_back(MAD2); dqs.push_back(MCM1); dqs.push_back(NET1T); dqs.push_back(PE); @@ -740,6 +737,9 @@ void Chen2004SbmlOdeSystem::Initialise(double time) SetDefaultInitialCondition(34, SWI5P); SetDefaultInitialCondition(35, TEM1GTP); + mParameters.push_back(BUB2); + mParameters.push_back(LTE1); + mParameters.push_back(MAD2); mParameters.push_back(b0); mParameters.push_back(bub2h); mParameters.push_back(bub2l); @@ -1076,8 +1076,23 @@ double Chen2004SbmlOdeSystem::ProcessModelEvents(double time, const std::vector< if (mEventTriggered[1]) { [[maybe_unused]] double event_priority = std::numeric_limits::max(); - MAD2 = mad2h; - BUB2 = bub2h; + // MAD2 = mad2h + if (!mEventAdjustedParameters[2] + || event_priority <= mEventAdjustedParameterPriority[2]) + { + mEventAdjustedParameters[2] = true; + mEventAdjustedParameterValues[2] = mad2h; + mEventAdjustedParameterPriority[2] = event_priority; + } + + // BUB2 = bub2h + if (!mEventAdjustedParameters[0] + || event_priority <= mEventAdjustedParameterPriority[0]) + { + mEventAdjustedParameters[0] = true; + mEventAdjustedParameterValues[0] = bub2h; + mEventAdjustedParameterPriority[0] = event_priority; + } } } @@ -1160,9 +1175,32 @@ double Chen2004SbmlOdeSystem::ProcessModelEvents(double time, const std::vector< if (mEventTriggered[2]) { [[maybe_unused]] double event_priority = std::numeric_limits::max(); - MAD2 = mad2l; - LTE1 = lte1h; - BUB2 = bub2l; + // MAD2 = mad2l + if (!mEventAdjustedParameters[2] + || event_priority <= mEventAdjustedParameterPriority[2]) + { + mEventAdjustedParameters[2] = true; + mEventAdjustedParameterValues[2] = mad2l; + mEventAdjustedParameterPriority[2] = event_priority; + } + + // LTE1 = lte1h + if (!mEventAdjustedParameters[1] + || event_priority <= mEventAdjustedParameterPriority[1]) + { + mEventAdjustedParameters[1] = true; + mEventAdjustedParameterValues[1] = lte1h; + mEventAdjustedParameterPriority[1] = event_priority; + } + + // BUB2 = bub2l + if (!mEventAdjustedParameters[0] + || event_priority <= mEventAdjustedParameterPriority[0]) + { + mEventAdjustedParameters[0] = true; + mEventAdjustedParameterValues[0] = bub2l; + mEventAdjustedParameterPriority[0] = event_priority; + } } } @@ -1254,7 +1292,15 @@ double Chen2004SbmlOdeSystem::ProcessModelEvents(double time, const std::vector< mEventAdjustedStatePriority[22] = event_priority; } - LTE1 = lte1l; + // LTE1 = lte1l + if (!mEventAdjustedParameters[1] + || event_priority <= mEventAdjustedParameterPriority[1]) + { + mEventAdjustedParameters[1] = true; + mEventAdjustedParameterValues[1] = lte1l; + mEventAdjustedParameterPriority[1] = event_priority; + } + // BUD = 0.0 if (!mEventAdjustedStateVars[0] || event_priority <= mEventAdjustedStatePriority[0]) @@ -1317,149 +1363,152 @@ std::vector Chen2004SbmlOdeSystem::RunModelEquations(double time, const SWI5P = rStateVariables[34]; TEM1GTP = rStateVariables[35]; - b0 = GetParameter(0); - bub2h = GetParameter(1); - bub2l = GetParameter(2); - C0 = GetParameter(3); - CDC15T = GetParameter(4); - Dn3 = GetParameter(5); - ebudb5 = GetParameter(6); - ebudn2 = GetParameter(7); - ebudn3 = GetParameter(8); - ec1b2 = GetParameter(9); - ec1b5 = GetParameter(10); - ec1k2 = GetParameter(11); - ec1n2 = GetParameter(12); - ec1n3 = GetParameter(13); - ef6b2 = GetParameter(14); - ef6b5 = GetParameter(15); - ef6k2 = GetParameter(16); - ef6n2 = GetParameter(17); - ef6n3 = GetParameter(18); - eicdhb2 = GetParameter(19); - eicdhb5 = GetParameter(20); - eicdhn2 = GetParameter(21); - eicdhn3 = GetParameter(22); - eorib2 = GetParameter(23); - eorib5 = GetParameter(24); - esbfb5 = GetParameter(25); - esbfn2 = GetParameter(26); - esbfn3 = GetParameter(27); - ESP1T = GetParameter(28); - IET = GetParameter(29); - J20ppx = GetParameter(30); - Jacdh = GetParameter(31); - Jaiep = GetParameter(32); - Jamcm = GetParameter(33); - Jasbf = GetParameter(34); - Jatem = GetParameter(35); - Jd2c1 = GetParameter(36); - Jd2f6 = GetParameter(37); - Jicdh = GetParameter(38); - Jiiep = GetParameter(39); - Jimcm = GetParameter(40); - Jisbf = GetParameter(41); - Jitem = GetParameter(42); - Jn3 = GetParameter(43); - Jpds = GetParameter(44); - Jspn = GetParameter(45); - ka15_p = GetParameter(46); - ka15_p_p = GetParameter(47); - ka15p = GetParameter(48); - ka20_p = GetParameter(49); - ka20_p_p = GetParameter(50); - kacdh_p = GetParameter(51); - kacdh_p_p = GetParameter(52); - kaiep = GetParameter(53); - kamcm = GetParameter(54); - kasb2 = GetParameter(55); - kasb5 = GetParameter(56); - kasbf = GetParameter(57); - kasesp = GetParameter(58); - kasf2 = GetParameter(59); - kasf5 = GetParameter(60); - kasrent = GetParameter(61); - kasrentp = GetParameter(62); - kaswi = GetParameter(63); - kd14 = GetParameter(64); - kd1c1 = GetParameter(65); - kd1f6 = GetParameter(66); - kd1pds_p = GetParameter(67); - kd20 = GetParameter(68); - kd2c1 = GetParameter(69); - kd2f6 = GetParameter(70); - kd2pds_p_p = GetParameter(71); - kd3c1 = GetParameter(72); - kd3f6 = GetParameter(73); - kd3pds_p_p = GetParameter(74); - kdb2_p = GetParameter(75); - kdb2_p_p = GetParameter(76); - kdb2p = GetParameter(77); - kdb5_p = GetParameter(78); - kdb5_p_p = GetParameter(79); - kdbud = GetParameter(80); - kdcdh = GetParameter(81); - kdib2 = GetParameter(82); - kdib5 = GetParameter(83); - kdiesp = GetParameter(84); - kdif2 = GetParameter(85); - kdif5 = GetParameter(86); - kdirent = GetParameter(87); - kdirentp = GetParameter(88); - kdn2 = GetParameter(89); - kdnet = GetParameter(90); - kdori = GetParameter(91); - kdppx_p = GetParameter(92); - kdppx_p_p = GetParameter(93); - kdspn = GetParameter(94); - kdswi = GetParameter(95); - KEZ = GetParameter(96); - KEZ2 = GetParameter(97); - ki15 = GetParameter(98); - kicdh_p = GetParameter(99); - kicdh_p_p = GetParameter(100); - kiiep = GetParameter(101); - kimcm = GetParameter(102); - kisbf_p = GetParameter(103); - kisbf_p_p = GetParameter(104); - kiswi = GetParameter(105); - kkpnet_p = GetParameter(106); - kkpnet_p_p = GetParameter(107); - kppc1 = GetParameter(108); - kppf6 = GetParameter(109); - kppnet_p = GetParameter(110); - kppnet_p_p = GetParameter(111); - ks14 = GetParameter(112); - ks1pds_p_p = GetParameter(113); - ks20_p = GetParameter(114); - ks20_p_p = GetParameter(115); - ks2pds_p_p = GetParameter(116); - ksb2_p = GetParameter(117); - ksb2_p_p = GetParameter(118); - ksb5_p = GetParameter(119); - ksb5_p_p = GetParameter(120); - ksbud = GetParameter(121); - ksc1_p = GetParameter(122); - ksc1_p_p = GetParameter(123); - kscdh = GetParameter(124); - ksf6_p = GetParameter(125); - ksf6_p_p = GetParameter(126); - ksf6_p_p_p = GetParameter(127); - ksn2_p = GetParameter(128); - ksn2_p_p = GetParameter(129); - ksnet = GetParameter(130); - ksori = GetParameter(131); - kspds_p = GetParameter(132); - ksppx = GetParameter(133); - ksspn = GetParameter(134); - ksswi_p = GetParameter(135); - ksswi_p_p = GetParameter(136); - lte1h = GetParameter(137); - lte1l = GetParameter(138); - mad2h = GetParameter(139); - mad2l = GetParameter(140); - mdt = GetParameter(141); - TEM1T = GetParameter(142); + BUB2 = GetParameter(0); + LTE1 = GetParameter(1); + MAD2 = GetParameter(2); + b0 = GetParameter(3); + bub2h = GetParameter(4); + bub2l = GetParameter(5); + C0 = GetParameter(6); + CDC15T = GetParameter(7); + Dn3 = GetParameter(8); + ebudb5 = GetParameter(9); + ebudn2 = GetParameter(10); + ebudn3 = GetParameter(11); + ec1b2 = GetParameter(12); + ec1b5 = GetParameter(13); + ec1k2 = GetParameter(14); + ec1n2 = GetParameter(15); + ec1n3 = GetParameter(16); + ef6b2 = GetParameter(17); + ef6b5 = GetParameter(18); + ef6k2 = GetParameter(19); + ef6n2 = GetParameter(20); + ef6n3 = GetParameter(21); + eicdhb2 = GetParameter(22); + eicdhb5 = GetParameter(23); + eicdhn2 = GetParameter(24); + eicdhn3 = GetParameter(25); + eorib2 = GetParameter(26); + eorib5 = GetParameter(27); + esbfb5 = GetParameter(28); + esbfn2 = GetParameter(29); + esbfn3 = GetParameter(30); + ESP1T = GetParameter(31); + IET = GetParameter(32); + J20ppx = GetParameter(33); + Jacdh = GetParameter(34); + Jaiep = GetParameter(35); + Jamcm = GetParameter(36); + Jasbf = GetParameter(37); + Jatem = GetParameter(38); + Jd2c1 = GetParameter(39); + Jd2f6 = GetParameter(40); + Jicdh = GetParameter(41); + Jiiep = GetParameter(42); + Jimcm = GetParameter(43); + Jisbf = GetParameter(44); + Jitem = GetParameter(45); + Jn3 = GetParameter(46); + Jpds = GetParameter(47); + Jspn = GetParameter(48); + ka15_p = GetParameter(49); + ka15_p_p = GetParameter(50); + ka15p = GetParameter(51); + ka20_p = GetParameter(52); + ka20_p_p = GetParameter(53); + kacdh_p = GetParameter(54); + kacdh_p_p = GetParameter(55); + kaiep = GetParameter(56); + kamcm = GetParameter(57); + kasb2 = GetParameter(58); + kasb5 = GetParameter(59); + kasbf = GetParameter(60); + kasesp = GetParameter(61); + kasf2 = GetParameter(62); + kasf5 = GetParameter(63); + kasrent = GetParameter(64); + kasrentp = GetParameter(65); + kaswi = GetParameter(66); + kd14 = GetParameter(67); + kd1c1 = GetParameter(68); + kd1f6 = GetParameter(69); + kd1pds_p = GetParameter(70); + kd20 = GetParameter(71); + kd2c1 = GetParameter(72); + kd2f6 = GetParameter(73); + kd2pds_p_p = GetParameter(74); + kd3c1 = GetParameter(75); + kd3f6 = GetParameter(76); + kd3pds_p_p = GetParameter(77); + kdb2_p = GetParameter(78); + kdb2_p_p = GetParameter(79); + kdb2p = GetParameter(80); + kdb5_p = GetParameter(81); + kdb5_p_p = GetParameter(82); + kdbud = GetParameter(83); + kdcdh = GetParameter(84); + kdib2 = GetParameter(85); + kdib5 = GetParameter(86); + kdiesp = GetParameter(87); + kdif2 = GetParameter(88); + kdif5 = GetParameter(89); + kdirent = GetParameter(90); + kdirentp = GetParameter(91); + kdn2 = GetParameter(92); + kdnet = GetParameter(93); + kdori = GetParameter(94); + kdppx_p = GetParameter(95); + kdppx_p_p = GetParameter(96); + kdspn = GetParameter(97); + kdswi = GetParameter(98); + KEZ = GetParameter(99); + KEZ2 = GetParameter(100); + ki15 = GetParameter(101); + kicdh_p = GetParameter(102); + kicdh_p_p = GetParameter(103); + kiiep = GetParameter(104); + kimcm = GetParameter(105); + kisbf_p = GetParameter(106); + kisbf_p_p = GetParameter(107); + kiswi = GetParameter(108); + kkpnet_p = GetParameter(109); + kkpnet_p_p = GetParameter(110); + kppc1 = GetParameter(111); + kppf6 = GetParameter(112); + kppnet_p = GetParameter(113); + kppnet_p_p = GetParameter(114); + ks14 = GetParameter(115); + ks1pds_p_p = GetParameter(116); + ks20_p = GetParameter(117); + ks20_p_p = GetParameter(118); + ks2pds_p_p = GetParameter(119); + ksb2_p = GetParameter(120); + ksb2_p_p = GetParameter(121); + ksb5_p = GetParameter(122); + ksb5_p_p = GetParameter(123); + ksbud = GetParameter(124); + ksc1_p = GetParameter(125); + ksc1_p_p = GetParameter(126); + kscdh = GetParameter(127); + ksf6_p = GetParameter(128); + ksf6_p_p = GetParameter(129); + ksf6_p_p_p = GetParameter(130); + ksn2_p = GetParameter(131); + ksn2_p_p = GetParameter(132); + ksnet = GetParameter(133); + ksori = GetParameter(134); + kspds_p = GetParameter(135); + ksppx = GetParameter(136); + ksspn = GetParameter(137); + ksswi_p = GetParameter(138); + ksswi_p_p = GetParameter(139); + lte1h = GetParameter(140); + lte1l = GetParameter(141); + mad2h = GetParameter(142); + mad2l = GetParameter(143); + mdt = GetParameter(144); + TEM1T = GetParameter(145); BCK2 = b0 * MASS; // CDC14T = CDC14 + RENT + RENTP; // @@ -1855,9 +1904,6 @@ void CellwiseOdeSystemInformation::Initialise() this->mDerivedQuantityNames.push_back("BCK2"); this->mDerivedQuantityUnits.push_back("non-dim"); - this->mDerivedQuantityNames.push_back("BUB2"); - this->mDerivedQuantityUnits.push_back("non-dim"); - this->mDerivedQuantityNames.push_back("CDC14T"); this->mDerivedQuantityUnits.push_back("non-dim"); @@ -1882,12 +1928,6 @@ void CellwiseOdeSystemInformation::Initialise() this->mDerivedQuantityNames.push_back("IE"); this->mDerivedQuantityUnits.push_back("non-dim"); - this->mDerivedQuantityNames.push_back("LTE1"); - this->mDerivedQuantityUnits.push_back("non-dim"); - - this->mDerivedQuantityNames.push_back("MAD2"); - this->mDerivedQuantityUnits.push_back("non-dim"); - this->mDerivedQuantityNames.push_back("MCM1"); this->mDerivedQuantityUnits.push_back("non-dim"); @@ -2411,6 +2451,15 @@ void CellwiseOdeSystemInformation::Initialise() this->mDerivedQuantityUnits.push_back("non-dim"); // PARAMETERS + this->mParameterNames.push_back("BUB2"); + this->mParameterUnits.push_back("non-dim"); + + this->mParameterNames.push_back("LTE1"); + this->mParameterUnits.push_back("non-dim"); + + this->mParameterNames.push_back("MAD2"); + this->mParameterUnits.push_back("non-dim"); + this->mParameterNames.push_back("b0"); this->mParameterUnits.push_back("non-dim"); diff --git a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.hpp b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.hpp index 787accbe..47f083d9 100644 --- a/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.hpp +++ b/chaste_sbml/SbmlRefModels/src/reference/Chen2004/Chen2004SbmlOdeSystem.hpp @@ -27,6 +27,9 @@ class Chen2004SbmlOdeSystem : public AbstractSbmlOdeSystem } // PARAMETERS + double BUB2; // BUB2 + double LTE1; // LTE1 + double MAD2; // MAD2 double b0; // b0 double bub2h; // bub2h double bub2l; // bub2l @@ -249,7 +252,6 @@ class Chen2004SbmlOdeSystem : public AbstractSbmlOdeSystem // DERIVED QUANTITIES double cell; // cell double BCK2; // BCK2 - double BUB2; // BUB2 double CDC14T; // CDC14T double CDC15i; // CDC15i double CDC6T; // CDC6T @@ -258,8 +260,6 @@ class Chen2004SbmlOdeSystem : public AbstractSbmlOdeSystem double CLB5T; // CLB5T double CLN3; // CLN3 double IE; // IE - double LTE1; // LTE1 - double MAD2; // MAD2 double MCM1; // MCM1 double NET1T; // NET1T double PE; // PE diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 3a7c84f1..f92a7cfc 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -928,12 +928,12 @@ test,status,notes 00927,pass,ok 00928,pass,ok 00929,pass,ok -00930,fail_known,event priority/ordering -00931,fail_known,event priority/ordering +00930,pass,ok (event-modified species as parameter) +00931,pass,ok (event-modified species as parameter) 00932,unsupported,delay 00933,unsupported,delay -00934,fail_known,event priority/ordering -00935,fail_known,event priority/ordering +00934,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) +00935,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 00936,unsupported,delay 00937,unsupported,delay function 00938,unsupported,delay @@ -951,7 +951,7 @@ test,status,notes 00950,fail_known,non-finite expected results (INF/NaN) 00951,fail_known,non-finite expected results (INF/NaN) 00952,unsupported,random event execution -00953,fail_known,event priority/ordering +00953,unsupported,event execution semantics (competing simultaneous events) 00954,pass,ok (no-ode placeholder) 00955,pass,ok (no-ode placeholder) 00956,pass,ok (no-ode placeholder) @@ -961,11 +961,11 @@ test,status,notes 00960,pass,ok (no-ode placeholder) 00961,pass,ok (no-ode placeholder) 00962,unsupported,random event execution -00963,fail_known,event priority/ordering +00963,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 00964,unsupported,random event execution 00965,unsupported,random event execution 00966,unsupported,random event execution -00967,fail_known,event priority/ordering +00967,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 00968,pass,ok 00969,pass,ok 00970,pass,ok @@ -976,7 +976,7 @@ test,status,notes 00975,pass,ok 00976,pass,ok 00977,pass,ok -00978,fail_known,event priority/ordering +00978,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 00979,pass,ok (no-ode placeholder) 00980,unsupported,delay 00981,unsupported,delay @@ -1329,13 +1329,13 @@ test,status,notes 01328,unsupported,delay 01329,unsupported,delay 01330,pass,ok (no-ode placeholder) -01331,fail_known,event priority/ordering +01331,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 01332,pass,ok (no-ode placeholder) -01333,fail_known,event priority/ordering +01333,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) 01334,pass,ok (no-ode placeholder) 01335,unsupported,delay 01336,pass,ok (no-ode placeholder) -01337,fail_known,event priority/ordering +01337,unsupported,event execution semantics (persistent triggers) 01338,pass,ok 01339,pass,ok 01340,pass,ok diff --git a/chaste_sbml/SbmlRefModels/test/reference/TestChen2004SbmlOdeSystem.hpp b/chaste_sbml/SbmlRefModels/test/reference/TestChen2004SbmlOdeSystem.hpp index 3e932785..3d62b88f 100644 --- a/chaste_sbml/SbmlRefModels/test/reference/TestChen2004SbmlOdeSystem.hpp +++ b/chaste_sbml/SbmlRefModels/test/reference/TestChen2004SbmlOdeSystem.hpp @@ -66,7 +66,7 @@ class TestChen2004SbmlOdeSystem : public AbstractCellBasedTestSuite { private: const unsigned ODE_SIZE = 36u; - const unsigned NUM_DERIVED_QUANTITIES = 187u; + const unsigned NUM_DERIVED_QUANTITIES = 184u; std::vector default_initial_conditions = { 0.008473, // BUD @@ -470,7 +470,6 @@ class TestChen2004SbmlOdeSystem : public AbstractCellBasedTestSuite std::vector dq_names = { "cell", "BCK2", - "BUB2", "CDC14T", "CDC15i", "CDC6T", @@ -479,8 +478,6 @@ class TestChen2004SbmlOdeSystem : public AbstractCellBasedTestSuite "CLB5T", "CLN3", "IE", - "LTE1", - "MAD2", "MCM1", "NET1T", "PE", @@ -524,7 +521,6 @@ class TestChen2004SbmlOdeSystem : public AbstractCellBasedTestSuite std::vector dqs_expected = { 1.0, // cell 0.065125026, // BCK2 - 0.2, // BUB2 2.117884, // CDC14T 0.34346699999999997, // CDC15i 0.3866693, // CDC6T @@ -533,8 +529,6 @@ class TestChen2004SbmlOdeSystem : public AbstractCellBasedTestSuite 0.1289119, // CLB5T 0.06694509131879892, // CLN3 0.8985, // IE - 0.1, // LTE1 - 0.01, // MAD2 0.4690076182110798, // MCM1 2.638456, // NET1T 0.6986870000000001, // PE diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index f1429ab9..befaefd5 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -1108,8 +1108,11 @@ def _add_species_dynamics( self._add_equation( var=state_var.derivative_id, math=parseL3Formula(rhs), eq_type=EquationType.DERIVATIVE ) + elif species_id in self._event_assigned_ids: + # Otherwise constant but modified by an event: model as a (variable) parameter + self._add_parameter(species_id, label, initial_value, units) else: - # Truly constant: no reactions, no rules, constant compartment + # Truly constant: no reactions, no rules, no events, constant compartment self._add_derived_quantity(species_id, label, initial_value, units) def _formula_to_string(self, math: "ASTNode", local_parameters: Optional[list["LocalParameter"]] = None) -> str: @@ -1220,6 +1223,12 @@ def build(self) -> None: self._events = [] self._functions = [] + # Ids assigned by some event, so species classification can model an otherwise-constant + # event-modified species as a (variable) parameter. + self._event_assigned_ids = { + ea.getVariable() for event in self._sbml_events for ea in event.getListOfEventAssignments() + } + # Names already claimed by real SBML entities and the Chaste base classes. Synthetic # identifiers (derivatives, amount/concentration conversions, initial-assignment # intermediates) are allocated against this so they never collide. diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index a4b6f182..eaedbbf0 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -240,6 +240,44 @@ def test_build_no_ode_model_preserves_real_outputs(tmp_path): assert "S" in {p.id for p in builder._parameters} +def test_build_event_modified_constant_species_is_a_parameter(tmp_path): + """An otherwise-constant species modified by an event is a parameter, not a derived quantity. + + A derived quantity is recomputed from the current member each step, so an event's change would + be reported at every earlier time point; a parameter is time-resolved in the recorded solution. + """ + doc = libsbml.SBMLDocument(3, 2) + model = doc.createModel() + compartment = model.createCompartment() + compartment.setId("C") + compartment.setConstant(True) + compartment.setSize(1.0) + compartment.setSpatialDimensions(3) + species = model.createSpecies() + species.setId("S") + species.setCompartment("C") + species.setConstant(False) + species.setBoundaryCondition(False) + species.setHasOnlySubstanceUnits(False) + species.setInitialAmount(0.0) + event = model.createEvent() + event.setUseValuesFromTriggerTime(True) + trigger = event.createTrigger() + trigger.setMath(parseL3Formula("time >= 1")) + trigger.setInitialValue(True) + trigger.setPersistent(True) + assignment = event.createEventAssignment() + assignment.setVariable("S") + assignment.setMath(parseL3Formula("5")) + path = tmp_path / "event_species.xml" + libsbml.writeSBMLToFile(doc, str(path)) + + builder = _build(path) + + assert "S" in {p.id for p in builder._parameters} + assert "S" not in {d.id for d in builder._derived_quantities} + + def test_build_does_not_add_placeholder_when_model_has_odes(): """A model with genuine ODEs is left untouched: no placeholder state variable is added.""" builder = _build(REFERENCE / "Goldbeter1991" / "Goldbeter1991.xml") From c2e5aa6ef118fc8f993eae8267599366ad1bfd21 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 07:06:23 +0100 Subject: [PATCH 09/19] #42 Reject flux-balance (fbc package) models with a clear error Flux-balance-constraint models cannot be translated to an ODE system. Detect the SBML fbc package and raise NotImplementedError rather than silently generating a model that ignores the constraints. Covers 34/35 flux-balance test-suite cases (no supported case uses fbc). The remaining unsupported categories (event execution semantics, random event execution) cannot be distinguished from supported models by an SBML feature, so are left generating. Also pass -DChaste_UPDATE_PROVENANCE=OFF in infra/configure.sh. Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/infra/configure.sh | 4 ++-- chaste_sbml/_model_builder.py | 12 ++++++++++- chaste_sbml/tests/test_model_builder.py | 21 ++++++++++++++++++++ 3 files changed, 34 insertions(+), 3 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/infra/configure.sh b/chaste_sbml/SbmlRefModels/infra/configure.sh index 5f0dd878..485447ca 100755 --- a/chaste_sbml/SbmlRefModels/infra/configure.sh +++ b/chaste_sbml/SbmlRefModels/infra/configure.sh @@ -19,7 +19,7 @@ fi if [[ ! -d "${chaste_dir}" ]]; then echo "Error: Chaste source directory not found at '${chaste_dir}'." >&2 - echo "Set CHASTE_SOURCE_DIR to override the default sibling checkout path." >&2 + echo "Set CHASTE_SOURCE_DIR to override the default checkout path." >&2 exit 1 fi @@ -27,5 +27,5 @@ mkdir -p "${output_dir}" "${build_dir}" export CHASTE_TEST_OUTPUT="${output_dir}" cd "${build_dir}" -nice -n 9 cmake "${chaste_dir}" +nice -n 9 cmake -DChaste_UPDATE_PROVENANCE=OFF "${chaste_dir}" diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index befaefd5..35e145c4 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -575,6 +575,14 @@ def _format_events(self) -> None: priority, ) + def _reject_unsupported_packages(self) -> None: + """Raise if the model uses an SBML package the generator cannot translate. + + :raises NotImplementedError: if the model uses the flux-balance-constraints (fbc) package. + """ + if self._sbml_model.getPlugin("fbc") is not None: + raise NotImplementedError("Flux balance constraint models (SBML fbc package) are not supported.") + def _reject_unsupported_delay(self, event) -> None: """Raise if the event has a non-zero delay (delays are not supported).""" if not event.isSetDelay(): @@ -1205,6 +1213,8 @@ def group(dq: "DerivedQuantity") -> int: def build(self) -> None: """Process the SBML model to set up the formatted variables for templates.""" + self._reject_unsupported_packages() + self._variable_types = {} self._index_of = None # id -> index cache, built lazily by _get_variable_index self._odes = {} @@ -1224,7 +1234,7 @@ def build(self) -> None: self._functions = [] # Ids assigned by some event, so species classification can model an otherwise-constant - # event-modified species as a (variable) parameter. + # event-modified species as a (variable) parameterp. self._event_assigned_ids = { ea.getVariable() for event in self._sbml_events for ea in event.getListOfEventAssignments() } diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index eaedbbf0..c8d393b7 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -367,6 +367,27 @@ def test_reject_unsupported_delay_allows_missing_delay_math(): builder._reject_unsupported_delay(event) # should not raise +def test_reject_unsupported_packages_rejects_fbc(): + """A model using the flux-balance-constraints (fbc) package is rejected.""" + doc = libsbml.SBMLDocument(libsbml.SBMLNamespaces(3, 1, "fbc", 2)) + builder = _builder_without_init() + builder._sbml_model = doc.createModel() + builder._keep_doc = doc + + with pytest.raises(NotImplementedError, match="fbc"): + builder._reject_unsupported_packages() + + +def test_reject_unsupported_packages_allows_plain_model(): + """A model without an unsupported package is not rejected.""" + doc = libsbml.SBMLDocument(3, 2) + builder = _builder_without_init() + builder._sbml_model = doc.createModel() + builder._keep_doc = doc + + builder._reject_unsupported_packages() # should not raise + + def test_format_function_definitions_skips_missing_body(): """A function definition with no MathML body is skipped rather than crashing.""" builder = _event_builder() From c5e3031ea08d4559fee1cdd95a35feea31a64d72 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 09:35:51 +0100 Subject: [PATCH 10/19] #42 Support the implies operator and single-argument min/max Map SBML 'implies' to a new sm::implies helper (implies(a, b) == (not a) or b), and add single-argument overloads for sm::min/sm::max (the min/max of one value is that value). Both were unmapped and failed to compile. Adds C++ SbmlMath tests for each. Marks cases 01274/01275/01276/01279/01280/01281/01497 passing. Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/src/SbmlMath.hpp | 9 +++++++-- chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp | 10 ++++++++++ .../test/data/sbml_test_suite_status.csv | 14 +++++++------- chaste_sbml/_expressions.py | 1 + chaste_sbml/tests/test_expressions.py | 5 +++++ 5 files changed, 30 insertions(+), 9 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp b/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp index 35001701..60ab6272 100644 --- a/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp +++ b/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp @@ -95,6 +95,9 @@ constexpr bool or_(Args... args); // not_ inline bool not_(bool x) { return !x; } +// implies (a implies b == (not a) or b) +inline bool implies(double antecedent, double consequent) { return !antecedent || consequent; } + // xor_ template constexpr bool xor_(Args... args); @@ -149,11 +152,13 @@ inline double asech(double x) { return std::acosh(1.0 / x); } // factorial inline double factorial(double x) { return std::tgamma(x + 1.0); } -// max +// max (the max of a single value is that value) +constexpr double max(double only) { return only; } template constexpr double max(double first, double second, Args... rest); -// min +// min (the min of a single value is that value) +constexpr double min(double only) { return only; } template constexpr double min(double first, double second, Args... rest); diff --git a/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp b/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp index f592cb63..074ef9d2 100644 --- a/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp +++ b/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp @@ -154,6 +154,14 @@ class TestSbmlMath : public CxxTest::TestSuite TS_ASSERT_EQUALS(sm::not_(true), false); TS_ASSERT_EQUALS(sm::not_(false), true); } + void TestImplies() + { + // implies(a, b) == (not a) or b + TS_ASSERT_EQUALS(sm::implies(true, true), true); + TS_ASSERT_EQUALS(sm::implies(true, false), false); + TS_ASSERT_EQUALS(sm::implies(false, true), true); + TS_ASSERT_EQUALS(sm::implies(false, false), true); + } void TestXor() { @@ -361,6 +369,7 @@ class TestSbmlMath : public CxxTest::TestSuite void TestMax() { + TS_ASSERT_EQUALS(sm::max(2.0), 2.0); // single argument TS_ASSERT_EQUALS(sm::max(1.0, 2.0), 2.0); TS_ASSERT_EQUALS(sm::max(2.0, 1.0), 2.0); TS_ASSERT_EQUALS(sm::max(1.0, 2.0, 3.0), 3.0); @@ -374,6 +383,7 @@ class TestSbmlMath : public CxxTest::TestSuite void TestMin() { + TS_ASSERT_EQUALS(sm::min(2.0), 2.0); // single argument TS_ASSERT_EQUALS(sm::min(1.0, 2.0), 1.0); TS_ASSERT_EQUALS(sm::min(2.0, 1.0), 1.0); TS_ASSERT_EQUALS(sm::min(3.0, 2.0, 1.0), 1.0); diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index f92a7cfc..9dd7426a 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1272,14 +1272,14 @@ test,status,notes 01271,pass,ok (missing mathml treated as no-op) 01272,pass,ok (no-ode placeholder) 01273,pass,ok (no-ode placeholder) -01274,fail_known,codegen: implies operator -01275,fail_known,codegen: implies operator -01276,fail_known,codegen: implies operator +01274,pass,ok (implies operator) +01275,pass,ok (implies operator) +01276,pass,ok (implies operator) 01277,pass,ok (no-ode placeholder) 01278,pass,ok (no-ode placeholder) -01279,fail_known,codegen: single-arg min -01280,fail_known,codegen: implies operator -01281,fail_known,codegen: implies operator +01279,pass,ok (single-arg min/max) +01280,pass,ok (implies operator) +01281,pass,ok (implies operator) 01282,pass,ok (no-ode placeholder) 01283,pass,ok (no-ode placeholder) 01284,pass,ok (no-ode placeholder) @@ -1495,7 +1495,7 @@ test,status,notes 01494,pass,ok (nary relational via sm:: helpers) 01495,pass,ok (no-ode placeholder) 01496,pass,ok (no-ode placeholder) -01497,fail_known,codegen: implies operator +01497,pass,ok (implies operator) 01498,pass,ok 01499,unsupported,algebraic rule 01500,unsupported,algebraic rule diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index 72a830c5..a6aa45f5 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -313,6 +313,7 @@ def rewrite_power(node: Optional["ASTNode"]) -> None: "factorial": "factorial", "geq": "geq", "gt": "gt", + "implies": "implies", "leq": "leq", "log": "log", "lt": "lt", diff --git a/chaste_sbml/tests/test_expressions.py b/chaste_sbml/tests/test_expressions.py index 623bcd62..44f9cd14 100644 --- a/chaste_sbml/tests/test_expressions.py +++ b/chaste_sbml/tests/test_expressions.py @@ -77,6 +77,11 @@ def test_formula_to_string_functions_constants_and_power(): assert out == "k * std::pow(S, 2.0) + M_PI" +def test_formula_to_string_implies_maps_to_sbmlmath(): + """The MathML 'implies' logical operator maps to the sm::implies helper.""" + assert formula_to_string(parseL3Formula("implies(a, b)"), {}, state_variables=[]) == "sm::implies(a, b)" + + @pytest.mark.parametrize( ("formula", "expected"), [ From a24d29a12eac27804446d2c8f5ac3bc503186961 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 09:35:51 +0100 Subject: [PATCH 11/19] #42 Derive variable types from records; use a DerivedQuantity kind enum - Drop the parallel _variable_types dict: a variable's VarType is implied by which record collection it belongs to, so derive it (_variable_type_map) instead of maintaining a duplicate updated in five places. - Replace the DerivedQuantity is_conversion/is_reaction boolean pair with a single DerivedQuantityKind enum (NORMAL/REACTION/CONVERSION) -- mutually exclusive by construction, template-friendly, self-documenting. No behaviour change (reference generation output is byte-identical). Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/_config.py | 13 +++++ chaste_sbml/_model_builder.py | 70 +++++++++++++++---------- chaste_sbml/_records.py | 5 +- chaste_sbml/_rendering.py | 2 + chaste_sbml/templates/ode/ode.hpp | 2 +- chaste_sbml/tests/test_model_builder.py | 29 ++++++---- 6 files changed, 78 insertions(+), 43 deletions(-) diff --git a/chaste_sbml/_config.py b/chaste_sbml/_config.py index c3c16ec0..cc37b5fb 100644 --- a/chaste_sbml/_config.py +++ b/chaste_sbml/_config.py @@ -51,6 +51,19 @@ class VarType(Enum): UNKNOWN = 5 +class DerivedQuantityKind(Enum): + """Kind of derived quantity, which determines how it is declared and computed. + + NORMAL is declared as a member and computed from an equation; REACTION is a reaction flux + (declared and computed by the reaction machinery, only exposed here as an output); CONVERSION + is an amount/concentration conversion computed in ComputeDerivedQuantities. + """ + + NORMAL = 0 + REACTION = 1 + CONVERSION = 2 + + class ModelType(Enum): """Enumeration of model types for code generation.""" diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index 35e145c4..f5900691 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -35,6 +35,7 @@ NON_DIM_UNITS, PLACEHOLDER_STATE_ID, PREFIX_SEP, + DerivedQuantityKind, EquationType, EventType, VarType, @@ -115,7 +116,7 @@ def _add_amount(self, species: "Species") -> None: amt_rhs = f"{species_id} * {compartment_id}" amt_math = parseL3Formula(amt_rhs) - self._add_derived_quantity(amt_id, amt_label, None, amt_units, is_conversion=True) + self._add_derived_quantity(amt_id, amt_label, None, amt_units, kind=DerivedQuantityKind.CONVERSION) self._add_equation(var=amt_id, math=amt_math, eq_type=EquationType.CONVERSION) def _add_concentration(self, species: "Species") -> None: @@ -136,7 +137,7 @@ def _add_concentration(self, species: "Species") -> None: conc_rhs = f"{species_id} / {compartment_id}" conc_math = parseL3Formula(conc_rhs) - self._add_derived_quantity(conc_id, conc_label, None, conc_units, is_conversion=True) + self._add_derived_quantity(conc_id, conc_label, None, conc_units, kind=DerivedQuantityKind.CONVERSION) self._add_equation(var=conc_id, math=conc_math, eq_type=EquationType.CONVERSION) def _add_assignment_rule(self, id_: str, label: str, var: str, math: "ASTNode") -> None: @@ -155,8 +156,7 @@ def _add_derived_quantity( label: str, initial_value: Optional[float], units: str = NON_DIM_UNITS, - is_conversion: bool = False, - is_reaction: bool = False, + kind: DerivedQuantityKind = DerivedQuantityKind.NORMAL, ) -> None: """Add a derived quantity to the template variables. @@ -164,11 +164,10 @@ def _add_derived_quantity( :param label: The variable description. :param initial_value: The variable initial value. :param units: The variable units. - :param is_conversion: True if the derived quantity is an amount/concentration - conversion (computed in ComputeDerivedQuantities rather than stored as a member). - :param is_reaction: True if the derived quantity is a reaction flux exposed as an - output. Its member is already declared by the reaction loop, so it is not - re-declared, and its variable type stays VarType.REACTION. + :param kind: The kind of derived quantity. CONVERSION is an amount/concentration conversion + (computed in ComputeDerivedQuantities rather than stored as a member); REACTION is a + reaction flux exposed as an output -- its member is already declared by the reaction + loop, so it is not re-declared and its variable type stays VarType.REACTION. """ self._derived_quantities.append( DerivedQuantity( @@ -176,14 +175,11 @@ def _add_derived_quantity( label=label, index=len(self._derived_quantities), initial_value=initial_value, - is_conversion=is_conversion, - is_reaction=is_reaction, units=units, + kind=kind, name=self._names.sbml_name(id_), ) ) - if not is_reaction: - self._variable_types[id_] = VarType.DERIVED_QUANTITY def _add_equation( self, @@ -245,7 +241,6 @@ def _add_function(self, id_: str, label: str, args: str, body: str) -> None: :param body: The function body. """ self._functions.append(Function(id=id_, label=label, index=len(self._functions), args=args, body=body)) - self._variable_types[id_] = VarType.FUNCTION def _add_initial_assignment(self, id_: str, label: str, var: str, math: Optional["ASTNode"] = None) -> None: """Add an initial assignment to the template variables. @@ -284,7 +279,6 @@ def _add_parameter( name=self._names.sbml_name(id_), ) ) - self._variable_types[id_] = VarType.PARAMETER def _add_rate_rule(self, id_: str, label: str, var: str, math: "ASTNode") -> None: """Add a rate rule to the template variables. @@ -303,7 +297,6 @@ def _add_reaction(self, id_: str, label: str) -> None: :param label: The variable description. """ self._reactions.append(Reaction(index=len(self._reactions), id=id_, label=label)) - self._variable_types[id_] = VarType.REACTION def _add_state_variable( self, id_: str, label: str, initial_value: Optional[float], units: str = NON_DIM_UNITS @@ -327,7 +320,6 @@ def _add_state_variable( ) self._state_variables.append(state_var) - self._variable_types[id_] = VarType.STATE_VARIABLE return state_var @@ -865,8 +857,8 @@ def _format_reactions(self) -> None: # A reaction's ID denotes its rate (flux), an observable that rules, events and # test outputs can read. Expose it as a derived quantity. The flux member is # already declared and computed by the reaction machinery, so it keeps its - # VarType.REACTION type and is not re-declared (is_reaction=True). - self._add_derived_quantity(id_, label, None, NON_DIM_UNITS, is_reaction=True) + # VarType.REACTION type and is not re-declared (kind=REACTION). + self._add_derived_quantity(id_, label, None, NON_DIM_UNITS, kind=DerivedQuantityKind.REACTION) kinetic_law = reaction.getKineticLaw() if kinetic_law is None: @@ -1133,7 +1125,7 @@ def _formula_to_string(self, math: "ASTNode", local_parameters: Optional[list["L :param local_parameters: Local parameters in scope for the formula. :return: The equivalent C++ string. """ - return formula_to_string(math, self._variable_types, self._state_variables, local_parameters) + return formula_to_string(math, self._variable_type_map(), self._state_variables, local_parameters) def _get_timescale_multiplier(self) -> float: """Get the timescale multiplier. @@ -1178,17 +1170,38 @@ def _build_variable_index(self) -> dict[str, int]: for record in collection: index_of[record.id] = record.index for derived_quantity in self._derived_quantities: - if not derived_quantity.is_reaction: + if derived_quantity.kind != DerivedQuantityKind.REACTION: index_of[derived_quantity.id] = derived_quantity.index return index_of + def _variable_type_map(self) -> dict[str, VarType]: + """Map each variable id to its VarType, derived from the built record collections. + + The type of a variable is implied by which collection its record belongs to, so no separate + mapping is kept. A reaction flux is both a Reaction (typed VarType.REACTION) and an + REACTION-kind DerivedQuantity; the latter is skipped so the reaction id keeps its type. + """ + types: dict[str, VarType] = {} + for derived_quantity in self._derived_quantities: + if derived_quantity.kind != DerivedQuantityKind.REACTION: + types[derived_quantity.id] = VarType.DERIVED_QUANTITY + for function in self._functions: + types[function.id] = VarType.FUNCTION + for parameter in self._parameters: + types[parameter.id] = VarType.PARAMETER + for reaction in self._reactions: + types[reaction.id] = VarType.REACTION + for state_variable in self._state_variables: + types[state_variable.id] = VarType.STATE_VARIABLE + return types + def _get_variable_type(self, var_id: str) -> VarType: """Get the type of a variable based on its ID. :param var_id: The variable ID. :return: The variable's VarType, or VarType.UNKNOWN if it is not a known variable. """ - return self._variable_types.get(var_id, VarType.UNKNOWN) + return self._variable_type_map().get(var_id, VarType.UNKNOWN) def _order_derived_quantities(self) -> None: """Order derived quantities: normal quantities, then reactions, then conversions. @@ -1201,9 +1214,9 @@ def _order_derived_quantities(self) -> None: """ def group(dq: "DerivedQuantity") -> int: - if dq.is_conversion: + if dq.kind == DerivedQuantityKind.CONVERSION: return 2 - if dq.is_reaction: + if dq.kind == DerivedQuantityKind.REACTION: return 1 return 0 @@ -1215,7 +1228,6 @@ def build(self) -> None: """Process the SBML model to set up the formatted variables for templates.""" self._reject_unsupported_packages() - self._variable_types = {} self._index_of = None # id -> index cache, built lazily by _get_variable_index self._odes = {} @@ -1295,10 +1307,12 @@ def _check_name_conflicts(self) -> None: identifiers.append((var.id, "state variable")) identifiers.append((var.derivative_id, "state-variable derivative")) for dq in self._derived_quantities: - if dq.is_reaction: + if dq.kind == DerivedQuantityKind.REACTION: continue # Declared under reactions; counted there. - kind = "amount/concentration conversion" if dq.is_conversion else "derived quantity" - identifiers.append((dq.id, kind)) + kind_label = ( + "amount/concentration conversion" if dq.kind == DerivedQuantityKind.CONVERSION else "derived quantity" + ) + identifiers.append((dq.id, kind_label)) identifiers += [(s.id, "stoichiometry variable") for s in self._stoichiometry_variables] identifiers += [(r.id, "reaction") for r in self._reactions] identifiers += [(f.id, "function") for f in self._functions] diff --git a/chaste_sbml/_records.py b/chaste_sbml/_records.py index d98b4050..4abba602 100644 --- a/chaste_sbml/_records.py +++ b/chaste_sbml/_records.py @@ -9,7 +9,7 @@ from dataclasses import dataclass from typing import TYPE_CHECKING, Optional -from ._config import EquationType, EventType, VarType +from ._config import DerivedQuantityKind, EquationType, EventType, VarType if TYPE_CHECKING: from libsbml import ASTNode @@ -71,9 +71,8 @@ class DerivedQuantity: label: str index: int initial_value: Optional[float] - is_conversion: bool - is_reaction: bool units: str + kind: DerivedQuantityKind = DerivedQuantityKind.NORMAL # SBML id reported to Chaste; equals ``id`` unless the id was escaped for C++ (see NameManager). name: str = "" diff --git a/chaste_sbml/_rendering.py b/chaste_sbml/_rendering.py index fff44b51..29aba91e 100644 --- a/chaste_sbml/_rendering.py +++ b/chaste_sbml/_rendering.py @@ -18,6 +18,7 @@ CONCENTRATION_PREFIX, PREFIX_SEP, ROOT_DIR, + DerivedQuantityKind, EquationType, EventType, VarType, @@ -35,6 +36,7 @@ class CodeRenderer: ) _env.globals["AMOUNT_PREFIX"] = AMOUNT_PREFIX _env.globals["CONCENTRATION_PREFIX"] = CONCENTRATION_PREFIX + _env.globals["DerivedQuantityKind"] = DerivedQuantityKind _env.globals["EquationType"] = EquationType _env.globals["EventType"] = EventType _env.globals["PREFIX_SEP"] = PREFIX_SEP diff --git a/chaste_sbml/templates/ode/ode.hpp b/chaste_sbml/templates/ode/ode.hpp index 76d03632..a84752ae 100644 --- a/chaste_sbml/templates/ode/ode.hpp +++ b/chaste_sbml/templates/ode/ode.hpp @@ -42,7 +42,7 @@ class {{ ode_class_name }} : public AbstractSbmlOdeSystem // DERIVED QUANTITIES {% for dq in derived_quantities %} -{% if dq["is_conversion"] is false() and dq["is_reaction"] is false() %} +{% if dq["kind"] == DerivedQuantityKind.NORMAL %} double {{ dq["id"] }}; // {{ dq["label"] }} {% endif %} {% endfor %} diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index c8d393b7..e4234af4 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -9,6 +9,7 @@ PLACEHOLDER_STATE_ID, PREFIX_SEP, ROOT_DIR, + DerivedQuantityKind, EquationType, VarType, ) @@ -101,10 +102,8 @@ def _make_state_variable(*, id, derivative_id): return StateVariable(index=0, id=id, derivative_id=derivative_id, label="", initial_value=None, units="") -def _make_derived_quantity(*, id, is_conversion=False, is_reaction=False): - return DerivedQuantity( - id=id, label="", index=0, initial_value=None, is_conversion=is_conversion, is_reaction=is_reaction, units="" - ) +def _make_derived_quantity(*, id, kind=DerivedQuantityKind.NORMAL): + return DerivedQuantity(id=id, label="", index=0, initial_value=None, units="", kind=kind) def test_check_name_conflicts_passes_for_clean_model(): @@ -117,7 +116,7 @@ def test_check_name_conflicts_flags_synthetic_vs_real_collision(): """A real id equal to a synthesised amount-conversion name is caught.""" builder = _builder_with_names( state_variables=["X"], - derived_quantities=[{"id": "amt__X", "is_conversion": True, "is_reaction": False}], + derived_quantities=[{"id": "amt__X", "kind": DerivedQuantityKind.CONVERSION}], parameters=["amt__X"], ) with pytest.raises(NameConflictError, match="amt__X"): @@ -143,7 +142,7 @@ def test_check_name_conflicts_ignores_reaction_flux_duplicate(): builder = _builder_with_names( state_variables=["C"], reactions=["reaction1"], - derived_quantities=[{"id": "reaction1", "is_conversion": False, "is_reaction": True}], + derived_quantities=[{"id": "reaction1", "kind": DerivedQuantityKind.REACTION}], ) builder._check_name_conflicts() # should not raise @@ -301,7 +300,7 @@ def test_build_extracts_parameters_and_reactions(): def test_build_adds_amount_conversions_for_amount_species(): """Amount species get an amt__ conversion derived quantity.""" builder = _build(REFERENCE / "Goldbeter1991" / "Goldbeter1991.xml") - conversions = {d.id for d in builder._derived_quantities if d.is_conversion} + conversions = {d.id for d in builder._derived_quantities if d.kind == DerivedQuantityKind.CONVERSION} assert conversions == {"amt__C", "amt__M", "amt__X"} @@ -329,8 +328,10 @@ def _event_builder(): builder._keep_doc = doc # keep the owning document (and thus the model) alive builder._sbml_events = builder._sbml_model.getListOfEvents() builder._sbml_function_definitions = builder._sbml_model.getListOfFunctionDefinitions() - builder._variable_types = {} builder._state_variables = [] + builder._parameters = [] + builder._derived_quantities = [] + builder._reactions = [] builder._events = [] builder._functions = [] return builder @@ -404,8 +405,11 @@ def _distance_builder() -> ModelBuilder: doc = libsbml.SBMLDocument(3, 2) builder._sbml_model = doc.createModel() builder._keep_doc = doc # keep the owning document (and thus the model) alive - builder._variable_types = {} builder._state_variables = [] + builder._parameters = [] + builder._derived_quantities = [] + builder._reactions = [] + builder._functions = [] return builder @@ -456,8 +460,11 @@ def _resize_builder(*, assignment_rules=()) -> ModelBuilder: builder._sbml_compartments = model.getListOfCompartments() builder._sbml_species = model.getListOfSpecies() builder._assignment_rules = list(assignment_rules) - builder._variable_types = {"S": VarType.STATE_VARIABLE} - builder._state_variables = [] + builder._state_variables = [_make_state_variable(id="S", derivative_id="d_S_dt")] + builder._parameters = [] + builder._derived_quantities = [] + builder._reactions = [] + builder._functions = [] builder._index_of = {"S": 0} return builder From 86212a80b61fd5b00e284c08ebdaf6361a74db01 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 10:00:23 +0100 Subject: [PATCH 12/19] #42 Fix xor on doubles, piecewise without otherwise, and map arcsec/arcsech - sm::xor_ used bitwise ^ directly on its operands, invalid for the double (numeric-boolean) arguments the generated code passes; cast each to bool first. - Add a 2-argument sm::piecewise overload (a piece with no otherwise clause): the value when the condition holds, NaN otherwise -- also the recursion base case for an even argument count. - Map the arcsec/arcsech MathML functions (missing, unlike the other arc* inverses). Adds C++ SbmlMath tests. Marks 01116 and 01486 passing; 01490/01491 compile now but hit a CVODE convergence failure on their degenerate no-ODE models (re-noted). Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/src/SbmlMath.hpp | 13 ++++++++++++- chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp | 11 +++++++++++ .../test/data/sbml_test_suite_status.csv | 8 ++++---- chaste_sbml/_expressions.py | 2 ++ 4 files changed, 29 insertions(+), 5 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp b/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp index 60ab6272..51dc4494 100644 --- a/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp +++ b/chaste_sbml/SbmlRefModels/src/SbmlMath.hpp @@ -165,6 +165,10 @@ constexpr double min(double first, double second, Args... rest); // piecewise constexpr double piecewise(double otherwise); +// A piecewise with no otherwise clause: undefined (NaN) when the condition does not hold. Also the +// recursion base case for an even number of arguments (pieces without a trailing otherwise). +constexpr double piecewise(double value, bool condition); + constexpr double piecewise(double value, bool condition, double otherwise); template @@ -211,7 +215,9 @@ constexpr bool sbmlmath::or_(Args... args) template constexpr bool sbmlmath::xor_(Args... args) { - return (false ^ ... ^ args); + // Cast each operand to bool first: the arguments are numeric booleans (0/non-0), and ^ is + // bitwise (undefined for double), so xor them as bools -- true when an odd number are true. + return (false ^ ... ^ static_cast(args)); } // Relational (variadic templates) ============== @@ -301,6 +307,11 @@ constexpr double sbmlmath::piecewise(double otherwise) return otherwise; } +constexpr double sbmlmath::piecewise(double value, bool condition) +{ + return condition ? value : std::numeric_limits::quiet_NaN(); +} + constexpr double sbmlmath::piecewise(double value, bool condition, double otherwise) { return condition ? value : otherwise; diff --git a/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp b/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp index 074ef9d2..daa40883 100644 --- a/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp +++ b/chaste_sbml/SbmlRefModels/test/TestSbmlMath.hpp @@ -174,6 +174,11 @@ class TestSbmlMath : public CxxTest::TestSuite TS_ASSERT_EQUALS(sm::xor_(false, false), false); TS_ASSERT_EQUALS(sm::xor_(false, true, false), true); TS_ASSERT_EQUALS(sm::xor_(true, false, true), false); + + // Numeric-boolean (double) operands: cast to bool before xoring. + TS_ASSERT_EQUALS(sm::xor_(1.0, 0.0), true); + TS_ASSERT_EQUALS(sm::xor_(1.0, 1.0), false); + TS_ASSERT_EQUALS(sm::xor_(1.0, 0.0, 1.0), false); } // Relational ================================= @@ -404,6 +409,12 @@ class TestSbmlMath : public CxxTest::TestSuite TS_ASSERT_EQUALS(sm::piecewise(1, true, 2), 1); TS_ASSERT_EQUALS(sm::piecewise(1, false, 2), 2); + + // No otherwise clause: the value when the condition holds, NaN when it does not. + TS_ASSERT_EQUALS(sm::piecewise(1.0, true), 1.0); + TS_ASSERT(std::isnan(sm::piecewise(1.0, false))); + TS_ASSERT_EQUALS(sm::piecewise(1.0, true, 2.0, false), 1.0); + TS_ASSERT(std::isnan(sm::piecewise(1.0, false, 2.0, false))); } void TestQuotient() diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 9dd7426a..31d0f1fa 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1114,7 +1114,7 @@ test,status,notes 01113,pass,ok (no-ode placeholder) 01114,pass,ok (no-ode placeholder) 01115,pass,ok (no-ode placeholder) -01116,fail_known,codegen: piecewise arity +01116,pass,ok (piecewise without otherwise) 01117,pass,ok 01118,pass,ok 01119,unsupported,delay @@ -1484,12 +1484,12 @@ test,status,notes 01483,unsupported,algebraic rule 01484,unsupported,algebraic rule 01485,pass,ok (no-ode placeholder) -01486,fail_known,codegen: unconverted ^ operator +01486,pass,ok (xor / arcsec / uncommon mathml) 01487,pass,ok (no-ode placeholder) 01488,pass,ok (no-ode placeholder) 01489,pass,ok (no-ode placeholder) -01490,fail_known,codegen: unconverted ^ operator -01491,fail_known,codegen: unconverted ^ operator +01490,fail_known,CVODE convergence failure (degenerate no-ode model) +01491,fail_known,CVODE convergence failure (degenerate no-ode model) 01492,pass,ok 01493,pass,ok 01494,pass,ok (nary relational via sm:: helpers) diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index a6aa45f5..591b4c6f 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -305,6 +305,8 @@ def rewrite_power(node: Optional["ASTNode"]) -> None: "arccoth": "acoth", "arccsc": "acsc", "arccsch": "acsch", + "arcsec": "asec", + "arcsech": "asech", "cot": "cot", "coth": "coth", "csc": "csc", From 6ef1b0935c19f4096ee892185db914ea460ae689 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 10:33:29 +0100 Subject: [PATCH 13/19] #42 Rename SBML ids that collide with C++ macros An SBML id spelled like a C++ macro (e.g. a parameter named NAN or INFINITY) cannot be emitted verbatim: `double NAN;` is rewritten by the preprocessor into the macro's expansion (which contains a string literal). Treat these macros like keywords/reserved names in NameManager, so such an id is escaped for the emitted code (NAN -> NAN_) while still reported to Chaste under its original id. Marks 01812 passing; 01813 stays fail_known -- its reference output for a initial assignment is genuinely NaN, which the finite-tolerance test cannot validate (non-finite expected results). Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 4 ++-- chaste_sbml/_names.py | 11 +++++++--- chaste_sbml/tests/test_names.py | 21 +++++++++++++++++++ 3 files changed, 31 insertions(+), 5 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 31d0f1fa..0fe31e65 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1810,8 +1810,8 @@ test,status,notes 01809,pass,ok 01810,pass,ok (no-ode placeholder) 01811,fail_known,non-finite expected results (INF/NaN) -01812,fail_known,codegen: string literal in math -01813,fail_known,codegen: string literal in math +01812,pass,ok (id spelled like a C++ macro) +01813,fail_known,non-finite expected results (INF/NaN) 01814,pass,ok (no-ode placeholder) 01815,pass,ok (no-ode placeholder) 01816,pass,ok (no-ode placeholder) diff --git a/chaste_sbml/_names.py b/chaste_sbml/_names.py index bbc35dce..3a6a6dff 100644 --- a/chaste_sbml/_names.py +++ b/chaste_sbml/_names.py @@ -156,8 +156,13 @@ } ) +# C++ macros (from /) that expand when used as a bare identifier, so an SBML id +# spelled like one cannot be emitted verbatim (e.g. a parameter NAN becomes `double NAN;`, which +# the preprocessor rewrites into the macro's expansion). +CPP_MACROS = frozenset({"NAN", "INFINITY", "HUGE_VAL", "EOF", "NULL"}) + # Every name the generator may not use for one of its own identifiers. -RESERVED_NAMES = CPP_KEYWORDS | CHASTE_RESERVED_NAMES +RESERVED_NAMES = CPP_KEYWORDS | CHASTE_RESERVED_NAMES | CPP_MACROS _IDENTIFIER_RE = re.compile(r"^[A-Za-z_][A-Za-z0-9_]*$") @@ -239,7 +244,7 @@ def resolve_cpp_name(base: str, taken) -> str: :return: A C++-safe identifier not present in ``taken``. """ candidate = base - while candidate in CPP_KEYWORDS or candidate in CHASTE_RESERVED_NAMES: + while candidate in CPP_KEYWORDS or candidate in CHASTE_RESERVED_NAMES or candidate in CPP_MACROS: candidate += "_" return unique_name(candidate, taken) @@ -293,7 +298,7 @@ def resolve_real_id_conflicts(self) -> None: if not element.isSetId(): continue old_id = element.getId() - if old_id not in CPP_KEYWORDS and old_id not in CHASTE_RESERVED_NAMES: + if old_id not in CPP_KEYWORDS and old_id not in CHASTE_RESERVED_NAMES and old_id not in CPP_MACROS: continue new_id = resolve_cpp_name(old_id, taken) element.setId(new_id) diff --git a/chaste_sbml/tests/test_names.py b/chaste_sbml/tests/test_names.py index cb8be04e..d2a57cf8 100644 --- a/chaste_sbml/tests/test_names.py +++ b/chaste_sbml/tests/test_names.py @@ -187,6 +187,27 @@ def test_name_manager_sbml_name_recovers_original_id(): assert manager.sbml_name("k") == "k" # an un-renamed id is returned unchanged +def test_name_manager_resolve_renames_cpp_macro(): + """An id spelling a C++ macro (e.g. NAN) is renamed but still reported under its original id. + + ``double NAN;`` would be rewritten by the preprocessor into the macro's expansion, so the id + must be escaped for the emitted code while keeping ``NAN`` as its reported name. + """ + doc = libsbml.SBMLDocument(3, 2) + sbml_model = doc.createModel() + parameter = sbml_model.createParameter() + parameter.setId("NAN") # the NAN macro from + parameter.setConstant(True) + parameter.setValue(0.1) + + manager = NameManager(sbml_model) + manager.resolve_real_id_conflicts() + + assert sbml_model.getParameter("NAN") is None + assert sbml_model.getParameter("NAN_") is not None + assert manager.sbml_name("NAN_") == "NAN" + + def test_generate_header_guard(): """Test header guard generation.""" test_cases = [ From 764c88afb97272d11bd297312166445a244b0d33 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 11:00:02 +0100 Subject: [PATCH 14/19] #42 Rewrite function-call refs when renaming a keyword-clashing function libsbml's renameSIdRefs updates variable references but not function-call names (a call f(x) is an AST_FUNCTION node, not an AST_NAME), so renaming a function definition to escape a keyword/reserved/macro clash left its callers invoking the old, now-invalid name. Walk the model's math and rewrite those calls explicitly. Fixes case 01312 (nested function definitions with a function named `this`). Co-Authored-By: Claude Opus 4.8 --- .../test/data/sbml_test_suite_status.csv | 2 +- chaste_sbml/_names.py | 27 +++++++++++++++++++ chaste_sbml/tests/test_names.py | 24 +++++++++++++++++ 3 files changed, 52 insertions(+), 1 deletion(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 0fe31e65..5a2a4ee3 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1310,7 +1310,7 @@ test,status,notes 01309,pass,ok 01310,pass,ok (no-ode placeholder) 01311,pass,ok (no-ode placeholder) -01312,fail_known,codegen: expression called as function +01312,pass,ok (nested function definitions; a function named after a C++ keyword) 01313,pass,ok (no-ode placeholder) 01314,pass,ok (no-ode placeholder) 01315,pass,ok (no-ode placeholder) diff --git a/chaste_sbml/_names.py b/chaste_sbml/_names.py index 3a6a6dff..1c940376 100644 --- a/chaste_sbml/_names.py +++ b/chaste_sbml/_names.py @@ -23,6 +23,8 @@ import re +import libsbml + # C++ keywords (and alternative operator spellings); none may be used as an identifier. CPP_KEYWORDS = frozenset( { @@ -304,9 +306,34 @@ def resolve_real_id_conflicts(self) -> None: element.setId(new_id) for referrer in model.getListOfAllElements(): referrer.renameSIdRefs(old_id, new_id) + # ``renameSIdRefs`` rewrites variable references but not function-call names (a call + # ``f(x)`` is an ``AST_FUNCTION`` node, not an ``AST_NAME``), so a renamed function + # definition would leave its callers dangling. Rewrite those calls explicitly. + self._rename_function_call_refs(old_id, new_id) taken.add(new_id) self._sbml_names[new_id] = old_id + def _rename_function_call_refs(self, old_id: str, new_id: str) -> None: + """Rename every ``old_id(...)`` function call in the model's math to ``new_id(...)``. + + Function calls reference a function definition by name through an ``AST_FUNCTION`` node, + which ``renameSIdRefs`` leaves untouched; without this a function definition renamed to + avoid a keyword/reserved/macro clash (e.g. one called ``this``) keeps being called under + its old, now-invalid name. + """ + + def rename(node: "libsbml.ASTNode") -> None: + if node is None: + return + if node.getType() == libsbml.AST_FUNCTION and node.getName() == old_id: + node.setName(new_id) + for i in range(node.getNumChildren()): + rename(node.getChild(i)) + + for element in self._sbml_model.getListOfAllElements(): + if hasattr(element, "isSetMath") and element.isSetMath(): + rename(element.getMath()) + def sbml_name(self, cpp_id: str) -> str: """Return the original SBML id for a (possibly renamed) C++ identifier. diff --git a/chaste_sbml/tests/test_names.py b/chaste_sbml/tests/test_names.py index d2a57cf8..ca4bbf28 100644 --- a/chaste_sbml/tests/test_names.py +++ b/chaste_sbml/tests/test_names.py @@ -208,6 +208,30 @@ def test_name_manager_resolve_renames_cpp_macro(): assert manager.sbml_name("NAN_") == "NAN" +def test_name_manager_resolve_renames_function_calls(): + """A function definition renamed to dodge a keyword also has its callers rewritten. + + A call ``this(x)`` is an ``AST_FUNCTION`` node, which ``renameSIdRefs`` leaves untouched, so + renaming the ``this`` function definition must explicitly rewrite the call inside ``caller`` -- + otherwise the emitted C++ would call a function that no longer exists (and ``this`` is a keyword). + """ + doc = libsbml.SBMLDocument(3, 2) + sbml_model = doc.createModel() + callee = sbml_model.createFunctionDefinition() + callee.setId("this") # 'this' is a C++ keyword + callee.setMath(libsbml.parseL3Formula("lambda(x, x)")) + caller = sbml_model.createFunctionDefinition() + caller.setId("caller") + caller.setMath(libsbml.parseL3Formula("lambda(x, this(x))")) + + NameManager(sbml_model).resolve_real_id_conflicts() + + assert sbml_model.getFunctionDefinition("this") is None + assert sbml_model.getFunctionDefinition("this_") is not None + # The body of 'caller' now calls the renamed function. + assert "this_(x)" in libsbml.formulaToL3String(sbml_model.getFunctionDefinition("caller").getMath()) + + def test_generate_header_guard(): """Test header guard generation.""" test_cases = [ From 7a7c7afbffef1111daf852cbbb0a3747e4b3b850 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 11:16:47 +0100 Subject: [PATCH 15/19] #42 Handle non-finite values in results and MathML constants Reference results and models can be non-finite in two ways, both mishandled: - Expected-results columns may be INF/-INF/NaN (a division by zero, 0/0). These were emitted verbatim into the C++ initialiser list (not valid tokens), and a delta/tolerance check cannot validate them anyway (|INF-INF| is NaN). Emit the INFINITY/NAN macros, and compare non-finite expected values by kind and sign instead of by tolerance. - The MathML constants and are stored by libsbml as ordinary real nodes that render as "INF"/"NaN" -- textually identical to a reference to a parameter so named -- so after stringification they collided with such a parameter and resolved to its value. Detect them by their non-finite value on the AST and rename to placeholders, exactly as pi/avogadro are, so the constant and a same-named parameter stay distinct. Fixes cases 00950, 00951, 01811, 01813. Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/generate_cases.py | 24 +++++++++++-- .../test/data/sbml_test_suite_status.csv | 8 ++--- chaste_sbml/_expressions.py | 36 ++++++++++++------- chaste_sbml/templates/cases/semantic.hpp | 24 ++++++++++--- chaste_sbml/tests/test_expressions.py | 22 ++++++++++-- 5 files changed, 90 insertions(+), 24 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/generate_cases.py b/chaste_sbml/SbmlRefModels/generate_cases.py index 044c052b..e20e57e0 100644 --- a/chaste_sbml/SbmlRefModels/generate_cases.py +++ b/chaste_sbml/SbmlRefModels/generate_cases.py @@ -9,7 +9,7 @@ from typing import TYPE_CHECKING from chaste_sbml import ChasteSbmlModel -from chaste_sbml._config import ROOT_DIR, ModelType +from chaste_sbml._config import ModelType from chaste_sbml._names import generate_header_guard if TYPE_CHECKING: @@ -20,6 +20,24 @@ sbml_versions = ["l2v5", "l3v2", "l3v1"] +def _cpp_double_literal(value: str) -> str: + """Render a reference-CSV cell as a C++ double literal. + + Reference results can be non-finite -- a division by zero gives INF and 0/0 gives NaN -- and + the SBML Test Suite writes these as ``INF``/``-INF``/``NaN``, which are not valid C++ tokens. + Map them to the ```` macros so the expected-results initialiser list compiles; finite + numbers are already valid literals and pass through unchanged. + """ + token = value.strip() + magnitude = token.lstrip("+-").lower() + sign = "-" if token.startswith("-") else "" + if magnitude == "inf": + return f"{sign}INFINITY" + if magnitude == "nan": + return "NAN" + return token + + class TestType(Enum): """Enumeration of test types for code generation.""" @@ -50,7 +68,9 @@ def __init__(self, sbml_file: str, sbml_version: str, test_params: dict[str, "An self._test_hpp_filename = f"Test{self._model_name}.hpp" test_result_columns = ", ".join(f'"{col}"' for col in test_params["results"][0]) - test_result_data = ["{ " + ", ".join(row) + " }" for row in test_params["results"][1:]] + test_result_data = [ + "{ " + ", ".join(_cpp_double_literal(cell) for cell in row) + " }" for row in test_params["results"][1:] + ] test_result_data = ",\n".join(test_result_data) test_amounts = test_params["settings"]["amount"].split(",") diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index 5a2a4ee3..bd5f2b9e 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -948,8 +948,8 @@ test,status,notes 00947,pass,ok 00948,pass,ok 00949,pass,ok (no-ode placeholder) -00950,fail_known,non-finite expected results (INF/NaN) -00951,fail_known,non-finite expected results (INF/NaN) +00950,pass,ok (non-finite results: INF/-INF/NaN) +00951,pass,ok (non-finite results: INF/-INF/NaN) 00952,unsupported,random event execution 00953,unsupported,event execution semantics (competing simultaneous events) 00954,pass,ok (no-ode placeholder) @@ -1809,9 +1809,9 @@ test,status,notes 01808,pass,ok 01809,pass,ok 01810,pass,ok (no-ode placeholder) -01811,fail_known,non-finite expected results (INF/NaN) +01811,pass,ok (parameters named INF/Inf/inf vs the infinity constant) 01812,pass,ok (id spelled like a C++ macro) -01813,fail_known,non-finite expected results (INF/NaN) +01813,pass,ok (parameters named NAN/NaN/nan vs the notanumber constant) 01814,pass,ok (no-ode placeholder) 01815,pass,ok (no-ode placeholder) 01816,pass,ok (no-ode placeholder) diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index 591b4c6f..11a0cb79 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -6,6 +6,7 @@ into the model. """ +import math import re from typing import TYPE_CHECKING, Optional @@ -81,6 +82,9 @@ def substitute_ast_names(node: "ASTNode", replacements: dict) -> "ASTNode": PI_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}pi" E_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}exponentiale" TIME_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}time" +INFINITY_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}infinity" +NEG_INFINITY_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}neg_infinity" +NAN_PLACEHOLDER = f"{CHASTE_PREFIX}{PREFIX_SEP}nan" _SYMBOL_PLACEHOLDERS = { AST_NAME_AVOGADRO: AVOGADRO_PLACEHOLDER, @@ -91,12 +95,15 @@ def substitute_ast_names(node: "ASTNode", replacements: dict) -> "ASTNode": def replace_constant_symbols(node: "ASTNode") -> None: - """Recursively rename SBML math-symbol nodes (avogadro, pi, exponentiale, time) to placeholders. + """Recursively rename SBML math-symbol nodes (avogadro, pi, exponentiale, time, inf, NaN). Each of these symbols renders identically to a same-named parameter, so it is renamed to a unique placeholder that maps to its C++ value, keeping it distinct from a parameter of the same spelling. Handling this at the AST level (rather than string-matching the rendered name) means a - real variable named ``t``, ``pi`` etc. is never mistaken for the symbol. + real variable named ``t``, ``pi`` etc. is never mistaken for the symbol. The MathML constants + ```` and ```` are especially treacherous: libsbml stores them as ordinary + real nodes that render as ``INF``/``NaN``, indistinguishable from a reference to a parameter so + named (see test-suite cases 1811/1813), so they are matched here by their non-finite value. :param node: The root AST node. """ @@ -104,6 +111,14 @@ def replace_constant_symbols(node: "ASTNode") -> None: if placeholder is not None: node.setType(AST_NAME) node.setName(placeholder) + elif node.getType() == AST_REAL: + value = node.getReal() + if math.isnan(value): + node.setType(AST_NAME) + node.setName(NAN_PLACEHOLDER) + elif math.isinf(value): + node.setType(AST_NAME) + node.setName(NEG_INFINITY_PLACEHOLDER if value < 0 else INFINITY_PLACEHOLDER) for i in range(node.getNumChildren()): replace_constant_symbols(node.getChild(i)) @@ -239,21 +254,18 @@ def rewrite_power(node: Optional["ASTNode"]) -> None: rewrite_power(node.getChild(i)) -# SBML symbolic constants replaced with their C++ equivalents, keyed by the -# placeholder each symbol is renamed to (see replace_constant_symbols) plus the -# infinity/NaN literals formulaToL3String emits. (avogadro/pi/exponentiale/time -# are handled via placeholders; true/false deliberately absent.) +# SBML symbolic constants replaced with their C++ equivalents, keyed by the placeholder each symbol +# is renamed to (see replace_constant_symbols). Every one is matched at the AST level and renamed to +# its placeholder, so a same-named parameter (e.g. one literally called pi, inf or NaN) is never +# mistaken for the constant. (true/false deliberately absent.) _CONSTANTS = { AVOGADRO_PLACEHOLDER: "sm::AVOGADRO", PI_PLACEHOLDER: "M_PI", E_PLACEHOLDER: "M_E", TIME_PLACEHOLDER: "time", - "inf": "std::numeric_limits::infinity()", - "INF": "std::numeric_limits::infinity()", - "infinity": "std::numeric_limits::infinity()", - "nan": "NAN", - "NaN": "NAN", - "notanumber": "NAN", + INFINITY_PLACEHOLDER: "std::numeric_limits::infinity()", + NEG_INFINITY_PLACEHOLDER: "-std::numeric_limits::infinity()", + NAN_PLACEHOLDER: "std::numeric_limits::quiet_NaN()", } # SBML functions whose name matches the C++ (std::) equivalent. diff --git a/chaste_sbml/templates/cases/semantic.hpp b/chaste_sbml/templates/cases/semantic.hpp index fedd2f0e..36299f24 100644 --- a/chaste_sbml/templates/cases/semantic.hpp +++ b/chaste_sbml/templates/cases/semantic.hpp @@ -160,11 +160,27 @@ class Test{{ model_name }} : public CxxTest::TestSuite TS_ASSERT_EQUALS(values.size(), expected_result_data.size()); for (unsigned i = 0; i < expected_result_data.size(); i++) { - double delta = std::abs(expected_result_data[i][j] - values[i]); - double tol = tol_absolute + tol_relative * std::abs(expected_result_data[i][j]); - std::string msg(sth::ToString(values[i]) + " vs " + sth::ToString(expected_result_data[i][j]) + double expected = expected_result_data[i][j]; + double actual = values[i]; + std::string msg(sth::ToString(actual) + " vs " + sth::ToString(expected) + " at " + sth::ToString(ode_solution.rGetTimes()[i], 3) + " for " + var_name); - TSM_ASSERT_LESS_THAN_EQUALS(msg.c_str(), delta, tol); + // Non-finite expected values (e.g. a division by zero gives +/-INF, 0/0 gives + // NaN) cannot be checked with a tolerance: |INF - INF| is NaN, and any + // comparison against NaN is false. Match them by kind (and sign) instead. + if (std::isnan(expected)) + { + TSM_ASSERT(msg.c_str(), std::isnan(actual)); + } + else if (std::isinf(expected)) + { + TSM_ASSERT(msg.c_str(), std::isinf(actual) && std::signbit(actual) == std::signbit(expected)); + } + else + { + double delta = std::abs(expected - actual); + double tol = tol_absolute + tol_relative * std::abs(expected); + TSM_ASSERT_LESS_THAN_EQUALS(msg.c_str(), delta, tol); + } } } diff --git a/chaste_sbml/tests/test_expressions.py b/chaste_sbml/tests/test_expressions.py index 44f9cd14..70d57838 100644 --- a/chaste_sbml/tests/test_expressions.py +++ b/chaste_sbml/tests/test_expressions.py @@ -189,15 +189,33 @@ def test_replace_constant_symbols_leaves_other_names(): ("", "M_PI"), ("", "M_E"), ("", "std::numeric_limits::infinity()"), - ("", "NAN"), + ("", "std::numeric_limits::quiet_NaN()"), ], ) def test_formula_to_string_math_constants(mathml_symbol, expected): - """SBML math constants (incl. the capitalised INF/NaN libsbml emits) map to their C++ form.""" + """SBML math constants map to their C++ form.""" math = libsbml.readMathMLFromString(f'{mathml_symbol}') assert formula_to_string(math, {}, state_variables=[]) == expected +@pytest.mark.parametrize( + ("mathml_symbol", "param", "expected"), + [ + ("", "INF", "std::numeric_limits::infinity()"), + ("", "NaN", "std::numeric_limits::quiet_NaN()"), + ], +) +def test_formula_to_string_infinity_nan_distinct_from_parameter(mathml_symbol, param, expected): + """The / constants stay distinct from a parameter spelled INF/NaN. + + libsbml stores these constants as real nodes that render as ``INF``/``NaN`` -- the same text as + a reference to a parameter so named -- so the constant must be resolved from the AST, not its + rendered spelling (test-suite cases 1811/1813). + """ + math = libsbml.readMathMLFromString(f'{mathml_symbol}') + assert formula_to_string(math, {param: VarType.PARAMETER}, state_variables=[]) == expected + + def test_formula_to_string_pi_constant_and_parameter_are_distinct(): """The constant maps to M_PI while a same-named parameter keeps its own name.""" math = libsbml.readMathMLFromString( From fdd5e88fdb4e2d1c89d80efa89e98e1bf589ef85 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 11:37:21 +0100 Subject: [PATCH 16/19] #42 Initialise the no-ODE placeholder state variable Initialise() reads each state variable's member to set the solver's initial condition via SetDefaultInitialCondition. Real state variables are assigned their value first (by an INITIAL_VALUE equation), but the synthesised no-ODE placeholder had only a derivative equation -- its member was never assigned, so its initial condition was read from uninitialised memory. When that memory happened to hold a NaN, CVODE integrated from a NaN state and failed with CV_CONV_FAILURE (undefined behaviour that only surfaced on the no-ODE cases that integrate over time rather than solving for a steady state). Add an INITIAL_VALUE equation so the placeholder member is deterministically zeroed before it is used. Fixes cases 01490, 01491. Co-Authored-By: Claude Opus 4.8 --- .../SbmlRefModels/test/data/sbml_test_suite_status.csv | 4 ++-- chaste_sbml/_model_builder.py | 2 ++ chaste_sbml/tests/test_model_builder.py | 6 ++++++ 3 files changed, 10 insertions(+), 2 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index bd5f2b9e..a67e99bc 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -1488,8 +1488,8 @@ test,status,notes 01487,pass,ok (no-ode placeholder) 01488,pass,ok (no-ode placeholder) 01489,pass,ok (no-ode placeholder) -01490,fail_known,CVODE convergence failure (degenerate no-ode model) -01491,fail_known,CVODE convergence failure (degenerate no-ode model) +01490,pass,ok (no-ODE model integrated over time; placeholder initial condition) +01491,pass,ok (no-ODE model integrated over time; placeholder initial condition) 01492,pass,ok 01493,pass,ok 01494,pass,ok (nary relational via sm:: helpers) diff --git a/chaste_sbml/_model_builder.py b/chaste_sbml/_model_builder.py index f5900691..864eb5a0 100644 --- a/chaste_sbml/_model_builder.py +++ b/chaste_sbml/_model_builder.py @@ -334,6 +334,8 @@ def _add_placeholder_state_variable(self) -> None: """ placeholder_id = self._names.reserve(PREFIX_SEP.join([CHASTE_PREFIX, PLACEHOLDER_STATE_ID])) state_var = self._add_state_variable(placeholder_id, "Placeholder state variable (model has no ODEs)", 0.0) + # Initialise the member to zero + self._add_equation(var=placeholder_id, math=parseL3Formula("0.0"), eq_type=EquationType.INITIAL_VALUE) self._add_equation(var=state_var.derivative_id, math=parseL3Formula("0.0"), eq_type=EquationType.DERIVATIVE) def _add_stoichiometry_variable(self, species_reference: "SpeciesReference") -> None: diff --git a/chaste_sbml/tests/test_model_builder.py b/chaste_sbml/tests/test_model_builder.py index e4234af4..ea76aff2 100644 --- a/chaste_sbml/tests/test_model_builder.py +++ b/chaste_sbml/tests/test_model_builder.py @@ -228,6 +228,12 @@ def test_build_synthesises_placeholder_state_variable_for_no_ode_model(tmp_path) assert deriv_equations[0].var == placeholder.derivative_id assert deriv_equations[0].rhs == "0.0" + # The placeholder's member must be initialised to zero: Initialise() reads it to set the + # solver's initial condition, so without this equation the state would be uninitialised memory. + init_equations = [e for e in builder._equations if e.type == EquationType.INITIAL_VALUE] + assert placeholder.id in [e.var for e in init_equations] + assert next(e for e in init_equations if e.var == placeholder.id).rhs == "0.0" + def test_build_no_ode_model_preserves_real_outputs(tmp_path): """The placeholder does not displace the model's real outputs (assignment rule, species).""" From 1648c5d76649623f47876f608aeb33a43138ab15 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 12:15:39 +0100 Subject: [PATCH 17/19] Tidy SBML test-suite status notes - Collapse all passing-case notes to "ok"; once a case is green the specific reason it passes is not worth tracking. - Regroup the unsupported notes so each names one unsupported SBML feature and cases are grouped consistently: * merge "delay function" into "delay" -- the event element and the delay() csymbol are unsupported for the same reason (no delayed/historical values) and the split was already inconsistent; * consolidate the deterministic event cases under "event execution semantics" (the previous sub-notes misfiled persistence and t0-firing cases as priority-ordering), keeping the stochastic "random event execution" cases separate. - Fix case 01197: it is a trivial 7-dimension-compartment model that passes, not a flux-balance case -> pass/ok. - Document in src/cases/README.md (matching test/cases) where case support is tracked. Co-Authored-By: Claude Opus 4.8 --- chaste_sbml/SbmlRefModels/src/cases/README.md | 7 +- .../test/data/sbml_test_suite_status.csv | 432 +++++++++--------- 2 files changed, 221 insertions(+), 218 deletions(-) diff --git a/chaste_sbml/SbmlRefModels/src/cases/README.md b/chaste_sbml/SbmlRefModels/src/cases/README.md index f23f57cb..951602bf 100644 --- a/chaste_sbml/SbmlRefModels/src/cases/README.md +++ b/chaste_sbml/SbmlRefModels/src/cases/README.md @@ -1,7 +1,10 @@ # SBML Test Suite Cases ## Semantic - -## Stochastic +Case support tracked in test/data/sbml_test_suite_status.csv ## Syntactic +Covered by libSBML + +## Stochastic +Not supported diff --git a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv index a67e99bc..ab3c1b76 100644 --- a/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv +++ b/chaste_sbml/SbmlRefModels/test/data/sbml_test_suite_status.csv @@ -7,11 +7,11 @@ test,status,notes 00006,pass,ok 00007,pass,ok 00008,pass,ok -00009,pass,ok (no-ode placeholder) +00009,pass,ok 00010,pass,ok 00011,pass,ok 00012,pass,ok -00013,pass,ok (no-ode placeholder) +00013,pass,ok 00014,pass,ok 00015,pass,ok 00016,pass,ok @@ -27,8 +27,8 @@ test,status,notes 00026,pass,ok 00027,pass,ok 00028,pass,ok -00029,pass,ok (no-ode placeholder) -00030,pass,ok (no-ode placeholder) +00029,pass,ok +00030,pass,ok 00031,pass,ok 00032,pass,ok 00033,pass,ok @@ -116,7 +116,7 @@ test,status,notes 00115,pass,ok 00116,pass,ok 00117,pass,ok -00118,pass,ok (no-ode placeholder) +00118,pass,ok 00119,pass,ok 00120,pass,ok 00121,pass,ok @@ -172,7 +172,7 @@ test,status,notes 00171,pass,ok 00172,pass,ok 00173,pass,ok -00174,pass,ok (no-ode placeholder) +00174,pass,ok 00175,pass,ok 00176,pass,ok 00177,pass,ok @@ -211,7 +211,7 @@ test,status,notes 00210,pass,ok 00211,pass,ok 00212,pass,ok -00213,pass,ok (no-ode placeholder) +00213,pass,ok 00214,pass,ok 00215,pass,ok 00216,pass,ok @@ -226,7 +226,7 @@ test,status,notes 00225,pass,ok 00226,pass,ok 00227,pass,ok -00228,pass,ok (no-ode placeholder) +00228,pass,ok 00229,pass,ok 00230,pass,ok 00231,pass,ok @@ -577,7 +577,7 @@ test,status,notes 00576,unsupported,algebraic rule 00577,pass,ok 00578,pass,ok -00579,pass,ok (no-ode placeholder) +00579,pass,ok 00580,pass,ok 00581,pass,ok 00582,pass,ok @@ -918,54 +918,54 @@ test,status,notes 00917,pass,ok 00918,pass,ok 00919,pass,ok -00920,pass,ok (no-ode placeholder) -00921,pass,ok (no-ode placeholder) -00922,pass,ok (no-ode placeholder) -00923,pass,ok (no-ode placeholder) -00924,pass,ok (no-ode placeholder) -00925,pass,ok (no-ode placeholder) +00920,pass,ok +00921,pass,ok +00922,pass,ok +00923,pass,ok +00924,pass,ok +00925,pass,ok 00926,pass,ok 00927,pass,ok 00928,pass,ok 00929,pass,ok -00930,pass,ok (event-modified species as parameter) -00931,pass,ok (event-modified species as parameter) +00930,pass,ok +00931,pass,ok 00932,unsupported,delay 00933,unsupported,delay -00934,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) -00935,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) +00934,unsupported,event execution semantics +00935,unsupported,event execution semantics 00936,unsupported,delay -00937,unsupported,delay function +00937,unsupported,delay 00938,unsupported,delay 00939,unsupported,delay -00940,unsupported,delay function -00941,unsupported,delay function -00942,unsupported,delay function -00943,unsupported,delay function +00940,unsupported,delay +00941,unsupported,delay +00942,unsupported,delay +00943,unsupported,delay 00944,pass,ok 00945,pass,ok 00946,pass,ok 00947,pass,ok 00948,pass,ok -00949,pass,ok (no-ode placeholder) -00950,pass,ok (non-finite results: INF/-INF/NaN) -00951,pass,ok (non-finite results: INF/-INF/NaN) +00949,pass,ok +00950,pass,ok +00951,pass,ok 00952,unsupported,random event execution -00953,unsupported,event execution semantics (competing simultaneous events) -00954,pass,ok (no-ode placeholder) -00955,pass,ok (no-ode placeholder) -00956,pass,ok (no-ode placeholder) -00957,pass,ok (no-ode placeholder) -00958,pass,ok (no-ode placeholder) -00959,pass,ok (no-ode placeholder) -00960,pass,ok (no-ode placeholder) -00961,pass,ok (no-ode placeholder) +00953,unsupported,event execution semantics +00954,pass,ok +00955,pass,ok +00956,pass,ok +00957,pass,ok +00958,pass,ok +00959,pass,ok +00960,pass,ok +00961,pass,ok 00962,unsupported,random event execution -00963,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) +00963,unsupported,event execution semantics 00964,unsupported,random event execution 00965,unsupported,random event execution 00966,unsupported,random event execution -00967,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) +00967,unsupported,event execution semantics 00968,pass,ok 00969,pass,ok 00970,pass,ok @@ -976,8 +976,8 @@ test,status,notes 00975,pass,ok 00976,pass,ok 00977,pass,ok -00978,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) -00979,pass,ok (no-ode placeholder) +00978,unsupported,event execution semantics +00979,pass,ok 00980,unsupported,delay 00981,unsupported,delay 00982,unsupported,delay @@ -993,9 +993,9 @@ test,status,notes 00992,pass,ok 00993,unsupported,algebraic rule 00994,pass,ok -00995,pass,ok (no-ode placeholder) +00995,pass,ok 00996,pass,ok -00997,pass,ok (no-ode placeholder) +00997,pass,ok 00998,pass,ok 00999,pass,ok 01000,unsupported,delay @@ -1110,11 +1110,11 @@ test,status,notes 01109,pass,ok 01110,pass,ok 01111,pass,ok -01112,pass,ok (no-ode placeholder) -01113,pass,ok (no-ode placeholder) -01114,pass,ok (no-ode placeholder) -01115,pass,ok (no-ode placeholder) -01116,pass,ok (piecewise without otherwise) +01112,pass,ok +01113,pass,ok +01114,pass,ok +01115,pass,ok +01116,pass,ok 01117,pass,ok 01118,pass,ok 01119,unsupported,delay @@ -1122,21 +1122,21 @@ test,status,notes 01121,pass,ok 01122,pass,ok 01123,pass,ok -01124,pass,ok (no-ode placeholder) -01125,pass,ok (no-ode placeholder) -01126,pass,ok (no-ode placeholder) +01124,pass,ok +01125,pass,ok +01126,pass,ok 01127,pass,ok -01128,pass,ok (no-ode placeholder) -01129,pass,ok (no-ode placeholder) -01130,pass,ok (no-ode placeholder) -01131,pass,ok (no-ode placeholder) -01132,pass,ok (no-ode placeholder) -01133,pass,ok (no-ode placeholder) -01134,pass,ok (no-ode placeholder) -01135,pass,ok (no-ode placeholder) -01136,pass,ok (no-ode placeholder) -01137,pass,ok (no-ode placeholder) -01138,pass,ok (no-ode placeholder) +01128,pass,ok +01129,pass,ok +01130,pass,ok +01131,pass,ok +01132,pass,ok +01133,pass,ok +01134,pass,ok +01135,pass,ok +01136,pass,ok +01137,pass,ok +01138,pass,ok 01139,pass,ok 01140,pass,ok 01141,pass,ok @@ -1144,23 +1144,23 @@ test,status,notes 01143,pass,ok 01144,pass,ok 01145,pass,ok -01146,pass,ok (no-ode placeholder) -01147,pass,ok (no-ode placeholder) +01146,pass,ok +01147,pass,ok 01148,pass,ok -01149,pass,ok (no-ode placeholder) +01149,pass,ok 01150,pass,ok 01151,pass,ok -01152,pass,ok (no-ode placeholder) -01153,pass,ok (no-ode placeholder) -01154,pass,ok (no-ode placeholder) -01155,pass,ok (no-ode placeholder) +01152,pass,ok +01153,pass,ok +01154,pass,ok +01155,pass,ok 01156,unsupported,delay -01157,pass,ok (no-ode placeholder) -01158,pass,ok (no-ode placeholder) +01157,pass,ok +01158,pass,ok 01159,pass,ok 01160,pass,ok 01161,pass,ok -01162,pass,ok (no-ode placeholder) +01162,pass,ok 01163,pass,ok 01164,pass,ok 01165,pass,ok @@ -1169,9 +1169,9 @@ test,status,notes 01168,pass,ok 01169,pass,ok 01170,pass,ok -01171,pass,ok (no-ode placeholder) +01171,pass,ok 01172,pass,ok -01173,unsupported,delay function +01173,unsupported,delay 01174,unsupported,algebraic rule 01175,pass,ok 01176,unsupported,delay @@ -1195,7 +1195,7 @@ test,status,notes 01194,unsupported,flux balance 01195,unsupported,flux balance 01196,unsupported,flux balance -01197,unsupported,flux balance +01197,pass,ok 01198,pass,ok 01199,pass,ok 01200,pass,ok @@ -1207,16 +1207,16 @@ test,status,notes 01206,pass,ok 01207,pass,ok 01208,pass,ok -01209,pass,ok (no-ode placeholder) -01210,pass,ok (no-ode placeholder) -01211,pass,ok (no-ode placeholder) -01212,pass,ok (nary relational via sm:: helpers) +01209,pass,ok +01210,pass,ok +01211,pass,ok +01212,pass,ok 01213,unsupported,delay -01214,pass,ok (no-ode placeholder) +01214,pass,ok 01215,pass,ok -01216,pass,ok (nary relational via sm:: helpers) -01217,pass,ok (no-ode placeholder) -01218,pass,ok (no-ode placeholder) +01216,pass,ok +01217,pass,ok +01218,pass,ok 01219,pass,ok 01220,pass,ok 01221,pass,ok @@ -1232,22 +1232,22 @@ test,status,notes 01231,pass,ok 01232,pass,ok 01233,pass,ok -01234,pass,ok (no-ode placeholder) -01235,pass,ok (no-ode placeholder) -01236,pass,ok (no-ode placeholder) -01237,pass,ok (no-ode placeholder) -01238,pass,ok (missing mathml treated as no-op) -01239,pass,ok (missing mathml treated as no-op) -01240,pass,ok (no-ode placeholder) -01241,pass,ok (missing mathml treated as no-op) -01242,pass,ok (no-ode placeholder) -01243,pass,ok (no-ode placeholder) +01234,pass,ok +01235,pass,ok +01236,pass,ok +01237,pass,ok +01238,pass,ok +01239,pass,ok +01240,pass,ok +01241,pass,ok +01242,pass,ok +01243,pass,ok 01244,unsupported,algebraic rule -01245,pass,ok (no-ode placeholder) +01245,pass,ok 01246,pass,ok -01247,pass,ok (no-ode placeholder) -01248,pass,ok (no-ode placeholder) -01249,pass,ok (no-ode placeholder) +01247,pass,ok +01248,pass,ok +01249,pass,ok 01250,pass,ok 01251,pass,ok 01252,pass,ok @@ -1269,27 +1269,27 @@ test,status,notes 01268,unsupported,delay 01269,unsupported,delay 01270,unsupported,delay -01271,pass,ok (missing mathml treated as no-op) -01272,pass,ok (no-ode placeholder) -01273,pass,ok (no-ode placeholder) -01274,pass,ok (implies operator) -01275,pass,ok (implies operator) -01276,pass,ok (implies operator) -01277,pass,ok (no-ode placeholder) -01278,pass,ok (no-ode placeholder) -01279,pass,ok (single-arg min/max) -01280,pass,ok (implies operator) -01281,pass,ok (implies operator) -01282,pass,ok (no-ode placeholder) -01283,pass,ok (no-ode placeholder) -01284,pass,ok (no-ode placeholder) -01285,pass,ok (no-ode placeholder) -01286,pass,ok (no-ode placeholder) +01271,pass,ok +01272,pass,ok +01273,pass,ok +01274,pass,ok +01275,pass,ok +01276,pass,ok +01277,pass,ok +01278,pass,ok +01279,pass,ok +01280,pass,ok +01281,pass,ok +01282,pass,ok +01283,pass,ok +01284,pass,ok +01285,pass,ok +01286,pass,ok 01287,unsupported,delay 01288,pass,ok -01289,pass,ok (no-ode placeholder) +01289,pass,ok 01290,pass,ok -01291,pass,ok (no-ode placeholder) +01291,pass,ok 01292,unsupported,algebraic rule 01293,pass,ok 01294,pass,ok @@ -1298,50 +1298,50 @@ test,status,notes 01297,pass,ok 01298,pass,ok 01299,unsupported,delay -01300,pass,ok (no-ode placeholder) -01301,pass,ok (no-ode placeholder) +01300,pass,ok +01301,pass,ok 01302,pass,ok -01303,pass,ok (no-ode placeholder) -01304,pass,ok (no-ode placeholder) +01303,pass,ok +01304,pass,ok 01305,unsupported,delay -01306,pass,ok (no-ode placeholder) +01306,pass,ok 01307,pass,ok 01308,pass,ok 01309,pass,ok -01310,pass,ok (no-ode placeholder) -01311,pass,ok (no-ode placeholder) -01312,pass,ok (nested function definitions; a function named after a C++ keyword) -01313,pass,ok (no-ode placeholder) -01314,pass,ok (no-ode placeholder) -01315,pass,ok (no-ode placeholder) -01316,pass,ok (no-ode placeholder) -01317,pass,ok (no-ode placeholder) -01318,unsupported,delay function -01319,unsupported,delay function +01310,pass,ok +01311,pass,ok +01312,pass,ok +01313,pass,ok +01314,pass,ok +01315,pass,ok +01316,pass,ok +01317,pass,ok +01318,unsupported,delay +01319,unsupported,delay 01320,unsupported,delay 01321,pass,ok 01322,pass,ok -01323,pass,ok (no-ode placeholder) +01323,pass,ok 01324,unsupported,delay 01325,unsupported,delay 01326,unsupported,delay 01327,unsupported,delay 01328,unsupported,delay 01329,unsupported,delay -01330,pass,ok (no-ode placeholder) -01331,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) -01332,pass,ok (no-ode placeholder) -01333,unsupported,event execution semantics (priority ordering / useValuesFromTriggerTime) -01334,pass,ok (no-ode placeholder) +01330,pass,ok +01331,unsupported,event execution semantics +01332,pass,ok +01333,unsupported,event execution semantics +01334,pass,ok 01335,unsupported,delay -01336,pass,ok (no-ode placeholder) -01337,unsupported,event execution semantics (persistent triggers) +01336,pass,ok +01337,unsupported,event execution semantics 01338,pass,ok 01339,pass,ok 01340,pass,ok -01341,pass,ok (no-ode placeholder) -01342,pass,ok (no-ode placeholder) -01343,pass,ok (no-ode placeholder) +01341,pass,ok +01342,pass,ok +01343,pass,ok 01344,pass,ok 01345,pass,ok 01346,pass,ok @@ -1351,29 +1351,29 @@ test,status,notes 01350,unsupported,algebraic rule 01351,pass,ok 01352,pass,ok -01353,pass,ok (no-ode placeholder) -01354,pass,ok (no-ode placeholder) -01355,pass,ok (no-ode placeholder) -01356,pass,ok (no-ode placeholder) -01357,pass,ok (no-ode placeholder) +01353,pass,ok +01354,pass,ok +01355,pass,ok +01356,pass,ok +01357,pass,ok 01358,unsupported,delay 01359,unsupported,algebraic rule 01360,pass,ok 01361,pass,ok -01362,pass,ok (no-ode placeholder) -01363,pass,ok (no-ode placeholder) -01364,pass,ok (no-ode placeholder) -01365,pass,ok (no-ode placeholder) -01366,pass,ok (no-ode placeholder) +01362,pass,ok +01363,pass,ok +01364,pass,ok +01365,pass,ok +01366,pass,ok 01367,unsupported,delay 01368,unsupported,algebraic rule 01369,pass,ok 01370,pass,ok -01371,pass,ok (no-ode placeholder) -01372,pass,ok (no-ode placeholder) -01373,pass,ok (no-ode placeholder) -01374,pass,ok (no-ode placeholder) -01375,pass,ok (no-ode placeholder) +01371,pass,ok +01372,pass,ok +01373,pass,ok +01374,pass,ok +01375,pass,ok 01376,unsupported,delay 01377,unsupported,algebraic rule 01378,pass,ok @@ -1388,11 +1388,11 @@ test,status,notes 01387,pass,ok 01388,pass,ok 01389,pass,ok -01390,pass,ok (no-ode placeholder) -01391,pass,ok (no-ode placeholder) -01392,pass,ok (no-ode placeholder) -01393,pass,ok (no-ode placeholder) -01394,pass,ok (no-ode placeholder) +01390,pass,ok +01391,pass,ok +01392,pass,ok +01393,pass,ok +01394,pass,ok 01395,pass,ok 01396,unsupported,fast reaction 01397,unsupported,fast reaction @@ -1414,7 +1414,7 @@ test,status,notes 01413,unsupported,delay 01414,unsupported,delay 01415,unsupported,delay -01416,unsupported,delay function +01416,unsupported,delay 01417,unsupported,delay 01418,unsupported,delay 01419,unsupported,delay @@ -1459,13 +1459,13 @@ test,status,notes 01458,pass,ok 01459,pass,ok 01460,pass,ok -01461,pass,ok (no-ode placeholder) +01461,pass,ok 01462,pass,ok 01463,pass,ok 01464,pass,ok 01465,pass,ok -01466,pass,ok (no-ode placeholder) -01467,pass,ok (no-ode placeholder) +01466,pass,ok +01467,pass,ok 01468,pass,ok 01469,pass,ok 01470,pass,ok @@ -1483,19 +1483,19 @@ test,status,notes 01482,unsupported,algebraic rule 01483,unsupported,algebraic rule 01484,unsupported,algebraic rule -01485,pass,ok (no-ode placeholder) -01486,pass,ok (xor / arcsec / uncommon mathml) -01487,pass,ok (no-ode placeholder) -01488,pass,ok (no-ode placeholder) -01489,pass,ok (no-ode placeholder) -01490,pass,ok (no-ODE model integrated over time; placeholder initial condition) -01491,pass,ok (no-ODE model integrated over time; placeholder initial condition) +01485,pass,ok +01486,pass,ok +01487,pass,ok +01488,pass,ok +01489,pass,ok +01490,pass,ok +01491,pass,ok 01492,pass,ok 01493,pass,ok -01494,pass,ok (nary relational via sm:: helpers) -01495,pass,ok (no-ode placeholder) -01496,pass,ok (no-ode placeholder) -01497,pass,ok (implies operator) +01494,pass,ok +01495,pass,ok +01496,pass,ok +01497,pass,ok 01498,pass,ok 01499,unsupported,algebraic rule 01500,unsupported,algebraic rule @@ -1520,7 +1520,7 @@ test,status,notes 01519,unsupported,delay 01520,unsupported,delay 01521,unsupported,delay -01522,unsupported,delay function +01522,unsupported,delay 01523,unsupported,delay 01524,unsupported,delay 01525,unsupported,delay @@ -1528,10 +1528,10 @@ test,status,notes 01527,pass,ok 01528,unsupported,delay 01529,unsupported,delay -01530,pass,ok (no-ode placeholder) -01531,partly_supported,clustered events +01530,pass,ok +01531,unsupported,event execution semantics 01532,unsupported,delay -01533,pass,ok (no-ode placeholder) +01533,pass,ok 01534,unsupported,delay 01535,unsupported,delay 01536,unsupported,delay @@ -1594,7 +1594,7 @@ test,status,notes 01593,unsupported,delay 01594,unsupported,delay 01595,unsupported,delay -01596,pass,ok (no-ode placeholder) +01596,pass,ok 01597,unsupported,delay 01598,unsupported,delay 01599,unsupported,random event execution @@ -1633,13 +1633,13 @@ test,status,notes 01632,pass,ok 01633,pass,ok 01634,pass,ok -01635,pass,ok (no-ode placeholder) +01635,pass,ok 01636,pass,ok 01637,pass,ok -01638,pass,ok (no-ode placeholder) +01638,pass,ok 01639,pass,ok 01640,pass,ok -01641,pass,ok (no-ode placeholder) +01641,pass,ok 01642,pass,ok 01643,pass,ok 01644,pass,ok @@ -1656,13 +1656,13 @@ test,status,notes 01655,pass,ok 01656,pass,ok 01657,pass,ok -01658,pass,ok (no-ode placeholder) +01658,pass,ok 01659,unsupported,delay 01660,unsupported,delay 01661,unsupported,delay -01662,pass,ok (no-ode placeholder) -01663,pass,ok (no-ode placeholder) -01664,pass,ok (no-ode placeholder) +01662,pass,ok +01663,pass,ok +01664,pass,ok 01665,pass,ok 01666,pass,ok 01667,pass,ok @@ -1691,14 +1691,14 @@ test,status,notes 01690,unsupported,delay 01691,unsupported,delay 01692,unsupported,delay -01693,pass,ok (no-ode placeholder) -01694,pass,ok (no-ode placeholder) -01695,pass,ok (no-ode placeholder) -01696,pass,ok (no-ode placeholder) -01697,pass,ok (no-ode placeholder) -01698,pass,ok (no-ode placeholder) -01699,pass,ok (no-ode placeholder) -01700,pass,ok (no-ode placeholder) +01693,pass,ok +01694,pass,ok +01695,pass,ok +01696,pass,ok +01697,pass,ok +01698,pass,ok +01699,pass,ok +01700,pass,ok 01701,unsupported,delay 01702,unsupported,delay 01703,unsupported,delay @@ -1777,10 +1777,10 @@ test,status,notes 01776,pass,ok 01777,pass,ok 01778,pass,ok -01779,pass,ok (no-ode placeholder) -01780,pass,ok (no-ode placeholder) -01781,pass,ok (nary relational via sm:: helpers) -01782,pass,ok (no-ode placeholder) +01779,pass,ok +01780,pass,ok +01781,pass,ok +01782,pass,ok 01783,pass,ok 01784,pass,ok 01785,unsupported,algebraic rule @@ -1794,29 +1794,29 @@ test,status,notes 01793,unsupported,algebraic rule 01794,unsupported,algebraic rule 01795,pass,ok -01796,pass,ok (no-ode placeholder) +01796,pass,ok 01797,pass,ok 01798,unsupported,delay 01799,pass,ok 01800,pass,ok 01801,pass,ok 01802,pass,ok -01803,pass,ok (no-ode placeholder) -01804,pass,ok (no-ode placeholder) -01805,pass,ok (no-ode placeholder) -01806,pass,ok (no-ode placeholder) +01803,pass,ok +01804,pass,ok +01805,pass,ok +01806,pass,ok 01807,pass,ok 01808,pass,ok 01809,pass,ok -01810,pass,ok (no-ode placeholder) -01811,pass,ok (parameters named INF/Inf/inf vs the infinity constant) -01812,pass,ok (id spelled like a C++ macro) -01813,pass,ok (parameters named NAN/NaN/nan vs the notanumber constant) -01814,pass,ok (no-ode placeholder) -01815,pass,ok (no-ode placeholder) -01816,pass,ok (no-ode placeholder) -01817,pass,ok (no-ode placeholder) -01818,pass,ok (no-ode placeholder) -01819,pass,ok (pi constant vs parameter distinguished) -01820,pass,ok (no-ode placeholder) -01821,pass,ok (no-ode placeholder) +01810,pass,ok +01811,pass,ok +01812,pass,ok +01813,pass,ok +01814,pass,ok +01815,pass,ok +01816,pass,ok +01817,pass,ok +01818,pass,ok +01819,pass,ok +01820,pass,ok +01821,pass,ok From 5b7c58168c1533c1169c88959c87203fcab7ad87 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 12:15:56 +0100 Subject: [PATCH 18/19] Rebalance SBML test-suite CI matrix into 5 equal batches Re-split the semantic-test matrix so each of the 5 batches runs exactly 269 passing cases (was ~230 and unbalanced after retriage). The ranges are contiguous case-id spans covering the whole suite: 1-276, 277-678, 679-1030, 1031-1364, 1365-1821. Co-Authored-By: Claude Opus 4.8 --- .github/workflows/sbml_test_suite.yml | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/.github/workflows/sbml_test_suite.yml b/.github/workflows/sbml_test_suite.yml index c3d0a7a4..f34c5b19 100644 --- a/.github/workflows/sbml_test_suite.yml +++ b/.github/workflows/sbml_test_suite.yml @@ -20,18 +20,18 @@ jobs: strategy: fail-fast: false - # matrix is split by number of passing test cases (~230 per batch) + # matrix is split by number of passing test cases (269 per batch) matrix: include: - first_case: "1" - last_case: "245" - - first_case: "246" - last_case: "578" - - first_case: "579" - last_case: "890" - - first_case: "891" - last_case: "1267" - - first_case: "1268" + last_case: "276" + - first_case: "277" + last_case: "678" + - first_case: "679" + last_case: "1030" + - first_case: "1031" + last_case: "1364" + - first_case: "1365" last_case: "1821" steps: From defe36a8e39f82e4dfd16fbdc7b166bd9a445143 Mon Sep 17 00:00:00 2001 From: Kwabena Amponsah Date: Mon, 6 Jul 2026 12:59:49 +0100 Subject: [PATCH 19/19] #42 Add test cases README --- chaste_sbml/SbmlRefModels/test/cases/.gitignore | 3 +++ chaste_sbml/SbmlRefModels/test/cases/README.md | 10 ++++++++++ chaste_sbml/_expressions.py | 2 +- chaste_sbml/_model_builder.py | 2 +- 4 files changed, 15 insertions(+), 2 deletions(-) create mode 100644 chaste_sbml/SbmlRefModels/test/cases/.gitignore create mode 100644 chaste_sbml/SbmlRefModels/test/cases/README.md diff --git a/chaste_sbml/SbmlRefModels/test/cases/.gitignore b/chaste_sbml/SbmlRefModels/test/cases/.gitignore new file mode 100644 index 00000000..855b75d4 --- /dev/null +++ b/chaste_sbml/SbmlRefModels/test/cases/.gitignore @@ -0,0 +1,3 @@ +semantic/ +stochastic/ +syntactic/ diff --git a/chaste_sbml/SbmlRefModels/test/cases/README.md b/chaste_sbml/SbmlRefModels/test/cases/README.md new file mode 100644 index 00000000..1731ab44 --- /dev/null +++ b/chaste_sbml/SbmlRefModels/test/cases/README.md @@ -0,0 +1,10 @@ +# SBML Test Suite Cases + +## Semantic +Case support tracked in test/data/sbml_test_suite_status.csv + +## Syntactic +Covered by libSBML + +## Stochastic +Not supported \ No newline at end of file diff --git a/chaste_sbml/_expressions.py b/chaste_sbml/_expressions.py index 11a0cb79..06faecad 100644 --- a/chaste_sbml/_expressions.py +++ b/chaste_sbml/_expressions.py @@ -391,7 +391,7 @@ def formula_to_string( rewrite_power(math) formula = formulaToL3String(math) - # Convert all integer literals to doubles to fix integer division. Thelookbehinds + # Convert all integer literals to doubles to fix integer division. The lookbehinds # keep the digits of a scientific-notation exponent (e.g. 1e-5, 1e+5) intact. formula = re.sub(r"(? None: self._functions = [] # Ids assigned by some event, so species classification can model an otherwise-constant - # event-modified species as a (variable) parameterp. + # event-modified species as a (variable) parameter. self._event_assigned_ids = { ea.getVariable() for event in self._sbml_events for ea in event.getListOfEventAssignments() }