状態FB用LMI(H∞制御)

状態FB用LMI(H∞制御)…Homework

[1] n次系\dot{x}=Ax+Bu,\ y=Cx+DuH_\inftyノルムが\gammaより小となるためのLMI条件は、次の通りでした。

\displaystyle{(15)\quad \begin{array}{lll} && \sup_{\omega\in{\rm\bf R}}||\hat{G}(j\omega)||_2<\gamma \nonumber\\ &\Leftrightarrow& \exists X>0:\ \left[\begin{array}{ccc} A^TX+XA & XB & C^T \\ B^TX & -\gamma I & D^T \\ C & D & -\gamma I \end{array}\right]<0 \nonumber\\ &\Leftrightarrow& \exists Y>0:\ \left[\begin{array}{ccc} YA^T+AY & B & YC^T \\ B^T & -\gamma I & D^T \\ CY & D & -\gamma I \end{array}\right]<0\nonumber \end{array} }

ただし

\displaystyle{(16)\quad \hat{G}(s)=C(sI-A)^{-1}B+D }


図1 2ポートシステム

いま、次の2ポート表現されたn次系を考えます。

\displaystyle{(17)\quad P: \left\{\begin{array}{l} \dot{x}=Ax+B_1u_1+B_2u_2 \\ \underbrace{ \left[\begin{array}{c} y\\ u \end{array}\right] }_{y_1} = \underbrace{ \left[\begin{array}{c} C\\ 0 \end{array}\right] }_{C_1} x+ \underbrace{ \left[\begin{array}{c} 0\\ 0 \end{array}\right] }_{D_{11}}u_1+ \underbrace{ \left[\begin{array}{c} 0\\ I \end{array}\right] }_{D_{12}} u_2 \\ y_2=x \end{array}\right. }

これに対する状態フィードバック

\displaystyle{(18)\quad K: u_2=-Fy_2 }

による閉ループ系

\displaystyle{()\quad P_{CL}: \left\{\begin{array}{l} \dot{x}=(A-B_2F)x+B_1u_1 \\ y_1= \underbrace{ \left[\begin{array}{c} C\\ -F \end{array}\right] }_{C_1-D_{12}F} x \end{array}\right. }

において、そのH_\inftyノルムが\gammaより小、すなわち

\displaystyle{(19)\quad \sup_{\omega\in{\rm\bf R}} \underbrace{||\hat{P}_{CL}(j\omega)||_2}_{\sigma_1(\hat{P}_{CL}(j\omega))} <\gamma }

ただし

\displaystyle{(20)\quad \hat{P}_{CL}(s)= \left[\begin{array}{c} \hat{G}_y(s)\\ \hat{G}_u(s) \end{array}\right]= \left[\begin{array}{c} C(sI-A-B_2F)^{-1}B_1\\ -F(sI-A-B_2F)^{-1}B_1 \end{array}\right] }

となるように状態フィードバックゲインFを求める問題を考えます。

この問題のLMI条件は、次のようになります。

\displaystyle{(21)\quad \exists Y>0:\ }
\displaystyle{ \left[\begin{array}{ccc} Y(A-B_2F)^T+(A-B_2F)Y & B_1 & Y(C_1-D_{12}F)^T \\ B_1^T & -\gamma I & D_{11}^T \\ (C_1-D_{12}F)Y & D_{11} & -\gamma I \end{array}\right]<0 }

ここで、変数変換Z=FYを行うと、次のようなYZに関するLMIとなります。

\displaystyle{(22)\quad \exists Y>0:\ }
\displaystyle{ \left[\begin{array}{ccc} (AY-B_2Z)^T+AY-B_2Z & B_1 & (C_1Y-D_{12}Z)^T \\ B_1^T & -\gamma I & D_{11}^T \\ C_1Y-D_{12}Z & D_{11} & -\gamma I \end{array}\right]<0 }

これを解いてYZを求め、次式によって状態フィードバックゲインを決定します。

\displaystyle{(23)\quad F=ZY^{-1} }

演習B33…Flipped Classroom
1^\circ 次のコードを参考にして、H_\infty制御(状態FB)を求める関数を作成せよ。

MATLAB
%sf_syn_lmi5.m 
%----- 
 A=[0 1;-1 -2*0.01]; B1=[0;1]; B2=B1;
 C1=[1 0;0 0]; D11=[0;0]; D12=[0;1]; 
 [n,m]=size(B2);  
%----- 
 setlmis([]);
 gam=lmivar(1,[1 0]); 
 Y=lmivar(1,[n 1]); 
 Z=lmivar(2,[m n]); 
%----- 
 lmi1=newlmi;
 lmiterm([lmi1,1,1,Y],A,1,'s');      %#1:A*Y+Y*A' 
 lmiterm([lmi1,1,1,Z],-B2,1,'s');    %#1:-(B2*Z+Z*B2') 
 lmiterm([lmi1,1,2,0],B1);           %#1:B1 
 lmiterm([lmi1,2,2,gam],-1,1);       %#1:-gam 
 lmiterm([lmi1,3,1,Y],C1,1);         %#1:C1*Y 
 lmiterm([lmi1,3,1,Z],-D12,1);       %#1:D12*Z 
 lmiterm([lmi1,3,2,0],D11);          %#1:D11 
 lmiterm([lmi1,3,3,gam],-1,1);       %#1:-gam 
%----- 
 lmi2=newlmi;
 lmiterm([-lmi2,1,1,Y],1,1);         %#2:Y   
%----- 
 lmi3=newlmi;
 lmiterm([lmi3,1,1,gam],1,1);        %#3:gam 
 lmiterm([-lmi3,1,1,0],1e3);         %#3:1000    
%----- 
 LMIs=getlmis; 
 cobj=zeros(1,decnbr(LMIs)); 
 cobj(1)=1; 
 [cost,xopt]=mincx(LMIs,cobj); 
 Y=dec2mat(LMIs,xopt,Y);  
 Z=dec2mat(LMIs,xopt,Z);  
 F=Z/Y;  
%----- 
 pl=eig(A-B2*F)
 close all,figure(1) 
 dregion(0,0,20,pi/4,20*[-1,1,-1,1]) 
 plot(real(pl),imag(pl),'*') 
 figure(2) 
 ol=ltisys(A,B1,C1,D11); 
 cl=ltisys(A-B2*F,B1,C1-D12*F,D11); 
 om=logspace(-2,2,200); 
 splot(ssub(cl,1,1),'sv',om),hold on 
 splot(ssub(ol,1,1),'sv',om),grid 
%----- 
%eof

2^\circ 上の2次振動系のH_\infty制御について、閉ループ系の固有値に制約を付けた場合の効果を確かめよ。

MATLAB

%sf_syn_lmi6.m 
%----- 
 clear all, close all
 A=[0 1;-1 -2*0.01]; B1=[0;1]; B2=B1;
 C1=[1 0;0 0]; D11=[0;0]; D12=[0;1]; 
 [n,m]=size(B2);  
%----- 
 setlmis([]);
 gam=lmivar(1,[1 0]); 
 Y=lmivar(1,[n 1]); 
 Z=lmivar(2,[m n]); 
%----- 
 lmi1=newlmi;
 lmiterm([lmi1,1,1,Y],A,1,'s');      %#1:A*Y+Y*A' 
 lmiterm([lmi1,1,1,Z],-B2,1,'s');    %#1:-(B2*Z+Z*B2') 
 lmiterm([lmi1,1,2,0],B1);           %#1:B1 
 lmiterm([lmi1,2,2,gam],-1,1);       %#1:-gam 
 lmiterm([lmi1,3,1,Y],C1,1);         %#1:C1*Y 
 lmiterm([lmi1,3,1,Z],-D12,1);       %#1:D12*Z 
 lmiterm([lmi1,3,2,0],D11);          %#1:D11 
 lmiterm([lmi1,3,3,gam],-1,1);       %#1:-gam 
%----- 
 alpha=0.1;
 lmi2=newlmi;  
 lmiterm([lmi2,1,1,Y],A,1,'s');      %#2:A*Y+Y*A'   
 lmiterm([lmi2,1,1,Z],-B2,1,'s');    %#2:-(B2*Z+Z'*B2')   
 lmiterm([lmi2,1,1,Y],2*alpha,1);    %#2:2*alpha*Y   
%----- 
 r=2;
 lmi3=newlmi;  
 lmiterm([lmi3,1,1,Y],-r,1);         %#3:-r*Y  
 lmiterm([lmi3,1,2,Y],A,1);          %#3:A*Y   
 lmiterm([lmi3,1,2,Z],-B2,1);        %#3:-B2*Z  
 lmiterm([lmi3,2,2,Y],-r,1);         %#3:-r*Y  
%----- 
 theta=pi/4; sth=sin(theta); cth=cos(theta);
 lmi4=newlmi;  
 lmiterm([lmi4,1,1,Y],sth*A,1,'s');  %#4:sth*(A*Y+Y*A')  
 lmiterm([lmi4,1,1,Z],-sth*B2,1,'s');%#4:-sth*(B2*Z+Z'*B2')  
 lmiterm([lmi4,1,2,Y],cth*A,1);      %#4:cth*A*Y  
 lmiterm([lmi4,1,2,Y],1,-cth*A');    %#4:-cth*Y*A'  
 lmiterm([lmi4,1,2,Z],-cth*B2,1);    %#4:-cth*B2*Z  
 lmiterm([lmi4,1,2,-Z],1,cth*B2');   %#4:cth*Z'*B2'  
 lmiterm([lmi4,2,2,Y],sth*A,1,'s');  %#4:sth*(A*Y+Y*A')  
 lmiterm([lmi4,2,2,Z],-sth*B2,1,'s');%#4:-sth*(B2*Z+Z'*B2')  
%----- 
 lmi5=newlmi;
 lmiterm([-lmi5,1,1,Y],1,1);         %#5:Y   
%----- 
 lmi6=newlmi;
 lmiterm([lmi6,1,1,gam],1,1);        %#6:gam 
 lmiterm([-lmi6,1,1,0],1e3);         %#6:1000    
%----- 
 LMIs=getlmis; 
 cobj=zeros(1,decnbr(LMIs)); 
 cobj(1)=1; 
 [cost,xopt]=mincx(LMIs,cobj); 
 gopt=dec2mat(LMIs,xopt,gam);  
 Y=dec2mat(LMIs,xopt,Y);  
 Z=dec2mat(LMIs,xopt,Z);  
 F=Z/Y;  
%----- 
 pl=eig(A-B2*F)
 close all,figure(1) 
 dregion(alpha,0,r,theta,r*[-1,1,-1,1]) 
 plot(real(pl),imag(pl),'*') 
 figure(2) 
 ol=ltisys(A,B1,C1,D11); 
 cl=ltisys(A-B2*F,B1,C1-D12*F,D11); 
 om=logspace(-2,2,200); 
 splot(ssub(cl,1,1),'sv',om),hold on 
 splot(ssub(ol,1,1),'sv',om),grid 
%----- 
%eof