# -*- coding: utf-8 -*-
"""
Created on Fri Aug 14 11:28:18 2026

@author: KAPOTA
"""

# 1- importation des bibliothèques
import numpy as np
import matplotlib.pyplot as plt
import numpy.random as rd

# 2- conditions initiales et variables
C0=0.15 #mol/L
Cper = 4.0E-2 #mol/L
V0 = 10 #mL 
D_V0=0.1 #mL
D_Veq=0.05 #mL
u_Cper=0.01*Cper #mol/L
VVeq= 8 #mL pour le calcul d'incertitude à Veq=8mL
N=100000
D_V0=0.1 #mL
D_Veq=0.05 #mL
u_Cper=0.01*Cper #mol/L
VVeq= 8 #mL pour le calcul d'incertitude à Veq=8mL
N=100000

# 2- Données expérimentales définies comme des listes mais si vous préférez des array, à vous de les modifier
t = [50,100,150,200,250,300,350,400,450,500]
Veq = [13.6,12.1,11.1,10.1,8.9,8.2,7.5,6.6,6,5.4]

# 3- Calcul de C
C=[]
for i in range(len(Veq)) :
    C.append(5*Cper*Veq[i]/(2*V0))

#3- Calcul de lnC puisqu'on suppose un ordre 1
lnC=[]
for i in range(len(C)):
    lnC.append(np.log(C[i]))

# 3- Calcul de la regression linéaire attendue et exploitation
reglin=np.polyfit(t,lnC,1)
k=-reglin[0]
lnC_theo=np.polyval(reglin,t)
residu=[]
for i in range(len(lnC)):
    residu.append(lnC[i]-lnC_theo[i])

#3- Calcul de N valeur de C8, pour cela on va générer N valeurs de Cper, de Veq8 et de V0
Cper_MC=rd.normal(Cper, u_Cper,N)
V0_MC=V0+D_V0*rd.uniform(-1,1,N)
Veq_MC=[]
C_MC=[]
u_lnC=[]
for i in range(len(Veq)):
    Veq_MC.append(Veq[i]+D_Veq*rd.uniform(-1,1,N))
    C_MC=5*Cper_MC*Veq_MC[i]/(2*V0_MC)
    lnC_MC=np.log(C_MC)
    u_lnC.append(np.std((lnC_MC),ddof=1))


#3- calcul des barres d'erreurs
barre_lnC=[]
for i in range(len(t)):
    barre_lnC.append(2*u_lnC[i])

#4- affichage de k
print('u_lnC=',u_lnC)

# 4- affichage de k
print('k={:.2e} s-1'. format(k))

#tracé de C=f(t)
plt.plot(t,C,'+b')
plt.title('C=f(t)')
plt.xlabel('t (s)')
plt.ylabel('C (mol/L)')
plt.grid() #parce que c'est plus joli
plt.show()

# tracé de lnC et sa regression linéaire
plt.plot(t, lnC_theo,':r')
plt.errorbar(t,lnC,yerr=barre_lnC,fmt='o',capsize=5)
plt.title('vérification ordre 1')
plt.xlabel('t (s)')
plt.ylabel('lnC')
plt.grid() #parce que c'est plus joli
plt.show()

# tracé des résidus
plt.errorbar(t,residu,yerr=barre_lnC,fmt='o',capsize=5)
plt.title('résidu=f(t)')
plt.xlabel('t (s)')
plt.ylabel('résidu')
plt.axhline()
plt.grid() #parce que c'est plus joli
plt.show()