反馈控制器的程序实现问题

zyj_030 · 发表于 2009-9-9 15:01

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

您需要登录才可以下载或查看，没有账号？我要加入

x

对于矩阵不等式：

采用锥补线性化方法等价于求解下述的优化问题：

求该优化问题的算法：

zyj_030 · 发表于 2009-9-9 15:10

对于

在论文中，对于h2=0.1，经过212次迭代后，获得可使系统稳定的最大时滞界h1=0.65，相应的控制器为

K=Y*L-1 =[-0.3543 -0.2372
0.0515
-0.2752]

我编的程序如下，但得不到论文的结果，可能是程序有误，不知错在何处，请各位大侠帮忙

clc;

clear all;

close all;

A0=[0 0 0 0;0 0.5 0 0;-0.5 0 0.3 0;0 0 0 1];

A1=[-2 -0.5 0 0;-0.2 -1 0 0; 0.5 0 -2 -0.5;0 0 0 -1];

B1=[1 1 1 0]';

B2=[0 1 1 1]';

U=zeros(4);

h2=0.1; %固定h2

h1=0.05; %初始化h1

k=0;

a=1;

setlmis([]);

L=lmivar(1,[4 1]);

Q1=lmivar(1,[4 1]);

Q2=lmivar(1,[4 1]);

S1=lmivar(1,[4 1]);

S2=lmivar(1,[4 1]);

R1=lmivar(1,[4 1]);

R2=lmivar(1,[4 1]);

T1=lmivar(1,[4 1]);

T2=lmivar(1,[4 1]);

J=lmivar(1,[4 1]);

Z1=lmivar(1,[4 1]);

Z2=lmivar(1,[4 1]);

Y=lmivar(2,[1 4]);

N11=lmivar(2,[4 4]);

N12=lmivar(2,[4 4]);

N21=lmivar(2,[4 4]);

N22=lmivar(2,[4 4]);

lmiterm([-1 1 1 L],1,1);

lmiterm([-2 1 1 Q1],1,1);

lmiterm([-3 1 1 Q2],1,1);

lmiterm([-4 1 1 S1],1,1);

lmiterm([-5 1 1 S2],1,1);

lmiterm([-6 1 1 R1],1,1);

lmiterm([-7 1 1 R2],1,1);

lmiterm([-8 1 1 T1],1,1);

lmiterm([-9 1 1 T2],1,1);

lmiterm([-10 1 1 J],1,1);

lmiterm([-11 1 1 Z1],1,1);

lmiterm([-12 1 1 Z2],1,1);

lmiterm([13 1 1 L],A0,1,'s');

lmiterm([13 1 1 Y],B1,1,'s');

lmiterm([13 1 1 N11],1,1,'s');

lmiterm([13 1 1 N12],1,1,'s');

lmiterm([13 1 1 Q1],1,1);

lmiterm([13 1 1 Q2],1,1);

lmiterm([13 1 2 L],A1,1);

lmiterm([13 1 2 -N11],-1,1);

lmiterm([13 1 2 N21],1,1);

lmiterm([13 1 3 Y],B2,1);

lmiterm([13 1 3 -N12],-1,1);

lmiterm([13 1 3 N22],1,1);

lmiterm([13 1 4 -L],h1,A0');

lmiterm([13 1 4 -Y],h1,B1');

lmiterm([13 1 5 -L],h2,A0');

lmiterm([13 1 5 -Y],h2,B1');

lmiterm([13 1 6 -N11],h1,1);

lmiterm([13 1 7 -N12],h2,1);

lmiterm([13 2 2 Q1],-1,1);

lmiterm([13 2 2 N21],-1,1,'s');

lmiterm([13 2 4 -L],h1,A1');

lmiterm([13 2 5 -L],h2,A1');

lmiterm([13 2 6 -N21],h1,1);

lmiterm([13 3 3 Q2],-1,1);

lmiterm([13 3 3 N22],-1,1,'s');

lmiterm([13 3 4 -Y],h1,B2');

lmiterm([13 3 5 -Y],h2,B2');

lmiterm([13 3 7 -N22],h2,1);

lmiterm([13 4 4 R1],-h1,1);

lmiterm([13 5 5 R2],-h2,1);

lmiterm([13 6 6 S1],-h1,1);

lmiterm([13 7 7 S2],-h2,1);

lmiterm([-14 1 1 T1],1,1);

lmiterm([-14 1 2 J],1,1);

lmiterm([-14 2 2 Z1],1,1);

lmiterm([-15 1 1 T2],1,1);

lmiterm([-15 1 2 J],1,1);

lmiterm([-15 2 2 Z2],1,1);

lmiterm([-16 1 1 S1],1,1);

lmiterm([-16 1 2 0],1);

lmiterm([-16 2 2 T1],1,1);

lmiterm([-17 1 1 S2],1,1);

lmiterm([-17 1 2 0],1);

lmiterm([-17 2 2 T2],1,1);

lmiterm([-18 1 1 L],1,1);

lmiterm([-18 1 2 0],1);

lmiterm([-18 2 2 J],1,1);

lmiterm([-19 1 1 R1],1,1);

lmiterm([-19 1 2 0],1);

lmiterm([-19 2 2 Z1],1,1);

lmiterm([-20 1 1 R2],1,1);

lmiterm([-20 1 2 0],1);

lmiterm([-20 2 2 Z2],1,1);

lmis=getlmis;

[t,xopt]=feasp(lmis);

L0=dec2mat(lmis,xopt,L);

Q10=dec2mat(lmis,xopt,Q1);

Q20=dec2mat(lmis,xopt,Q2);

S10=dec2mat(lmis,xopt,S1);

S20=dec2mat(lmis,xopt,S2);

R10=dec2mat(lmis,xopt,R1);

R20=dec2mat(lmis,xopt,R2);

T10=dec2mat(lmis,xopt,T1);

T20=dec2mat(lmis,xopt,T2);

J0=dec2mat(lmis,xopt,J);

Z10=dec2mat(lmis,xopt,Z1);

Z20=dec2mat(lmis,xopt,Z2);

Y0=dec2mat(lmis,xopt,Y);

N110=dec2mat(lmis,xopt,N11);

N120=dec2mat(lmis,xopt,N12);

N210=dec2mat(lmis,xopt,N21);

N220=dec2mat(lmis,xopt,N22);

zyj_030 · 发表于 2009-9-9 15:16

while a>0.001,
%循环条件

n=decnbr(lmis);

c=zeros(n,1);

for j=1:n,

[Lj,Q1j,Q2j,S1j,S2j,R1j,R2j,T1j,T2j,Jj,Z1j,Z2j,Yj,N11j,N12j,N21j,N22j]=defcx(lmis,j,L,Q1,Q2,S1,S2,R1,R2,T1,T2,J,Z1,Z2,Y,N11,N12,N21,N22);

c(j)=trace c(j)=trace(L0*Jj+J0*Lj+S10*T1j+T10*S1j+S20*T2j+T20*S2j+R10*Z1j+Z10*R1j+R20*Z2j+Z20*R2j);

end

options=[1e-5,0,0,0,0];

[copt,xopt]=mincx(lmis,c,options);

L1=dec2mat(lmis,xopt,L);

Q11=dec2mat(lmis,xopt,Q1);

Q21=dec2mat(lmis,xopt,Q2);

S11=dec2mat(lmis,xopt,S1);

S21=dec2mat(lmis,xopt,S2);

R11=dec2mat(lmis,xopt,R1);

R21=dec2mat(lmis,xopt,R2);

T11=dec2mat(lmis,xopt,T1);

T21=dec2mat(lmis,xopt,T2);

J1=dec2mat(lmis,xopt,J);

Z11=dec2mat(lmis,xopt,Z1);

Z21=dec2mat(lmis,xopt,Z2);

Y1=dec2mat(lmis,xopt,Y);

N111=dec2mat(lmis,xopt,N11);

N121=dec2mat(lmis,xopt,N12);

N211=dec2mat(lmis,xopt,N21);

N221=dec2mat(lmis,xopt,N22);

Z=[A0*L1+L1*A0'+B1*Y1+Y1'*B1'+N111'+N111+N121'+N121+Q11+Q21
A1*L1-N111'+N211
B2*Y1-N121'+N221
h1*(L1'*A0'+Y1'*B1')
h2*(L1'*A0'+Y1'*B1')
h1*N111'
h2*N121';

(A1*L1-N111'+N211)'
-Q11-N211'-N211
U
h1*L1'*A1'
h2*L1'*A1'
h1*N211'
U
;

(B2*Y1-N121'+N221)'
U
-Q21-N221'-N221
h1*Y1'*B2'
h2*Y1'*B2'
U
h2*N221';

h1*(L1'*A0'+Y1'*B1')'
h1*(L1'*A1')'
h1*(Y1'*B2')'
-h1*R11
U
U
U
;

h2*(L1'*A0'+Y1'*B1')'
h2*(L1'*A1')'
h2*(Y1'*B2')'
U
-h2*R21
U
U
;

h1*N111
h1*N211
U
U
U
-h1*L1*inv(R11)*L1
U
;

h2*N121
U
h2*N221
U
U
U
-h2*L1*inv(R21)*L1];
%Z就是公式（4.25）的左端

M=eig(Z);% 计算Z的特征值，进而判断（4.25）是否满足

m=0;

for i=1:length(M),

if M(i)<0,

m=m+1;

end

if m==length(M), %（4.25）满足

h1=h1+0.03;
%适当增加h1，返回到step2