UNIVERSITE DE DOUALA ENSET DE DOUALA
2024/2025 DEPARTEMENT/FILIERE ET NIVEAU : ESB/IM5
CC NUMERO 2 DE CONTROL OPTIMAL :
ENSEIGNANT : Pr. BOUM ALEXANDRE
NOMS ET PRENOMS : MEZATIO TIAGHO VANEL
MATRICULE : 23NIML04A
CONSIGNE : En utilisant la méthode du Tir direct (Direct single shooting method) résolvons les
problèmes ci-dessous avec Matlab.
EXERCICE 1 :
Soit le problème de commande optimale suivant :
1
min J (u)= 1 ∫ ( 3 x 2 ( t ) +u 2(t))dt,
u 2 0
sous contrainte : ẋ (t) = - x(t) + u(t), x(0) = 0, x(1) = 2
RESOLUTON
CODE MATLAB
% Paramètres du problème
tf = 1; % Temps final
N = 100; % Nombre de points de discrétisation
t = linspace(0, tf, N); % Discrétisation du temps
dt = t(2) - t(1); % Pas de temps
% Fonction coût
J = @(u, x) 0.5 * sum(3 * x.^2 + u.^2) * dt;
% Fonction dynamique
dynamics = @(x, u) -x + u;
% Conditions initiales et finales
x0 = 0; % x(0)
xf = 2; % x(1)
% Tir initial pour u
u = zeros(1, N); % Initialisation de u
% Optimisation avec fmincon
options = optimoptions('fmincon', 'Display', 'iter', 'Algorithm', 'sqp');
u_opt = fmincon(@(u) optimize(u, dynamics, x0, xf, N, dt, J), u, [], [], [], [], [], [], [], options);
% Résultat final
x_opt = zeros(1, N);
x_opt(1) = x0;
for i = 1:N-1
x_opt(i+1) = x_opt(i) + dt * dynamics(x_opt(i), u_opt(i));
end
% Affichage des résultats
figure;
subplot(2,1,1);
plot(t, x_opt, 'b-', 'LineWidth', 1.5); hold on;
xlabel('Temps t'); ylabel('x(t)');
title('État x(t)');
subplot(2,1,2);
plot(t, u_opt, 'r-', 'LineWidth', 1.5); hold on;
xlabel('Temps t'); ylabel('u(t)');
title('Contrôle u(t)');
% Fonction pour l'optimisation
function J_total = optimize(u, dynamics, x0, xf, N, dt, J)
x = zeros(1, N);
x(1) = x0;
for i = 1:N-1
x(i+1) = x(i) + dt * dynamics(x(i), u(i));
end
% Pénalité pour contrainte finale
penalty = 1e6 * (x(end) - xf)^2;
J_total = J(u, x) + penalty;
end
RESULTAT DE LA SIMULATION DU CODE
EXERCICE2
Soit le problème de commande optimale suivant :
1
min J(u(t)) = ∫ ( x 1 + x 2+ 0.005 u (t) ) dt ,
2 2 2
Sous contraintes : x˙ 1(t) = x2(t), x1(0) = 0, x˙ 2(t) = −x2(t) + u(t), x2(0) = −1
RESOLUTION
CODE MATLAB
% Méthode de tir direct pour l'Exercice 2
clc; clear; close all;
% Paramètres du problème
tf = 1; % Temps final
N = 100; % Nombre de points de discrétisation
t = linspace(0, tf, N); % Discrétisation du temps
dt = t(2) - t(1); % Pas de temps
% Conditions initiales
x0 = [0; -1]; % Conditions initiales : [x1(0); x2(0)]
% Initialisation du contrôle
u_init = zeros(1, N); % Contrôle initial (vecteur)
% Options pour l'optimisation
options = optimoptions('fmincon', 'Display', 'iter', 'Algorithm', 'sqp', ...
'MaxFunctionEvaluations', 1e5, 'MaxIterations', 1e3);
% Optimisation du contrôle avec fmincon
u_opt = fmincon(@(u) cost_function(u, x0, t, dt), u_init, [], [], [], [], [], [], ...
@(u) constraints(u, x0, t, dt), options);
% Résolution finale des équations d'état avec u_opt
[x1_opt, x2_opt] = integrate_state(x0, u_opt, t, dt);
% Affichage des résultats
figure;
subplot(3,1,1);
plot(t, x1_opt, 'b-', 'LineWidth', 1.5);
xlabel('Temps t'); ylabel('x_1(t)');
title('État x_1(t) optimisé');
grid on;
subplot(3,1,2);
plot(t, x2_opt, 'g-', 'LineWidth', 1.5);
xlabel('Temps t'); ylabel('x_2(t)');
title('État x_2(t) optimisé');
grid on;
subplot(3,1,3);
plot(t, u_opt, 'r-', 'LineWidth', 1.5);
xlabel('Temps t'); ylabel('u(t)');
title('Contrôle u(t) optimisé');
grid on;
%% Fonction coût pour la méthode de tir direct
function J = cost_function(u, x0, t, dt)
% Intégration des états avec le contrôle donné
[x1, x2] = integrate_state(x0, u, t, dt);
% Calcul du coût J(u)
J = 0.5 * sum(x1.^2 + x2.^2 + 0.005 * u.^2) * dt;
end
%% Intégration des équations d'état
function [x1, x2] = integrate_state(x0, u, t, dt)
N = length(t);
x1 = zeros(1, N);
x2 = zeros(1, N);
x1(1) = x0(1); % Condition initiale pour x1
x2(1) = x0(2); % Condition initiale pour x2
% Intégration des équations d'état avec schéma d'Euler explicite
for i = 1:N-1
x1(i+1) = x1(i) + dt * x2(i); % Dynamique : x1_dot = x2
x2(i+1) = x2(i) + dt * (-x2(i) + u(i)); % Dynamique : x2_dot = -x2 + u
end
end
%% Contraintes du problème
function [c, ceq] = constraints(u, x0, t, dt)
% Intégration des états avec le contrôle donné
[~, x2] = integrate_state(x0, u, t, dt);
% Contrainte d'inégalité sur x2
c = x2 - (8 * (t - 0.5).^2 + 0.5); % x2 <= 8(t - 0.5)^2 + 0.5
% Pas de contrainte d'égalité (conditions aux frontières sont implicites)
ceq = [];
end
RESULTAT DE LA SILULATION DU CODE
S
