TP 5 : Résolution numérique des systèmes différentiels
Exercice 1 :
Soit l’équation différentielle suivante :
𝑢̇ (t) − t 2 − 𝑢2 − 1 = 0
{
𝑢(0) = 0
1- Ecrire l’équation sous forme de Cauchy.
2- Ecrire un programme qui donne une approximation de 𝑢(0.2) avec un pas ℎ = 0.1𝑠 par
les algorithmes de : Euler explicite, Euler Implicite et Runge Kutta d’orde 2
%%%%%===============================================================
% Schemas d'integration temporelle
%%%%%===============================================================
clear
clc
% Donnees initiales
u0=0;
t0=0;
h=0.1;
% Fonction de Cauchy
du=@(u,t) t^2+u^2+1;
uex=zeros(1, 4); % Solutions d'Euler explicite
uim=zeros(1, 4); % Solutions d'Euler implicite
urk=zeros(1, 4); % Solutions de Runge-Kutta 2
t=zeros(1, 4); % Vecteur temps
% Initialisation des vecteurs solutions
uex(1)=u0; uim(1)=u0; urk(1)=u0;
t(1)=t0;
% 2- Calcul
%pour n = 1; t1 = 0.1
t(2)=t(1)+h;
% Euler explicite
uex(2)=uex(1)+h*du(t(1),uex(1));
% Euler implicite
uim(2)=uim(1)+h*du(t(2),uex(2)); % Remarquez qu'on a utilisé
uex(2) à droite car uim(2) est une inconnue
% Rk2
k1=du(t(1),urk(1));
k2=du(t(1)+h,urk(1)+h.*k1);
urk(2)=urk(1)+(h/2).*(k1+k2);
%pour n = 2; t1 = 0.2
t(3)=t(2)+h;
% Euler explicite
uex(3)=uex(2)+h*du(t(1),uex(2));
% Euler implicite
1 2022/2023 Pr O. Bourihane
uim(3)=uim(2)+h*du(t(2),uex(3));
% Rk2
k1=du(t(2),urk(2));
k2=du(t(2)+h,urk(2)+h.*k1);
urk(3)=urk(2)+(h/2).*(k1+k2);
fprintf('A l''instant t = %3.2f s\n', t(3))
fprintf('Solution d''Euler explicite: u(0.2s) = %6.5f\n', uex(3))
fprintf('Solution d''Euler implicite: u(0.2s) = %6.5f\n', uim(3))
fprintf('Solution de Runge-Kutta 2: u(0.2s) = %6.5f\n', urk(3))
Exercice 2 :
Considérant l’équation différentielle suivante :
𝑢̇ − 𝑡 × 𝑠𝑖𝑛(𝑢) = 1
{
𝑢(0) = 1
1- Ecrivez des fonctions qui permettent la résolution d’un système différentiel dans un
intervalle de temps quelconque, par les méthodes citées dans l’exercice 1.
%%%%%===============================================================
Euler Explicite
%%%%%===============================================================
function [t,u]=EulerExp(t0,tf,h,u0,du)
% Resolution par methode Euler Explicite
n=ceil((tf-t0)/h);
t=zeros(1, n+1);
u=zeros(1, n+1);
u(1)=u0;
t(1)=t0;
%%%%%
for i=1:n
t(i+1)=t(i)+h;
u(i+1)=u(i)+h*du(t(i),u(i));
end
end
%%%%%===============================================================
Euler Implicite
%%%%%===============================================================
n=ceil((tf-t0)/h);
t=zeros(1, n+1);
u=zeros(1, n+1);
u(1)=u0;
t(1)=t0;
%%%%%
for i=1:n
t(i+1)=t(i)+h;
uex = u(i)+h*du(t(i),u(i));
2 2022/2023 Pr O. Bourihane
u(i+1)=u(i)+h*du(t(i+1),uex); % De meme u(i+1) est
inconnue, on l'approxime par uex
end
end
%%%%%===============================================================
Runge-kutta d’ordre 2
%%%%%===============================================================
function [t,u]=RungeKutta2(t0,tf,h,u0,du)
% Resolution par methode de Runge-Kutta
n=ceil((tf-t0)/h);
t=zeros(1, n);
u=zeros(1, n);
u(1)=u0;
t(1)=t0;
%%%%%
for i=1:n
t(i+1)=t(i)+h;
k1 = du(t(i),u(i));
k2 = du(t(i)+h,u(i)+h*k1);
u(i+1)=u(i) + h/2*(k1+k2);
end
end
2- Comparez les résultats des trois méthodes dans l‘intervalle [0 10𝑠] en considérant les
pas de temps suivants : ℎ = 0.1; 0.05; 0.01𝑠.
%%%%%=========================================================
======
% Comparaison des Schemas
%%%%%=========================================================
======
% Main
clear
clc
% Donnees initiales
u0=0;
t0=0;
tf = 12;
h1=0.1; h2=0.05; h3=0.01;
% Fonction de Cauchy
du=@(u,t) t*sin(u)+1;
% Calcul
% h = 0.1
[~,uex1]=EulerExp(t0,tf,h1,u0,du);
[~,uim1]=EulerImp(t0,tf,h1,u0,du);
[t1,urk1]=RungeKutta2(t0,tf,h1,u0,du);
% h = 0.05
[~,uex2]=EulerExp(t0,tf,h2,u0,du);
3 2022/2023 Pr O. Bourihane
[~,uim2]=EulerImp(t0,tf,h2,u0,du);
[t2,urk2]=RungeKutta2(t0,tf,h2,u0,du);
% h = 0.1
[~,uex3]=EulerExp(t0,tf,h3,u0,du);
[~,uim3]=EulerImp(t0,tf,h3,u0,du);
[t3,urk3]=RungeKutta2(t0,tf,h3,u0,du);
% Traçage
figure(1)
hold on
plot(t1,uex1,'b',t1,uim1,'r',t1,urk1,'k')
plot(t2,uex2,'b',t2,uim2,'r',t2,urk2,'k')
plot(t3,uex3,'b',t3,uim3,'r',t3,urk3,'k')
hold off
Exercice 3 : Homework
Résolvez l’équation différentielle ci-dessous par les méthodes citées avant :
𝑚𝑥̈ + 𝑐𝑥̇ + 𝑘𝑥 = 𝐹
{ 𝑥(𝑡0 ) = 𝑥0
𝑥̇ (𝑡0 ) = 𝑣0
Etudiez l’effet des paramètres 𝑚, 𝑐, 𝑘, 𝐹. 𝑥̇
4 2022/2023 Pr O. Bourihane