#Taux d'avancement en fonction de la température
#Exemple : Procédé Deacon

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import bisect



#Données
R=8.314
DrH0= -112400         # noter la valeur de l'enthalpie standard de réaction à 298K en J/mol
DrS0= -129.1          #entrer la valeur de l'entropie standard de réaction à 298K en J/K/mol
P = 1                #entrer la valeur de la pression totale en bar

#Valeur de la constante d'équilibre à 298K
DrG0=DrH0-298*DrS0
K298= np.exp(-DrG0/R/298)
print ( "la valeur de la constante d'équilibre à 298K est : ",K298)


# introduction de T comme variable
T=np.linspace( 298, 1000,100)

# Valeurs de la constante d'équilibre en fonction de T  et graphe K°(T)
A=np.log(K298)+ DrH0/R/298
def F1(x):
    return (np.exp(A-DrH0/R/x))
Y1=F1(T)
plt.figure(0)
plt.plot(T,Y1,'+b',label="constante d'équilibre")
plt.grid()
plt.title('Variations de K° avec la température')
plt.xlabel("T(K)")
plt.ylabel("K°")
plt.show()


# valeur du taux d'avancement  pour T = 298 K
def F(x):
    return (K298-(x**4*(5-x)/16/(1-x)**5/P))
tau0=bisect(F,0.001,0.999)
print ( "la valeur du taux d'avancement à 298K est : ",tau0)



# Valeurs du taux d'avancement en fonction de la température et graphe tau (T)
Y2=np.zeros(100)
def F2(x,t):
    return (F1(t)-(x**4*(5-x)/16/(1-x)**5/P))
for i in range (len(T)):
    Y2[i]=bisect(F2,0.001,0.999,args=T[i])

plt.figure(1)
plt.plot(T,Y2,'+g')
plt.grid()
plt.title("variations du taux d'avancement avec la température")
plt.xlabel("T(K)")
plt.ylabel("tau")
plt.show ()


#intervalle de températures pour lequel le taux d'avancement est supérieur à 0,5
Y=np.zeros(100)
for i in range (len(T)):
    Y[i]=0.5
plt.figure(1)
plt.plot(T,Y2,'+g')
plt.plot(T,Y,'+r')
plt.grid()
plt.xlabel("T(K)")
plt.ylabel("tau")
plt.show ()

#Influence de la pression
Y4=np.zeros(100)
Y5=np.zeros(100)

def F3(x,t,p):
    return (F1(t)-(x**4*(5-x)/16/(1-x)**5/p))
for i in range (len(T)):
    Y4[i]=bisect(F3,0.001,0.999,args=(T[i],1))
    Y5[i]=bisect(F3,0.001,0.999,args=(T[i],10))

plt.figure(2)
plt.plot(T,Y4,'+g',label='P=1 bar')
plt.plot(T,Y5,'xb',label='P=10 bars')
plt.grid()
plt.legend()
plt.title('Influence de la pression')
plt.xlabel("T(K)")
plt.ylabel("tau")
plt.show ()
