TP XXI
Méthode d’Euler
1 Méthode d’Euler explicite
Programme XXI.1 – Euler explicite
1 def Euler ( phi , y0 , T ) :
2 """ Entree : une fonction f de 2 variables
3 un reel y0 , la condition initiale en t0
4 une liste de r é els , monotone , dont le premier
terme est t0
5 Sortie : une liste de points d é finissant
6 une solution approchee de y ’ = phi (y , t )
7 avec les conditions initiales ( t0 , y0 )
8 y ( T [ i ]) est approch é par Y [ i ] """
9 n = len ( T ) # n valeurs à calculer
10 Y = [0]* n # cr é ation de la liste
11 Y [0] = y0 # ordonnee initiale
12 for k in range (n -1) : # il reste (n -1) points à trouver
13 pas = T [ k +1] - T [ k ]
14 pente = phi ( Y [ k ] , T [ k ])
15 Y [ k +1] = Y [ k ] + pas * pente # On applique la formule
16 return Y
1.1 Peu de points
On commence par illustrer la méthode d’Euler sur un exemple pour lequel on connaît la solution.
Pour mettre en évidence la convergence des calculs on ne considérera pas un grand nombre de
points. Cette situation est l’opposée de l’utilisation usuelle de la méthode : si on veut approcher
une solution qu’on ne sait pas calculer on choisit un grand nombre de points.
On considère l’équation y � = y + t avec y(0) = 0.5 dont la solution est t �→ 12 et − t − 1.
Exercice XXI.1
Écrire la fonction phi1(y, t) qui renvoie y + t
et la fonction f (t) qui renvoie 12 et − t − 1
XXI. Méthode d’Euler
Exercice XXI.2
Compléter les instructions suivantes pour qu’elles calculent les solutions calculées sur [0; 3] de
l’équation pour des listes de n points avec n ∈ {3, 6, 10, 20, 40, 100} et les affichent sur un même
graphe puis calcule les listes avec 10000 points qui permettent de tracer le graphe de la solution,
sur le même graphe.
Programme XXI.2 – Tracés multiples
liste_n = [3 , ............]
p = len ( liste_n )
for n in liste_n :
T = np . linspace (......)
Y = ........... # solution approch é e avec Euler
plt . plot (T , Y , label = " Euler avec {} points " . format ( n )
T = ...........
Y = sol ( T ) # solution exacte
plt . plot (T , Y , label = " Solution exacte " )
plt . legend ()
plt . show ()
1.2 Étude d’un filtre actif passe bas
R1
On peut réaliser un filtre actif passe-
ue −
bas avec un A.O. (Amplificateur
opérationnel). us
Le fonctionnement de l’amplificateur +
opérationnel implique que la tension R3
de sortie, us , est solution de C
R2 + R 3
R1 Cy � (t) + y(t) = ue (t) R2
R2
où ue est la tension d’entrée.
Les valeurs numériques sont les suivantes :R2 = R3 = 1000Ω et C = 1µF.
Différentes valeurs de R1 seront testées afin de mettre en évidence les différents comportements.
Le signal d’entrée ue (t) est un signal rectangulaire de fréquence 500 Hz.
Exercice XXI.3
Écrire une fonction creneau(t,T) qui calcule la valeur d’un signal rectangulaire de période T : la
valeur devra être 1 sur [0; T2 [ et −1 sur [ T2 ; T [.
On pourra utiliser l’expression t%T qui renvoie le réel x appartenant à [0; T [ tel que t − x est un
multiple entier de T .
XXI-2
XXI.2-Méthode d’Euler implicite
Exercice XXI.4
En utilisant la méthode d’Euler déterminer alors les solutions de l’équation : on prendra une liste
de temps de 1000 points entre 0 et 0,02s avec une condition initiale nulle.
On calculera les solutions pour R1 = 100Ω, R1 = 200Ω, R1 = 500Ω, R1 = 1000Ω et R1 = 2000Ω.
On pourra tracer les solutions en même temps que le signal d’entrée.
2 Méthode d’Euler implicite
La méthode d’Euler implicite utilise la formule
yk+1 − yk = (tk+1 − tk )ϕ(yk+1 , tk+1 )
Dans cette équation les temps, tk et tk+1 , sont connus ainsi que la valeur (approchée) de la solution
en tk , yk . Elle demande de résoudre une équation en yk+1 .
Elle permet parfois une meilleure stabilité de la solution.
2.1 Exemple simple
�
y � = −20y
On peut illustrer la stabilité avec l’équation .
y(0) = y0
Exercice XXI.5
Déterminer et tracer la solution approchée de l’équation ci-dessus par la méthode d’Euler explicite
sur [0; 10] avec 100 points pour y0 = 1.
On pourra aussi regarder le comportement pour 101 points puis 102.
On ne retrouve pas le résultat attendu, la solution est t �→ e−20t .
Exercice XXI.6
Résoudre (mathématiquement) l’équation d’inconnue yk+1 ,
yk+1 = yk + (tk+1 − tk )ϕ(yk+1 , tk+1 ) = yk − 20(tk+1 − tk )yk+1
En déduire une fonction solution1(y0, T) de résolution de l’équation en fonction de y0 et de la
liste des temps puis tracer la solution pour y0 = 1 avec 100 points entre 0 et 1.
2.2 Autre exemple
�
y � = ty − y 2
On considère l’équation .
y(0) = y0
Exercice XXI.7
Montrer que la méthode d’Euler implicite donne le calcul
�
htk+1 − 1 ± (htk+1 − 1)2 + 4hyk
yk+1 =
2h
où h est le pas de temps tk+1 − tk . �
htk+1 − 1 + (htk+1 − 1)2 + 4hyk
Pourquoi est-il judicieux de choisir la solution ?
2h
Exercice XXI.8
Déterminer une fonction de résolution de l’équation (E2 ) en fonction de y0 et de la liste des temps :
solution2(y0, T).
Comparer graphiquement les solutions par les deux méthodes (implicite et explicite) sur [0; 4] avec
50 puis 500 points pour y0 = 10.
XXI-3
XXI. Méthode d’Euler
3 Autres méthodes de résolution
3.1 Bibliothèque scipy
Dans la bibliothèque scipy, on dispose d’une fonction odeint
from scipy . integrate import odeint
La fonction odeint est alors substituable, sans modification, à la fonction Euler.
Elle donnera des résultats plus précis.
3.2 Schéma de Heun
La méthode de Heun améliore le schéma d’Euler en approchant l’intégrale par la méthode des
� b
g(a) + g(b)
trapèzes : g(u)du � (b − a) . On obtient
a 2
� � � �
ϕ y(tk ), tk + ϕ y(tk+1 ), tk+1
y(tk+1 ) − y(tk ) � (tk+1 − tk )
2
Cependant cela devient une équation implicite en y(tk+1 ) qu’on ne souhaite pas résoudre.
On va donc procéder en approchant la valeur de y(tk+1 ) dans le second membre par celle que
calcule la méthode d’Euler.
• On calcule zk = yk + ϕ(yk , tk )(tk+1 − tk ) comme valeur approchée de y(tk+1 )
ϕ(yk , tk ) + ϕ(zk , tk+1 )
• On calcule la valeur approchée yk+1 = yk + (tk+1 − tk ).
2
On remarque que zk est la valeur approchée par la méthode d’Euler.
Exercice XXI.9
Écrire une fonction Heun(phi, y0, T) qui approche, aux points de T, une solution de l’équation
y � = ϕ(y, t) avec la condition initiale y(t0 ) = y0 où t0 est la première valeur de la liste T. Le
docstring est le même que celui de la fonction Euler.
3.3 Méthode de Runge-Kutta
La méthode de Runge-Kutta s’inspire de la méthode de Simpson :
� b � �
g(a) + 4g a+b2 + g(b)
g(u)du � (b − a)
a 6
tk +tk+1
On obtient, pour mk = 2 et h = tk+1 − tk ,
� � � � � �
ϕ y(tk ), tk + 4ϕ y(mk ), mk + ϕ y(tk+1 ), tk+1
y(tk+1 ) � y(tk ) + h
6
On va, ici encore, approcher provisoirement les valeurs de y(mk ) et y(tk+1 ).
On va en fait utiliser
� deux � approximations de y(mk ).
• a = yk + h2 ϕ yk , tk est une première valeur approchée de y(mk )
� �
• b = yk + h2 ϕ a, tk + h2 est aussi une valeur approchée de y(mk )
� �
• c = yk + hϕ b, tk + h2 est une valeur approchée de y(tk+1 ) obtenue en prenant comme valeur
moyenne la valeur au point milieu.
• On pose alors � � � � � � � �
ϕ yk , tk + 2ϕ a, tk + h2 + 2ϕ b, tk + h2 + ϕ c, tk + h
yk+1 = yk + h
6
Exercice XXI.10
Écrire une fonction RK(phi, y0, T) avec le même docstring que celui de la fonction Euler.
XXI-4
XXI.4-Solutions
4 Solutions
Solution de l’exercice XXI.1 -
def phi1 (y , t ) :
return y + t
def sol ( t ) :
return exp ( t ) /2 - t - 1
Solution de l’exercice XXI.2 -
liste_n = [3 , 6 , 10 , 20 , 40 , 100]
for n in liste_n :
T = np . linspace (0 , 3 , n )
Y = Euler ( phi1 , 0.5 , T )
plt . plot (T , Y , label = " Euler avec {} points " . format ( n ) )
T = np . linspace (0 , 3 , 10000)
Y = sol ( T )
plt . plot (T , Y , label = " Solution exacte " )
plt . legend ()
plt . show ()
Solution de l’exercice XXI.3 -
def creneau (t , T ) :
x = t%T
if x < T /2:
return 1
else :
return -1
XXI-5
XXI. Méthode d’Euler
Solution de l’exercice XXI.4 -
def filtre (y , t ) :
return (( R2 + R3 ) / R2 * creneau (t , per ) - y ) / R1 / C
R2 = 1000
R3 = 1000
f = 500
per = 1/ f
C = 1e -6
R = [100 , 200 , 500 , 1000 , 2000]
T = liste_abs (0 , 0.02 , 1000)
n = len ( R )
for i in range ( n ) :
R1 = R [ i ]
Y = Euler ( filtre , 0 , T )
plt . plot (T , Y , label = " $R_1$ = {} " . format ( R1 ) )
def cr ( t ) :
return creneau (t , per )
Y = appliquer ( cr , T )
plt . plot (T , Y , label = " Entr é e " )
plt . legend ()
plt . show ()
XXI-6
XXI.4-Solutions
Solution de l’exercice XXI.5 -
def phi1 (y , t ) :
return -20* y
for n in range (100 , 103) :
T = liste_abs (0 , 10 , n )
Y = Euler ( phi1 , 1 , T )
plt . plot (T , Y , label = " Euler explicite avec {} points " .
format ( n ) )
Y = [ exp ( -20* t ) for t in T ]
plt . plot (T , Y , label = " Solution exacte " )
plt . legend ()
plt . show ()
yk
Solution de l’exercice XXI.6 - L’équation donne yk+1 = .
1 + 20(tk+1 − tk )
def solution1 ( y0 , T ) :
n = len ( T )
Y = [0]* n
Y [0] = y0
for k in range (n -1) :
pas = T [ k +1] - T [ k ]
Y [ k +1] = Y [ k ]/(1+20* pas )
return Y
T = liste_abs (0 , 10 , 100)
Yi = solution1 (1 , T )
plt . plot (T , Yi )
plt . show ()
Solution de l’exercice XXI.7 - � �
2
L’équation devient yk+1 − yk = h. tk+1 yk+1 − yk+1 donc yk+1 est solution de hY 2 − (htk+1 −
1)Y − yk = 0 d’où la formule.
XXI-7
XXI. Méthode d’Euler
Si on fait un développement
� limité à l’ordre 1� en h on trouve
htk+1 − 1 ± 1 + 2hyk − htk+1 + o(h)
yk+1 = .
2h
Seul le signe + donne une limite finie (yk ) quand h tend vers 0.
Solution de l’exercice XXI.8 -
def solution2 ( y0 , T ) :
n = len ( T )
Y = [0]* n
Y [0] = y0
for k in range (n -1) :
pas = T [ k +1] - T [ k ]
delta = ( pas * T [ k +1] -1) **2 + 4* pas * Y [ k ]
Y [ k +1] = ( pas * T [ k +1] - 1 + delta **0.5) /(2* pas )
return Y
def phi2 (y , t ) :
return t * y - y **2
for n in [50 , 500]:
T = liste_abs (0 , 4 , n )
Ye = Euler ( phi2 , 10 , T )
Yi = solution2 (10 , T )
plt . plot (T , Ye , label = " M é thode explicite avec {} points "
. format ( n ) )
plt . plot (T , Yi , label = " M é thode implicite avec {} points "
. format ( n ) )
plt . legend ()
XXI-8
XXI.4-Solutions
Solution de l’exercice XXI.9 -
def Heun ( phi , y0 , T ) :
n = len ( T )
Y = [0]* n
Y [0] = y0
for k in range (n -1) :
pas = T [ k +1] - T [ k ]
pente1 = phi ( Y [ k ] , T [ k ])
z = Y [ k ] + pas * pente1
pente2 = phi (z , T [ k +1])
pente = ( pente1 + pente2 ) /2
Y [ k +1] = Y [ k ] + pas * pente
return Y
Solution de l’exercice XXI.10 -
def RK (f , y0 , T ) :
n = len ( T )
Y = [0]* n
Y [0] = y0
for i in range (n -1) :
y = Y[i]
t = T[i]
t1 = T [ i +1]
m = ( t + t1 ) /2
pas = t1 - t
pente1 = f (y , t )
a = y + pente1 * pas /2 # 1 è re valeur au milieu
pente2 = f (a , m ) # 1 è re pente au milieu
b = y + pente2 * pas /2 # 2 è me valeur au milieu
pente3 = f (b , m ) # 2 è me pente au milieu
c = y + pente3 * pas # valeur à droite
pente4 = f (c , t1 ) # pente à droite
pente = ( pente1 + 2* pente2 + 2* pente3 + pente4 ) /6
Y [ i +1] = y + pas * pente
return Y
XXI-9