# équation de diddusion en dimension 1
import numpy as np
import matplotlib.pyplot as plt
from math import *

#données numériques
L=1# longueur de la barre
tau=4000 #durée totale d'évolution
D=1.2*10**(-4)#diffusivité thermique
nx=100 #nombre de tronçons découpés dans la barre#nt=40000 #nombre d'intervalles de temps
Deltax=L/nx # pas discrétisation spatiale
Deltat=0.1*(Deltax**2)/(2*D)#pas de discrétisation spatial compatible avec la condition CFL
nt=int(tau/Deltat)#nombre de points axe des temps
T1=25.0 #température initiale de la barre
T2=100

#construction axe des abscisses
x=np.linspace(0,L,nx)

#condition initiale
T=np.linspace(T1,T1,nx)# température initiale uniforme

#température en régime stationnaire
def Tsta(x):
    return (T1-T2)*x/L+T2


#schéma numérique de calcul des accroissements de température et des températures en chaque points bouclage sur le temps
T[0]=T2
for n in range(nt): #boucle évolution du temps pas à pas
    for m in range(1,nx-1):#boucle calcul des accroissements de température à un instant donné aux différents points
        T[m]=T[m]+((Deltat*D)/(Deltax**2))*(T[m+1]-2*T[m]+T[m-1])
    if (n%10000==0):
        plotlabel="t = %1.0f s"%(n*Deltat)#plotlabel="t = %1.2f s"%(n*Deltat) le 1.2f permet de choisir l'écriture de t nombre de chiffre avant et après la virgule avec le type float
        plt.plot(x,T,label=plotlabel)

#tracé évolution températureen régime stationnaire
plt.plot(x,Tsta(x),'+')
plt.grid()
plt.xlabel(r'x(m)',fontsize=26)
plt.ylabel(r'T(K)',fontsize=26)
plt.title("Evolution de la température dans la barre au cours du temps après de passage de T(0) de T1 à T2")
plt.axis([0,L,0,105])
plt.legend()
plt.show()