KM02_1.mw

Méthode d'Euler

> restart;

>

> suite_aet suite_t sont des variables de stockage
an et tn sont des variables de calcul (toute nouvelle valeur écrase la précédente)
A chaque tour de Boucle (i-ième tour), on calcule la nouvelle valeur approchée an (correspondant à l'indice i) idem pour tn
Puis on empile cette nouvelle valeur dans la variable de stockage
le Résultat de la procedure est une liste contenant la liste des instants et celles des valeurs approchées

> Euler:=proc(f,t_0,y_0,N,delta)
local i,z,tn,an,suite_a,suite_t;
suite_a:=y_0;suite_t:=t_0;an:=y_0;tn:=t_0;             
for i from 1 to N do                                        
   an:=an+delta*f(tn,an); tn:=tn+delta;                
  suite_a:=suite_a,an; suite_t:=suite_t,tn;
od:        
[suite_t],[suite_a];                                        
end:

>

Testons la procédure

> N:=10: delta:=1:t_0:=0;y_0:=1;f:=(t,y)->1:         

f est l'expression définisant l'equation différentielle y'=f(t,y)... la condition initiale étant y(t_0)=y_0

t_0 := 0

y_0 := 1

> EDO:=diff(y(t),t)=f(t,y(t));                         

Voici l'Equation Différentielle de référence

EDO := diff(y(t), t) = 1

> Eul:=Euler(f,t_0,y_0,N,delta):

Voici les listes correspondants aux intstants caculés et aux valeurs approchées obtenues pour ces instants

> suite_des_temps:=Eul[1];suite_des_valeurs_approchées:=Eul[2];

suite_des_temps := [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

suite_des_valeurs_approchées := [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11]

Testons la procédure

> N:=10: delta:=1:t_0:=0;y_0:=1;f:=(t,y)->y:         

t_0 := 0

y_0 := 1

> EDO:=diff(y(t),t)=f(t,y(t));                       

EDO := diff(y(t), t) = y(t)

> Eul:=Euler(f,t_0,y_0,N,delta):

> suite_des_temps:=Eul[1];suite_des_valeurs_approchées:=Eul[2];

suite_des_temps := [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

suite_des_valeurs_approchées := [1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024]

On constate qu'il s'agit là d'une suite géométrique de raison 2

Calcul approché de exp(-1/2)

> N:=10: delta:=1/N:t_0:=0:y_0:=1:f:=(t,y)->-t*y:      

> EDO:=diff(y(t),t)=f(t,y(t));                       

EDO := diff(y(t), t) = -t*y(t)

> Eul:=Euler(f,t_0,y_0,N,delta):

> suite_des_temps:=Eul[1];suite_des_valeurs_approchées:=Eul[2];

suite_des_temps := [0, 1/10, 1/5, 3/10, 2/5, 1/2, 3/5, 7/10, 4/5, 9/10, 1]

suite_des_valeurs_approchées := [1, 1, 99/100, 4851/5000, 470547/500000, 1411641/1562500, 26821179/31250000, 1260595413/1562500000, 117235373409/156250000000, 2696413588407/3906250000000, 245373636545...

Augmentons le nombre de valeurs approchées calculées sur [0,1] ce qui affine la précision de la valeur approchée de y(1)=exp(-1/2)

> N:=1000: delta:=1/N:

> Eul:=Euler(f,t_0,y_0,N,delta):

> evalf(Eul[2][N+1]-exp(-1/2));

0.2023117e-3

En calculant la méthode d'Euler avec un pas de 1/1000 on obtient une précision de l'ordre de 10^(-3) pour la valeur de y(1)

>

> point_suite:=NULL:
for N from 10 to 1000 by 50 do
if N<>0 then delta:=1/N: else delta:=0; fi;
Eul:=Euler(f,t_0,y_0,N,delta):
point_suite:=point_suite,[N,evalf(Eul[2][N+1])];
od:

On trace ici les valeurs de la suite de valeurs approchées calculées à partir de la méthode d'Euler

> plot([point_suite]);

[Plot]

>