Skip to content

Commit

Permalink
fix #2642 (#2643)
Browse files Browse the repository at this point in the history 
* fix #2642

* fixup

* fixup

* add issue as regression test

* fix path
  • Loading branch information
@FFroehlich
FFroehlich authored Feb 5, 2025
1 parent 64de63a commit c9f698d
Show file tree
Hide file tree
Showing 3 changed files with 176 additions and 3 deletions.
26 changes: 23 additions & 3 deletions python/sdist/amici/sbml_import.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2292,9 +2292,29 @@ def _make_initial(

sym_math, rateof_to_dummy = _rateof_to_dummy(sym_math)

for species_id, species in self.symbols[SymbolId.SPECIES].items():
if "init" in species:
sym_math = smart_subs(sym_math, species_id, species["init"])
# we can't rely on anything else being properly initialized at this point, so we need to
# compute all initial values from scratch, recursively
for var in sym_math.free_symbols:
element_id = str(var)
# already recursive since _get_element_initial_assignment calls _make_initial
if (
ia := self._get_element_initial_assignment(element_id)
) is not None:
sym_math = sym_math.subs(var, ia)
elif (species := self.sbml.getSpecies(element_id)) is not None:
# recursive!
init = self._make_initial(get_species_initial(species))
sym_math = sym_math.subs(var, init)
elif var in self.symbols[SymbolId.SPECIES]:
sym_math = sym_math.subs(
var, self.symbols[SymbolId.SPECIES][var]["init"]
)
elif (
element := self.sbml.getElementBySId(element_id)
) and self.is_rate_rule_target(element):
# no need to recurse here, as value is numeric
init = sp.Float(element.getValue())
sym_math = sym_math.subs(var, init)

sym_math = smart_subs(sym_math, sbml_time_symbol, sp.Float(0))

Expand Down
130 changes: 130 additions & 0 deletions python/tests/sbml_models/regression_2642.xml
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,130 @@
<?xml version="1.0" encoding="UTF-8"?>
<sbml xmlns="http://www.sbml.org/sbml/level3/version2/core" level="3" version="2">
<model metaid="debug" id="debug" substanceUnits="substance" timeUnits="time_unit" volumeUnits="volume" areaUnits="area" lengthUnits="length">
<listOfUnitDefinitions>
<unitDefinition id="length">
<listOfUnits>
<unit kind="metre" exponent="1" scale="0" multiplier="1"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="area">
<listOfUnits>
<unit kind="metre" exponent="2" scale="0" multiplier="1"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="volume">
<listOfUnits>
<unit kind="litre" exponent="1" scale="0" multiplier="1"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="time_unit" name="time">
<listOfUnits>
<unit kind="second" exponent="1" scale="1" multiplier="6"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="substance">
<listOfUnits>
<unit kind="mole" exponent="1" scale="-9" multiplier="1"/>
</listOfUnits>
</unitDefinition>
</listOfUnitDefinitions>
<listOfCompartments>
<compartment id="cytosol" spatialDimensions="3" size="0.0005" constant="true"/>
</listOfCompartments>
<listOfSpecies>
<species id="MR_i" compartment="cytosol" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
<species id="MR_m" compartment="cytosol" hasOnlySubstanceUnits="false" boundaryCondition="false" constant="false"/>
</listOfSpecies>
<listOfParameters>
<parameter id="synthesis_MR" constant="true"/>
<parameter id="recycling" value="0.0000045" constant="true"/>
<parameter id="binding" constant="true"/>
<parameter id="MR_i_t0" value="1400000" constant="true"/>
<parameter id="MR_m_t0" value="600000" constant="true"/>
</listOfParameters>
<listOfInitialAssignments>
<initialAssignment symbol="MR_i">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> MR_i_t0 </ci>
</math>
</initialAssignment>
<initialAssignment symbol="MR_m">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> MR_m_t0 </ci>
</math>
</initialAssignment>
<initialAssignment symbol="synthesis_MR">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<divide/>
<cn> 0.00585177319092733 </cn>
<cn> 0.00606 </cn>
</apply>
</math>
</initialAssignment>
<initialAssignment symbol="binding">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<divide/>
<apply>
<minus/>
<apply>
<times/>
<ci> recycling </ci>
<ci> MR_i </ci>
</apply>
<ci> synthesis_MR </ci>
</apply>
<ci> MR_m </ci>
</apply>
</math>
</initialAssignment>
</listOfInitialAssignments>
<listOfReactions>
<reaction id="_J0" reversible="true">
<listOfProducts>
<speciesReference species="MR_i" stoichiometry="1" constant="true"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> synthesis_MR </ci>
</math>
</kineticLaw>
</reaction>
<reaction id="_J1" reversible="true">
<listOfReactants>
<speciesReference species="MR_i" stoichiometry="1" constant="true"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="MR_m" stoichiometry="1" constant="true"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> recycling </ci>
<ci> MR_i </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="_J2" reversible="true">
<listOfReactants>
<speciesReference species="MR_m" stoichiometry="1" constant="true"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="MR_i" stoichiometry="1" constant="true"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> binding </ci>
<ci> MR_m </ci>
</apply>
</math>
</kineticLaw>
</reaction>
</listOfReactions>
</model>
</sbml>
23 changes: 23 additions & 0 deletions python/tests/test_sbml_import.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -852,3 +852,26 @@ def test_import_same_model_name():

assert model_module_1c.get_model().getParameters()[0] == 1.0
assert model_module_1c.get_model().module is model_module_1c


@skip_on_valgrind
def test_regression_2642():
sbml_file = Path(__file__).parent / "sbml_models" / "regression_2642.xml"
sbml_importer = amici.SbmlImporter(sbml_file)
model_name = "regression_2642"
with TemporaryDirectory(prefix="regression_2642") as outdir:
sbml_importer.sbml2amici(
model_name=model_name,
output_dir=outdir,
)
module = amici.import_model_module(
module_name=model_name, module_path=outdir
)
model = module.getModel()
solver = model.getSolver()
model.setTimepoints(np.linspace(0, 1, 3))
r = amici.runAmiciSimulation(model, solver)
assert (
len(np.unique(r.w[:, model.getExpressionIds().index("binding")]))
== 1
)

0 comments on commit c9f698d

Please sign in to comment.