#!/usr/bin/env python
# coding: utf-8

# ## Comparing the dynamics of the reaction `A <-> B` , with and without an enzyme
# 
# LAST REVISED: June 23, 2024 (using v. 1.0 beta36)

# In[1]:


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


# In[2]:


from life123 import ChemData
from life123 import UniformCompartment


# # 1. WITHOUT ENZYME
# ### `A` <-> `B`

# In[3]:


# Initialize the system
chem_data = ChemData(names=["A", "B"])

# Reaction A <-> B , with 1st-order kinetics, favorable thermodynamics in the forward direction, 
# and a forward rate that is much slower than it would be with the enzyme - as seen in part 2, below
chem_data.add_reaction(reactants="A", products="B",
                       forward_rate=1., delta_G=-3989.73)

chem_data.describe_reactions()


# ### Set the initial concentrations of all the chemicals

# In[4]:


dynamics = UniformCompartment(chem_data=chem_data, preset="fast")
dynamics.set_conc(conc={"A": 20.},
                  snapshot=True)
dynamics.describe_state()


# ### Take the initial system to equilibrium

# In[5]:


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

dynamics.single_compartment_react(duration=3.0,
                                  initial_step=0.1, variable_steps=True)


# In[6]:


#dynamics.explain_time_advance()


# In[7]:


dynamics.plot_history(colors=['darkorange', 'green'], show_intervals=True, title_prefix="WITHOUT enzyme")


# #### Note how the time steps get automatically adjusted, as needed by the amount of change - including a complete step abort/redo at time=0

# In[8]:


dynamics.curve_intersect("A", "B", t_start=0, t_end=1.0)


# In[9]:


# Verify that the reaction has reached equilibrium
dynamics.is_in_equilibrium()


# In[ ]:





# # 2. WITH ENZYME `E`
# ### `A` + `E` <-> `B` + `E`

# ### Note: for the sake of the demo, we'll completely ignore the concomitant (much slower) direct reaction A <-> B
# This in an approximation that we'll drop in later experiments

# In[10]:


# Initialize the system
chem_data = ChemData(names=["A", "B", "E"])

# Reaction A + E <-> B + E , with 1st-order kinetics, and a forward rate that is faster than it was without the enzyme
# Thermodynamically, there's no change from the reaction without the enzyme
chem_data.add_reaction(reactants=["A", "E"], products=["B", "E"],
                       forward_rate=10., delta_G=-3989.73)

chem_data.describe_reactions()  # Notice how the enzyme `E` is noted in the printout below


# ### Notice how, while the ratio kF/kR is the same as it was without the enzyme (since it's dictated by the energy difference), the individual values of kF and kR now are each 10 times bigger than before

# ### Set the initial concentrations of all the chemicals (to what they originally were)

# In[11]:


dynamics = UniformCompartment(chem_data=chem_data, preset="mid")
dynamics.set_conc(conc={"A": 20., "B": 0., "E": 30.},
                  snapshot=True)      # Plenty of enzyme `E`
dynamics.describe_state()


# ### Take the initial system to equilibrium

# In[12]:


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

dynamics.single_compartment_react(duration=0.1,
                                  initial_step=0.1, variable_steps=True)


# #### Note how the (proposed) initial step - kept the same as the previous run - is now _extravagantly large_, given the much-faster reaction dynamics.  However, the variable-step engine intercepts and automatically corrects the problem!

# In[13]:


#dynamics.explain_time_advance()


# In[14]:


dynamics.plot_history(colors=['darkorange', 'green', 'violet'], show_intervals=True, title_prefix="WITH enzyme")


# In[15]:


dynamics.curve_intersect("A", "B", t_start=0, t_end=0.02)


# In[16]:


# Verify that the reaction has reached equilibrium
dynamics.is_in_equilibrium()


# ## Thanks to the (abundant) enzyme, the reaction reaches equilibrium roughly around t=0.02, far sooner than the roughly t=3.5 without enzyme
# The concentrations of `A` and `B` now become equal (cross-over) at t=0.00246 , rather than t=0.740

# In[17]:


dynamics.get_history()


# In[ ]: