基于混合灵敏度的H∞鲁棒控制器的设计仿真

clc

close all

K=31.31

b0=50.88

b1=57.81

Ag=[0 1;-b0 -b1]

Bg=[0;1]

Cg=[K 0]

Dg=[0];

[num,den]=ss2tf(Ag,Bg,Cg,Dg) Gs=tf(num,den)

[z,p,k]=tf2zp(num,den)

Gs1=zpk(z,p,k)

%W1S=tf([0.2 3],[1 0.001])

%W2S=tf([0.002],[1])

%W3S=tf([0.002 0.00001],[1]) %W1S=tf([0.2 4],[1 0.002])

%W2S=tf([0.003],[1])

%W3S=tf([0.001 0.00001],[1]) %加权函数

W1S=tf([0.2 10],[1 0.001])

W2S=tf([0.001],[1])

W3S=tf([0.001 0.00001],[1]) W3NUM=W3S.num{1,1};

W3DEN=W3S.den{1,1};

i=0;

j=0;

a=length(W3NUM);

b=length(W3DEN);

for c=1:1:a

if W3NUM(c-i)~=0

break;

else

W3NUM(c-i)=[];

i=i+1;

end

end

for d=1:1:b

if W3DEN(d-j)~=0

break;

else

W3DEN(d-j)=[];

j=j+1;

end

end

Gg=ss(Ag,Bg,Cg,Dg);

[Aw1,Bw1,Cw1,Dw1]=tf2ss(W1S.num{1,1},W1S.den{1,1});

[Aw2,Bw2,Cw2,Dw2]=tf2ss(W2S.num{1,1},W2S.den{1,1});

if Polyorder(W3NUM)>Polyorder(W3DEN)

[Q,R]=deconv(W3NUM,W3DEN)

W3poly=Q

[Aw3,Bw3,Cw3,Dw3]=tf2ss(R,W3DEN)

else

[Aw3,Bw3,Cw3,Dw3]=tf2ss(W3NUM,W3DEN)

W3poly=[]

end

W1=[Aw1,Bw1;Cw1,Dw1]

W2=[Aw2,Bw2;Cw2,Dw2]

W3=[Aw3,Bw3;Cw3,Dw3]

[A,B1,B2,C1,C2,D11,D12,D21,D22]=augss(Ag,Bg,Cg,Dg,Aw1,Bw1,Cw1,Dw1,Aw2,Bw2,Cw2,Dw2,A w3,Bw3,Cw3,Dw3,W3poly)

G=[A,B1,B2;C1,D11,D12;C2,D21,D22];

[Acp,Bcp,Ccp,Dcp,Ac1,Bc1,Cc1,Dc1]=hinf(A,B1,B2,C1,C2,D11,D12,D21,D22)

[knum,kden]=ss2tf(Acp,Bcp,Ccp,Dcp) %控制器的状态空间形式转换成传递函数形式

K=tf(knum,kden) %控制器的传递函数

[z,p,k]=tf2zp(knum,kden) %控制器的传递函数转换成零极点增益形式

Kzpk=zpk(z,p,k)

%{

[gnum,gden]=ss2tf(Ag,Bg,Cg,Dg)

Gs=tf(gnum,gden) %标称系统的传递函数

P=nd2sys(gnum,gden,1) %标称系统的系统矩阵形式

K1=nd2sys(knum,kden,1) %鲁棒控制器的系统矩阵形式

[type,out,in,n]=minfo(P)

I=eye(out)

S=minv(madd(I,mmult(P,K1))) %灵敏度函数(1+GK)^-1 ，其中P和K1都是系统矩阵的形式T=msub(I,S)

W1S1=nd2sys(W1S.num{1,1},W1S.den{1,1},1)

INVW1S1=minv(W1S1) %加权函数W1的逆W1^-1

INVW3S1=nd2sys(W3S.den{1,1},W3S.num{1,1},1) %加权函数W3的逆W3^-1

w=logspace(-5,5,500)

Sw=vsvd(frsp(S,w)) %灵敏度函数的频域响应的奇异值

Tw=vsvd(frsp(T,w)) %补灵敏度函数的频域响应的奇异值

INVW1S1w=vsvd(frsp(INVW1S1,w))

INVW3S1w=vsvd(frsp(INVW3S1,w))

figure(1)

vplot('liv,lm',Sw,'r-',INVW1S1w,'b--')

title('Singular values of sensivity function S and W1^{-1}')

set(gca,'color','w')

xlabel('Frequency(rad/sec)')

ylabel('Amplitude')

grid on

figure(2)

vplot('liv,lm',Tw,'r-',INVW3S1w,'b--')

title('Singular values of complementary sensivity function T and W3^{-1}')

set(gca,'color','w')

xlabel('Frequency(rad/sec)')

ylabel('Amplitude')

grid on

%}

w=logspace(-5,5,500);

figure(3)

bode(W1S,'-b',W3S,'-r',w)

legend('W1','W3')

title('Bode Diagram Of Weighed Function W1 and W3')

grid on

[Acg,Bcg,Ccg,Dcg]=series(Acp,Bcp,Ccp,Dcp,Ag,Bg,Cg,Dg)

[As,Bs,Cs,Ds]=feedbk(Acg,Bcg,Ccg,Dcg,1) %灵敏度函数的状态空间(1+GK)^-1 [At,Bt,Ct,Dt]=feedbk(Acg,Bcg,Ccg,Dcg,2) %补灵敏度函数的状态空间GK(1+GK)^-1

svs=sigma(As,Bs,Cs,Ds,w);

svs=20*log10(svs);

[Aw1i,Bw1i,Cw1i,Dw1i]=unpck(minv(pck(Aw1,Bw1,Cw1,Dw1)));

svw1i=sigma(Aw1i,Bw1i,Cw1i,Dw1i,w);

svw1i= 20*log10(svw1i);

figure(4)

semilogx(w,svw1i,'b--',w,svs,'r-')

title('Singular values of sensivity function S and W1^{-1}')

xlabel('Frequency(rad/sec)')

ylabel('Amplitude(db)')

legend('W1^{-1}','S');

grid on

[Aw3i,Bw3i,Cw3i,Dw3i]=tf2ss(W3S.den{1,1},W3S.num{1,1})

svw3i=sigma(Aw3i,Bw3i,Cw3i,Dw3i,w);

svw3i=20*log10(svw3i);

svt=sigma(At,Bt,Ct,Dt,w);

svt=20*log10(svt);

figure(5)

semilogx(w,svw3i,'b--',w,svt,'r-')

title('Singular values of complementary sensivity function T and W3^{-1}')

xlabel('Frequency(rad/sec)')