import numpy as np
import matplotlib.pyplot as plt

L=2e-2 # plateau de la table de 2 cm d'épaisseur
N=100
dx=L/N
duree=3*3600

lamb=0.18 # bois
D=0.15e-7

h=1

T1=30 # température cutanée moyenne
T0=20


T=np.zeros(N)
T2=np.zeros(N)
T3=np.zeros((N,10))

# conditions aux limites initiale
T[1:N-1]=T0
T[0]=T1

# convergence si D*dt/(dx**2)<0.5
P=int(3*D*duree/dx**2)
dt=duree/P
print("nombres d'iterations : "+str(P))

j=0
for k in range(P):
    T2[0]=T[0]
    for i in range(1,N-1):
        T2[i]=T[i]+D*dt/(dx**2)*(T[i+1]+T[i-1]-2*T[i])

    T2[N-1]=20 #(h*T0+lamb/dx*T[N-2])/(h+lamb/dx) condition convection
    T=T2.copy()

    if k%int(P/10)==0 and k>0:
        T3[:,j]=T.copy()
        j=j+1


x=np.linspace(0,L*1000,N) # en mm

plt.figure(figsize=(10,8))
plt.xlabel("en mm")
plt.ylabel("en °C")
plt.title("Profil de température au cours du temps")

for i in range(9):
    plt.plot(x,T3[:,i],label="t = "+str(round(duree/10*i,4))+" s")
    plt.legend()

plt.show()