h8.m

来自「h2 and h-infinite control code」· M 代码 · 共 61 行

clear;
clc;

A=[0 1;0 -1];
B=[0 1]';
C1=[1 0];
D1=[0 1];
C2=[0 1];
D2=0;
E=[0 0;1 0];
ep=0.01;
gama=inf;
D1tutor=[D1 0 0 ep];
Etutor=[E ep*eye(2) zeros(2,1)];
C2tutor=[C2;ep*eye(2);zeros(1,2)];
D2tutor=[D2;zeros(2,1);ep];

P = h8care(A,B,C2tutor,D2tutor,E,gama);
EigP=eig(P);
Q = h8care(A',C1',E',D1',C2tutor',gama);
EigQ=eig(Q);

if norm(P*Q)<gama^2
    disp('gama=');
    disp(gama);
    disp('ep=');
    disp(ep);
    disp('Eig(P)=');
    disp(EigP);
    disp('Eig(Q)=');
    disp(EigQ);
    
    F=-inv(D2tutor'*D2tutor)*(D2tutor'*C2tutor+B'*P);
    K=-(Q*C1'+E*D1')*inv(D1*D1');

    Ac=A+E*E'*P/(gama*gama)+B*F+inv(eye(2)-Q*P/(gama*gama))*K*(C1+D1*E'*P/(gama*gama));
    Bc=-inv(eye(2)-Q*P/(gama*gama))*K;
    Cc=F;
    Dc=0;
    
    Acl=[A+B*Dc*C1 B*Cc;Bc*C1 Ac];
    Bcl=[E+B*Dc*D1;Bc*D1];
    Ccl=[C2+D2*Dc*C1 D2*Cc];
    Dcl=D2*Dc*D1;
    cllsys=ss(Acl,Bcl,Ccl,Dcl);
    w=-2:0.1:2;
    for i=1:length(w)
        w(i)=10^w(i);
    end
    for i=1:length(w)
        sigma(i)=norm(Ccl*inv(w(i)*(0+1i)*eye(4)-Acl)*Bcl);
    end
    plot(w,sigma);
    xlabel('Frequency (rad/sec)');
    ylabel('Singular value of Hinf closed-loop system');
    grid;
else
    disp('Unfeasible gama');
end