Bonjour tous le monde j’étais entrain de faire une estimation de paramètres par la méthode de Marquardt d'un procédé pilote, le modèle linaire est données par:
d2y/dt2=-teta1*dy/dt-teta2*y+teta3*u avec u=0.1 amplitude échelon j'ai estimé les paramètres les programmes sont les suivants
1)
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
clear all 
close all
clc
% méthode de marquardt
load donneesbo.dat
global u tspan h3
tspan=(donneesbo(:,1));
h3=donneesbo(:,3)-donneesbo(1,3);
u=0.1;
teta0=[1 1 1]'/1000;
lamda=1
esp=0.001
%calcul de J0,G0 et H0
[J0,H0,G0]=feval('tp',teta0);
 
while(1)
 
    teta=teta0-inv(H0+lamda*eye(3))*G0;%eye(3)=identité
    [J,H,G]=feval('tp',teta);
    if(abs(J-J0)<eps)
        break;
    end
    if(J>=J0)
        lamda=10*lamda;
    else
        lamda=lamda/10;
        teta0=teta;
        H0=H;
        J0=J;
        G0=G;
    end
    i=i+1;
end
2)
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
function [J,H,G]=tp(teta)
global  teta1 teta2 teta3 tspan h3
%boucle de simulation 
teta1=teta(1);
teta2=teta(2);
teta3=teta(3);
x0=zeros(8,1);
[t,x]=ode23s('modele',tspan,x0);
plot(t,x(:,1),t,h3)
xlabel('temps en (s)');
% Create ylabel
ylabel('sortie en (cm)');
% Create title
title('comparaison modèle et mesures');
legend('mesure','modele')
%plot(h3)
pause(1)
%calcul de J0,H0,G0
M=[x(:,1)-h3 x(:,3) x(:,4) x(:,5)];
P=M'*M;
J=P(1,1);%l'erreur 
G=P(2:4,1);%gradiant
H=P(2:4,2:4);%Hessien
end
3)
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
function [dx] = modele(t,x)
global teta1 teta2 teta3 u
dx1=x(2);%y=x1= hauteur h3
dx2=-teta1*x(2)-teta2*x(1)+teta3*u;%u=0,1
dx3=x(6);
dx4=x(7);
dx5=x(8);
dx6=-x(2)-teta1*x(6)-teta2*x(3);
dx7=-teta1*x(7)-x(1)-teta2*x(4);
dx8=u-teta1*x(8)-teta2*x(5);
dx=[dx1;dx2;dx3;dx4;dx5;dx6;dx7;dx8];
end
ce modele etant validé je voulais identifier les parametres du regulateur PID c'est a dire que je vais remplacer u par : u=Kp*epsilon+KI*epsilon*dt+Kd*depsilon/dt avec Kp=teta1 KI=teta2 Kd=teta3 et espilon=consigne -X1
mais je ne sais pas par ou commencer est ce que quelqu'un pourrait m'expliquer la procedure
Merci d'avance