线性系统结构分析与分解及标准型

1、广西大学实验报告纸实验内容：线性系统结构分析与分解及标准型【实验时间】2014年5月28日【实验地点】宿舍(课外)【实验目的】掌握线性系统状态空间标准型、解及其模型转换。【实验设备与软件】1、 MATLAB数值分析软件【实验原理】1、 标准型变换、矩阵Jordan型变换阵、特征值。其格式分别为csys=canon（sys，type）；V J=Jordan（A）；V J=eig（A）、cv=eig（A）。2、 状态模型的相似变换：命令格式 sysb=ss2ss（sys，T）、传递函数模型与状态空间模型之间的相互转换： 命令格式 A,B,C,D=tf2ss（num，den） num,den=ss2

2、tf（A,B,C,D,iu）、zpk模型与状态空间模型之间的互相转换： 命令格式 A,B,C,D=zp2ss（z,p,k） z,p,k=ss2zp（A,B,C,D,iu）3、 线性定常系统的可控制性与可观性及结构分解（1） 可控性和可观性一般采用能控性或能观性矩阵类别状态可观性判别子函数代码如下：Function str=pdobsv(A,C) Qo=obsv(A,C); r=rank(Qo); 1=size(A,1); if r=1 Str=系统是状态完全可观的！ else str=系统不是状态完全可观的！ end（2） 可控性和可观性Gram矩阵可由下面函数求得W=gram（sys，typ

3、e）-sys是系统的状态空间模型，type可以是c和o。（3） 结构分解a. 在matlab中调用ctrbf（）函数对系统按能控性分解Abar,Bbar,Cbar,T,K=ctrbf(A,B,C)Abar,Bbar,Cbar,T,K=ctrbf(A,B,C,TOL)b. 在matlab中调用ctrbf（）函数对系统按能观性分解Abar,Bbar,Cbar,T,K=obsvf(A,B,C)Abar,Bbar,Cbar,T,K=obsvf(A,B,C,TOL)c. 按能控、能观性分解-Kalman分解4、 定常线性系统的标准型（1） Jordan标准型a. 互异根的情况下，代码可以采用如下的形式，

4、出可以有采用Jordan指令。num=.;den=.;G=tfnum,den;r,p,k=residue(num,den);A=diag(p);B=ones(length(r),1);C=k;Gss=SS(A,B,C,D)b.有重根的情况num=.;den=.;Gtf=tf(num,den)Gs=ss(G)V J=jordan(Gs.a)Gss=ss2ss(Gs,inv(V)（2） 能控标准型-限于SISO系统转换成能控标准型代码：A=.;B=.;C=.;D=.;Gs=ss(A,B,C,D);T=ctrb(Gs.a,Gs.b)Abar=inv(T)*A*T；Bar=inv(T)*B;Cbar=

5、C*T,Dar=D;Gss=ss(Abar,Bbar,Cbar,Dbar)（3） 能观标准型-限于SISO系统转换成能观标准型代码：A=.;B=.;C=.;D=.;Gs=ss(A,B,C,D);Tinv=obsv(Gs.a,Gs.b);T=inv(Tinv);Abar=inv(T)*A*T;Bar=inv(T)*B;Cbar=C*T,Dar=D;Gss=ss(Abar,Bbar,Cbar,Dbar)【实验内容、方法、过程与分析】1. 实验内容1. 已知线性系统（1） 判断其状态可控性、可观性和传递函数的关系，并加以说明分析。程序代码如下：子函数如下：function str=pdctrb(a,

6、b)Qc=ctrb(a,b);r=rank(Qc);l=length(a);if r=l str=系统是状态完全可控的！;else str=系统是状态不完全可控的！;endfunction str=pdobsv(a,c)Qo=obsv(a,c);r=rank(Qo);l=size(a,1);if r=l str=系统是状态完全可观的！ ;else str=系统是状态不完全可观的！;end主函数为：a=-6 -0.625 0.75;8 0 0;0 2 0;b=1 0 0;c=1 -0.25 0.0625;d=1;G=ss(a,b,c,d);Gz=zpk(G)y,t,x=step(G,36);pl

7、ot(t,x,r,t,y,b);grid on;axis(0 36 0 10);title(系统状态和单位阶跃输出曲线);xlabel(t (s);ylabel(y(t)/x(t);str=pdctrb(a,b)str=pdobsv(a,c)运行结果：Zero/pole/gain:(s+6.236) (s+1.764) (s-1)- (s+4) (s+3) (s-1)str =系统是状态完全可控的！str =系统是状态不完全可观的！说明分析：由运行结果表明，传递函数发生了s=1的零极点对消。因此，系统是不完全可观的。（2） 对系统分别按能控性分解、能观性分解以及能控能观性分解。能控性分解： a

8、bar,bbar,cbar,T,K=ctrbf(a,b,c)abar = 0 2.0000 0 0 0 -8.0000 -0.7500 0.6250 -6.0000bbar = 0 0 1cbar = -0.0625 0.2500 1.0000T = 0 0 -1 0 -1 0 1 0 0K = 1 1 1能观性分解： abar,bbar,cbar,T,K=obsvf(a,b,c)abar = 1.0000 2.2772 -2.7134 0.0000 0.3407 -5.5238 -0.0000 2.6251 -7.3407bbar = 0.0558 0.2432 -0.9684cbar =

9、0.0000 0.0000 -1.0327T = 0.0558 0.4465 0.8930 0.2432 0.8614 -0.4459 -0.9684 0.2421 -0.0605K = 1 1 0能控能观性分解：function Gk,T,K=kalmdec(G)G=ss(G);A=G.a;B=G.b;C=G.c;Ac,Bc,Cc,Tc,Kc=ctrbf(A,B,C);nc=rank(ctrb(A,B);n=length(A);ic=n-nc+1:n;Ao1,Bo1,Co1,To1,Ko1=obsv(Ac(ic,ic),Bc(ic),Cc(ic);if ncn inc=1:n-nc; Ao2

10、,Bo2,Co2,To2,Ko2=obsvf(Ac(inc,inc),Bc(inc),Cc(inc);endm1,n1=size(To1);m2,n2=size(To2);To=To2,zeros(m2,n1);zeros(m1,n2),To1;T=To*Tc;n1=rank(obsv(Ac(ic,ic),Cc(ic);n2=rank(obsv(Ac(inc,inc),Cc(inc);K=zeros(1,n-nc-n2),ones(1,n2),.,2*ones(1,nc-n1),3*ones(1,n1);Ak=T*A*inv(T);Bk=T*B;Ck=C*inv(T);Gk=ss(Ak,Bk,

11、Ck,G.d);2. 在MATLAB建立并运行如下的.m代码，回答下面的问题。num=1 2 3;den=conv(1 6 25,1 12 35);G=tf(num,den)Gs=ss(G)V J=jordan(Gs.a)Gss=ss2ss(Gs,inv(V)Gsm=canon(Gs,model)Gsf=canon(Gs,companion)运行结果：Transfer function: s2 + 2 s + 3-s4 + 18 s3 + 132 s2 + 510 s + 875 a = x1 x2 x3 x4 x1 -18 -8.25 -3.984 -3.418 x2 16 0 0 0 x3

12、 0 8 0 0 x4 0 0 2 0b = u1 x1 0.25 x2 0 x3 0 x4 0c = x1 x2 x3 x4 y1 0 0.25 0.0625 0.04688d = u1 y1 0Continuous-time model.V = -0.4883 -1.3398 0.4570 + 0.1719i 0.4570 - 0.1719i 1.5625 3.0625 -0.4375 - 1.5000i -0.4375 + 1.5000i -2.5000 -3.5000 -1.5000 + 2.0000i -1.5000 - 2.0000i 1.0000 1.0000 1.0000 1.

13、0000 J = -5.0000 0 0 0 0 - 7.0000 0 0 0 0 -3.0000 + 4.0000i 0 0 0 0 -3.0000 - 4.0000ia = x1 x2 x3 x4 x1 -5+8.4e-015i -1.07e-014+8.55e-015i 1.9e-014+8.38e-015i 2.19e-014+1.19e-014i x2 -3.55e-015-4.3e-015i -7-4.43e-015i -5.7e-015-6.06e-015i -4.07e-015+1.57e-016i x3 4.97e-016-3.55e-015i 1.81e-015-3.55e

14、-015i -3+4i 4.44e-016-4.22e-015i x4 3.91e-016-1.78e-015i 4.41e-015-1.78e-015i -4.88e-015-3.11e-015i -3-4i b = u1 x1 1.6+2.22e-016i x2 -1-1.67e-016i x3 -0.3+0.1i x4 -0.3-0.1ic = x1 x2 x3 x4 y1 0.281+1.25e-016i 0.594+1.94e-016i -0.156-0.25i -0.156+0.25id = u1 y1 0Continuous-time model.a = x1 x2 x3 x4

15、x1 -3 4 0 0 x2 -4 -3 0 0 x3 0 0 -7 0 x4 0 0 0 -5b = u1 x1 0.9492 x2 -0.1931 x3 4.942 x4 5.042c = x1 x2 x3 x4 y1 0.1699 0.09054 -0.1201 0.08925d = u1 y1 0Continuous-time model.a = x1 x2 x3 x4 x1 0 0 0 -875 x2 1 0 0 -510 x3 0 1 0 -132 x4 0 0 1 -18b = u1 x1 1 x2 0 x3 0 x4 0 c = x1 x2 x3 x4 y1 0 1 -16 1

16、59d = u1 y1 0（2） canon得到的对角型系统为：num=2 1;den=1 7 14 8;G=tf(num,den);Gs=ss(G);Gsm=canon(Gs,model)运行结果：a = x1 x2 x3 x1 -4 0 0 x2 0 -2 0 x3 0 0 -1b = u1 x1 -3.83 x2 -4.899 x3 1.915c = x1 x2 x3 y1 0.3046 -0.3062 -0.1741d = u1 y1 0说明：传递函数中只含有各个相异的实数极点，可以转换为对角型系统。（3） 将原理中给出的能控标准型与能观标准型转换代码写成子函数的形式，并通过调用你所编

17、写的子函数将.m文件中给出的模型变换成能控标准、型和能观标准、型，并从结果说明能控与能观标准型间关系。因为系统可控不能观，故只能能转化成能控型。转换成能控标准型代码：A=-6 -0.625 0.75;8 0 0;0 2 0;B=1 0 0;C=1 -0.25 0.0625;D=1;Gs=ss(A,B,C,D);T=ctrb(Gs.a,Gs.b)Abar=inv(T)*A*T;Bbar=inv(T)*B;Cbar=C*T,Dbar=D;Gss=ss(Abar,Bbar,Cbar,Dbar)运行结果：T = 1 -6 31 0 8 -48 0 0 16Cbar = 1 -8 44a = x1 x2

18、 x3 x1 0 0 12 x2 1 0 -5 x3 0 1 -6b = u1 x1 1 x2 0 x3 0c = x1 x2 x3 y1 1 -8 44d = u1 y1 1转换成能观标准型代码： A=-6 -0.625 0.75;8 0 0;0 2 0;B=1 0 0;C=1 -0.25 0.0625;D=1;Gs=ss(A,B,C,D);Tinv=obsv(Gs.a,Gs.c);T=inv(Tinv);Abar=inv(T)*A*T;Bbar=inv(T)*B;Cbar=C*T,Dbar=D;Gss=ss(Abar,Bbar,Cbar,Dbar)运行结果：Cbar = 0.7500 0 -0.0313a = x1 x2 x3 x1 4 6 0.5 x2