/////////////////////////////////////////////
/// Equation des ondes 1D, CL Dirichlet, schéma différences finies, tracé de l'erreur
/////////////////////////////////////////////
//
alpha=evstr(x_dialog('CFL (<1)','0.5'));
Tf=evstr(x_dialog('Temps final','1'));
//
deff('y=U0(x)','y=sin(%pi*x)');// condition initiale  sur u (u'=0)
deff('y=Uex(t,x)','y=cos(%pi*t)*sin(%pi*x)')// solution exacte
//
//
err=[];
Ntot=[];
//
///////////////// début boucle sur N (nombre de pts intérieurs)//////////
//
for N=10:20:200
//
dx=1/(N+1);
dt=alpha*dx;
//
M=int(Tf/dt);//nombre d'itérations en temps
//
A1=toeplitz([-2,1,zeros(1,N-2)]);
//
A=2*eye(N,N)+(dt^2)/(dx^2)*A1;
u0=zeros(N,1);
ufin=zeros(N,1);
//
//
for i=1:N
   u0(i)=U0(i*dx);  
   ufin(i)=Uex(dt*M,i*dx);
end
//
//besoin  des 2 premier temps de u pour demarrer le schema (n=0 et n=1)
//
//
u11=u0;// approche u au temps n=1 (u'=0)
//
Uint2=zeros(N,1);
for nt=1:M
    Uint2=u11;//sauvegarde du temps u au temps n
    u11=A*u11-u0; // mise a jour de u 
    u0=Uint2; 
    //
end
//
err=[err,max(abs(u0-ufin))];
Ntot=[Ntot,N];
//
end
//
clf();
plot2d1('gll',Ntot,err)
plot2d1('gll',Ntot,ones(Ntot)./(Ntot.^2),2);
xtitle('decroissance de l''erreur à l''instant T (noir). En bleu: pente -2')
