"Simulation d''une EDP (équation de la chaleur) sur le tore R²/Z²"

"1-Saisie des paramètres"

T = 1; //Durée de la simulation (en s)
J = 100; //Nombre de subdivisions en espace
D = 0.2; //Coefficient de diffusion thermique (en m²/s)
N = ceil(4*D*T*J^2); //Nombre de subdivisions en temps

"2-Pas de subdivision et intervalles de résolution"

ht = T/N; //Pas de subdivision en temps
hs = 1/J; //Pas de subdivision en espace
TT = 0:ht:T; //Intervalle de résolution en temps
XX = 0:hs:1; //Intervalle de résolution en abscisse
YY = 0:hs:1; //Intervalle de résolution enordonnée


"3-Condition initiale"

U0 = exp(-100*(XX-0.5).^2)'*exp(-100*(YY-0.5).^2);


"4-Création des matrices"

A = sparse([(1:J+1)',(1:J+1)'],-2*ones(J+1,1)) + sparse([(1:J+1)',[(2:J+1)';1]],ones(J+1,1)) + sparse([[(2:J+1)';1],(1:J+1)'],ones(J+1,1)); //Matrice de discrétisation du laplacien sur le tore
A1 = speye(J+1,J+1) - (D*ht/hs^2)*A;
A2 = speye(J+1,J+1) + (D*ht/hs^2)*A;
U = zeros(J+1,(N+1)*(J+1)); //Matrice de la solution approchée de l'équation
V = zeros(J+1,J+1);
V = zeros(J+1,J+1);


"5-Calcul de la solution"
U(:,1:J+1) = U0;
for n=1:N
    V = U(:,(n-1)*(J+1)+1:n*(J+1));
    W = A1\(V*A2 + 0.005*V.^2.*(1-V)); //Schéma semi-explicite-semi-implicite
    U(:,n*(J+1)+1:(n+1)*(J+1)) = W;
end


"6-Tracé de la solution"

f = scf(1);
cmap = jetcolormap(512);
f.color_map = cmap;
Min = min(U);
Max = max(U);
UU = U(:,1:(J+1));
plot3d1(XX,YY,UU,270,0,"x@y@T",[-6,1,1]);
title('Température sur le tore (en °C) à t = 0 s');
xlabel('x (en m)');
ylabel('y (en m)');
for n = 2:50:N+1
    UU = U(:,n*(J+1)+1:(n+1)*(J+1));
    cmap2 = cmap(1:ceil(max(UU)*512/Max),:);
    drawlater();
    f=gcf();
    f.color_map=cmap2;
    title('Température sur la plaque (en °C) à t = '+string(round(1000*n*ht)/1000)+' s');
    plot3d1(XX,YY,UU,270,0,"x@y@T",[-6,3,1],[0,1,0,1,0,100]);
    drawnow();
end
