Source code for sbmlsim.examples.example_units

"""
Example for handling units in simulations and results.
"""

import numpy as np
from matplotlib import pyplot as plt
from sbmlutils.console import console

from sbmlsim.resources import DEMO_SBML
from sbmlsim.simulation import Dimension, ScanSim, Timecourse, TimecourseSim
from sbmlsim.simulator.rr_simulator_serial import SimulatorSerialRR
from sbmlsim.units import Quantity, UnitsInformation
from sbmlsim.xresult import XResult



[docs]def run_demo_example(): """Run various timecourses.""" simulator = SimulatorSerialRR.from_sbml(DEMO_SBML) # units information uinfo = UnitsInformation.from_sbml(DEMO_SBML) Q_ = uinfo.Q_ # 1. simple timecourse simulation print("*** setting concentrations and amounts ***") # FIXME: units of changes are not used for conversion tc_scan = ScanSim( simulation=TimecourseSim( [ Timecourse( start=0, end=10, steps=100, changes={ "[e__A]": Q_(10, "mM"), "[e__B]": Q_(1, "mmole/litre"), "[e__C]": Q_(1, "mole/m**3"), "c__A": Q_(1e-5, "mole"), "c__B": Q_(10, "µmole"), "Vmax_bA": Q_(300.0, "mole/min"), }, ) ] ), dimensions=[ Dimension( "dim1", index=np.arange(20), changes={"[e__A]": Q_(np.linspace(5, 15, num=20), "mM")}, ) ], ) # print(tc_sim) xres: XResult = simulator.run_scan(tc_scan) xres.set_units(udict=uinfo.udict) console.log(xres) # create figure fig, (ax1, ax2, ax3, ax4) = plt.subplots(nrows=1, ncols=4, figsize=(20, 5)) fig.subplots_adjust(wspace=0.3, hspace=0.3) axes = (ax1, ax2, ax3, ax4) axes_units = [ {"xunit": "s", "yunit": "mM"}, {"xunit": "ms", "yunit": "mole/litre"}, {"xunit": "min", "yunit": "nM"}, {"xunit": "hr", "yunit": "pmole/cm**3"}, ] ax: plt.Axes for ax, ax_units in dict(zip(axes, axes_units)).items(): xunit = ax_units["xunit"] yunit = ax_units["yunit"] for key in ["[e__A]", "[e__B]", "[e__C]", "[c__A]", "[c__B]", "[c__C]"]: # => How to better handle units !!! # FIXME: correct handling of units; and subsequent conversion # Mapping of model to units needed # FIXME: correct reduction of additional dimensions! ax.plot( Q_(xres["time"].values, xres.units["time"]).to(xunit).m, Q_(xres[key].values, xres.units[key]).to(yunit).m, label=f"{key} [{yunit}]", ) ax.legend() ax.set_xlabel(f"time [{xunit}]") ax.set_ylabel(f"concentration [{yunit}]") plt.show()
if __name__ == "__main__": run_demo_example()