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

# ## The attenuation over time of a single-frequency in the initial concentration
# 
# ### The initial system state is a sine wave of frequency 2 (i.e. 2 cycles across the system's length), of amplitude 10, with a baseline (bias) of 30
# ### Afterward, the process is restarted and repeated with a frequency 5 times larger  
#   
#   
# 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 BioSim1D

import plotly.express as px
import plotly.graph_objects as go

from life123 import ChemData as chem


# In[3]:


# Initialize the system.  We use a RELATIVELY LARGE NUMBER OF BINS, 
# to captures the finer changes in the frequency components of the concentration function
chem_data = chem(names=["A"], diffusion_rates=[0.5])
bio = BioSim1D(n_bins=500, chem_data=chem_data)


# ## Initial Preparation -
# ### Start with a sinusoidal concentration, with exactly 2 cycles over the length of the system
# #### (of amplitude 10 and bias value of 30)

# In[4]:


bio.inject_sine_conc(species_name="A", amplitude=10, bias=30, frequency=2)


# In[5]:


bio.show_system_snapshot()


# In[6]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= "Initial System State",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[7]:


# Show as heatmap
fig = px.imshow(bio.system_snapshot().T, 
                title= "Initial System State (as a heatmap)", 
                labels=dict(x="Bin number", y="Chem. species", color="Concentration"),
                text_auto=False, color_continuous_scale="gray_r") 

fig.data[0].xgap=0
fig.data[0].ygap=0

fig.show()


# In[8]:


# Take a look at the frequency domain of the concentration values
frequency_data = bio.frequency_analysis(species_name="A")
frequency_data


# In[9]:


ratio = frequency_data.loc[1, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
ratio


# In[10]:


bio.add_snapshot(data_snapshot={"ratio": ratio})
bio.get_history()


# # Start the diffusion - to time t=10

# In[11]:


bio.diffuse(total_duration=10, n_steps=100)


# In[12]:


# Take a look at the frequency domain of the concentration values
frequency_data = bio.frequency_analysis(species_name="A")
frequency_data


# In[13]:


ratio = frequency_data.loc[2, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
ratio


# In[14]:


bio.add_snapshot(data_snapshot={"ratio": ratio})
bio.get_history()


# ## Do 49 more rounds of diffusion - to time t = 500

# In[15]:


for i in range(49):
    bio.diffuse(total_duration=10, n_steps=100)
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[2, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[16]:


bio.get_history()


# In[17]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= f"System State at time t={bio.system_time}",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# #### Note how the curve is beginning to flatten; the peaks are no longer between 20 and 40

# ## Do 150 more rounds of diffusion - to time t = 2000

# In[18]:


for i in range(150):
    bio.diffuse(total_duration=10, n_steps=75)    # Notice the gradual decreas of the number of intermediate steps, given the smaller gradient
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[2, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[19]:


bio.get_history()


# In[20]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
          title= f"Diffusion. System snapshot at time t={bio.system_time}",
          color_discrete_sequence = ['red'],
          labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# #### Note how the curve is flatter still - and even beginning to lose shape at the boundary 

# ## Do 800 more rounds of diffusion - to time t = 10000

# In[21]:


for i in range(800):
    bio.diffuse(total_duration=10, n_steps=20)   # Note how we're gradually increasing the time steps, because the gradiants are now smaller
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[2, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[22]:


bio.get_history()


# In[23]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
          title= f"Diffusion. System snapshot at time t={bio.system_time}",
          color_discrete_sequence = ['red'],
          labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# #### Getting decisively flatter, as it approaches equilibrium

# ## Now, let's look how the amplitude of the sine signal, relative to the constant bias, changed over time

# In[24]:


fig = px.line(data_frame=bio.get_history(), y=["ratio"], 
          title= "Component of frequency=2, relative to amplitude of constant bias",
          color_discrete_sequence = ['purple'],
          labels={"value":"ratio", "variable":"Freq. 2", "index":"Time (groups of 10 units)"})
fig.show()


# In[25]:


fig = px.line(data_frame=bio.get_history(), y=["ratio"], 
              title= "Component of frequency=2, relative to amplitude of constant bias (log scale in y-axis)",
              color_discrete_sequence = ['purple'],
              labels={"value":"ratio", "variable":"Freq. 2", "index":"Time (groups of 10 units)"},
              log_y=True)
fig.show()


# In[26]:


saved_freq_2_df = bio.get_history().copy(deep=True)   # Save the whole dataframe, for later use


# # Start a NEW system
# #### Everything same as before EXCEPT the frequency is now 5 times bigger (10 cycles over the length of the system)

# In[27]:


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


# In[28]:


bio.inject_sine_conc(species_name="A", amplitude=10, bias=30, frequency=10)   # x5 higher frequency than before


# In[29]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= "Initial System State",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[30]:


# Show as heatmap
fig = px.imshow(bio.system_snapshot().T, 
                title= "Initial System State (as a heatmap)", 
                labels=dict(x="Bin number", y="Chem. species", color="Concentration"),
                text_auto=False, color_continuous_scale="gray_r") 

fig.data[0].xgap=0
fig.data[0].ygap=0

fig.show()


# In[31]:


# Take a look at the frequency domain of the concentration values
frequency_data = bio.frequency_analysis(species_name="A")
frequency_data


# In[32]:


ratio = frequency_data.loc[1, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
ratio


# In[33]:


bio.add_snapshot(data_snapshot={"ratio": ratio})
bio.get_history()


# In[34]:


bio.diffuse(total_duration=10, n_steps=100)


# In[35]:


# Take a look at the frequency domain of the concentration values
frequency_data = bio.frequency_analysis(species_name="A")
frequency_data.loc[10]


# In[36]:


ratio = frequency_data.loc[10, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
ratio


# In[37]:


bio.add_snapshot(data_snapshot={"ratio": ratio})
bio.get_history()


# In[38]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= "Initial System State",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[39]:


for i in range(4):
    bio.diffuse(total_duration=10, n_steps=100)
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[10, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[40]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= f"System State at time t={bio.system_time}",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[41]:


bio.get_history()


# In[42]:


for i in range(10):
    bio.diffuse(total_duration=10, n_steps=100)
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[10, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[43]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= f"System State at time t={bio.system_time}",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[44]:


bio.get_history()


# In[45]:


for i in range(35):
    bio.diffuse(total_duration=10, n_steps=100)
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[10, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[46]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= f"System State at time t={bio.system_time}",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# In[47]:


bio.get_history()


# In[48]:


for i in range(50):
    bio.diffuse(total_duration=10, n_steps=100)
    frequency_data = bio.frequency_analysis(species_name="A")
    ratio = frequency_data.loc[10, "Relative Amplitude"] / frequency_data.loc[0, "Relative Amplitude"]
    bio.add_snapshot(data_snapshot={"ratio": ratio})


# In[49]:


fig = px.line(data_frame=bio.system_snapshot(), y=["A"], 
              title= f"System State at time t={bio.system_time}",
              color_discrete_sequence = ['red'],
              labels={"value":"concentration", "variable":"Chemical", "index":"Bin number"})
fig.show()


# **Most of the concentration graph is now flat as a board!**

# In[50]:


bio.get_history()


# ## Like done earlier for the smaller frequency, let's look how the amplitude of the sine signal, relative to the constant bias, changed over time

# In[51]:


fig = px.line(data_frame=bio.get_history(), y=["ratio"], 
          title= "Component of frequency=10, relative to amplitude of constant bias",
          color_discrete_sequence = ['green'],
          labels={"value":"ratio", "variable":"Freq. 10", "index":"Time (groups of 10 units)"})
fig.show()


# In[52]:


fig = px.line(data_frame=bio.get_history(), y=["ratio"], 
              title= "Component of frequency=10, relative to amplitude of constant bias (log scale in y-axis)",
              color_discrete_sequence = ['green'],
              labels={"value":"ratio", "variable":"Freq. 10", "index":"Time (groups of 10 units)"},
              log_y=True)
fig.show()


# ### Resurrect the earlier plot for frequency=2 (using a saved dataset)
# ### and then superpose them

# In[53]:


fig_2 = px.line(data_frame=saved_freq_2_df, y=["ratio"], 
              title= "Component of frequency=2, relative to amplitude of constant bias (log scale in y-axis)",
              color_discrete_sequence = ['purple'],
              labels={"value":"ratio", "variable":"Freq. 2", "index":"Time (groups of 10 units)"},
              log_y=True)
fig_2.show()


# In[54]:


# Combine the last 2 plots
all_fig = go.Figure(data=fig.data + fig_2.data)    # Note that the + is concatenating lists
all_fig.update_layout(title="Superposed plots for frequency 2 (purple) and 10 (green)")
all_fig.show()


# ## Note how more more substantial the attenuation of the higher frequency (green) is!

# In[ ]: