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

# # Reaction-Diffusion in 1-D: `A + B <-> C` in 1-D, taken to equilibrium
# #### Mostly forward reaction, with 1st-order kinetics for each species
# 
# Initial concentrations of `A` and `B` are spatially separated to the opposite ends of the system;
# as a result, no `C` is being generated.
# 
# But, as soon as `A` and `B`, from their respective distant originating points at the edges, 
# diffuse into the middle - and into each other - the reaction starts,
# consuming both `A` and `B`,
# until an equilibrium is reached in both diffusion and reactions.

# ### TAGS :  "reactions 1D", "diffusion 1D", "quick-start"

# In[1]:


LAST_REVISED = "May 19, 2025"
LIFE123_VERSION = "1.0.0rc3"        # Library version this experiment is based on


# In[2]:


#import set_path                    # Using MyBinder?  Uncomment this before running the next cell!


# In[3]:


#import sys
#sys.path.append("C:/some_path/my_env_or_install")   # CHANGE to the folder containing your venv or libraries installation!
# NOTE: If any of the imports below can't find a module, uncomment the lines above, or try:  import set_path   

from life123 import check_version, BioSim1D, ChemData, UniformCompartment


# In[4]:


check_version(LIFE123_VERSION)


# In[ ]:





# In[5]:


# Initialize the system
chem_data = ChemData(names=["A", "B", "C"], diffusion_rates=[50., 50., 1.],   # `A` and `B` diffuse fast; `C` diffuses slowly
                     plot_colors=["red", "blue", "purple"])                   # Color choice is a reminder that red + blue = purple


# In[6]:


bio = BioSim1D(n_bins=7, chem_data=chem_data)


# In[7]:


reactions = bio.get_reactions()

# Reaction A + B <-> C , with 1st-order kinetics for each species; note that it's mostly in the forward direction
reactions.add_reaction(reactants=["A", "B"], products="C", forward_rate=20., reverse_rate=2.)
reactions.describe_reactions()


# In[ ]:





# # TIME 0 : Inject initial concentrations of `A` and `B` at opposite ends of the system

# In[8]:


bio.set_bin_conc(bin_address=0, chem_label="A", conc=20.)
bio.set_bin_conc(bin_address=6, chem_label="B", conc=20.)

bio.describe_state()


# In[9]:


bio.show_system_snapshot()     # A more streamlined alternate view


# In[10]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[11]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ## Enable History

# In[12]:


# Let's take a peek at the current concentrations of all chemicals in the bins with the initial concentration injections,
# as well as at the bin in the very center
bio.selected_concentrations(bins=[0, 6, 3])


# In[13]:


# Let's enable history for those same 3 bins
bio.enable_history(bins=[0, 6, 3], frequency=2, take_snapshot=True)     # Taking a snapshot to include the current initial state in the history


# In[ ]:





# In[ ]:





# ### <a name="sec_2_first_step"></a>First step : advance to time t=0.002

# In[14]:


delta_t = 0.002   # This will be our time "quantum" for this experiment


# In[15]:


# First step
bio.react_diffuse(time_step=delta_t, n_steps=1)
bio.describe_state()


# In[16]:


bio.show_system_snapshot()


# In[17]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[18]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_2"></a>Several more steps : advance to time t=0.016

# In[19]:


# Continue with several delta_t steps
for _ in range(7):
    bio.react_diffuse(time_step=delta_t, n_steps=1)
    bio.describe_state(concise=True)


# In[20]:


bio.show_system_snapshot()


# In[21]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# `A` is continuing to diffuse from the left.  
# `B` is continuing to diffuse from the right.  
# They're finally beginning to overlap in the middle bins!

# In[22]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_3"></a>Several groups of longer runs : advance to time t=0.096

# In[23]:


# Now, do several group of longer runs
for _ in range(4):
    print("\n\n+ 10 steps later:")
    bio.react_diffuse(time_step=delta_t, n_steps=10)
    bio.describe_state(concise=True)


# In[24]:


bio.show_system_snapshot()


# In[25]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# `A` is continuing to diffuse from the left.  
# `B` is continuing to diffuse from the right.  
# The overlap of `A` and `B`, especially in the central bins, is by now extensive, and the reaction is proceeding in earnest.  
# Notice the continue symmetry about the center of the system

# In[26]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_4"></a>Advance to time t=0.336

# In[27]:


# Continue the simulation
for _ in range(4):
    print("\n\n+++ 30 steps later:")
    bio.react_diffuse(time_step=delta_t, n_steps=30)
    bio.describe_state(concise=True)


# In[28]:


bio.show_system_snapshot()


# In[29]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[30]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_5"></a>Advance to time t=0.736

# In[31]:


# Continue the simulation
for _ in range(4):
    print("\n+++++ 50 steps later:")
    bio.react_diffuse(time_step=delta_t, n_steps=50)
    bio.describe_state(concise=True)


# In[32]:


bio.show_system_snapshot()


# In[33]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[34]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_6"></a>Advance to time t=1.936

# In[35]:


# Continue the simulation
for _ in range(4):
    print("\n+++++++++++++++ 150 steps later:")
    bio.react_diffuse(time_step=delta_t, n_steps=150)
    bio.describe_state(concise=True)


# In[36]:


bio.show_system_snapshot()


# In[37]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[38]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# ### <a name="sec_7"></a>Advance to time t=5.936

# In[39]:


# Continue the simulation
for _ in range(2):
    print("\n++++++++++ ... ++++++++++ 1,000 steps later:")
    bio.react_diffuse(time_step=delta_t, n_steps=1000)
    bio.describe_state(concise=True)


# In[40]:


bio.show_system_snapshot()


# In[41]:


bio.visualize_system(title_prefix="Reaction-Diffusion A + B <-> C")   # Line curve view


# In[42]:


bio.system_heatmaps(title_prefix="Reaction-Diffusion A + B <-> C")


# In[ ]:





# In[ ]:





# ## Equilibrium in both diffusion and reaction  
# All bins now have essentially uniform concentration, for each of the chemicals

# In each bin at equilibrium, [A] = 0.49, [B] = 0.49, [C] = 2.37   
# Let's verify that it's consistent with our (mostly-forward) equation `A + B <-> C`:

# In[43]:


bio.reaction_in_equilibrium(bin_address=0, rxn_index=0, explain=True)  # Choice of bin is immaterial now, because they have all equilibrated


# In[ ]:





# **Mass conservation**: The initial total quantities of `A` and `B` were 20 each, and zero `C`.   Here's how they have changed:

# In[44]:


20. - bio.chem_quantity(chem_label="A")  # Tot. decrease in `A`


# In[45]:


20. - bio.chem_quantity(chem_label="B")  # Tot. decrease in `B`


# In[46]:


bio.chem_quantity(chem_label="C")  # Tot. increase in `C`


# All consistent with the stoichiometry of our reaction `A + B <-> C`

# In[ ]:





# In[ ]:





# # Plots of changes of concentration with time

# In[47]:


bio.plot_history_single_bin(bin_address=0)


# In[48]:


bio.plot_history_single_bin(bin_address=6)


# In[49]:


bio.plot_history_single_bin(bin_address=3)


# `A` and `B` overlap on the plot, due to the symmetry of the system.  
# Initially, in the middle bin, neither `A` nor `B` are present; over time they diffuse there... but then they react and get consumed (producing `C`), to an equilibrium value.  
# Meanwhile, `C` gradually diffuses to uniformity.

# In[ ]: