% TP : Commande basée sur un observateur d'état

clear; clc;

% === 1. Définition des matrices du système ===

A = [0 -1; 1 -1]; % Matrice dynamique

B = [1; 0]; % Matrice d'entrée

C = [0 1]; % Matrice de sortie

% === 2. Gain de l'observateur ===

% Les pôles souhaités pour l'observateur sont : -1 ± j*sqrt(3)

poles_obs = [-1 + 1i*sqrt(3), -1 - 1i*sqrt(3)];

L = place(A', C', poles_obs)'; % Calcul du gain de l'observateur

% Affichage du gain de l'observateur

disp('Gain de l''observateur :');

disp(L);

% === 3. Fonction de transfert du système avec observateur ===

% Modification de la matrice A pour inclure l'observateur

A_obs = A - L * C;

% Déterminant du système

syms p;

determinant = det(p * eye(2) - A_obs);

disp('Dénominateur de la fonction de transfert (polynôme caractéristique) :');

disp(determinant);

% Fonction de transfert H(p) = C * (pI - A_obs)^(-1) * B

H = simplify(C * inv(p * eye(2) - A_obs) * B);

disp(H);

% === 4. Simulation de la réponse indicielle ===

% Conversion de H(p) en une fonction de transfert MATLAB

num = [1]; % Numérateur : 1

den = [1 2 4]; % Dénominateur obtenu dans les calculs : p^2 + 2p + 4

sys = tf(num, den); % Fonction de transfert

% Réponse indicielle

figure;

step(sys);

title('Réponse indicielle du système');

xlabel('Temps (s)');

ylabel('Amplitude');

% === 5. Gain de retour d'état ===

% Les pôles souhaités pour le système bouclé sont : -0.5 ± j*sqrt(3)/2

poles_ctrl = [-0.5 + 1i*sqrt(3)/2, -0.5 - 1i*sqrt(3)/2];

K = place(A, B, poles_ctrl); % Calcul du gain de retour d'état

% Affichage du gain de retour d'état

disp('Gain de retour d''état :');

disp(K);

% === 6. Commande basée sur l'observateur ===

% Matrice augmentée pour la commande

A_cl = A - B * K - L * C;

% Fonction de transfert pour le système commandé

disp('Fonction de transfert avec commande basée sur observateur H_aug(p) :');

disp(H_aug);

% Réponse indicielle du système observateur-commandé

sys_obs = tf(num, den); % Vous pouvez recalculer le numérateur/dénominateur pour ce système

figure;

step(sys_obs);

title('Réponse indicielle du système observateur-commandé');

xlabel('Temps (s)');

ylabel('Amplitude');

% === 7. Dépassement et amortissement ===

% Calcul des coefficients d'amortissement et pulsation naturelle

D_percent = 10; % Dépassement maximal en %

zeta = -log(D_percent/100) / sqrt(pi^2 + (log(D_percent/100))^2); % Coefficient d'amortissement

wn = sqrt(4) / zeta; % Pulsation propre du système

disp('Coefficient d''amortissement :');

disp(zeta);

disp('Pulsation propre du système :');

disp(wn);
% Script MATLAB - Atelier 2

% 1. Création du vecteur de temps Tsim

t0 = 0; t1 = 5; t2 = 30; t3 = 31; tmax = 80;

n_sim = 0.1; % Pas d'échantillonnage

Tsim = [t0:n_sim:t1, t1+n_sim:n_sim:t2, t2+n_sim:n_sim:t3, t3+n_sim:n_sim:tmax]';

% 2. Création du vecteur de la consigne d'entrée U_def

U = 10;
tol = 1e-6; % Tolérance pour comparaison

U_def = [U * ones(1, find(abs(Tsim - t1) < tol)), ... % Valeur U jusqu'à t1

0 * ones(1, find(abs(Tsim - t2) < tol) - find(abs(Tsim - t1) < tol)), ... % Défaut d'actionneur entre
t1 et t2

U * ones(1, find(abs(Tsim - t3) < tol) - find(abs(Tsim - t2) < tol)), ... % Valeur U restaurée jusqu'à
t3

0 * ones(1, length(Tsim) - find(abs(Tsim - t3) < tol))]'; % Défaut au-delà de t3

% 3. Création de la matrice de données Data

Data = [Tsim, U_def];

disp('Matrice Data :');

disp(Data);

% 4. Création de la matrice X0 et de la matrice Nsim

X0 = [0; 0; 0; 0];

Nsim = [0:n_sim:tmax]';

% 5. Utilisation de ode45 pour résoudre une équation différentielle

% Exemple d'une équation différentielle x' = A*x + B*u

global A B C L U_ob Data

A = [0 1 0 0; 0 0 1 0; 0 0 0 1; -1 -2 -3 -4]; % Exemple de matrice A

B = [0; 0; 0; 1]; % Exemple de matrice B

C = [1 0 0 0]; % Exemple de matrice C

L = 1; % Exemple de gain L

U_ob = 10; % Exemple d'entrée observée

% Résolution avec ode45

tspan = [0 tmax];

[t, y] = ode45(@(t, x) adfunc(t, x), tspan, X0);

% 6. Affichage des résultats

figure;
subplot(2, 1, 1);

plot(Tsim, U_def, 'r', 'LineWidth', 1.5);

xlabel('Temps (s)');

ylabel('Entrée U\_def');

title('Consigne d''entrée');

grid on;

subplot(2, 1, 2);

plot(t, y, 'b', 'LineWidth', 1.5);

xlabel('Temps (s)');

ylabel('État du système');

title('Réponse du système');

grid on;

% Fonction utilisée par ode45

function dx = adfunc(t, x)

global A B U_ob Data

u = interp1(Data(:, 1), Data(:, 2), t); % Interpolation de l'entrée à l'instant t

dx = A * x + B * u; % Équation différentielle

end
% Initialisation des paramètres

t0 = 0; t1 = 5; t2 = 30; t3 = 31; tmax = 80; n_sim = 0.1;

Tsim = [t0; t1; t1 + n_sim; t2; t2 + n_sim; t3; t3 + n_sim; tmax]; % Colonne

u = 10; fd = 0; u_def = u + fd; u_obs = 10;

U_def = [u; u; u; u; u_def; u_def; u_def; u_def]; % Colonne

Data = [Tsim, U_def];

x0 = [0; 0; 0; 0]; % Conditions initiales

Nsim = 0:n_sim:tmax; % Plage de simulation

global A B C L u_obs Data

% Définition des matrices système

B = [1; 0];

C = [0, 1];

L = [3; 0];

% Résolution du système différentiel

[T, X] = ode45(@Sys_mode, Nsim, x0);

% Extraction des variables d'état

x1 = X(:, 1);

x2 = X(:, 2);

x1_obs = X(:, 3);

x2_obs = X(:, 4);

% Calcul des résidus

r1 = x1 - x1_obs;

r2 = x2 - x2_obs;

% Affichage des résultats

figure;

plot(T, x1, 'b-', T, x1_obs, 'g--');

xlabel('Temps (s)');

title('Evolution de x1 et x1_{obs}');

legend('x1', 'x1_{obs}');

figure;

plot(T, x2, 'r-', T, x2_obs, 'y--');

xlabel('Temps (s)');

title('Evolution de x2 et x2_{obs}');

legend('x2', 'x2_{obs}');
figure;

plot(T, r1, 'm-', T, r2, 'k--');

xlabel('Temps (s)');

title('Evolution des residus r1 et r2');

legend('r1', 'r2');

% Simulation des sorties

Cdef = 0.2;

fd1 = 0;

C_def = [Cdef, Cdef + fd1];

y_sys = zeros(size(T)); % Préallocation

y_obs = zeros(size(T)); % Préallocation

for i = 1:length(T)

if T(i) < max(T)/2

y_sys(i) = C * [X(i, 1); X(i, 2)];

else

y_sys(i) = C_def * [X(i, 1); X(i, 2)];

end

y_obs(i) = C * [X(i, 3); X(i, 4)];

end

% Affichage des sorties et des résidus

figure;

plot(T, y_sys, 'b-', T, y_obs, 'g--');

title('y_{sys} et y_{obs}');

xlabel('Temps (s)');

legend('y_{sys}', 'y_{obs}');

rs = y_sys - y_obs;

figure;
plot(T, rs, 'b-');

title('Residus rs');

xlabel('Temps (s)');

% Fonction de définition du système

function dx = Sys_mode(t, x)

global A B C L u_obs Data

% Lecture des données de commande

u = interp1(Data(:,1), Data(:,2), t, 'linear', 'extrap');

% Etat observateur

x_obs = x(3:4);

x_real = x(1:2);

% Modèle système avec observateur

dx_real = A * x_real + B * u;

dx_obs = A * x_obs + B * u + L * (C * x_real - C * x_obs);

% Concaténation

dx = [dx_real; dx_obs];

end