Université A.
Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
TD n°2
Résolution Numérique des Equations Différentielles
Exercice 1
Trouver par les méthodes (Euler, point milieu, Euler modifié et Runge-kutta d’ordre 4) la solution de l’équation
suivante:
y (t ) = − y (t ) + t + 1
on prend h = 0.1
y (0) = 1
• Sur un tableau, écrire les résultats obtenus pour chaque méthode (6 itérations) en calculant l’erreur relative
sachant que la solution analytique est : y(t ) = e−t + t .
• Donner l’algorithme de chaque méthode.
Exercice 2
• Trouver la solution de l’équation différentielle suivante :
y (t ) = 2ty (t ) + 2 y (t )
y (0) = 1, y (0) = 0
• Proposer l’algorithme de résolution.
Exercice 3
Soit le système électrique schématisé par la figure (fig.1).
Fig.1
• Trouver le courant électrique aux instants 0.1s et 0.2s et 0.3s, sachant que la condition initiale au
problème est i(o)=0A.
Mme ABDELLI.R Page 1
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
RAPPEL
( t0 , y0 ) donné
• Méthode d’Euler :
yn +1 = yn + hf (tn , yn )
( t0 , y0 ) donné
h
Méthode du point milieu : yn +1 = yn + 2 f ( tn , yn )
*
•
h
yn +1 = yn + hf tn + , yn* +1
2
( t0 , y0 ) donné
• Méthode d’Euler-Cauchy (Euler modifié):
yn +1 = yn + hf ( tn , yn )
*
h
(
yn +1 = yn + 2 f (tn , yn ) + f (tn +1 , yn +1 )
*
)
k1 = hf (tn , yn )
k
h
k2 = hf (tn + , yn + 1 )
2 2
h k2
k3 = hf (tn + , yn + )
• Méthode de Runge-Kutta d’ordre4: 2 2
k4 = hf (tn + h, yn + k3 )
y = y + 1 ( k + 2k + 2k + k )
n +1 n
6
1 2 3 4
tn +1 = tn + h
Solution
Exo 1
L’objectif de cet exercice est de comparer entre les différentes méthodes de résolution des équations différentielles
les plus utilisées commençant par la méthode d’Euler qui est une méthode classique est basique mais qui reste
moins précise par rapport aux deux autres méthodes d’Euler avancées à savoir la méthode du point milieu et
d’Euler modifié, ces méthodes sont précises mais leur principal inconvénient est la dépendance vis-à-vis le pas du
calcul, à l’inverse de la méthode de Runge Kutta qui est plus précise, ce qui nous permettra de choisir un pas assez
grand pour R-K-4 par rapport aux 2 autres méthodes d’Euler.
Méthode d’Euler
y(0.1) y(0) + 0.1(− y(0) + 0 + 1) 1
y(0.2) y(0.1) + 0.1(− y(0.1) + 0.1 + 1) 1.01
Méthode d’Euler modifié
Mme ABDELLI.R Page 2
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
y* (0.1) y (0) + 0.1(− y (0) + 0 + 1) 1
0.1
y (0.1) y (0) + (− y (0) + 0 + 1) + (− y* (0.1) + 0.1 + 1 1.005
2
y (0.2) y (0.1) + 0.1(− y (0.1) + 0.1 + 1) 1.0145
*
0.1
y (0.2) = y (0.1) + (− y (0.1) + 0.1 + 1) + (− y* (0.2) + 0.2 + 1 1.019025
2
Méthode du point milieu
y* (0.1) y (0) + 0.05(− y (0) + 0 + 1) 1
0.1
y (0.1) y (0) + 0.1 − y* (0.1) + + 1 1.005
2
y (0.2) y (0.1) + 0.05(− y (0.1) + 0.1 + 1) 1.00975
*
y (0.2) y (0.1) + 0.1 − y* (0.2) + 0.15 + 1 1.019025
Méthode de Runge Kutta d’ordre 4
Pour y(0.1)
k1 = 0.1(− y (0) + 0 + 1) = 0
0
k 2 = 0.1 −( y (0) + ) + 0.05 + 1 = 0.005
2
0.005
k 3 = 0.1 −( y (0) + ) + 0.05 + 1 = 0.00475
2
k 4 = 0.1 −( y (0) + 0.00475) + 0.1 + 1 = 0.009525
1
y (0.1) y (0) + 0 + 2(0.005 + 0.00475) + 0.009525 = 1.0048375
6
Pour y(0.2)
k1 = 0.1(− y (0.1) + 0.1 + 1) = 0.00951625
k1
k 2 = 0.1 −( y (0.1) + ) + 0.15 + 1 = 0.014040437
2
k2
k 3 = 0.1 −( y (0.1) + ) + 0.15 + 1 = 0.013814228
2
k 4 = 0.1 −( y (0.1) + k 3) + 0.2 + 1 = 0.018134827
1
y (0.2) y (0.1) + 0 + 2(k 2 + k 3) + k 4 = 1.018730901
6
On continue les calculs et on trace le tableau :
Euler Euler modifié Point milieu R-K-4 Solution
analytique
Y(0.1) 1 1.005 1.005 1.0048375 1.004837418
Y(0.2) 1.01 1.019025 1.019025 1.018730901 1.018730753
Y(0.3)
Mme ABDELLI.R Page 3
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
On peut remarquer que la méthode de R-K-4 est la plus précise, les deux méthodes d’Euler
modifié et du point milieu donnent des résultats très satisfaisants mais la méthode d’Euler n’est
pas du tout précise.
Les algorithmes
• Méthode d’Euler
a=input('donner a=') : On donne les bornes de l’intervalle de te temps
b=input('donner b=')
n=input('donner n=') : on précise le nombre de points
h=(b-a)/n; : on calcule le pas
y(1)=input('donner la condition initiale')
t(1)=a;
for i=2:n+1 : la boucle de remplissage des vecteurs
y(i)=y(i-1)+h*f(y(i-1),t(i-1));
t(i)=t(i-1)+h;
end
• Méthode d’Euler modifié
a=input('donner a=')
b=input('donner b=')
n=input('donner n=')
h=(b-a)/n;
y(1)=input('donner la condition initiale')
t(1)=a;
for i=2:n+1
y0=y(i-1)+h*f(t(i-1),y(i-1));
y(i)=y(i-1)+(h/2)*(f(t(i-1),y(i-1))+f(t(i)+h,y0);
t(i)=t(i-1)+h;
end
• Méthode du point milieu
a=input('donner a=')
b=input('donner b=')
n=input('donner n=')
h=(b-a)/n;
y(1)=input('donner la condition initiale')
t(1)=a;
for i=2:n+1
y0=y(i-1)+(h/2)*f(t(i-1),y(i-1));
y(i)=y(i-1)+h*(f(t(i)+(h/2),y0);
t(i)=t(i-1)+h;
end
• Méthode de Runge-Kutta
a=input('donner a=')
b=input('donner b=')
n=input('donner n=')
h=(b-a)/n;
y(1)=input('donner la condition initiale')
t(1)=a;
for i=2:n+1
k1=h*f(t(i-1),y(i-1));
k2=h*f((t(i-1)+h/2),(y(i-1)+k1/2));
k3=h*f((t(i-1)+h/2),(y(i-1)+k2/2));
k4=h*f((t(i-1)+h),(y(i-1)+k3));
y(i)=y(i-1)+(k1+2*k2+2*k3+k4)/6;
t(i)=t(i-1)+h;
Mme ABDELLI.R Page 4
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
end
Exo3
Dans cet exercice, on étudiera un simple circuit électrique dont la recherche du courant
électrique est demandée.
L’écriture de l’équation régissant le fonctionnement de ce circuit nous donne le problème
différentiel suivant :
di (t )
= 100 − 10i (t )
dt
i (0) = 0
−10t
La solution analytique de cette équation est connue i(t ) = 10(1 − e ) , l’objectif de cet exercice
est de résoudre le problème avec les méthodes numériques et de comparer pour tirer des
conclusions et approfondir la compréhension.
Pour calculer la valeur du courant à l’instant t=0.1s donc i(0.1), on peut choisir entre les trois
méthodes étudiées (point milieu ou euler modifié ou R-K-4) mais il faut faire très attention aux
pas choisis, on sait bien que les méthodes d’Euler dépendent fortement du pas de calcul donc
on peut utiliser une des deux méthodes en prenant un pas de calcul de 0.01s par exemple (ce
pas reste très grand, on choisit de tels pas car on ne peut pas encore réduire puisque les méthodes
ne sont pas programmées et le calcul se fait à la main) mais si on choisit la méthode R-K-4,cette
méthode puisque elle ne dépend pas trop du pas de calcul, on peut prendre par exemple un pas
de 0.05 ou 0.025 pour des systèmes qui nécessitent trop de précision.
Pour cette solution, je prends un pas h=0.05 car le calcul du courant dans ce circuit peut tolérer
une petite erreur.
Avant de calculer i(0.1), je passe d’abords par le calcul de i(0.05).
Pour t=0.05 :
k1 = 0.05 100 − 10.i (0) = 5
k 2 = 0.05 100 − 10 ( i(0) + 2.5 ) = 3.75
3.75
k 3 = 0.05 100 − 10 i(0) + = 4.0625
2
k 4 = 0.05 100 − 10 ( i(0) + 4.0625 ) = 2.96875
1
i (0.05) = i (0) + 5 + 2(3.75 + 4.0625) + 2.96875 = 3.932291667 A
6
Pour t=0.1
k1 = 3.033854167
k 2 = 2.275390625
k 3 = 2.46500651
k 4 = 1.801350911
i (0.1) 6.318291561A
Mme ABDELLI.R Page 5
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
La solution exacte est : 6.321205588A
Exo 2
Le problème à résoudre est :
y ''(t ) = 2ty '(t ) + 2 y (t )
y (0) = 1, y '(0) = 0
u1 (t ) = y (t ) u1 '(t ) = y '(t ) = u2 (t )
On prend :
u2 (t ) = y '(t )
Donc le nouveau problème à résoudre est donné par :
u1 '(t ) = u2 (t )
u2 '(t ) = 2t.u2 (t ) + 2u1 (t )
u (0) = 1, u (0) = 0
1 2
On utilise bien sur la méthode de Runge-Kutta avec un pas h quelconque :
Pour calculer y(h), on aura deux étapes :
u2 ( h ) = ?
k1 = h. ( 2*0* u2 (0) + 2u1 (0) )
u1 (h) = ?
k = h.u (0) = 0 k2 = h. 2* 0 + * u2 (0) + 1 + 2u1 (0)
h k
1 2
2 2
k2 = h.u2 (0) = 0
h k2
k3 = h.u2 (0) = 0 k3 = h. 2* 0 + * u2 (0) + + 2u1 (0)
k = h.u (0) = 0 2 2
4 2 k = h. ( 2* ( 0 + h ) * ( u (0) + k ) + 2u (0) )
u1 (h) = 1 4 2 3 1
1
u2 (h) = u2 (0) + ( k1 + 2k2 + 2k3 + k4 )
6
Exemple : h=0.1
Mme ABDELLI.R Page 6
Université A.Mira de BEJAIA
Faculté de Technologie, Département GE
Module UEF 1122
Année 2022/2023
y (0.1) = u2 (0.1) = ?
k1 = 0.1. ( 2*0* u2 (0) + 2u1 (0) ) = 0.2
u1 (0.1) = ?
0.2 k = 0.1.u (0) = 0
k2 = 0.1. 2* ( 0 + 0.05 ) * 0 + + 2 = 0.201
2 1 2
k2 = 0.1.u2 (0) = 0
0.201
k3 = 0.1. 2* ( 0 + 0.05 ) * 0 + + 2 = 0.201005 k3 = 0.1.u2 (0) = 0
2 k = 0.1.u (0) = 0
k = 0.1. ( 2* ( 0 + 0.1) * ( 0 + 0.201005 ) + 2 ) 4 2
4
u1 (0.1) = 1
1
y (0.1) = u2 (0.1) = 0 + ( k1 + 2k2 + 2k3 + k4 ) = 0.20133835
6
y (0.2) = u2 (0.2) = ?
k1 = 0.1. ( 2*0.1* u2 (0.1) + 2u1 (0.1) )
k2 = 0.1. 2* ( 0.1 + 0.05 ) * u2 (0.1) + 1 + 2u1 (0.1)
k
2
k2
k3 = 0.1. 2* ( 0.1 + 0.05 ) * u2 (0.1) + + 2u1 (0.1)
2
k = 0.1. ( 2* ( 0.1 + 0.1) * ( u (0.1) + k ) + 2u (0.1) )
4 2 3 1
1
y (0.2) = u2 (0.2) = u2 (0.1) + ( k1 + 2k2 + 2k3 + k4 )
6
u1 (0.2) = ?
k1 = 0.1.u2 (0.1) = 0.020133835
k2 = 0.1.u2 (0.1) = 0.020133835
k3 = 0.1.u2 (0.1) = 0.020133835
k = 0.1.u (0.1) = 0.020133835
4 2
1
u1 (0.2) = 1 + (4*0.020133835) = 1.013422556667
6
Mme ABDELLI.R Page 7
