# -*- coding: utf-8 -*-
"""
Created on Tue Nov 21 15:38:23 2023

@author: KAPOTA
"""

# 1- importation des bibliothèques
import numpy as np
import scipy.optimize as op
import matplotlib.pyplot as plt

# 2- Définition des paramètres du problème
DrH0 = (2*-241.8-4*-92.8)*10**3 #J.mol-1
DrS0 = (2*223.1+2*188.7-4*186.9-205.1) #J/K/mol
T0 = 298.00 #K
TF= 2000 #K
R= 8.314 #J/K/mol
P0=1 #bar
P=1 #bar
N=100

# 3- Fonction DrG0 dans l'approx d'Ellingham
def DrG0(T):
    return DrH0 - T*DrS0

# 3- Fonction K dans l'approx d'Ellingham
def K(T) :
    return np.exp(-DrG0(T)/(R*T))

# 3- Fonction à annuler si Q=num/denom=K alors f(T,tau)=num - K*denom
def f(T,P,tau):
    return tau**4*(4 - 2*tau)*P0**2-K(T)*4**2*(1-tau)**5*P

# 3-définir la fonction dichotomie 
def dichoto(min,max,ff,precision):
    while (max-min) > precision:
        m=(max+min)/2
        if ff(min)*ff(m) < 0 :
            max=m
        else :
            min=m
    return m

# 3-liste des abscisses :
T_abs=np.linspace(T0,TF,N)
Ptot=(1,10,20,50,100)

# 3- calcul de la liste de tau par la fonction dichoto
# utilisation de op.bisect(g,min,max) avec g une fonction à 1 valeur
# on va stocker les valeurs dans la liste tau_dichoto
for i in range(len(Ptot)):
    def fp(T,tau):
        return f(T,Ptot[i],tau)
    tau_dichoto=[]
    for j in range(len(T_abs)):
        def g(tau):
            return fp(T_abs[j],tau)
        tau_dichoto.append(dichoto(0,1,g,0.0001))
#4-tracé de chaque courbe
    titre = 'P = {} bar'.format(Ptot[i])
    plt.plot(T_abs,tau_dichoto,label=titre)
plt.xlabel('T (K)')
plt.ylabel('taux davancement (mol)' )
plt.legend()
plt.title('évolution de tau avec T pour l équilibre de Deacon méthode Dichoto')
plt.grid()
plt.show()