From finer to coarser resolution in advancing the two coupled reactions:¶

`2 S <-> U` and `S <-> X` (both mostly forward)¶

1st-order kinetics throughout.

Notes:

for an exploration of instabilities, see the experiment negative_concentrations_1
for an accurate longer run of the same reactions (with the same initial conditions), see the experiment "up_regulate_3"

LAST REVISED: Feb. 11, 2023 THIS IS AN ARCHIVED EXPERIMENT

IMPORTANT: DO NOT ATTEMPT TO RUN THIS NOTEBOOK!¶

This is a "frozen run" that depends on an old version of Life123, for demonstration purposes¶

(newer versions tend to recover from instability more gracefully.)
If you bypass the execution exit in the first cell, and run the other cells, you WILL NOT REPLICATE the results below!

In [1]:

# To stop the current and subsequent cells: USED TO PREVENT ACCIDENTAL RUNS OF THIS NOTEBOOK!

class StopExecution(Exception):
    def _render_traceback_(self):
        return []

raise StopExecution     # See: https://stackoverflow.com/a/56953105/5478830

In [ ]:

In [1]:

import set_path      # Importing this module will add the project's home directory to sys.path

Added 'D:\Docs\- MY CODE\BioSimulations\life123-Win7' to sys.path

In [2]:

from experiments.get_notebook_info import get_notebook_basename

from src.modules.reactions.reaction_data import ReactionData as chem
from src.modules.reactions.reaction_dynamics import ReactionDynamics

import numpy as np
import plotly.express as px
from src.modules.visualization.graphic_log import GraphicLog

In [3]:

# Initialize the HTML logging
log_file = get_notebook_basename() + ".log.htm"    # Use the notebook base filename for the log file

# Set up the use of some specified graphic (Vue) components
GraphicLog.config(filename=log_file,
                  components=["vue_cytoscape_1"],
                  extra_js="https://cdnjs.cloudflare.com/ajax/libs/cytoscape/3.21.2/cytoscape.umd.js")

-> Output will be LOGGED into the file 'large_time_steps_2.log.htm'

Initialize the system¶

In [4]:

# Initialize the system
chem_data = chem(names=["U", "X", "S"])

# Reaction 2 S <-> U , with 1st-order kinetics for all species (mostly forward)
chem_data.add_reaction(reactants=[(2, "S")], products="U",
                       forward_rate=8., reverse_rate=2.)

# Reaction S <-> X , with 1st-order kinetics for all species (mostly forward)
chem_data.add_reaction(reactants="S", products="X",
                       forward_rate=6., reverse_rate=3.)

chem_data.describe_reactions()

# Send the plot of the reaction network to the HTML log file
graph_data = chem_data.prepare_graph_network()
GraphicLog.export_plot(graph_data, "vue_cytoscape_1")

Number of reactions: 2 (at temp. 25 C)
0: 2 S <-> U  (kF = 8 / kR = 2 / Delta_G = -3,436.56 / K = 4) | 1st order in all reactants & products
1: S <-> X  (kF = 6 / kR = 3 / Delta_G = -1,718.28 / K = 2) | 1st order in all reactants & products
[GRAPHIC ELEMENT SENT TO LOG FILE `large_time_steps_2.log.htm`]

Run 1 : extremely small fixed time steps (no substeps)¶

In [5]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.001, stop_time=0.8)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

800 total step(s) taken
From time 0 to 0.8, in 800 FULL steps of 0.001
(for a grand total of the equivalent of 800 FULL steps)

In [6]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Run 2 : very small time steps, with dynamic substeps¶

In [7]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.01, stop_time=0.8,
                                  dynamic_substeps=4, rel_fast_threshold=50)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

single_compartment_react(): setting abs_fast_threshold to 50.0
80 total step(s) taken
From time 0 to 0.44, in 176 substeps of 0.0025 (each 1/4 of full step)
From time 0.44 to 0.8, in 36 FULL steps of 0.01
(for a grand total of the equivalent of 80 FULL steps)

In [8]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Run 3 : small-ish time steps, with dynamic substeps¶

In [9]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.08, stop_time=0.8,
                                  dynamic_substeps=4, rel_fast_threshold=250)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

single_compartment_react(): setting abs_fast_threshold to 31.25
11 total step(s) taken
From time 0 to 0.64, in 32 substeps of 0.02 (each 1/4 of full step)
From time 0.64 to 0.88, in 3 FULL steps of 0.08
(for a grand total of the equivalent of 11 FULL steps)

In [10]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Run 4 : same as previous run, but fewer dynamic substeps¶

In [11]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.08, stop_time=0.8,
                                  dynamic_substeps=2, rel_fast_threshold=150)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

single_compartment_react(): setting abs_fast_threshold to 18.75
11 total step(s) taken
From time 0 to 0.8, in 20 substeps of 0.04 (each 1/2 of full step)
From time 0.8 to 0.88, in 1 FULL step of 0.08
(for a grand total of the equivalent of 11 FULL steps)

In [12]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Run 5 : same as previous run, but slightly larger primary step¶

In [13]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.1, stop_time=0.8,
                                  dynamic_substeps=2, abs_fast_threshold=80.0)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

single_compartment_react(): setting rel_fast_threshold to 800.0
9 total step(s) taken
From time 0 to 0.3, in 6 substeps of 0.05 (each 1/2 of full step)
From time 0.3 to 0.9, in 6 FULL steps of 0.1
(for a grand total of the equivalent of 9 FULL steps)

In [14]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Note: the above run is identical to the last run in the experiment `negative_concentrations_1`¶

Run 6 : same as previous run, but no substeps¶

In [15]:

dynamics = ReactionDynamics(reaction_data=chem_data)
dynamics.set_conc(conc={"U": 50., "X": 100., "S": 0.})
#dynamics.describe_state()

dynamics.set_diagnostics()       # To save diagnostic information about the call to single_compartment_react()

dynamics.single_compartment_react(time_step=0.1, stop_time=0.8)

df = dynamics.get_history()
#df
dynamics.explain_time_advance()

******** CAUTION: negative concentration in chemical `S` (resulting from reaction 2 S <-> U)
         upon advancing from system time t=0.1 [Baseline value: 50 ; delta conc: -64]
         It will be AUTOMATICALLY CORRECTED with a reduction in the time step size
******** CAUTION: negative concentration resulting from the combined effect of multiple reactions, upon advancing reactions from system time t=0.45
         It will be AUTOMATICALLY CORRECTED with a reduction in the time step size
The computation took 2 extra step(s) - automatically added to prevent negative concentrations
10 total step(s) taken
From time 0 to 0.1, in 1 FULL step of 0.1
From time 0.1 to 0.15, in 1 substep of 0.05 (1/2 of full step)
From time 0.15 to 0.45, in 3 FULL steps of 0.1
From time 0.45 to 0.5, in 1 substep of 0.05 (1/2 of full step)
From time 0.5 to 0.9, in 4 FULL steps of 0.1
(for a grand total of the equivalent of 9 FULL steps)

Notice the automated detection - and correction - of negative concentrations arising from the excessively-large time steps¶

In [16]:

dynamics.plot_curves(colors=['green', 'orange', 'blue'])

Note: for an exploration of instabilities, see the experiment `negative_concentrations_1`¶