140
9
•
Méthodes multipas
9.4 Programmation
eqdm1
✞
☎
N=5;
nomf=’f1’;
t0=0;
T=5;
eta=1;
[t,u]=multip(nomf,t0,T,N,eta)
tt=t0:T/500:t0+T;
yex= exp(tt);
plot(t,u,’rx’,tt,yex,’b’)
✝
✆
✞
☎
function[t,u]=multip(nomf,t0,T,N,eta)
h=T/N;
t=[t0:h:t0+T];
u(1)=eta;
x=u(1)+h*feval(nomf,t0,u(1));
u(3)=u(1)+2*h*feval(nomf,t(2),x);
u(2)=u(1)/4+x/2+u(3)/4;
c1=9*h/4;
c2=3*h/4;
for n=3:length(t)-1
u(n+1)=u(n-2)+c1*feval(nomf,t(n),u(n))+c2
* feval(nomf,t(n-2),u(n-2));
end
✝
✆
Dans le cas où y
( t ) = ( − y ( t ) t
2
− y
2
( t )+2t ) / (1 − t
3
)avec y (0) = 1, la solution estdonnée par
y ( t ) = ( t
2
+1) / ( t +1)auvoisinage de 0. Apriori l’équation ne peut être vérifiée en 1, mais la
solution, elle peut être prolongée. En choisissantune valeur T > 1etun N qui évite la valeur
t i = 1, Matlab effectue des calculsaudelà de 1mais les résultats n’ont plus rien àvoir avec la
solution...
eqdm2.m,étude d’erreur
✞
☎
tabN=10:10:200;
for k=1:length(tabN)
[t,u]=multip(’f1’,0,5,tabN(k),1);
yex= exp(t);
taberr(k)=norm(u-yex,inf);
end
plot( log(tabN),log(taberr))
xlabel(’log(N)’,’FontSize’,16)
ylabel(’log(err)’,’FontSize’,16)
a=polyfit( log(tabN),log(taberr),1)
title([’pentedeeladroitedeeregression:’...
, num2str(a(1))],’FontSize’,16)
✝
✆