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

# In[1]:


import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import linalg
from scipy.integrate import odeint
sns.set()
get_ipython().run_line_magic('matplotlib', 'inline')


# In[2]:


t = 0;dt = 0.001; s4 = 0;
ci= 0 
s5=0.7; factive=0.2; s7=0.01; s6=1-factive-s7-s4;
tnow =[]
ixnow = []
F1=[]
F2=[]
F3=[]
F4=[]
for i in range(1,250000+1):
    t = t + dt
    ni = 0
    if t>30 and t<40:
        ni = 100
        ci =0.1
    if t>50 and t<60:
        ni = 100
        ci = 0.3
    if t>70 and t<80:
        ni = 100
        ci = 5
    if t>90 and t<100:
        ni = 100
        ci = 20
    if t>110 and t<120:
        ni = 100
        ci = 30
    if t>130 and t<140:
        ni = 100
        ci = 40
    if t>150 and t<160:
        ni = 100
        ci = 50
    if t>170 and t<180:
        ni = 100
        ci = 60
    if t>190 and t<200:
        ni = 100
        ci = 70
    #if t>310 and t<340:
        #ni = 100
        #ci = 70
    #if t>340 and t<370:
        #ni = 100
        #ci = 80
    #if t>370 and t<400:
        #ni = 100
        #ci =90
    #if t>400 and t<430:
        #ni = 100
        #ci = 95
    #if t>430 and t<460:
        #ni = 100
        #ci = 100
        
        
    f3n=ni**2.5/(ni**2.5+17**2.5)
    kcon1=0.1
    kcoff1=0.05
    kcon2=20 
    kcoff2=0.3
    kinact=0.05
    #kinact=1
    #kinact=0.001 
    #kinact=1
    #kinact=0.2
    
    s5=s5+ (s6*ci*kcon1-s5*kcoff1+factive*f3n*kinact-s5*0.3) *dt
    factive=factive+(s7*ci*kcon2-factive*kcoff2+s5*0.15-factive*f3n*kinact)*dt
    s7=s7+(factive*kcoff2+s6*0.1-s7*ci*kcon2-s7*f3n*kinact*25)*dt
    s6=1-s5-factive-s7
    
    F1.append(s5)
    F2.append(factive)
    F3.append(s7)
    F4.append(s6)
    
    incx=factive*f3n
    tnow.append(t)
    ixnow.append(incx)


# In[3]:


plt.figure(figsize = [10,5])
plt.plot(tnow,ixnow)
plt.xlabel("Time(t)")
plt.ylabel("Current")
plt.savefig("plot/current.png")
plt.savefig("plot/current.pdf")
plt.show()


# In[4]:


plt.figure(figsize = [10,5])
plt.plot(tnow,F1,label = "F1")
plt.plot(tnow,F2,label = "F2")
plt.plot(tnow,F3,label="F3")
plt.plot(tnow,F4,label="F4")
plt.xlabel("Time(t)",fontsize=15)
plt.ylabel("states",fontsize=15)
#plt.xticks([i*30 for i in range(10)])
plt.legend(fontsize = 8)
plt.savefig("plot/states.png")
plt.savefig("plot/states.pdf")
plt.show()


# In[ ]:





# In[ ]:





# In[ ]: