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

# ### Exploring the data structures of MEMBRANES, and reactions in them 
# #### - with NO DIFFUSION
# 
# LAST REVISED: June 23, 2024 (using v. 1.0 beta34.1)

# In[1]:


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


# In[2]:


from life123 import ChemData as chem

from life123 import BioSim1D


# In[3]:


chem_data = chem(names=["A", "B", "C"])     # NOTE: Diffusion not done
bio = BioSim1D(n_bins=5, chem_data=chem_data)

bio.set_membranes(membrane_pos=[1])   # A single membrane, passing thru bin 1

bio.set_all_uniform_concentrations(conc_list=[4., 8., 12.])

bio.describe_state()


# In[4]:


bio.set_bin_conc(bin_address=1, species_name="A", conc=10.)
bio.set_bin_conc(bin_address=1, species_name="A", conc=55., across_membrane=True)
bio.set_bin_conc(bin_address=1, species_name="B", conc=20.)
bio.set_bin_conc(bin_address=1, species_name="C", conc=30., both_sides=True)

bio.describe_state()


# In[5]:


# Make the last bin match all the concentrations of the "post-membrane" section of bin 1
bio.set_bin_conc(bin_address=4, species_name="A", conc=55.)
bio.set_bin_conc(bin_address=4, species_name="C", conc=30.)

bio.describe_state()


# In[6]:


# Reaction A + B <-> C , with 1st-order kinetics in both directions, mostly forward
bio.chem_data.add_reaction(reactants=["A", "B"], products=["C"], forward_rate=8., reverse_rate=2.)

bio.chem_data.describe_reactions()


# In[7]:


bio.react(time_step=0.002, n_steps=1)
bio.describe_state()


# ### Note how (in the absence of diffusion, which we're neglecting) the concentrations on the "post-membrane side" of bin 1 continue to match those of bin 5

# ## Now continue to reaction equilibrium

# In[8]:


bio.react(time_step=0.002, n_steps=100)
bio.describe_state()


# ### The system has now reached equilibrium
# ### in individual bins, which remain separate because we're NOT doing diffusion in this experiment

# Verify the equilibrium in each of the active bins

# In[9]:


bio.reaction_dynamics.is_in_equilibrium(rxn_index=0, conc={"A": 0.7932534523016195, "B": 4.793253452301622, "C": 15.206746547698396})
# A was largely the limiting reagent


# In[10]:


bio.reaction_dynamics.is_in_equilibrium(rxn_index=0, conc={"A": 0.897094737056602, "B": 10.897094737056594, "C": 39.10290526294337})
# A was largely the limiting reagent


# In[11]:


bio.reaction_dynamics.is_in_equilibrium(rxn_index=0, conc={"A": 47.20020986266437, "B": 0.20020986266437912, "C": 37.79979013733563})
# This time, with ample [A], the limiting reagent was B


# In[ ]: