%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%MAO 2004/2005
%   Animation de l'equation de la chaleur
%
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%function chaleur(M)
K=1:M;
%coefficients de la donnee initiale
ak0=-2*(sin(K*pi)-2*sin(1/2*K*pi))./(pi^2*K.^2);
%coefficients de la source f(x)=10*exp(-x)
fk=-10*(K*pi*exp(-1).*cos(K*pi)-K*pi)./(1+K.^2*pi^2);
%coeff de f(x)=exp(-x)*sin(2*Pi*x);
fk=10*2*pi*(-2*K*pi*exp(-1).*cos(K*pi)+2*K*pi)./((1+4*pi^2-4*K*pi^2+K.^2*pi^2).*(1+4*pi^2+4*K*pi^2+K.^2*pi^2));
% paramètres d'échantillonnage en temps
tmin=0;tmax=0.5;         % limites en temps
nt=40;dt=(tmax-tmin)/nt; % pas d'échantillonnage en temps
%
% paramètres d'échantillonnage en espace
xmin=0;xmax=1;        % limites suivant x
ymin=-1;ymax=+1;      % limites suivant y
dx=(xmax-xmin)/1000;  % pas d'échantillonnage suivant x
x=xmin:dx:xmax;       % échantillonnage suivant x
%
% initialisation du graphique
S = 0*[1:length(x)];
plot(x,S);
axis([xmin,xmax,ymin,ymax]);
%
for it=1:nt
%
    t = tmin+(it-1)*dt
    S = 0*[1:length(x)]; 
    %
    
        for n=1:M
        S = S+(ak0(n)*sin(pi*n*x))*exp(-n^2*pi^2*t) +...
        ((1-exp(-n^2*pi^2*t))/(pi^2*n^2))*fk(n)*sin(pi*n*x);
        end
        %
        
    %représentation graphique
    plot(x,S,'r','Erasemode','xor');
    axis([xmin,xmax,ymin,ymax]);drawnow;pause(0.1)
%
end

%
% Solution exacte
%
uex=(-4*pi*exp(-x).*cos(2*pi*x)+4*pi^2*exp(-x).*sin(2*pi*x)-exp(-x).*sin(2*pi*x)+4*(exp(-1)-1)*pi*x/(1+8*pi^2+16*pi^4)+32*(exp(-1)-1)*pi^3*x/(1+8*pi^2+16*pi^4)+64*(exp(-1)-1)*pi^5*x/(1+8*pi^2+16*pi^4)+4*pi/(1+8*pi^2+16*pi^4)+32*pi^3/(1+8*pi^2+16*pi^4)+64*pi^5/(1+8*pi^2+16*pi^4))/(1+8*pi^2+16*pi^4);