第3章 线性控制系统

1、第3章 线性控制系统的数学模型系统数学模型的重要性 系统仿真分析必须已知数学模型 系统设计必须已知数学模型 本课程数学模型是基础 系统数学模型的获取 建模方法：从已知的物理规律出发，用数学推导的方式建立起系统的数学模型 辨识方法：由实验数据拟合系统的数学模型 本章主要内容线性连续系统的数学模型与MATLABMATLAB表示线性离散时间系统的数学模型方框图描述系统的化简系统模型的相互转换线性系统的模型降阶线性系统的模型辨识本章要点简介 回顾：传递函数 零初始条件下，线性系统响应（即输出）量的拉普拉斯变换与激励（即输入）量的拉普拉斯变换之比。记作G（s）Y（s）U（s），其中Y（s）、U（s)分别

2、为输出量和输入量的拉普拉斯变换。传递函数是描述线性系统动态特性的基本数学工具之一，经典控制理论的主要研究方法，频率响应法和根轨迹法，都是建立在传递函数的基础之上。系统的传递函数与描述其运动规律的微分方程是对应的。可根据组成系统各单元的传递函数和它们之间的联结关系导出整体系统的传递函数，并用它分析系统的动态特性、稳定性，或根据给定要求综合控制系统，设计满意的控制器。以传递函数为工具分析和综合控制系统的方法称为频域法。它不但是经典控制理论的基础，而且在以时域方法为基础的现代控制理论发展过程中，也不断发展形成了多变量频域控制理论，成为研究多变量控制系统的有力工具。传递函数中的复变量s在实部为零、虚部

3、为角频率时就是频率响应。3.1 连续线性系统的数学模型与MATLAB表示-传递函数模型 (1)线性系统的微分方程Laplace变换MATLAB：num=b1, b2, , bm , bm+1; den=a1, a2, , an , an+1;G=tf (num, den);11211121( )( )( )( )( )( )( )( )nnnnnnmmmmmmd y tdy tdy taaaay tdtdtdtd u tdu tdu tbbbbu tdtdtdt11211121( )mmmmnnnnb sb sb sbG sa sa sa sa3.1 连续线性系统的数学模型与MATLAB表示-

4、传递函数模型 (2)例31 传递函数模型num=12 24 12 20; den=2 4 6 2 2; G=tf(num,den);还可以： s=tf(s); G=(12*s3+24*s2+12*s+20)/(2*s4+4*s3+6*s2+2*s+2);3243212241220( )24622sssG sssss3.1 连续线性系统的数学模型与MATLAB表示-传递函数模型 (3)例32 传递函数s=tf(s); G=3*(s2+3)/(s+2)3/(s2+2*s+1)/(s2+5)例33 传递函数 s=tf(s); G=(s3+2*s2+3*s+4)/(s3*(s+2)*(s+5)2+5)

5、应该根据给出传递函数形式选择输入方法 23223(3)( )(2) (21)(5)sG sssss3232234( )(2)(5)5sssG ssss3.1 连续线性系统的数学模型与MATLAB表示-传递函数模型 (4)MATLAB的tf对象还允许携带其他属性，set(tf)命令列出例33-4 延迟传递函数 G.ioDelay=3 或者 set(G,ioDelay,3)假设复域变量为p,则 G.Variable=p 或者 set(G,Variable,p)传递函数参数提取：num, den=tfdata(G,v) 或者num=G.numi,j; den=G.deni,j;%第i输入对第j输出间

6、的传函3( )sG s e3.1 连续线性系统的数学模型与MATLAB表示-线性系统的状态方程模型 (1)状态方程模型p路输入、q路输出、n个状态量、fi(.)和gi(.)可以是线性或非线性函数时变模型： 线性时不变模型：121121(,),1,(,),1,iinpiinpxf x xx uuinyg x xx uuiq ( )( ) ( )( ) ( )( )( ) ( )( ) ( )x tA t x tB t u ty tC t x tD t u t( )( )( )( )( )( )x tAx tBu ty tCx tDu t3.1 连续线性系统的数学模型与MATLAB表示-线性系统的

7、状态方程模型 (2)线性时不变模型的MATLAB描述G=ss(A,B,C,D); A矩阵为nn方阵，B为np矩阵，C为qn矩阵，D为qp矩阵 例3-51217.216.811.91.50.268.68.4610.3( )( )( )68.78.46215.98.68.3600.520.500.8( )( )0.30.30.21x tx tu ty tx t3.1 连续线性系统的数学模型与MATLAB表示-线性系统的状态方程模型 (3) A=-12, -17.2, -16.8, -11.9; 6.8, 6, 8.4, 6; 6, 8.7, 8.4, 6; -5.9, -8.6, -8.3, -6

8、 ; B=1.5, 0.2; 1, 0.3; 2, 1; 0, 0.5; C=2, 0.5, 0, 0.8; 0.3, 0.3, 0.2, 1; D=zeros(2, 2); G=ss (A, B, C, D);3.1 连续线性系统的数学模型与 MATLAB 表示-线性系统的状态方程模型 (4)带时间延迟的状态方程 数学模型MATLAB输入语句：( )( )()( )( )(),( )()iiox tAx tBu tz tCx tDu ty tz tio=ss( , , , , InputDelay,OutputDelay, ) GA B C D3.1 连续线性系统的数学模型与 MATLAB

9、表示-线性系统的零极点模型(1)零极点模型 是因式型传递函数模型 MATLAB表示:z=z1;z2;zm;p=p1;p2;pn;G=zpk(z,p,k);1212()()()( )()()()mnszszszG sKspspsp3.1 连续线性系统的数学模型与 MATLAB 表示-线性系统的零极点模型(2)例3-5 零极点模型MATLAB输入方法P=-1;-2;-3;-4; Z=-5;-2+2i;-2+2i;G=zpk(Z,P,6);另一种输入方法s=zpk(s,); G=6*(s+5)*(s+2+2i)*(s+2-2i)/ (s+1)*(s+2)*(s+3)*(s+4)6(5)(22 )(2

10、2 )( )(4)(3)(2)(1)ssj sjG sssss+-=+3.1 连续线性系统的数学模型与MATLAB表示-多变量系统传递函数矩阵模型(1)传递函数矩阵 为第i输出对第j输入的传递函数可以先定义子传递函数，再由矩阵定义G(s)111212122212( )( )( )( )( )( )( )( )( )( )ppqqqpgsgsgsgsgsgsG sgsgsgs轾犏犏犏=犏犏犏犏臌LLMMOML( )ijgs3.1 连续线性系统的数学模型与MATLAB表示-多变量系统传递函数矩阵模型(2)例37 多变量模型方法1：g11=tf(0.1134,1.78 4.48 1,ioDelay,

11、0.72);g12=tf(0.924, 2.07 1); g21=tf(0.3378, 0.361 1.09 1,ioDelay, 0.3); g22=tf(-0.318,2.93 1,ioDealy, 1.29);G=g11, g12; g21,g22;方法2：G=tf(0.1134,1.78 4.48 1), tf(0.924, 2.07 1); tf(0.3378, 0.361 1.09 1), tf(-0.318,2.93 1); G.ioDelay=0.72 0; 0.3 1.29;0.7220.31.2920.11340.9242.0711.784.481( )0.33780.31

12、82.9310.3611.091sssesssG seesssZ变换（Z-transformation），是对离散序列进行的一种数学变换。常用以求线性时不变差分方程的解。它在离散时间系统中的地位，如同拉普拉斯变换在连续时间系统中的地位。离散时间信号的Z变换已成为分析线性时不变离散时间系统问题的重要工具。在数字信号处理、计算机控制系统等领域有广泛的应用。3.2 线性离散时间系统的数学模型-离散传递函数模型(1)一般单变量离散系统 差分方程：Z变换代替Laplace变换：MATLAB表示：num=b0,b1, ,bn-1,bn;(采样周期T) den=a1,a2, ,an,an+1; H=tf(n

13、um, den, Ts,T);算子输入方法：z=tf(z, T)121011()(1) (1) () ()(1) (1) () nnnna y kTa y kTa y knTay kn Tb u kTbu kTbu knTb u kn T1011121( )nnnnnnb zb zbH za za za3.2 线性离散时间系统的数学模型-离散传递函数模型(2)例 38 离散传递函数，采样周期T=0.1秒输入方法1：z=tf(z,0.1); H=(6*z2-0.6*z-0.12)/(z4-z3+0.25*z-0.125);输入方法2：num=6 -0.6 -0.12; den=1 -1 0.25

14、 0.25 -0.125; H=tf(num, den, Ts,0.1);243260.60.12( )0.250.250.125zzH zzzzz3.2 线性离散时间系统的数学模型-离散传递函数模型(3)离散系统的时间延迟模型 延迟是采样周期的整数倍 MATLAB输入方法： H.ioDelay=m 或者 set(H, ioDelay,m)1011121( )nnmnnnnb zb zbH zza za za3.2 线性离散时间系统的数学模型-离散传递函数模型(4)滤波器型离散模型 记 ，则MATLAB表示方法：num=bn, bn-1, , b1, b0;den=an+1, an, an-1

15、 , , a1; H=tf(num, den, Ts,T,Variable, q); 11111011121()nnnnnnnnb zbzb zbH zaza za za 1qz1111011121( )nnnnnnnnb qbqb zbH qaza za za 3.2 线性离散时间系统的数学模型-离散传递函数模型(5)例 39 离线系统的零极点模型 采样周期T=0.1秒 z=1/2;1/2+j/2;1/2-j/2; p=-1/2;-1/3;-1/4;-1/5; H=zpk(z,p,1/20,Ts,0.1) (1/2)(1/2/2)(1/2/2)( )120(1/2)(1/3)(1/4)(1/

16、5)zzjzjH zzzzz3.2 线性离散时间系统的数学模型-离散状态方程模型一般的离散状态方程MATLAB表示： H=ss(F, G, C, D, Ts,T);带有延迟的离散系统状态方程MATLAB表示：H=ss(F, G, C, D, Ts, T, ioDelay, m);(1) ()()()()()x kTFx kTGu kTy kTCx kTDu kT(1) ()() ()()() x kTFx kTGu km Ty kTCx kTDu km T3.3 方框图描述系统的化简浙江大学电气学院系统系 包哲静前面介绍的传递函数、状态方程等都是单环节模型实际系统为若干个子模型互连的如何通过等

17、效变化进行化简主要内容控制系统的典型连接结构节点移动时的等效变换复杂系统模型的简化3.3 方框图描述系统的化简-控制系统的典型连接结构(1)系统串、并联 串联传递函数：并联传递函数：21( )( )( )G sG s G sG1(s)G2(s)u(t)y(t)G1(s)G2(s)u(t)y(t)(a) 串联结构(b) 并联结构12( )( )( )G sG sG s3.3 方框图描述系统的化简-控制系统的典型连接结构(2)串联系统的状态方程：并联系统的状态方程：1111221222112122120 xAxBuxB CAxB DxyD CCD Dux1111222211212200()xAxB

18、uxAxBxyCCDD ux3.3 方框图描述系统的化简-控制系统的典型连接结构(3)若一个模型为传递函数、另一个为状态方程，如何处理？【将二者变换成同样结构再计算！】基于MATLAB的计算方法串联 G=G2*G1 (多变量系统注意次序) 并联 G=G2+G1 3.3 方框图描述系统的化简-控制系统的典型连接结构(4)反馈连接正反馈 负反馈112( )( )1( )( )G sG sG s G s112( )( )1( )( )G sG sG s G sG1(s)G2(s)u(t)y(t)G1(s)G2(s)u(t)y(t)(a) 正反馈结构(a) 负反馈结构3.3 方框图描述系统的化简-控制

19、系统的典型连接结构(5)负反馈连接的数学模型：其中，若 ，则简化为1112112112212212221111212,()xAB ZD CB ZCxB ZuxB ZCAB D ZCxB D ZxyZCD ZCD Z ux112()ZID D120DD111211221221120, 0 xAB CxBuxB CAxxyCx 3.3 方框图描述系统的化简-控制系统的典型连接结构(6)反馈连接的MATLAB求解： G=feedback(G1,G2); %负反馈 G=feedback(G1,G2,1); %正反馈feedback()函数仅能用于G1,G2为具体参数指定的模型，通过适当扩展，可以处理符

20、号运算(置于sym目录下)：function H=feedback(G1,G2,key) if nargin=2 key=-1 end H=G1/(sym(1)-key*G1*G2); H=simple(H);3.3 方框图描述系统的化简-控制系统的典型连接结构(7)例3-10 s=tf(s); G=(12*s3+24*s2+12*s+20)/(2*s4+4*s3+6*s2+2*s+2); Gc=(5*s+3)/s; H=1000/(s+1000); GG=feedback(G*Gc,H)3243212241220( )24622sssG sssss53( )csG ssGc(s)H(s)R(

21、s)Y(s)G(s)1000( )1000H ss3.3 方框图描述系统的化简-控制系统的典型连接结构(8)例3-11 受控对象状态方程控制器为对角阵A=-12,-17.2,-16.8,-11.9; 6,8.6,8.4,6;6,8.7,8.4,6;-5.9 -8.6 -8.3 -6; B=-1.5,0.2;1,0.3;2,1;0,0.5; C=2,0.5,0,0.8;0.3,0.3,0.2,1; D=zeros(2,2);G=ss(A,B,C,D); s=tf(s); g11=(2*s+1)/s; g22=(5*s+2)/s; Gc=g11 0;0 g22; H=eye(2); GG=feed

22、back(G*Gc,H)1217.216.811.91.5268.68.4610.3( )( )( )68.78.46215.98.68.3600.520.500.8( )( )030.30.21x tx tu ty tx t(21)/0( )0(52)/cssG sss3.3 方框图描述系统的化简-结点移动的等效变换(1)考虑模型需要将A点等效移至输出端Y(s)G1(s)H(s)R(s)Y(s)G2(s)G3(s)G4(s)H2(s)H3(s)A3.3 方框图描述系统的化简-结点移动的等效变换(2)节点移动等效变换G1(s)G2(s)/G1(s)BAG1(s)G2(s) G1(s)BAG1(

23、s)G2(s) BAG1(s)G2(s)BA节点前移节点后移3.3 方框图描述系统的化简-复杂系统模型的简化(1)例3-12 原系统等效变换为syms G1 G2 G3 G4 H1 H2 H3 %定义各个子模块为符号变量 c1=feedback(G4*G3,H3); %最内层闭环 c2=feedback(c1*G2, H2/G4); %第二层闭环 G=feedback(c2*G1,H1); pretty(G) %总系统模型G1(s)H(s)R(s)Y(s)G2(s)G3(s)G4(s)H2(s)/ G4(s)H3(s)3.3 方框图描述系统的化简-复杂系统模型的简化(2)例3-13 电机拖动模

24、型输入r(t)单独激励syms Ka Kr c1 c2 c Ra T1 T2 Km Kb s Ga=feedback(1/Ra/(T1*s+1)*Km*1/c/(T2*s+1),Kb); G1=c1*feedback(Ka*Kr*Ga/s,c2); G1=collect(G1,s);a11/1RT s21/1cT skrc2r(t)n(t)1/skbc1kakmM(t)3.3 方框图描述系统的化简-复杂系统模型的简化(3) M(t)单独输入G2=-feedback(1/c/(T2*s+1)/s,Km/Ra/(T1*s+1)*(Kb*s+c2*Ka*Kr); G2=collect(simplif

25、y(G2),s);G=G1 G2;krn(t)kaM(t)1/s21/1cT sskbkma11/1RT sc23.4 系统模型的相互转换前面介绍的各种模型之间的相互等效变换主要内容连续模型和离散模型的相互转换系统传递函数的获取控制系统的状态方程实现状态方程的最小实现传递函数与符号表达式的相互转换3.4 系统模型的相互转换-连续模型与离散模型的相互转换(1)连续模型离散化连续状态方程的解析解：采样周期T，选择得到离散模型记则离散状态方程模型00()()( )(0)( )tA t tA ttx texeBud0,(1)tkT tkT (1)(1)(1)()( )kTATAkTkTxkTex kT

26、eBud0(1)()()TATAxkTex kTedBu kT0,TATAFeGedB(1) ()()x kTFx kTGu kT3.4 系统模型的相互转换-连续模型与离散模型的相互转换(2)MATLAB函数直接求解：还可采用Tustin变换(双线性变换)：例314 双输入模型，T=0.12(1)/ (1)szT z1217.216.811.91.50.268.68.4610.3( )( )( )68.78.46215.98.68.3600.520.500.8( )0.30.30.21kx tx tu tyx t1c2d(,)GG TG:传递函数或状态方程模型，还可处理延迟环节3.4 系统模型

27、的相互转换-连续模型与离散模型的相互转换(3)A=-12, -17.2,-16.8,-11.9;6,8.6,8.4,6;6,8.7,8.4,6;-5.9,-8.6,-8.3,-6; B=1.5,0.2;1,0.3;2,1;0,0.5; C=2,0.5,0,0.8;0.3,0.3,0.2,1; D=zeros(2,2); G=ss(A,B,C,D); T=0.1; Gd=c2d(G,T); 得到模型：10.149961.64811.60761.140.184230.127230.573541.8820.80180.573540.266840.103620.576450.836151.80590.

28、576450.36790.174010.566450.82610.795870.42360.165650.023263kkxx20.500.80.30.30.21kkkuyx3.4 系统模型的相互转换-连续模型与离散模型的相互转换(4)例315 时间延迟系统的离散化(采样周期T=0.1秒)s=tf(s); G=1/(s+2)3; G.ioDelay=2;分别用零阶保持器和Tustin算法进行离散化G1=c2d(G,0.1) 零阶保持器变换G2=c2d(G,0.1,tustin) %Tustin变换231( )(2)sG ses220320.00014360.00049460.0001064(

29、)2.4562.0110.5488ZOHzzGzzzzz532520329.391*100.00028170.00028179.391*10( )2.4552.0080.5477ZOHzzzGzzzzz仅从数值结果无法判断两个模型好坏3.4 系统模型的相互转换-连续模型与离散模型的相互转换(5)离散模型连续化 对变换 ，进行逆变换 Tustin反变换： MATLAB求解：例316 对前面的连续状态方程模型离散化，对结果再连续化A=-12, -17.2,-16.8,-11.9;6,8.6,8.4,6;6,8.7,8.4,6;-5.9,-8.6,-8.3,-6; B=1.5,0.2;1,0.3;2

30、,1;0,0.5; C=2,0.5,0,0.8;0.3,0.3,0.2,1; D=zeros(2,2); G=ss(A,B,C,D); T=0.1; Gd=c2d(G,T); G1=d2c(Gd)0,TATAFeGe d B11ln,()AF BFIAGT(1/2)(1/2)zsTsT1d2c()GG3.4 系统模型的相互转换-系统传递函数的获取(1)已知状态方程两端Laplace变换则得到传递函数 难点是求取 ，基于Fadeev-Fadeeva算法能得到较为精确的解由零极点模型，直接展开分子分母，得到传递函数模型用MATLAB求解：( )( )( )( )( )( )x tAx tBu ty

31、 tCx tDu t( )( )( )( )( )( )sIX sAX sBU sY sCX sDU s1( )()( )X ssIABU s11( )( )( )()G sY s UsC sIABD1()sIA1tf ()GG3.4 系统模型的相互转换-系统传递函数的获取(2)例3-17 多变量模型，求传递函数矩阵A=-12, -17.2,-16.8,-11.9;6,8.6,8.4,6;6,8.7,8.4,6;-5.9,-8.6,-8.3,-6; B=1.5,0.2;1,0.3;2,1;0,0.5; C=2,0.5,0,0.8;0.3,0.3,0.2,1; D=zeros(2,2); G=s

32、s(A,B,C,D); G1=tf(G);得到传递函数矩阵32324324323232432433.5144.120.690.83720.9564.139.1610.3740.350.050.00240.350.050.0024( )1.1536.326.2250.13390.8515.712.6190.045590.350.050.00240.3ssssssssssssssG sssssssssssss250.050.0024ss3.4 系统模型的相互转换-控制系统的状态方程实现(1)由传递函数到状态方程的转换 不同状态变量的选择，结果不唯一 默认状态方程实现，MATLAB函数：G1=ss(

33、G) 适用于有时间延迟的、离散的、多变量的系统3.4 系统模型的相互转换-控制系统的状态方程实现(2)例3-18 连续多变量模型g11=tf(0.1134,1.78 4.48 1, ioDelay,0.72); g12=tf(0.924,2.07 1); g22=tf(-0.318,2.93 1, ioDelay, 1.29); g21=tf(0.3378,0.361 1.09 1, ioDelay, 0.3); G=g11, g12;g21, g22; G1=ss(G);0.7220.31.2920.11340.9242.0711.784.481( )0.33780.3182.9310.36

34、11.091sssesssG seesss3.4 系统模型的相互转换-控制系统的状态方程实现(3)得到状态方程模型ioDelay矩阵122.51690.280900000.25020000000(0.3)003.01940.69252000.250( )( )( )0040000000000.48309001000000.341300.2500.12742000.446380( )0000.935730u tx tx tu tz t12(0.42)( ),( )(1.29)0.43413z tx ty tz t0.42001.29T3.4 系统模型的相互转换-控制系统的状态方程实现(4)该模型

35、可以转换回传递函数矩阵得出的转换结果21tf()GG0.72220.31.2920.063710.44640.48312.5170.56180.93750.10850.34133.0192.77sssesssGeesss 3.4 系统模型的相互转换-控制系统的状态方程实现(5)均衡实现为了将各个状态变量在整个控制系统中的重要程度表示出来，需要进一步变换 用MATLAB求解： 得出均衡实现的模型 得出排序的Gram矩阵 b, , balreal( )Gg TG3.4 系统模型的相互转换-状态方程的最小实现(1)例319 观察单变量系统传递函数模型若不进行变换，则不能发现该模型有哪些特点求取零极点

36、模型：G=tf(5 50 155 150,1 11 41 61 30); zpk(G) %得到零极点模型零极点模型32432550155150( )11416130sssG sssss5(3)(2)(5)( )(5)(3)(2)(1)sssG sssss3.4 系统模型的相互转换-状态方程的最小实现(2)零极点对消后，得到一阶模型完全对消相同零极点后的模型，又称为最小实现问题：若系统模型为多变量模型，如何获得最小实现？MATLAB中最小实现的方法：5( )1rG ssmminreal( )GG3.4 系统模型的相互转换-状态方程的最小实现(3)例3-20多变量模型不能直接看出是否是最小实现61

37、.5249.56462.52512.555( )( )( )50.250.53.59.753410.5011.502211233120.750.51.52.75( )( )01.251.51.52.25x tx tu ty tx t3.4 系统模型的相互转换-状态方程的最小实现(4)A=-6,-1.5,2,4,9.5;-6,-2.5,2,5,12.5;-5,0.25,-0.5,3.5,9.75;-1,0.5,0,-1,1.5;-2,-1,1,2,3 ; B=6,4;5,5;3,4;0,2;3,1; C=2,0.75,-0.5,-1.5,-2.75;0,-1.25,1.5,1.5,2.25; D

38、=zeros(2,2); G=ss(A,B,C,D); G1=minreal(G)2.41251.17290.170226.48434.0942( )0.739460.123330.37256( )5.15173.7888( )0.650671.67661.71083.2275.55720.842350.0737980.048876( )( )0.250850.361290.46861x tx tu ty tx t 3.4 系统模型的相互转换-传递函数与符号表达式相互转换传递函数到符号表达式 function P=tf2sym(G) P=poly2sym(G.num1,s)/poly2sym(

39、G.den1,s); 表达式到传递函数 function G=sym2tf(P) n,d=numden(P); G=tf(sym2poly(n),sym2poly(d);置于sym目录下符号表达式必须是系数已知的有理函数形式3.5 线性系统模型降阶与最小实现不同用低阶模型近似高阶模型最早由Edward J. Davison提出(1966)本节主要内容Pad与Routh降阶算法时间延迟模型的Pad近似带有延迟的最优降阶算法状态空间的降阶算法3.5 线性系统模型降阶-Pad与Routh降阶算法(1)原始模型：寻求降阶模寻求降阶模型型： ，其中假设展开原模型其中时间矩量 可以递推求出 1121/11

40、21( )rrrr kkkkkssGssss11211121( )mmmmnnnnbsb sb sbG sa sa sa sa11kkn2012( )G scc sc sic101,110,1,2,ikikijnijjcbcbc ai 3.5 线性系统模型降阶-Pad与Routh降阶算法(2)若已知状态方程模型，则时间距量的MATLAB求解：function M=timmomt(G,k)G=ss(G); C=G.c; B=G.b; iA=inv(G.a);iA1=iA; M=zeros(1,k);for i=1:k, M(i)=-C*iA1*B; iA1=iA*iA1; endPad降阶思想：

41、保留前r+k+1时间距量!(1)!01( ),0,1,!iisd G scCABiids iiiiiisc ssiiiiiic sssG(s)=C(SI-A)-1B+D3.5 线性系统模型降阶-Pad与Routh降阶算法(3)对比系数，则1010111 010211 01211000rrkrkrk rrkrk rrkrk rk rkk rrrcccccccccccccccc 3.5 线性系统模型降阶-Pad与Routh降阶算法(4)这样可以得出11121121.rrrkrrrkk rk rrr kcccccccccc 011010110010rkrrrk rcccccc 3.5 线性系统模型降

42、阶-Pad与Routh降阶算法(5)Pad降阶求解函数function G_r=pademod(G_Sys,r,k)c=timmomt(G_Sys,r+k+1); G_r=pade_app(c,r,k);function Gr=pade_app(c,r,k)w=-c(r+2:r+k+1);vv=c(r+1:-1:1);zeros(k-1-r,1);W=rot90(hankel(c(r+k:-1:r+1),vv);V=rot90(hankel(c(r:-1:1);x=1 (Ww); dred=x(k+1:-1:1)/x(k+1);y=c(1) x(2:r+1)*V+c(2:r+1);nred=y

43、(r+1:-1:1)/x(k+1); Gr=tf(nred,dred);3.5 线性系统模型降阶-Pad与Routh降阶算法(6)例3-21原始模型Pad近似G=tf(1 7 11 5,1 7 21 37 30); Gr=pademod(G,1,2);得到结果324327115( )7213730sssG sssss20.85440.6957( )1.0914.174sG sss3.5 线性系统模型降阶-Pad与Routh降阶算法(7)Pad降阶算法并不能保持原系统的稳定性例3-22 原始模型求取零极点模型num=0.067 0.6 1.5 2.016 1.55 0.6; den=0.067

44、0.7 3 6.67 7.93 4.63 1; G=tf(num,den);zpk(G)得到稳定的零极点模型5432654320.0670.61.52.0161.550.6( )0.0670.736.677.934.631sssssG sssssss22(5.92)(1.221)(0.897)(0.91711.381)( )(2.805)(1.856)(1.025)(0.501)(4.2615.582)sssssG sssssss3.5 线性系统模型降阶-Pad与Routh降阶算法(8)Pad近似Gr=pademod(G,1,3); zpk(Gr)得到不稳定的降阶模型 Pad降阶方法不能保证降

45、阶模型的稳定性20.6328(0.7695)( )(2.598)(1.1080.3123)rsG ssss3.5 线性系统模型降阶-Pad与Routh降阶算法(9)Routh算法(Hutton提出的)（较烦琐，从略）function G_r=routhmod(G_Sys,nr)num=G_Sys.num1;den=G_Sys.den1;n0=length(den);n1=length(num);a1=den(end:-1:1);b1=num(end:-1:1) zeros(1,n0-n1-1);for k=1:n0-1, k1=k+2;alpha(k)=a1(k)/a1(k+1);beta(k

46、)=b1(k)/a1(k+1); for i=k1:2:n0-1 a1(i)=a1(i)-alpha(k)*a1(i+1);b1(i)=b1(i)-beta(k)*a1(i+1);end,endnn=;dd=1;nn1=beta(1);dd1=alpha(1),1;nred=nn1;nred=dd1;for i=2:nr, nred=alpha(i)*nn1,beta(i);dred=alpha(i)*dd1,0; n0=length(dd); n1=length(dred); nred=nred+zeros(1,n1-n0),nn; dred=dred+zeros(1,n1-n0),dd;

47、nn=nn1; dd=dd1;nn1=nred;dd1=dred;EndG_r=tf(nred(nr:-1:1),dred(end:-1:1);3.5 线性系统模型降阶-Pad与Routh降阶算法(10)Routh算法中由于利用了Routh因子的近似方法对于稳定系统总能得到渐进稳定的降阶模型例3-23 仍考虑稳定模型 num=0.067 0.6 1.5 2.016 1.55 0.6; den=0.067 0.7 3 6.67 7.93 4.63 1; G=tf(num,den);Gr=zpk(routhmod(G,3)得到稳定的降阶模型：5432654320.0670.61.52.0161.5

48、50.6( )0.0670.736.677.934.631sssssG sssssss2r20.37792(0.94720.3423)( )(0.4658)(1.150.463)ssG ssssRouth算法得出的降阶模型分子阶次总是比分母阶次少13.5 线性系统模型降阶-时间延迟模型的Pad近似(1)假设，纯时间延迟项 的k阶传递函数矩阵为MATLAB近似函数se23123,23231/2()()( 1)()( )1/2()()()nknkknspspspsPsspspsps , pade( , )n dk缺点：分子和分母的阶次相同3.5 线性系统模型降阶-时间延迟模型的Pad近似(2)纯时

49、间延迟项可用Maclaurin级数近似：编写MATLAB函数function n,d=paderm(tau,r,k)c(1)=1;for i=2:r+k+1, c(i)=-c(i-1)*tau/(i-1);endGr=pade_app(c,r,k); n=Gr.num1(k-r+1:end);d=Gr.den1;其中r/k任意选择223311111!2!3!sesss 3.5 线性系统模型降阶-时间延迟模型的Pad近似(3)例3-24 纯延迟模型MATLAB求解tau=1; n1,d1=pade(tau,3);G1=tf(n1,d1) n2,d2=paderm(tau,1,3);G2=tf(n

50、2,d2) 两种方法得到结果( )sG se321321260120( )1260120sssG ssss132624( )61824sG ssss3.5 线性系统模型降阶-时间延迟模型的Pad近似(4)例3-25 已知带有延迟的线性模型通过以下命令得到近似模型cd=1;tau=2; for i=1:5, cd(i+1)=-tau*cd(i)/i; end; cd G=tf(3 1,1,3,3,1); c=timmomt(G,5); c_hat=conv(c,cd); Gr=zpk(pade_app(c_hat,1,3);2331( )(1)ssG ses20.20122(0.04545)(

51、)(0.04546)(0.40270.2012)sG ssss3.5 线性系统模型降阶-带有时间延迟的次最优降阶算法(1)降阶模型的降阶效果误差定义ISE准则( )TsG s e/( )sr mGs e( )r t( )e t222200( )( )( )hh t dtw t e t dt3.5 线性系统模型降阶-带有时间延迟的次最优降阶算法(2)原模型降阶模型降阶误差定义111111( )nTsTsnnnnnnb sbsbG s eesa sasa11/111( )rssrrr mmmmmssGs eesss/( ) ( )( ) ( )Tssr mE sG s eGs eR s3.5 线性

52、系统模型降阶-带有时间延迟的次最优降阶算法(3)参数向量误差MATLAB实现（略）调用格式111(, )mr220min( )( , )Jw t e tdt0opt_app( , ,key,)rGG r mG若降阶模型或原模型中有延迟环节，要对延迟采用Pad3.5 线性系统模型降阶-带有时间延迟的次最优降阶算法(4)例3-26 对传递函数进行降阶num=68.6131,80.3787,67.087,29.9339,8.8818,1; den=0.0462,3.5338,16.5609,28.4472,21.7611,7.6194,1; Gr=zpk(opt_app(tf(num,den),2,

53、3,0)得出降阶模型为2345234561 8.881829.933967.08780.378768.6131( )1 7.619421.761128.447216.56093.53380.0462sssssG sssssss221523.6536(0.34920.2482)( )(74.85)(3.8715.052)rssG ssss3.5 线性系统模型降阶-带有时间延迟的次最优降阶算法(5)例3-27 已知高阶模型den=conv(conv(conv(conv(5,1,2,1),0.7,1),1,1),0.4,1); G=tf(432,den);Gr=zpk(opt_app(G,0,2,1

54、)得到降阶模型432( )(51)(21)(0.71)(1)(0.41)G ssssss1.5r31.4907( )(0.3283)(0.222)sG sess3.5 线性系统模型降阶-状态方程模型的降阶算法(1)均衡实现模型的降阶算法MATLAB求解函数 Gr=modred(G,elim)111121111222122222,xAAxBxuyCCDuxAAxBx 1111112222111122221112222112222()()()()xAA A AxBA A B uyCC A AxDC A B u销去状态变量x23.5 线性系统模型降阶-状态方程模型的降阶算法(2)例3-28G=tf(

55、1,7,24,24),1,10,35,50,24); G_b,g=balreal(ss(G) 得到Gram向量g=0.5179,0.0309,0.0124,0.0006T 销去第3，4状态变量G_r=modred(G_b,3，4)；zpk(G_r) 得到r0.025974(4.307)(22.36)( )(1.078)(2.319)ssG sss3.5 线性系统模型降阶-状态方程模型的降阶算法(2)基于Schur均衡实现模型的降阶算法MATLAB求解函数：Gr=schmr(G, 1, k)例3-29 高阶传递函数思路：先转换成状态方程，再降阶 num=68.6131,80.3787,67.08

56、7,29.9339,8.8818,1; den=0.0462,3.5338,16.5609,28.4472,21.7611,7.6194,1; G=ss(tf(num,den); Gh=zpk(schmr(G,1,3)得到Schur降阶模型：2345234561 8.881829.933967.08780.378768.6131( )1 7.619421.761128.447216.56093.53380.0462sssssG sssssss2r21485.3076(0.17890.2601)( )(71.64)(3.8814.188)ssG ssss3.5 线性系统模型降阶-状态方程模型的降

57、阶算法(3)最优Hankel范数的降阶模型近似MATLAB求解函数： Gr=ohklmr(G, 1, k)例3-30 仍采用前面模型 num=68.6131,80.3787,67.087,29.9339,8.8818,1; den=0.0462,3.5338,16.5609,28.4472,21.7611,7.6194,1; G=ss(tf(num,den); Gr=zpk(ohklmr(G,1,3)得到2r21527.8048(0.27640.2892)( )(73.93)(3.8554.585)ssG ssss3.5 线性系统模型降阶-降阶算法综述状态方程方法不能任意选择分母分子阶次，而许

58、多传递函数方法可以降阶效果比较，下章给出时域响应比较频域响应比较降阶模型的应用仿真应用（用途越来越小）控制器设计应用3.6 线性系统的模型辨识模型辨识由已知实测数据获得系统模型的方法实测数据时域响应数据、频率响应数据主要内容离散系统辨识方法辨识信号生成多变量系统辨识离散系统在线辨识3.6 线性系统的模型辨识-离散系统的模型辨识(1)离散传递函数模型对应的差分方程模型已知实测信号： 输入 输出 12121111( )mmdmmnnnnb zb zbzbG zzza zaza12121( )(1)(2)()()(1)(1)( )nmy ta y ta y ta y tnbu tdb u tdbu

59、tdmt (1), (2), ()Tuuuu M (1), (2), ()Tyyyy M3.6 线性系统的模型辨识-离散系统的模型辨识(2)由数据可以得出矩阵形式 ，其中1111111(1)(0)(1)(1)()(1)(2)(1)(2)(2)(1)(2)()(1)nmnmya ya ynbudbu mdya ya ynbudbumdy Ma y M 11()()(1)()nma y Mnbu Mdbu MdmMy(0)(1)(1)()(1)(2)(2)(1)(1)()()(1)yynudu m dyynudum dy My Mnu Mdu Mm d 3.6 线性系统的模型辨识-离散系统的模型辨

60、识(3)定义残差最小指标最小二乘解系统辨识工具箱求解：T=arx(y,u,m,n,d) 其中，T为结构体变量，T.a,T.b,tf(T)当然，由前面的公式也能直接求解121,Tnmaaa bb 21min( )Mii1y 3.6 线性系统的模型辨识-离散系统的模型辨识(4)例3-31实测数据如表3-1所示基于MATLAB求解 u=;y=; t1=arx(y,u,4,4,1) 可得系统模型为 还可以写成81234112344.83 1060.59990.1196()10.250.250.125qqqqG zqqqq8324324.83 1060.59990.1196( )0.250.250.12