#!/bin/env python3

import matplotlib.pyplot as plt
import numpy as np
import scipy.optimize

"""
Date : 2019
Auteur : Joachim Galiana


Ce programme vise à tracer les branches d'un diagramme binaire solide/liquide possédant un eutectique.
Les deux branches sont calculées et tracées et leur intersection (eutectique) est calcuée et pointée de même.

Le fichier peut être adapté pour d'autres couples de composés, en modifiant la partie DÉFINITION DES GRANDEURS et CONSTANTES PHYSIQUES
Les légendes sont à adapter en conséquence.
"""

#DÉFINITION DES GRANDEURS et CONSTANTES PHYSIQUES
Tfus_d=352.4#en K
Tfus_p=372.4#en K
Delta_Fus_d=21e3#en J/mol
Delta_Fus_p=16.46e3#en J/mol
M_d=134.2#en g/mol
M_p=178.2#en g/mol

R=8.314#en J/K/mol

#DÉFINITION DES VARIABLES
w=np.linspace(0.01,0.99,99)
y=1-w

#DÉFINITION DES FONCTIONS À TRACER
#Température en fonction de la fraction massique en Durène
def Td(frac_w):
    return ( Tfus_d * Delta_Fus_d ) / ( Delta_Fus_d - R *Tfus_d * np.log( frac_w / (frac_w + (M_d / M_p) * (1-frac_w))))
#Température en fonction de la fraction massique en Phénanthrène
def Tp(frac_w):
    return ( Tfus_p * Delta_Fus_p ) / ( Delta_Fus_p - R *Tfus_p * np.log( frac_w / (frac_w + (M_p / M_d) * (1-frac_w))))
#Soustraction des deux fonctions pour en trouver l'intersection
def Soustraction(frac_w):
    return Td(frac_w) - Tp(1-frac_w)

#CALCUL DE L'INTERSECTION ENTRE LES COURBES DE TEMPÉRATURE
E=scipy.optimize.fsolve(Soustraction,0)
print(E)
print(Td(E))
print(Tp(1-E))

#TRACÉ DU DIAGRAMME BINAIRE
#Limites en température
plt.ylim(300,400)
#Pour différencier avant et après eutectique !
ante_w=np.array([a for a in w if a < E])
post_w=np.array([a for a in w if a > E])

#Tracé en rouge de la courbe associée au Durène
plt.plot(ante_w,Td(ante_w),"red",linestyle='--')
plt.plot(post_w,Td(post_w),"red")
#Tracé en bleu de la courbe associée au Phénanthrène
plt.plot(ante_w,Tp(1-ante_w),"blue")
plt.plot(post_w,Tp(1-post_w),"blue",linestyle='--')

#Tracé de l'eutectique
plt.plot(E,Td(E),"+")
plt.hlines(y=Td(E),xmin=0.01,xmax=0.99)

#Points expérimentaux
#w_exp=np.array([0.2,0.3,0.33,0.47,0.5,0.6])
#T_i=np.array([110.2,85.5,66.3,55.5,73,97])
#T_E=np.array([52.4,54,53.5,55.5,55,60])
#plt.plot(w_exp,T_i+273,"+")
#plt.plot(w_exp,T_E+273,"+")

#Légendes
plt.xlabel("Fraction massique en durène, w")
plt.ylabel("Température, T (K)")
plt.title("Diagramme binaire solide/liquide durène-phénanthrène\n Eutectique à w = 0,47, T = 324 K")
plt.grid()
plt.show()