"""
Authors : Raphael Rullan et Timothée Chauvirée

On simule le titrage d'un acide faible par une base forte
Remarque : on s'intéresse ici à la réaction AH + B = A + H20 (avec B = HO-), on a acide fort base/forte
"""
import matplotlib.pyplot as plt
import numpy as np
import widgets

parameters = {'Cini_AH' : widgets.FloatSlider(value=0.1, description='concentration en acide (mol/L)', min=0.00001, max=2),
              'V_ini'  : widgets.IntSlider(value = 10, description='volume initial (mL)', min = 0.1, max =100),
              'Cini_B'  : widgets.FloatSlider(value = 0.1, description='concentration en titrant (mol/L)', min = 0.00001, max = 2),
              'pKa' : widgets.FloatSlider(value = 4.3, description='pKa du couple', min = 0, max = 11)}

Cini_A = 0




V=np.linspace(0,24,500)


  
# affichage des courbes

fig=plt.figure(figsize=(12,10))
plt.subplots_adjust(left=0.125,bottom=0.2,right=0.8,top=0.9,wspace=0.3,hspace=0.5)
ax1 = plt.subplot(211)
ax2 = plt.subplot(212)

plt.subplots_adjust(left=0.125,bottom=0.2,right=0.8,top=0.9,wspace=0.3,hspace=0.5)
fig.suptitle("Titrage d'un acide fort par une base forte")

def plot_data(Cini_AH, Cini_B, V_ini,pKa):
    VE = Cini_AH*V_ini/Cini_B
    C_AH, C_B,C_A, pH = np.zeros(np.size(V),),np.zeros(np.size(V),),np.zeros(np.size(V),),np.zeros(np.size(V),) 
    def calcul_quantites_avant_et_a_equivalence(Cini_AH,V_ini,Cini_A,V, pKa):
        C_AH = (Cini_AH*V_ini-Cini_B*V)/(V_ini+V)
        C_B = np.zeros(np.size(V))
        C_A = (Cini_B*V)/(V_ini+V)
        x = Cini_B*V/(Cini_AH*V_ini)
        pH = pKa + np.log10(x/(1-x))
        return pH,C_AH,C_B, C_A
    def calcul_quantites_v0(Cini_AH,V_ini,Cini_A,V, pKa):
        C_AH = Cini_AH
        C_B = 0
        C_A = 0
        pH = 0.5*(pKa - np.log10(Cini_AH))
        return np.array([pH]),np.array([C_AH]),np.array([C_B]), np.array([C_A])

    def calcul_quantites_equivalence(Cini_AH,V_ini,Cini_A,VE,V, pKa) :
        K = (10**-pKa)/(10**-14)
        C_AH = np.sqrt((Cini_AH*V_ini)/(V_ini+V)/K)
        C_B = np.sqrt((Cini_AH*V_ini)/(V_ini+V)/K)
        C_A = (Cini_AH*V_ini)/(V_ini+V)
        pH = pKa + np.log10(C_A/C_AH)
        return pH,C_AH,C_B, C_A

    def calcul_quantites_apres_equivalence(Cini_AH,V_ini,Cini_B,VE,V, pKa) :
        C_AH = np.zeros(np.size(V))
        C_B =(Cini_B*(V-VE))/((V+V_ini))
        C_A = (Cini_AH*V_ini)/(V_ini+ V)
        x = Cini_B*V/(Cini_AH*V_ini)    
        pH = 14 + np.log10((Cini_AH*V_ini)/(V_ini+V))+np.log10(x-1)
        return pH,C_AH,C_B, C_A
    
    def dpH(pH):
        dpH = np.gradient(pH)
        return dpH
 
    
    V0 = np.array(V[np.where(V==0)])
    pH1,C_AH1,C_B1,C_A1 = calcul_quantites_v0(Cini_AH,V_ini,Cini_A,V0, pKa)
    
   
    Vav = V[np.where(np.logical_and(V>0,V<=VE))]
    
    #appelle une fonction qui calcule les quantités avant et à l’équivalence
    pH2,C_AH2,C_B2,C_A2 = calcul_quantites_avant_et_a_equivalence(Cini_AH,V_ini,Cini_B,Vav, pKa)
    
    Veq = V[np.where(V==VE)]
    pH3,C_AH3,C_B3,C_A3 = calcul_quantites_equivalence(Cini_AH,V_ini,Cini_A,VE,Veq, pKa)
    
    Vap = V[np.where(V>VE)]
    #appelle une fonction qui calcule les quantités après l’équivalence
    pH4,C_AH4,C_B4,C_A4 = calcul_quantites_apres_equivalence(Cini_AH,V_ini,Cini_B,VE,Vap, pKa)
    
    
          
    pH = np.concatenate((pH1,pH2,pH4))
    C_AH = np.concatenate((C_AH1,C_AH2,C_AH4))
    C_B = np.concatenate((C_B1,C_B2,C_B4))
    C_A = np.concatenate((C_A1,C_A2,C_A4))
    
    lines['Courbe AH'].set_data(V,C_AH)
    lines['Courbe A'].set_data(V,C_A)
    lines['Courbe B'].set_data(V,C_B)
    a['Courbe pH'].set_data(V, pH)
    lines['Courbe derpH'].set_data(V,dpH(pH))
    texts.set_text('VE (mL) = {:.2e}'.format(VE))
    ax1.set_xlim(0,np.max(V))
    ax1.set_ylim(-0.00005,np.maximum(np.max(C_AH),np.max(C_B))) 
    ax2.set_xlim(0,np.max(V))
    ax2.set_ylim(0, np.max(pH)+1)
    
lines = {}
lines['Courbe AH'], = ax1.plot([], [],'b-', label="Concentration Acide AH" )
lines['Courbe A'], = ax1.plot([], [],'g-', label="Concentration Base conjuguee A")
lines['Courbe B'], = ax1.plot([], [],'r-', label="Concentration Base B")
lines['Courbe derpH'], = ax2.plot([], [],'r', label="dérivée du pH")
a = {}
a['Courbe pH'], = ax2.plot([], [],'b', label="Evolution pH")


texts= ax1.text(0.851, 0.9, '$VE=$', bbox=dict(facecolor='red', alpha=0.5), fontsize=10, transform=ax1.transAxes)





ax1.set_xlabel('Volume ajouté (mL)')
ax1.set_ylabel('concentration en espèce (mol/L)')
ax1.legend(bbox_to_anchor=(0.05,0.1, 0.8,1), loc = 2)



ax2.set_xlabel('Volume ajouté (mL)')
ax2.set_ylabel("pH")
ax2.legend()

param_widgets = widgets.make_param_widgets(parameters, plot_data, slider_box=[0.30, 0.05, 0.4, 0.075])
choose_widget = widgets.make_choose_plot(lines, box=[0.9,0.2,0.08, 0.2])
reset_button = widgets.make_reset_button(param_widgets,box=[0.01, 0.01, 0.08, 0.05])

if __name__=='__main__':
    plt.show()