Skip to content

Instantly share code, notes, and snippets.

@winkmal
Last active February 23, 2023 16:56
Show Gist options
  • Select an option

  • Save winkmal/fb7e2b5d1deea97118636ba879c4bffe to your computer and use it in GitHub Desktop.

Select an option

Save winkmal/fb7e2b5d1deea97118636ba879c4bffe to your computer and use it in GitHub Desktop.
Porting simple mass-based 2frak VS model to SBML/SymBiology/AMICI
<?xml version="1.0" encoding="UTF-8"?>
<sbml xmlns="http://www.sbml.org/sbml/level2/version4" level="2" version="4">
<annotation>
<SimBiology xmlns="http://www.mathworks.com">
<Version Major="6" Minor="4" Point="0"/>
</SimBiology>
</annotation>
<model id="mwf36e8484_4f16_442a_9bf6_3c1e5dadb30d" name="mass2frak">
<listOfUnitDefinitions>
<unitDefinition id="MWBUILTINUNIT_liter" name="liter">
<listOfUnits>
<unit kind="metre" exponent="3"/>
<unit kind="dimensionless" multiplier="0.001"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_liter__kilogram" name="liter/kilogram">
<listOfUnits>
<unit kind="gram" exponent="-1"/>
<unit kind="metre" exponent="3"/>
<unit kind="dimensionless" multiplier="1e-06"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_kilogram__liter" name="kilogram/liter">
<listOfUnits>
<unit kind="gram" exponent="1"/>
<unit kind="metre" exponent="-3"/>
<unit kind="dimensionless" multiplier="1000000"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_kilogram__meter___3" name="kilogram/meter^3">
<listOfUnits>
<unit kind="gram" exponent="1"/>
<unit kind="metre" exponent="-3"/>
<unit kind="dimensionless" multiplier="1000"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_1__day" name="1/day">
<listOfUnits>
<unit kind="second" exponent="-1"/>
<unit kind="dimensionless" multiplier="1.15740740740741e-05"/>
</listOfUnits>
</unitDefinition>
</listOfUnitDefinitions>
<listOfCompartments>
<compartment id="mwe60e1361_4f09_4623_8996_333e3da58738" name="V_liq" spatialDimensions="3" size="0.001" units="MWBUILTINUNIT_liter" constant="false"/>
</listOfCompartments>
<listOfSpecies>
<species id="mwf7741186_33e4_4e9d_89f8_fa51d4ce19c3" name="m_DVS1" compartment="mwe60e1361_4f09_4623_8996_333e3da58738" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mw5d523312_56e6_4abd_9e23_754443bb2050" name="m_DVS2" compartment="mwe60e1361_4f09_4623_8996_333e3da58738" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mwd0c00b97_9eec_4627_a9e1_6e4cc5e8dfb2" name="m_NDVS" compartment="mwe60e1361_4f09_4623_8996_333e3da58738" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="true"/>
<species id="mwdc6f0d5f_17ac_4b3f_9950_cccfb65aba80" name="m_ch4" compartment="mwe60e1361_4f09_4623_8996_333e3da58738" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
</listOfSpecies>
<listOfParameters>
<parameter id="mw15cca539_ce49_4318_b193_c7a0b7729d6b" name="m_VS0" value="0.0024567" units="kilogram" constant="true"/>
<parameter id="mwce95f56f_84d6_46ea_9d9d_8dfc205e8084" name="f1" value="0.448" units="dimensionless" constant="true"/>
<parameter id="mw65cc611a_3b93_4da8_a1fe_e2f2b57bfbcc" name="f2" value="0.446" units="dimensionless" constant="true"/>
<parameter id="mw3c7aa987_eb84_4ac7_9692_94e615299756" name="Y_ch4" value="420" units="MWDERIVEDUNIT_liter__kilogram" constant="true"/>
<parameter id="mw8beed3dd_8f8e_4d56_bab5_58f3e0d5093f" name="V_liq0" value="0.4" units="MWBUILTINUNIT_liter" constant="true"/>
<parameter id="mw7f38f112_8721_46fa_a07e_4d9f1f6ac434" name="rho_liq" value="1" units="MWDERIVEDUNIT_kilogram__liter" constant="true"/>
<parameter id="mwfca314df_35c4_42c1_9132_0f00ee1865b4" name="rho_ch4" value="0.7168" units="MWDERIVEDUNIT_kilogram__meter___3" constant="true"/>
<parameter id="mw71887a17_54c2_44d4_94df_704e22d539b6" name="k1" value="0.5" units="MWDERIVEDUNIT_1__day" constant="true"/>
<parameter id="mwc81785c8_cc03_4508_8bbb_e88e1a454c9b" name="k2" value="0.05" units="MWDERIVEDUNIT_1__day" constant="true"/>
</listOfParameters>
<listOfInitialAssignments>
<initialAssignment symbol="mwe60e1361_4f09_4623_8996_333e3da58738">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> mw8beed3dd_8f8e_4d56_bab5_58f3e0d5093f </ci>
</math>
</initialAssignment>
<initialAssignment symbol="mwf7741186_33e4_4e9d_89f8_fa51d4ce19c3">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mwce95f56f_84d6_46ea_9d9d_8dfc205e8084 </ci>
<ci> mw15cca539_ce49_4318_b193_c7a0b7729d6b </ci>
</apply>
</math>
</initialAssignment>
<initialAssignment symbol="mw5d523312_56e6_4abd_9e23_754443bb2050">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mw65cc611a_3b93_4da8_a1fe_e2f2b57bfbcc </ci>
<ci> mw15cca539_ce49_4318_b193_c7a0b7729d6b </ci>
</apply>
</math>
</initialAssignment>
<initialAssignment symbol="mwd0c00b97_9eec_4627_a9e1_6e4cc5e8dfb2">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<apply>
<minus/>
<apply>
<minus/>
<cn type="integer"> 1 </cn>
<ci> mwce95f56f_84d6_46ea_9d9d_8dfc205e8084 </ci>
</apply>
<ci> mw65cc611a_3b93_4da8_a1fe_e2f2b57bfbcc </ci>
</apply>
<ci> mw15cca539_ce49_4318_b193_c7a0b7729d6b </ci>
</apply>
</math>
</initialAssignment>
</listOfInitialAssignments>
<listOfRules>
<rateRule variable="mwe60e1361_4f09_4623_8996_333e3da58738">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<divide/>
<apply>
<minus/>
<apply>
<plus/>
<apply>
<times/>
<ci> mw71887a17_54c2_44d4_94df_704e22d539b6 </ci>
<ci> mwf7741186_33e4_4e9d_89f8_fa51d4ce19c3 </ci>
</apply>
<apply>
<times/>
<ci> mwc81785c8_cc03_4508_8bbb_e88e1a454c9b </ci>
<ci> mw5d523312_56e6_4abd_9e23_754443bb2050 </ci>
</apply>
</apply>
</apply>
<ci> mw7f38f112_8721_46fa_a07e_4d9f1f6ac434 </ci>
</apply>
</math>
</rateRule>
</listOfRules>
<listOfReactions>
<reaction id="mw5c83fecf_db49_461b_abfa_850b4cb1d4ca" name="Reaction_1" reversible="false">
<listOfReactants>
<speciesReference species="mwf7741186_33e4_4e9d_89f8_fa51d4ce19c3" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mwdc6f0d5f_17ac_4b3f_9950_cccfb65aba80" stoichiometry="0.3011"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mw71887a17_54c2_44d4_94df_704e22d539b6 </ci>
<ci> mwf7741186_33e4_4e9d_89f8_fa51d4ce19c3 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="mw0e002e5f_68fb_41fd_aa4a_1a1651a8d323" name="Reaction_2" reversible="false">
<listOfReactants>
<speciesReference species="mw5d523312_56e6_4abd_9e23_754443bb2050" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mwdc6f0d5f_17ac_4b3f_9950_cccfb65aba80" stoichiometry="0.3011"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mwc81785c8_cc03_4508_8bbb_e88e1a454c9b </ci>
<ci> mw5d523312_56e6_4abd_9e23_754443bb2050 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
</listOfReactions>
</model>
</sbml>
<?xml version="1.0" encoding="UTF-8"?>
<sbml xmlns="http://www.sbml.org/sbml/level2/version4" level="2" version="4">
<annotation>
<SimBiology xmlns="http://www.mathworks.com">
<Version Major="6" Minor="4" Point="0"/>
</SimBiology>
</annotation>
<model id="mwe1e9b936_85f9_447c_9f44_6fd3235fedb0" name="mass2frak">
<listOfUnitDefinitions>
<unitDefinition id="MWBUILTINUNIT_liter" name="liter">
<listOfUnits>
<unit kind="metre" exponent="3"/>
<unit kind="dimensionless" multiplier="0.001"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_liter__kilogram" name="liter/kilogram">
<listOfUnits>
<unit kind="gram" exponent="-1"/>
<unit kind="metre" exponent="3"/>
<unit kind="dimensionless" multiplier="1e-06"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_kilogram__liter" name="kilogram/liter">
<listOfUnits>
<unit kind="gram" exponent="1"/>
<unit kind="metre" exponent="-3"/>
<unit kind="dimensionless" multiplier="1000000"/>
</listOfUnits>
</unitDefinition>
<unitDefinition id="MWDERIVEDUNIT_1__day" name="1/day">
<listOfUnits>
<unit kind="second" exponent="-1"/>
<unit kind="dimensionless" multiplier="1.15740740740741e-05"/>
</listOfUnits>
</unitDefinition>
</listOfUnitDefinitions>
<listOfCompartments>
<compartment id="mwd174a751_031f_4733_988a_5d2a58f842f0" name="V_liq" spatialDimensions="3" size="0.001" units="MWBUILTINUNIT_liter" constant="false"/>
</listOfCompartments>
<listOfSpecies>
<species id="mw6fe849d9_96e9_4fd0_9b94_46b102b37245" name="m_VS_feed" compartment="mwd174a751_031f_4733_988a_5d2a58f842f0" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f" name="m_DVS1" compartment="mwd174a751_031f_4733_988a_5d2a58f842f0" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mw2c4e3997_d475_4c13_b386_17d8aefeb632" name="m_DVS2" compartment="mwd174a751_031f_4733_988a_5d2a58f842f0" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mwa741a466_f218_4518_8e4e_df911c809a4a" name="m_NDVS" compartment="mwd174a751_031f_4733_988a_5d2a58f842f0" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
<species id="mw3e6298c9_2139_491d_9e8c_ba40279dac34" name="m_ch4" compartment="mwd174a751_031f_4733_988a_5d2a58f842f0" initialAmount="0" substanceUnits="kilogram" hasOnlySubstanceUnits="true" boundaryCondition="false" constant="false"/>
</listOfSpecies>
<listOfParameters>
<parameter id="mwe02ca726_d30b_4633_ad0d_a73f20615d17" name="m_VS0" value="0.0024567" units="kilogram" constant="true"/>
<parameter id="mwd21aa6fb_8550_478e_aaf9_29f30c5bd1fb" name="f1" value="0.48446" units="dimensionless" constant="true"/>
<parameter id="mw79899ef9_f66e_4da0_be8d_c70d176d52ee" name="f2" value="0.18212" units="dimensionless" constant="true"/>
<parameter id="mwc28a65d6_7e4a_4c1f_847e_186cfed9f576" name="Y_ch4" value="420" units="MWDERIVEDUNIT_liter__kilogram" constant="true"/>
<parameter id="mwc42cadec_8a86_466a_a818_4823058fdc5a" name="V_liq0" value="0.4" units="MWBUILTINUNIT_liter" constant="true"/>
<parameter id="mw4a7b0b4f_6842_464d_a837_19c29b798e0c" name="rho_liq" value="1" units="MWDERIVEDUNIT_kilogram__liter" constant="true"/>
<parameter id="mwc12c4f63_d57c_487e_84f7_1e7cfcf8ccfb" name="m_DVS10" value="1" units="kilogram" constant="true"/>
<parameter id="mwfd760d66_5e27_4919_ab62_909b24c704cf" name="m_DVS20" value="1" units="kilogram" constant="true"/>
<parameter id="mw42fffe24_7a65_4244_977c_fb43ef44eeae" name="m_NDVS0" value="1" units="kilogram" constant="true"/>
<parameter id="mwd928159d_a622_4546_b45c_494fc614dbde" name="kpseudo" value="10000" units="MWDERIVEDUNIT_1__day" constant="true"/>
<parameter id="mw6e40a5af_9692_49e3_a39c_ea84fb64637d" name="k1" value="0.5" units="MWDERIVEDUNIT_1__day" constant="true"/>
<parameter id="mwd9c7ab9e_2052_49f3_8f06_4bf2b66b8559" name="k2" value="0.05" units="MWDERIVEDUNIT_1__day" constant="true"/>
</listOfParameters>
<listOfInitialAssignments>
<initialAssignment symbol="mwd174a751_031f_4733_988a_5d2a58f842f0">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> mwc42cadec_8a86_466a_a818_4823058fdc5a </ci>
</math>
</initialAssignment>
<initialAssignment symbol="mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> mwc12c4f63_d57c_487e_84f7_1e7cfcf8ccfb </ci>
</math>
</initialAssignment>
<initialAssignment symbol="mw2c4e3997_d475_4c13_b386_17d8aefeb632">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> mwfd760d66_5e27_4919_ab62_909b24c704cf </ci>
</math>
</initialAssignment>
<initialAssignment symbol="mwa741a466_f218_4518_8e4e_df911c809a4a">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<ci> mw42fffe24_7a65_4244_977c_fb43ef44eeae </ci>
</math>
</initialAssignment>
</listOfInitialAssignments>
<listOfRules>
<rateRule variable="mwd174a751_031f_4733_988a_5d2a58f842f0">
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<divide/>
<apply>
<minus/>
<apply>
<plus/>
<apply>
<minus/>
<apply>
<times/>
<ci> mwd21aa6fb_8550_478e_aaf9_29f30c5bd1fb </ci>
<ci> mwd928159d_a622_4546_b45c_494fc614dbde </ci>
<ci> mw6fe849d9_96e9_4fd0_9b94_46b102b37245 </ci>
</apply>
<apply>
<times/>
<ci> mw6e40a5af_9692_49e3_a39c_ea84fb64637d </ci>
<ci> mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f </ci>
</apply>
</apply>
<apply>
<times/>
<ci> mw79899ef9_f66e_4da0_be8d_c70d176d52ee </ci>
<ci> mwd928159d_a622_4546_b45c_494fc614dbde </ci>
<ci> mw6fe849d9_96e9_4fd0_9b94_46b102b37245 </ci>
</apply>
</apply>
<apply>
<times/>
<ci> mwd9c7ab9e_2052_49f3_8f06_4bf2b66b8559 </ci>
<ci> mw2c4e3997_d475_4c13_b386_17d8aefeb632 </ci>
</apply>
</apply>
<ci> mw4a7b0b4f_6842_464d_a837_19c29b798e0c </ci>
</apply>
</math>
</rateRule>
</listOfRules>
<listOfReactions>
<reaction id="mw1e614f4b_6e55_41bf_8059_9751ffcdcfcf" name="Reaction_1" reversible="false">
<listOfReactants>
<speciesReference species="mw6fe849d9_96e9_4fd0_9b94_46b102b37245" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f" stoichiometry="1"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mwd21aa6fb_8550_478e_aaf9_29f30c5bd1fb </ci>
<ci> mwd928159d_a622_4546_b45c_494fc614dbde </ci>
<ci> mw6fe849d9_96e9_4fd0_9b94_46b102b37245 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="mw1530d3be_bad0_4a99_b007_c6bdbf7ea720" name="Reaction_2" reversible="false">
<listOfReactants>
<speciesReference species="mw6fe849d9_96e9_4fd0_9b94_46b102b37245" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mw2c4e3997_d475_4c13_b386_17d8aefeb632" stoichiometry="1"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mw79899ef9_f66e_4da0_be8d_c70d176d52ee </ci>
<ci> mwd928159d_a622_4546_b45c_494fc614dbde </ci>
<ci> mw6fe849d9_96e9_4fd0_9b94_46b102b37245 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="mwa52d878a_32ba_4f69_a447_0bfc7f5d637d" name="Reaction_3" reversible="false">
<listOfReactants>
<speciesReference species="mw6fe849d9_96e9_4fd0_9b94_46b102b37245" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mwa741a466_f218_4518_8e4e_df911c809a4a" stoichiometry="1"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<apply>
<minus/>
<apply>
<minus/>
<cn type="integer"> 1 </cn>
<ci> mwd21aa6fb_8550_478e_aaf9_29f30c5bd1fb </ci>
</apply>
<ci> mw79899ef9_f66e_4da0_be8d_c70d176d52ee </ci>
</apply>
<ci> mwd928159d_a622_4546_b45c_494fc614dbde </ci>
<ci> mw6fe849d9_96e9_4fd0_9b94_46b102b37245 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="mwf02945b4_04b7_4cc2_be4c_61fe304d8320" name="Reaction_4" reversible="false">
<listOfReactants>
<speciesReference species="mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mw3e6298c9_2139_491d_9e8c_ba40279dac34" stoichiometry="0.3011"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mw6e40a5af_9692_49e3_a39c_ea84fb64637d </ci>
<ci> mwdc320f45_e5ae_4a55_98ef_2f50f6b34b4f </ci>
</apply>
</math>
</kineticLaw>
</reaction>
<reaction id="mwe1256e69_68b4_403b_9ae6_f53563385fdc" name="Reaction_5" reversible="false">
<listOfReactants>
<speciesReference species="mw2c4e3997_d475_4c13_b386_17d8aefeb632" stoichiometry="1"/>
</listOfReactants>
<listOfProducts>
<speciesReference species="mw3e6298c9_2139_491d_9e8c_ba40279dac34" stoichiometry="0.3011"/>
</listOfProducts>
<kineticLaw>
<math xmlns="http://www.w3.org/1998/Math/MathML">
<apply>
<times/>
<ci> mwd9c7ab9e_2052_49f3_8f06_4bf2b66b8559 </ci>
<ci> mw2c4e3997_d475_4c13_b386_17d8aefeb632 </ci>
</apply>
</math>
</kineticLaw>
</reaction>
</listOfReactions>
</model>
</sbml>
sbml_importer = amici.SbmlImporter('../rawDataFiles/mass2frakVs_fromMatlabSimBiology.xml')
2023-02-14 21:16:18.366 - amici.sbml_import - ERROR - libSBML SBML component consistency (Error): The value of a <species>'s 'units' attribute is restricted.
Reference: L2V4 Section 4.8.5
The value of a <species>'s 'substanceUnits' attribute can only be one of the following: 'substance', 'mole', 'item', 'gram', 'kilogram', 'dimensionless', or the identifier of a <unitDefinition> derived from 'mole' (with an 'exponent' of '1'), 'item' (with an 'exponent' of '1'), 'gram' (with an 'exponent' of '1'), 'kilogram' (with an 'exponent' of '1'), or 'dimensionless'. The current value ('MWBUILTINUNIT_liter') is not allowed.
Traceback (most recent call last):
File "/tmp/ipykernel_2979/3590856256.py", line 1, in <module>
sbml_importer = amici.SbmlImporter('../rawDataFiles/mass2frakVs_fromMatlabSimBiology.xml')
File "/home/winkmal/.local/lib/python3.10/site-packages/amici/sbml_import.py", line 151, in __init__
self._process_document()
File "/home/winkmal/.local/lib/python3.10/site-packages/amici/sbml_import.py", line 184, in _process_document
_check_lib_sbml_errors(self.sbml_doc, self.show_sbml_warnings)
File "/home/winkmal/.local/lib/python3.10/site-packages/amici/sbml_import.py", line 2058, in _check_lib_sbml_errors
raise SBMLException(
SBMLException: SBML Document failed to load (see error messages above)
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Trying to use AMICI for a simple VS2CH4 model, to better understand syntax
@author: winkmal
"""
import amici
# Batch version of the model
sbml_importerB = amici.SbmlImporter('DVS_modelBatch_SimBio.xml')
sbml_importerB.sbml2amici('mass2frak', 'mass2frakVs')
# Conti version of the model
sbml_importerC = amici.SbmlImporter('DVS_modelConti_SimBio.xml')
# sbml_importerC.sbml2amici('mass2frak', 'mass2frakVs')
% Test (D)VS 2-frak model with SimBiology Toolbox, at first as pure batch model
% Code-Template taken from: https://mathworks.com/help/simbio/ug/rule-object.html#buerqv2
clc; clearvars
mass2frak = sbiomodel('mass2frak');
addcompartment(mass2frak,'liq1',1,Units='liter');
% VS = "volatile solids"
addspecies(mass2frak,'m_VS', Units='kilogram'); % 'Tag','m_{VS}'
addspecies(mass2frak,'m_DVS1', Units='kilogram');
addspecies(mass2frak,'m_DVS2', Units='kilogram');
% Set to const. right from the bat
addspecies(mass2frak,'m_NDVS', 'InitialAmountUnits','kilogram', 'Constant',true);
addspecies(mass2frak,'V_ch4', Units='liter', InitialAmount=0);
addparameter(mass2frak,'m_VS0', 2.9799e-3, Units='kilogram');
addparameter(mass2frak,'f1', 0.44145, Units='dimensionless');
addparameter(mass2frak,'f2', 0.27544, Units='dimensionless');
addparameter(mass2frak,'Y_ch4', 420, Units='liter/kilogram');
% Reactions, Rates, model parameters
addrule(mass2frak,'m_DVS1 = -k1 * m_DVS1', 'RuleType','rate'); % Reaction rate 1
addrule(mass2frak,'m_DVS2 = -k2 * m_DVS2', 'RuleType','rate'); % Reaction rate 2
addrule(mass2frak,'V_ch4 = Y_ch4*(k1*m_DVS1 + k2*m_DVS2)', 'RuleType','rate'); % Cumulative Gas production = volume
addparameter(mass2frak,'k1', 0.08608, Units='1/day');
addparameter(mass2frak,'k2', 0.89843, Units='1/day');
% Initial values
addrule(mass2frak,'m_DVS1 = f1*m_VS0', 'RuleType','initialAssignment');
addrule(mass2frak,'m_DVS2 = f2*m_VS0', 'RuleType','initialAssignment');
addrule(mass2frak,'m_NDVS = (1 - f1 - f2)*m_VS0', 'RuleType','initialAssignment');
% Algebraic Rules: Changed into repeatedAssignment
addrule(mass2frak,'m_VS = (m_DVS1 + m_DVS2 + m_NDVS)', 'RuleType','repeatedAssignment');
% Since q_ch4 and V_ch4_spec "only" algebraic output variables, use addobservable instead of addspecies
addobservable(mass2frak,'q_ch4', 'Y_ch4.*(k1.*m_DVS1 + k2.*m_DVS2)', Units='liter/day');
addobservable(mass2frak,'V_ch4_spec', 'V_ch4./m_VS0', Units='liter/kilogram');
% sbmlexport(mass2frak, 'mass2frakVs_fromMatlabSymBiology.xml')
%% Simulation (with createSimFunction instead of sbiosimulate)
cs1 = getconfigset(mass2frak); % configSet #1,
cs1.TimeUnits = 'day'; % units
cs1.StopTime = 15; % Set Final SimTime to 15 days
cs1.CompileOptions.UnitConversion = true;
inputs = {};
outputs = ["m_VS","m_DVS1","m_DVS2","q_ch4","V_ch4_spec"];
dose = [];
variant = [];
simfun = createSimFunction(mass2frak,inputs,outputs,dose,variant,AutoAccelerate=false);
stopTime = 15;
tic; sd = simfun([],stopTime); toc
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment