1、本章主要內(nèi)容線性系統(tǒng)定性分析線性系統(tǒng)時(shí)域響應(yīng)解析解法線性系統(tǒng)的數(shù)字仿真分析根軌跡分析線性系統(tǒng)頻域分析4.1 4.1 線性系統(tǒng)性質(zhì)分析-主要內(nèi)容 線性系統(tǒng)穩(wěn)定性分析線性反饋系統(tǒng)內(nèi)部穩(wěn)定性分析線性系統(tǒng)的相似變換線性系統(tǒng)可控性分析線性系統(tǒng)可觀測(cè)性分析Kalman分解系統(tǒng)狀態(tài)方程的標(biāo)準(zhǔn)型系統(tǒng)的范數(shù)測(cè)度及求解4.1 4.1 線性系統(tǒng)性質(zhì)分析-穩(wěn)定性分析(1)(1) ()()()()()x kTFx kTGu kTy kTCx kTDu kT00()()0( )( )( )tA t tA ttx tex teBud( )( )( )( )( )( )x tAx tBu ty tCx tDu t110()(

2、0)()kkk iix kTF xFGu iT 4.1 4.1 線性系統(tǒng)性質(zhì)分析-穩(wěn)定性分析(2)基于MATLAB的穩(wěn)定性直接判定狀態(tài)方程模型 由eig(A)可求出所有特征根 離散系統(tǒng)：abs(eig(A)傳遞函數(shù)模型：完全同樣方法基于MATLAB的圖解判定法連續(xù)系統(tǒng)：pzmap(G)離散系統(tǒng)：pzmap(G)，同時(shí)畫出單位圓4.1 4.1 線性系統(tǒng)性質(zhì)分析-穩(wěn)定性分析(3)例4-1 單位負(fù)反饋構(gòu)成的閉環(huán)系統(tǒng)穩(wěn)定性判定直接分析方法num=10 50 100 100 40; den=1 21 184 870 2384 3664 2496 0; G=tf(num,den);GG=feedback(

3、G,1); pzmap(GG) eig(GG)用zpk(GG)命令得出零極點(diǎn)模型，從而判斷穩(wěn)定性432765432105010010040( )21184870238436642496ssssG ssssssss4.1 4.1 線性系統(tǒng)性質(zhì)分析-穩(wěn)定性分析(4)例4-2 高階離散單位負(fù)反饋系統(tǒng)模型den=1 -1 0.25 0.25 -0.125; num=6 -0.6 -0.12; H=tf(num,den, Ts,0.1); z=tf(z, Ts,0.1); Gc=0.3*(z-0.6)/(z+0.8); GG=feedback(H*Gc,1); pzmap(GG), abs(eig(GG

4、)2432c60.60.12( ),250.6( )0.3,0.10.8zzH zzzzzzG zTz4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性反饋系統(tǒng)內(nèi)部穩(wěn)定性(1)輸入、輸出穩(wěn)定是不夠的，因?yàn)槿魞?nèi)部信號(hào)可能過大，對(duì)系統(tǒng)作硬件破壞應(yīng)該引入內(nèi)部穩(wěn)定性概念，保證內(nèi)部信號(hào)也是穩(wěn)定的。帶有擾動(dòng)的線性反饋控制系統(tǒng)如下圖所示c( )G sud1xr2x( )G sy( )H svn3x4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性反饋系統(tǒng)內(nèi)部穩(wěn)定性(2)從輸入信號(hào)(r,d,n)到內(nèi)部輸出信號(hào)(x1,x2,x3)都穩(wěn)定的系統(tǒng)稱為內(nèi)部穩(wěn)定系統(tǒng)傳遞函數(shù)矩陣 其中逐一判斷每個(gè)子傳遞函數(shù)的穩(wěn)定性很繁瑣1

5、2cc3c1( )( )( )1( )1( )( )( )( )( )( )1xG s H sH srxG sG s H sdM sxG s G sG sn c( )1( )( )( )M sG s G s H s 4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性反饋系統(tǒng)內(nèi)部穩(wěn)定性(3)閉環(huán)系統(tǒng)內(nèi)部穩(wěn)定的充要條件為 1+H(s)G(s)Gc(s)沒有不穩(wěn)定的零點(diǎn) H(s)G(s)Gc(s)沒有不穩(wěn)定零極點(diǎn)對(duì)消第一個(gè)條件等效于輸入輸出穩(wěn)定性判斷第2個(gè)條件即可可以編寫MATLAB函數(shù)判斷內(nèi)部穩(wěn)定性 key=intstable(G,Gc,H)4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性反饋系統(tǒng)內(nèi)部穩(wěn)定性(4)判定的

6、MATLAB函數(shù)function key=intstable(G,Gc,H) GG=minreal(feedback(G*Gc,H);G0=H*G*Gc;G01=minreal(G0);p=eig(GG);z0=eig(G0);z1=eig(G01);zz=setdiff(z0,z1);if(G.Ts1), %離散系統(tǒng)判定 key=any(abs(p)1); if key=0, key=2*any(abs(zz)1);endelse, %連續(xù)系統(tǒng)判斷 key=any(real(p)0); if key=0, key=2*any(real(zz)0);endend1:輸入輸出不穩(wěn)定；2：輸入輸出

7、穩(wěn)定、但內(nèi)部不穩(wěn)定；0：內(nèi)部穩(wěn)定4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性相似變換(1)系統(tǒng)的狀態(tài)方程表示稱為系統(tǒng)實(shí)現(xiàn)不同狀態(tài)選擇下，狀態(tài)方程不唯一相似變換非奇異矩陣T狀態(tài)變換z=T-1x新狀態(tài)方程模型1( )( )( ),(0)(0)( )( )( )ttttz tAz tBu tzTxy tC z tDu t4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性相似變換(2)狀態(tài)變換公式MATLAB求解方法G1=ss2ss(G,T)11,ttttATATBT BCCTDD輸入輸出都要是狀態(tài)方程對(duì)象4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性相似變換(3)例4-3 已知系統(tǒng)和轉(zhuǎn)換矩陣A=0 1 0 0;0 0 1 0

8、;0 0 0 1;-24 -50 -35 -10; G1=ss(A,0;0;0;1,24 24 7 1,0); T=fliplr(eye(4); G2=ss2ss(G1,T)0100010010( )010001( )( )012450351011( )242471( )x tx tu tTy tx t 4.1 4.1 線性系統(tǒng)性質(zhì)分析-線性相似變換(4)變換結(jié)果：相似變換能改變系統(tǒng)的結(jié)構(gòu)引入相似變換矩陣可以將已知系統(tǒng)轉(zhuǎn)換成其他的形式10355024110000( )( )( )0100000100( )172424( )z tz tu ty tz t 4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性

9、分析(1)可控性定義假設(shè)系統(tǒng)由狀態(tài)方程(A,B,C,D)給出，對(duì)任意的初始時(shí)刻t0，如果狀態(tài)空間中任一狀態(tài)xi(t)可以從初始狀態(tài)xi(t0)處，由有界的輸入信號(hào)u(t)的驅(qū)動(dòng)下，在有限時(shí)間tf內(nèi)能夠到達(dá)任意預(yù)先指定的狀態(tài)xi(tf)，則稱此狀態(tài)是可控的；如果系統(tǒng)中所有的狀態(tài)都是可控的，則稱該系統(tǒng)為完全可控的系統(tǒng)；系統(tǒng)的可控性就是指系統(tǒng)內(nèi)部的狀態(tài)是不是可以由外部輸出信號(hào)控制的性質(zhì)；4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析(2)可控性判定矩陣方法 若矩陣Tc滿秩，則系統(tǒng)完全可控 基于MATLAB的判定方法： rank(T) 構(gòu)造可控性判定矩陣： Tc=ctrb(A,B)21 ,ncTB AB

10、 A BAB4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析(3)例4-4 離散狀態(tài)方程的可控性A=-2.2 -0.7 1.5 -1;0.2 -6.3 6 -1.5;0.6 -0.9 -2 -0.5;1.4 -0.1 -1 -3.5; B=6 9;4 6;4 4;8 4; Tc=ctrb(A,B); rank(Tc)或者Tc=B,A*B,A2*B,A3*B; %用直接方法建立可控性矩陣 rank(Tc)1690.26.361.546(1) ()()0.60.920.5484x kTx kTu kT4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析(4)由Gram矩

11、陣的非奇異判定可控性 引入可控Gram矩陣： 該矩陣滿足Lyapunov方程 MATLAB求解：Lc=lyap(A,B*B) 可控Gram矩陣還可以由Gc=gram(G,c)直接求出0TAtTA tcLeBB edtTTccALL ABB 4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析(5)例4-5 求Gram矩陣num=6 -0.6 -0.12;den=1 -1 0.25 0.25 -0.125; H=tf(num,den, Ts,0.1); Lc=gram(ss(H), c)243260.60.12( )25zzH zzzzz4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析

12、(6)可控性階梯分解 對(duì)于不完全可控的系統(tǒng)進(jìn)行階梯分解 階梯標(biāo)準(zhǔn)型 MATLAB函數(shù)調(diào)用 Ac,Bc,Cc,Tc=ctrbf(A,B,C) 若原系統(tǒng)狀態(tài)方程完全可控，則不必分解2100,cccccccAABCCCBAA4.1 4.1 線性系統(tǒng)性質(zhì)分析-可控性分析(7)例4-6 不完全可控系統(tǒng)A=-2.2 -0.7 1.5 -1;0.2 -6.3 6 -1.5; . 0.6 -0.9 -2 -0.5; 1.4 -0.1 -1 -3.5; B=6 9;4 6;4 4;8 4;C=1 2 3 4; Ac,Bc,Cc,Tc=ctrbf(A,B,C);1690.26.361.546(1

13、) ()()0.60.920.5484x kTx kTu kT4.1 4.1 線性系統(tǒng)性質(zhì)分析-可觀測(cè)性分析(1)可觀測(cè)性定義假設(shè)系統(tǒng)由狀態(tài)方程(A,B,C,D)給出,對(duì)任意的初始時(shí)刻t0,如果狀態(tài)空間中任一狀態(tài)xi(t)在任意有限時(shí)刻tf的狀態(tài)xi(tf)可以由輸出信號(hào)在這一時(shí)間區(qū)間內(nèi)tt0, tf的值精確地確定出來，則稱此狀態(tài)是可觀測(cè)的。如果系統(tǒng)中所有的狀態(tài)都是可觀測(cè)的，則稱該系統(tǒng)為完全可觀測(cè)的系統(tǒng)。系統(tǒng)的可觀測(cè)性就是指系統(tǒng)內(nèi)部的狀態(tài)是不是可以由系統(tǒng)輸出信號(hào)重建起來的性質(zhì)。4.1 4.1 線性系統(tǒng)性質(zhì)分析-可觀測(cè)性分析(2)可觀測(cè)性與可控性是對(duì)偶問題可觀測(cè)性判定矩陣方

14、法 系統(tǒng)(A,C)的可觀測(cè)性問題等同于(AT ,CT) 的可控性問題 MATLAB求解可觀測(cè)矩陣obsv(), obsvf()21onCCATCACA4.1 4.1 線性系統(tǒng)性質(zhì)分析-可觀測(cè)性分析(3)Gram矩陣判定可觀測(cè)性 MATLAB求解：gram(G, o) Gram矩陣滿足Lyapunov方程 0TA tTAtoLeC CedtTTooA LL AC C 4.1 4.1 線性系統(tǒng)性質(zhì)分析-kalman規(guī)范分解(1)Kalman規(guī)范分解,1,2,3,13,2,3,4,4,2,0000000( )( )( )00( )00( )c oc oc oc oc oc oc oc oAAAz t

15、z tu tBAAAABAAy tCCz t4.1 4.1 線性系統(tǒng)性質(zhì)分析-kalman規(guī)范分解(2)子空間 為既不可控又不可觀測(cè)的子空間， 為不可控但可觀測(cè)的子空間，和 分別為可控但不可觀測(cè)的子空間和既可控又可觀測(cè)的子空間。既可控又可觀的子空間就是前面提及的最小實(shí)現(xiàn)模型,(,0,0)c oA,(,0,)c oc oAC,(,0)c oc oAB,(,)c oc oc oABC4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(1)常用標(biāo)準(zhǔn)型單變量系統(tǒng)的標(biāo)準(zhǔn)型可控標(biāo)準(zhǔn)型實(shí)現(xiàn)可觀測(cè)標(biāo)準(zhǔn)型實(shí)現(xiàn)Jordan標(biāo)準(zhǔn)型實(shí)現(xiàn)多變量系統(tǒng)Leunberge標(biāo)準(zhǔn)型側(cè)重點(diǎn)：如何用MATLAB直

16、接獲取標(biāo)準(zhǔn)型4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(2)單變量系統(tǒng)的標(biāo)準(zhǔn)型 可控標(biāo)準(zhǔn)型 可觀測(cè)標(biāo)準(zhǔn)型ccccxA xB uyC xDu12120100000000101nnxxuaaaybbb x ooooxA xB uyC xD u1122330001000100010,0,1nnababxaxbuabyx4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(3)可控標(biāo)準(zhǔn)型和可觀測(cè)標(biāo)準(zhǔn)型互為對(duì)偶形式 可控可觀測(cè)標(biāo)準(zhǔn)型MATLAB實(shí)現(xiàn) function Gs=sscanform(G,type) switch type case ctrl G=t

17、f(G);Gs=; G.num1=G.num1/G.den1(1); G.den1=G.den1/G.den1(1); d=G.num1(1); G1=G; G1.ioDelay=0; G1=G1-d; num=G1.num1;den=G1.den1;n=length(G.den1)-1;,TTTTcocococoAABCCBDD4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(4) A=zeros(n-1,1) eye(n-1);-den(end:-1:2); B=zeros(n-1,1);1; C=num(end:-1:2); D=d; Gs=ss(A,B,C,D,T

18、s,G.Ts,ioDelay,G.ioDelay); case obsv Gc=sscanform(G,ctrl); Gs=ss(Gc.a,Gc.c,Gc.b,Gc.d,Ts,G.Ts,ioDelay,G.ioDelay); otherwise error(Only options ctrl and obsv are applicable.) end 4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(5)Jordan標(biāo)準(zhǔn)型 假設(shè)系統(tǒng)矩陣A的特征根為1 ,2 , , n ，第i個(gè)特征根i 對(duì)應(yīng)特征向量為v i ，則 矩陣A對(duì)應(yīng)的模態(tài)矩陣定義為 MATLAB變換： G1,T=

19、canon(G, modal) ,1,2,iiiAvvin 12kJJJz=Tx4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(6)多變量系統(tǒng)的可控標(biāo)準(zhǔn)型實(shí)現(xiàn)-Luenberger標(biāo)準(zhǔn)型 可控性判定矩陣 構(gòu)造矩陣1211111122,ppSb AbAb bAbAb11211TTTllLl提取此行提取此行 : 保證前面各列線性無關(guān)的最大指數(shù)值，即最大可控性指數(shù)is4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(7)得出Luenberger變換矩陣編寫Luenberger 實(shí)現(xiàn)的MATLAB函數(shù): T=luenberger(A,B)111212111

20、TTTll ATlA11ss,ATAT BT B4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(8)MATLAB函數(shù)清單function T=luenberger (A,B)n=size(A,1); p=size(B,2); S=; sigmas=; k=1;for i=1:p for j=0:n-1 S=S, Aj*B(:,i) ; if rank(S)=k, k=k+1; else, sigmas(i)=j-1; S=S(:,1:end-1); break; end, end if kn, break; endendk=k-1; 4.1 4.1 線性系統(tǒng)性質(zhì)分析-系

21、統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(9)if kn % 如果不是完全可控，則用隨機(jī)數(shù)補(bǔ)足 while rank(S)=n, S(:,k+1:n)=rank(n,n-k); endendL=inv(S); iT=;for i=1:p for j=0:sigmas(i) iT=iT; L(i+sum(sigmas(1:i),:)*Aj; end,endif knum=6 0 2 8 10;den=2 0 6 4 8; G=tf(num,den); Gs=sscanform(G,obsv)例4-8A=15,6,-12,9;4,14,8,-4;2,4,10,-2;9,6,-12,15; B=3 3;2

22、 2;-2 -2;3 9;T=luenberger(A,B) A1=inv(T)*A*T, B1=inv(T)*B424262810( )648sssG ssss156129334148422( )( )( )241022296121539x tx tu t4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(12)系統(tǒng)的范數(shù)測(cè)度及求解 連續(xù)系統(tǒng)的范數(shù)定義 H2范數(shù)： H范數(shù)： H范數(shù)是系統(tǒng)頻域響應(yīng)幅值的峰值2211( ) ()2pjijiG sG jdj ( )sup ()G sG j4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(13)離散系統(tǒng)的范

23、數(shù)定義范數(shù)的MATLAB求解 norm(G) 或 norm(G,2) norm(G,inf)221( ) ()pjiiG zG ed( )sup ()jG zG e4.1 4.1 線性系統(tǒng)性質(zhì)分析-系統(tǒng)狀態(tài)方程標(biāo)準(zhǔn)型的MATLAB求解(14)例4-9 已知離散系統(tǒng)模型A=-2.2 -0.7 1.5 -1;0.2 -6.3 6 -1.5;0.6 -0.9 -2 -0.5;1.4 -0.1 -1 -3.5; B=6 9;4 6;4 4;8 4; C=1 2 3 4; G=ss(A,B,C,0 0); norm(G,2), norm(G,inf)1690.26.361.546(1)

24、 ()()0.60.920.5484()1234 ()x kTx kTu kTy kTx kT 4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法給線性系統(tǒng)一個(gè)激勵(lì)信號(hào)，輸出是什么？有兩大類方法解析解方法 求解微分方程、差分方程解析解數(shù)值解方法主要內(nèi)容基于狀態(tài)方程的解析解方法基于傳遞函數(shù)部分方式展開的解析解方法二階系統(tǒng)的解析解方法4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(1)狀態(tài)方程模型解析解求解難點(diǎn)( )( )( )( )( )( )x tAx tBu ty tCx tDu t00()()0( )( )( )tA t tA ttx tex teBud0(

25、)( )tA tteBud4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(2)狀態(tài)增廣方法消除矩陣B，變成自治系統(tǒng)單位階躍信號(hào) ，若假設(shè)有另外一個(gè)狀態(tài)變量 ，則其導(dǎo)數(shù)為增廣狀態(tài)方程自治系統(tǒng)可以直接求解析解 ( )1( )u tt1( )( )nxtu t1( )0nxt11( )( )( )( )00nnttABxtxtxx( )(0)Att = exx4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(3)一般輸入信號(hào)的系統(tǒng)增廣一般輸入信號(hào)模型引入增廣狀態(tài)變量11224340( )( )( )cos()sin()md tiiiu tu tu tctedd

26、tdd t111424(1)3131cos()sin()( ),( )d tnd tnmnn mxed txed txu txut4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(4)增廣狀態(tài)方程模型其中，解析解( )(0)Att = exx231441000001000000,000Ad Bd BBddddA 12343( )( )( )( )( )( )( )nnnnn mtxtxttxtxtxt xx( )(0)Att =exx01(0)10(0)!mccc mxx4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(5)MATLAB實(shí)現(xiàn)函數(shù) functi

27、on Ga,Xa=ss_augment(G,cc,dd,X) G=ss(G); Aa=G.a; Ca=G.c; Xa=X; Ba=G.b; D=G.d; if (length(dd)0 & sum(abs(dd)1e-5), if(abs(dd(4)1e-5), Aa=Aa dd(2)*Ba, dd(3)*Ba; zeros(2,length(Aa), dd(1),-dd(4); dd(4),dd(1); Ca=Ca dd(2)*D dd(3)*D; Xa=Xa;1;0; Ba=Ba;0;0; else, Aa= Aa dd(2)*B; zeros(1,length(Aa) dd(1)

28、; Ca= Ca dd(2)*D; Xa=Xa; 1; Ba=B;0; end end4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(6) if(length(cc)0 & sum(abs(cc)1e-5), M=length(cc); Aa=Aa Ba zeros(length(Aa),M-1); zeros(M-1,length(Aa)+1) eye(M-1);. zeros(1,length(Aa)+M); Ca=Ca D zeros(1, M-1); Xa=Xa; cc(1); ii=1; for i=2:M, ii=ii*I; Xa(length(Aa)+

29、i)=cc(i)*ii; end, end Ga=ss(Aa,zeros(size(Ca),Ca,D);調(diào)用格式信號(hào)描述00,ss_augment( , , ,)aGG xc d x010124,mc ccd d ddc =d4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(7)例4-10 連續(xù)系統(tǒng)模型初值輸入信號(hào)求解析解191616191211617190( )( )( )201716201201616192( )2,1,0,0 ( )ttu ty tt xxxT(0)0,1,1,2x3( )22sin(2 )tu tet4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方

30、程的解析解方法(8)系統(tǒng)增廣cc=2; dd=-3,0,2,2; x0=0;1;1;2; A=-19,-16,-16,-19; 21,16,17,19;20,17,16,20;-20,-16,-16,-19; B=1;0;1;2; C=2 1 0 0; D=0; G=ss(A,B,C,D); Ga, xx0=ss_augment(G,cc,dd,x0); Ga.a, xx0增廣模型 -19-16-16-1902102116171900012017162002011( ) =( ),(0) =-20-16-16-1904220000-3-20100002-30000000002tt xxx4.2

31、 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于狀態(tài)方程的解析解方法(9)解析解求解syms t; y=Ga.c*expm(Ga.a*t)*xx0; latex(y);解析解的數(shù)學(xué)形式323312711913577( )54574cos(2 )sin(2 )4884tttttty tteeet eetet4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(1)連續(xù)系統(tǒng)的解析解法 輸入信號(hào)的Laplace變換U(s) 輸出信號(hào)的Laplace變換Y(s)=G(s) U(s) 無重根時(shí)部分分式展開 由Laplace反變換求解析解112112121( )mmmmnnnnnbsb sb sbG

32、 ssa sa sasa1212( )mmrrrY sspspsp12-112( ) ( )mp tp tp tmy tY srerer eL4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(2)有重根時(shí)相應(yīng)項(xiàng)的解析解為部分分式的MATLAB求解：r,p,K=residue(num,den)112()()jjj mmjjjrrrspspsp11111111!(1)!111!(1)!jjjjmp tp tp tjmjjp tmjjjmrtr erteemrrtrtem-+-+-+-+-輊犏=+犏-臌LL4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(3)例4

33、-11 系統(tǒng)的傳遞函數(shù)為 ，輸入信號(hào)為階躍信號(hào)R(s)=1/s，則輸出信號(hào)的Laplace變換為MATLAB求解：num=1 7 3 4; den=1 7 17 17 6; R,P,K=residue(num, den, 0); R,P解析解：32432734( )717176sssG sssss+=+325432734( )717176sssY ssssss+=+32( )2.583395.753.50.667tttty teeete-=-+-+4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(4)例4-12 帶有復(fù)數(shù)極點(diǎn)的系統(tǒng)階躍響應(yīng)解析解num=1 3; den=1 2

34、 11 18 18; r,p,k=residue(num,den,0);r,p解析解4323( )2111818sG sssss+=+33( 1)( 1)( )(0.0020.0255 )(0.0020.0255 )( 0.08530.0088 )( 0.08530.0088 )0.1667jtjtj tj ty tj ej ej ej e- +- -=+-+ -+ -+4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(5)解析解的進(jìn)一步化簡(jiǎn) 基于Euler公式的化簡(jiǎn) 其中， 新MATLAB函數(shù)： r,p,K=pfrac(num,den) 若p(i)為實(shí)數(shù)，則(r(i),p(

35、i)同residue() 若p(i)為復(fù)數(shù)( )，則r(i),r(i+1)對(duì)返回A和()()()()sin()jtjttabj eabj eAetswswswF+-+-=+222,arctan(/ )Aaba bF= -+=-( )p ijsw=4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(6)MATLAB函數(shù)function R,P,K=pfrac(num,den)R,P,K=residue(num,den);for i=1:length(R), if imag(P(i)eps a=real(R(i); b=imag(R(i); R(i)=-2*sqrt(a2+b2);

36、 R(i+1)=-atan2(a,b); elseif abs(imag(P(i)num=1 3; den=1 2 11 18 18; r,p,k=pfrac(num,den,0);r,p解析解4323( )2111818sG sssss+=+( )0.0511sin(30.0768)0.1715sin(1.4677)0.1677ty ttet-= -+4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(7)基于Laplace變換的求解參見附錄A步驟： 定義符號(hào)變量 描述原函數(shù)表達(dá)式 調(diào)用ilaplace()函數(shù)求解 結(jié)果化簡(jiǎn)，如simple()函數(shù)求解舉例4.2 4.2 線性

37、系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(8)例1 求系統(tǒng) 的階躍響應(yīng)解析解 MATLAB求解syms s; G=(s3+7*s2+3*s+4)/(s4+7*s3+17*s2+17*s+6); y=ilaplace(G/s) latex(y)解析解32432734( )717176sssG sssss+=+3223123( )7/293124ttty tetee-驏=+ -+-桫4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(9)例2 求系統(tǒng) 的階躍響應(yīng)解析解 MATLAB求解：syms s; G=(s+3)/(s4+2*s3+11*s2+18*s+18); ilap

38、lace(G/s), latex(ans)解析解4323( )2111818sG sssss+=+1131293( )cos(3 )sin(3 )cossin2552556170170tty tttetet-=-+-4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(10)離散系統(tǒng)的解析解法無重根時(shí) 部分分式展開 解析解-111nqqzpp p-輊驏犏= -犏-桫臌L1211112( )mmrrrY zzpzpzp-=+-L-1121122111( )( )nnnmmmrrry nY zpppppp驏驏驏鼢瓏鼢= -瓏鼢瓏鼢瓏桫桫桫LL對(duì)z-1進(jìn)行展開4.2 4.2 線性系統(tǒng)時(shí)

39、域響應(yīng)解析解法-基于部分分式展開方法求解(11)考慮采樣周期例4-14 系統(tǒng)的離散傳遞函數(shù)為 ，輸入為階躍信號(hào) D=conv(1 -1/2,conv(1 -1/4, conv(1 1/5, 1 -1); %分母 N=0 0 conv(1 -1/3, 1 0); %分子，因?yàn)橐嫘蚺帕?，前面補(bǔ)足零 N=N(end:-1:1); D=D(end:-1:1); %因?yàn)樾枰獙?duì)z-1進(jìn)行展開，所以逆序排列 R,P,K=residue(N,D); R,P,-R./P-1121122111()( )kTkTkTmmmrrry kTY zpppppp驏驏驏鼢瓏鼢= -瓏鼢瓏鼢瓏桫桫桫LL(1/3)( )(1/

40、2)(1/4)(1/5)zG zzzz-=-+4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(12)輸出信號(hào)解析解Z變換求解步驟定義符號(hào)變量調(diào)用iztrans()函數(shù)求解化簡(jiǎn)11117.05473.95063.80951.4815( )5421Y zzzzz-=+-111( )1.41090.98771.90481.4815542nnny n驏驏驏鼢瓏=-+鼢瓏鼢瓏桫桫桫-4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(13)利用符號(hào)運(yùn)算工具箱求解syms z G=(z-1/3)/(z-1/2)*(z-1/4)*(z+1/5) iztrans(G*z/(

41、z-1)求解結(jié)果408080040( )(1/2)(1/4)( 1/5)218156727nnny n = -+-+4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(14)有重根時(shí)，部分分式表達(dá)式的Z反變換例4-15 系統(tǒng)的離散傳函為 ，求階躍響應(yīng)的解析解部分分式展開D=conv(1 -1/2,conv(1 -1/2,conv(1, -1/2,conv(1,-1/3,1,-1); N=0,0,0,5,-2,0; R, P=residue(N(end:-1:1), D(end:-1:1); R P -11( 1)(1)(2)(1)()(1)!()mmnmqqnnnmzpmp-+

42、輊-犏=+-犏-臌LL352( )(1/2) (1/3)zG zzz-=-4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(15)得到部分分式展開后Z變換表達(dá)式解析解11121313242409619236( )32(2)(2)1Y zzzzzz-=+-232324 1240 1961( )(1)3 322( 2)2192/2 1(1)(2)36( 2)211108( 126072)3632nnnnnny nnnnnn驏驏驏-鼢瓏=+鼢瓏鼢瓏桫桫桫-驏+桫-驏驏鼢瓏= -+-+鼢瓏鼢瓏桫桫4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(16)符號(hào)運(yùn)算求解s

43、yms z; G=(5*z-2)/(z-1/2)3/(z-1/3);iztrans(G*z/(z-1)解析解更直觀，不建議用前者求解，而直接采用Z反變換的符號(hào)運(yùn)算方法求解21111108726036123222nnnnnn驏驏驏驏鼢鼢瓏瓏-+-+-鼢鼢瓏瓏鼢鼢瓏瓏桫桫桫桫4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-基于部分分式展開方法求解(17)時(shí)間延遲系統(tǒng)的解析解法連續(xù)系統(tǒng)模型 : 求解G(s)的解析解，用t-L代替t即可離散系統(tǒng)傳遞函數(shù) : 求解H(z)的解析解，用n-k替代n即可例4-16無延遲解析解：有延遲解析解： ( )LsG s e-( )kH z z-55352( )(1/2) (

44、1/3)zG z zzzz-=-211( )108( 126072)3632nny nnn驏驏鼢瓏= -+-+鼢瓏鼢瓏桫桫55255211( )108 12(5)60(5)72)36 1(5)3211108( 126072)36 1(5)32nnnny nnnnnnn-驏驏鼢瓏= -+-+-鼢瓏鼢瓏桫桫驏驏鼢瓏= -+-+-鼢瓏鼢瓏桫桫4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-二階系統(tǒng)的階躍響應(yīng)及指標(biāo)(1)二階系統(tǒng)開環(huán)模型單位負(fù)反饋的閉環(huán)控制系統(tǒng)模型階躍響應(yīng)的解析解其中2nn( )( +2)oG ss swzw=2n22nn( )+2oG ssswzww=+ndnd()()2ndndnd( )

45、12tteey tzwwzwwwwzwwzww-+-驏=+-+-桫2dn1wz w=-4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-二階系統(tǒng)的階躍響應(yīng)及指標(biāo)(2)根據(jù) 的不同取值，分為如下幾種情況：u無阻尼振蕩 ，u欠阻尼振蕩 ，u臨界阻尼振蕩u過阻尼振蕩z0z=n( )1cos()y ttw=-01z22(1)(1)n222( )12111tteey tzzzzwzzzzz-+-驏=-+-桫4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-二階系統(tǒng)的階躍響應(yīng)及指標(biāo)(3)二階系統(tǒng)階躍響應(yīng)曲線wn=1; yy=;t=0:0.1:12;zet=0:0.1:1,2,3,5;for z=zet if z=0, y

46、=1-cos(wn*t); elseif (z0 & z1, dd=sqrt(z2-1); lam1=-z-dd; lam2=-z+dd; y=1-0.5*wn*(exp(lam1*t)/lam1-exp(lam2*t)/lam2)/dd; end yy=yy;y;endplot(t,yy)4.2 4.2 線性系統(tǒng)時(shí)域響應(yīng)解析解法-二階系統(tǒng)的階躍響應(yīng)及指標(biāo)(4)三維曲面繪制：i=find(zetG=tf(10 20,10 23 26 23 10,ioDelay,1); step(G,30);利用MATLAB 提供的功能，可以從曲線上得到更多的信息，如超調(diào)量等4321020( )1023

47、262310ssG sessss-+=+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖響應(yīng)(3)MATLAB求解解析解num=G.num1; den=G.den1; r,p=pfrac(num,den,0); r,p得到解析解數(shù)值解精度比較：y,t=step(tf(num,den); y0=r(1)*exp(p(1)*t)+r(2)*t.exp(p(2)*t)+r(5)*exp(p(5)*t)+ r(3)*exp(real(p(3)*t).*sin(imag(p(3)*t+r(4); norm(y-y0) 0.15( )1.76470.58821.248sin(0.98870

48、.1897)2ttty teteet-= -+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖響應(yīng)(4)例4-18 離散化采樣周期 T=0.01，0.1，0.5，1.2 G=tf(1,1 0.2 1,ioDelay,1); G1=c2d(G,0.01,zoh);G2=c2d(G,0.1); G3=c2d(G,0.5);G4=c2d(G,1.2); step(G,-,G2,-,G3,:,G4,-.,10)對(duì)得出的曲線進(jìn)行比較，發(fā)現(xiàn)若采樣周期選擇過大，離散化后可能丟失原來系統(tǒng)的信息21( )0.21sG sess-=+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖

49、響應(yīng)(5)例4-19 求多變量系統(tǒng)階躍響應(yīng)MATLAB實(shí)現(xiàn)g11=tf(0.1134,1.78 4.48 1,ioDelay,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,ioDelay,1.29);G=g11, g12;g21, g22; step(G)0.7220.31.2920.11340.9241.784.4812.071( )0.33780.3180.3611.0912.931sssesssG seesss-輊犏犏+犏=犏-犏犏+臌4.3 4.3 線性系統(tǒng)

50、的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖響應(yīng)(6)階躍響應(yīng)曲線是在兩路輸入單獨(dú)作用下得出的存在系統(tǒng)耦合4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖響應(yīng)(7)系統(tǒng)解耦引入前置補(bǔ)償矩陣 對(duì)延遲環(huán)節(jié)引入Pad近似 n1,d1=paderm(0.72,0,2); g11.ioDelay=0; g11=tf(n1,d1)*g11; n1,d1=paderm(0.3,0,2); g21.ioDelay=0; g21=tf(n1,d1)*g21; n1,d1=paderm(1.29,0,2); g22.ioDelay=0; g22=tf(n1,d1)*g22; G=g11,g12;g21

51、,g22;補(bǔ)償后系統(tǒng)的模型：Kp=0.1134,0.924;0.3378,-0.318; step(G*Kp);p0.11340.9240.33780.318K輊犏=犏-臌4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-線性系統(tǒng)的階躍響應(yīng)與脈沖響應(yīng)(8)解耦效果較好！若要使得多變量系統(tǒng)能直接設(shè)計(jì)，在設(shè)計(jì)前必須解耦系統(tǒng)的脈沖響應(yīng)曲線impluse()函數(shù)與step()函數(shù)調(diào)用結(jié)構(gòu)一致G=tf(10 20,10 23 26 23 10,ioDelay,1); impulse(G,30);4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-任意輸入下的系統(tǒng)響應(yīng)(1)可以利用step()和impulse()函數(shù)求解輸出信

52、號(hào)計(jì)算：Y(s)=G(s)R(s)如R(s)已知，則可以直接求解例4-20 斜坡響應(yīng)斜坡信號(hào)的Laplace變換為1/s2，故系統(tǒng)響應(yīng)可由兩種方法求出：1. G(s)/s系統(tǒng)的階躍響應(yīng)2. G(s)/s2系統(tǒng)的脈沖響應(yīng)4321020( )1023262310ssG sessss-+=+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-任意輸入下的系統(tǒng)響應(yīng)(2)MATLAB求解：G=tf(10 20,10 23 26 23 10,ioDelay,1); s=tf(s);step(G/s); %或者 impulse(G/s2)其他輸入的響應(yīng)可以由lsim()函數(shù)求?。簂sim(G,u,t)例4-21 多變量

53、系統(tǒng)輸入信號(hào)0.7220.31.2920.11340.9241.784.4812.071( )0.33780.3180.3611.0912.931sssesssG seesss-輊犏犏+犏=犏-犏犏+臌12( )1sin(31),( )sin( )cos(2)tu tetu ttt-=-+=+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-任意輸入下的系統(tǒng)響應(yīng)(3)MATLAB求解g11=tf(0.1134,1.78 4.48 1, ioDelay, 0.72);g12=tf(0.924,2.07 1);g21=tf(0.3378,0.361 1.09 1, ioDelay, 0.3);g22=tf(

54、-0.318,2.93 1, ioDelay, 1.29);G=g11,g12;g21,g22; t=0:0.1:15;u=1-exp(-t).*sin(3*t+1),sin(t).*cos(t+2);lsim(G,u,t);小節(jié)：多變量系統(tǒng)的時(shí)域響應(yīng)可以這樣求解比較容易理解曲線含義兩個(gè)信號(hào)共同作用下系統(tǒng)的時(shí)域響應(yīng)，與多變量系統(tǒng)階躍響應(yīng)的概念不同4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-降階模型的時(shí)域分析及比較(1)前面介紹了降階方法，但未比較效果例4-22MATLAB求解G=tf(1,7,11,5,1,7,21,37,30); Gr=pademod(G,1,2); Gr1=opt_app(G,

55、1,2,0) step(G,Gr,Gr1)通過階躍響應(yīng)曲線判斷降階模型的效果，可發(fā)現(xiàn)Pad近似算法的結(jié)果不是很理想，而次最優(yōu)降階算法效果較好。324327115( )7213730sssG sssss+=+4.3 4.3 線性系統(tǒng)的數(shù)字仿真分析-降階模型的時(shí)域分析及比較(2)例4-23 非最小相位系統(tǒng)MATLAB求解G=tf(10 -60 110 60,1 17 82 130 100); Gr=opt_app(G,1,2,1); Gr1=opt_app(G,1,2,0); %獲得帶有時(shí)間延遲和不帶時(shí)間延遲的最優(yōu)降階模型 step(G,Gr,Gr1)通過階躍響應(yīng)曲線判斷降階模型的效果324321

56、06011060( )1782130100sssG sssss-+=+4.4 根軌跡分析(1)單變量開環(huán)傳遞函數(shù)G(S)控制器增益K單位負(fù)反饋閉環(huán)系統(tǒng)特征方程1+kG(s)=0對(duì)k的不同取值，則可能繪制出每個(gè)特征根變化的曲線，這樣的曲線稱為系統(tǒng)的根軌跡。根軌跡用開環(huán)信息研究閉環(huán)特性4.4 根軌跡分析(2)MATLAB求解rlocus(G) rlocus(G,K) R,K=rlocus(G) rlocus(G1,-,G2,-.b,G3,:r)該函數(shù)可以用于單變量不含有時(shí)間延遲的連續(xù)、離散系統(tǒng)的根軌跡繪制，也可以用于帶有時(shí)間延遲的單變量離散系統(tǒng)的根軌跡繪制4.4 根軌跡分析(3)例4-24 開環(huán)系

57、統(tǒng)MATLAB求解num=1 4 8;den=1,18,120.3,357.5,478.5,306; G=tf(num,den); rlocus(G)兩個(gè)問題：1. 如何求解臨界增益？2. 閉環(huán)系統(tǒng)穩(wěn)定性如何變換2543248( )18120.3357.5478.5306ssG ssssss+=+4.4 根軌跡分析(4)例4-25根軌跡求解：s=tf(s); G=10/(s*(s+3)*(s2+3*s+4); rlocus(G), grid點(diǎn)擊根軌跡圖，求出阻尼比在=0.707處的增益求出臨界增益下的階躍響應(yīng)曲線：k=0.524;step(feedback(G*k,1)210( )(3)(34

58、)G ss sss=+浙江大學(xué)電氣工程學(xué)院系統(tǒng)系 包哲靜4.4 根軌跡分析(5)4.4 根軌跡分析(6)例4-26 繪制離散系統(tǒng)(T=0.1秒)根軌跡z=tf(z,Ts,0.1);G=-0.95*(z+0.51)*(z+0.68)*.(z+1.3)* (z2-0.84*z+0.196)/.(z+0.66)*(z+0.96)*.(z2-0.52*z+0.1117)*.(z2+1.36*z+0.7328);rlocus(G),grid2220.95(0.51)(0.68)(1.3)(0.840.196)( )(0.66)(0.96)(0.520.1117)(1.360.7328)zzzzzG zz

59、zzzzz-+-+=+-+4.4 根軌跡分析(7)例4-27 離散系統(tǒng)模型T=0.1秒z=tf(z,Ts,0.1);G=0.52*(z-0.49)*(z2+.1.28*z+0.4385)/.(z-0.78)*(z+0.29)*.(z2+0.7*z+0.1586);rlocus(G)根據(jù)根軌跡曲線求取系統(tǒng)的臨界增益220.52(0.49)(1.280.4385)( )(0.78)(0.29)(0.70.1586)zzzG zzzzz-+=-+4.4 根軌跡分析(8)系統(tǒng)假設(shè)有6步的純滯后，則z=tf(z,Ts,0.1); G=0.52*(z-0.49)*(z2+1.28*z+0.4385)/(z

60、-0.78)*(z+0.29)*(z2+0.7*z+0.1586); G.ioDelay=6; rlocus(G)根據(jù)根軌跡曲線求取系統(tǒng)的臨界增益，發(fā)現(xiàn)延遲后系統(tǒng)的臨界增益減小。4.4 根軌跡分析(9)例4-28 延遲狀態(tài)方程考慮無延遲系統(tǒng)根軌跡：A=-0.99,1.16,1.76,-0.16;-2.03, .-2.3,2.9,-2.45;-0.48,-3.96, .-2.05,-0.91;-0.43,1.23,2.26,-1.2; B=-1.3;-0.73;-0.57;0.62; C=-1.34,0.13,-1.11,0; G=ss(A,B,C,0); rlocus(G)0.991.161.732.