現(xiàn)代控制理論MATLAB上機(jī)實驗指導(dǎo)

1、精選優(yōu)質(zhì)文檔-傾情為你奉上現(xiàn)代控制理論MATLAB實踐指導(dǎo)書劉紅軍 編寫 現(xiàn)代控制理論MATLAB實踐指導(dǎo)書1 MATLAB概述MATLAB是MATrix LABoratory的縮寫，早期主要用于現(xiàn)代控制中復(fù)雜的矩陣、向量的各種運(yùn)算。由于 MATLAB提供了強(qiáng)大的矩陣處理和繪圖功能，很多專家因此在自己擅長的領(lǐng)域用它編寫了許多專門的MATLAB工具包（toolbox），如控制系統(tǒng)工具包（control systems toolbox）；系統(tǒng)辨識工具包（system identification toolbox）；信號處理工具包（signal processing toolbox）；魯棒控制工具包

2、（robust control toolbox）；最優(yōu)化工具包(optimization toolbox）等等。由于MATLAB功能的不斷擴(kuò)展，所以現(xiàn)在的MATLAB已不僅僅局限與現(xiàn)代控制系統(tǒng)分析和綜合應(yīng)用，它已是一種包羅眾多學(xué)科的功能強(qiáng)大的“技術(shù)計算語言（The Language of Technical Computing）”。MathWorks 公司于 1992 年推出了具有劃時代意義的 MATLAB 4.0 版本，并推出了交互式模型輸入與仿真系統(tǒng) SIMULINK，它使得控制系統(tǒng)的仿真與 CAD應(yīng)用更加方便、快捷，用戶可以方便地在計算機(jī)上建模和仿真實驗。1997 年MathWorks推

3、出的 MATLAB 5.0 版允許了更多的數(shù)據(jù)結(jié)構(gòu)， 1999 年初推出的 MATLAB 5.3 版在很多方面又進(jìn)一步改進(jìn)了 MATLAB 語言的功能。2000 年底推出的 MATLAB 6.0。最新版本是MATLAB7.0。MATLAB以矩陣作為基本編程單元，它提供了各種矩陣的運(yùn)算與操作，并有較強(qiáng)的繪圖功能。MATLAB集科學(xué)計算、圖像處理、聲音處理于一身，是一個高度的集成系統(tǒng)，有良好的用戶界面，并有良好的幫助功能。MATLAB不僅流行于控制界，在機(jī)械工程、生物工程、語音處理、圖像處理、信號分析、計算機(jī)技術(shù)等各行各業(yè)中都有極廣泛的應(yīng)用。2 如何獲得MATLAB幫助在MATLAB主窗口中鍵入h

4、elp，即可獲得第一層幫助：help 加重型內(nèi)容為用戶鍵入的內(nèi)容，其它為執(zhí)行后顯示的內(nèi)容。HELP topics:matlabgeneral - General purpose commands.matlabops - Operators and special characters.matlablang - Programming language constructs.matlabelmat - Elementary matrices and matrix manipulation.matlabelfun - Elementary math functions.matlabspecfun

5、- Specialized math functions.matlabmatfun - Matrix functions - numerical linear algebra.simulinksimulink - Simulink simulinkblocks - Simulink block library.simulinksimdemos - Simulink 3 demonstrations and samples.simulinkdee - Differential Equation EditorMATLABR11work - (No table of contents file)to

6、olboxlocal - Preferences.如果用戶對MATLAB的語言結(jié)構(gòu)lang感興趣，想進(jìn)一步了解，則鍵入：help langProgramming language constructs.Control flow.if - Conditionally execute statements.else - IF statement condition.elseif - IF statement condition.end - Terminate scope of FOR, WHILE, SWITCH, TRY and IF statements.for - Repeat statem

7、ents a specific number of times.while - Repeat statements an indefinite number of times. 如果想進(jìn)一步了解for語句，則鍵入：help forFOR Repeat statements a specific number of times.The general form of a FOR statement is:FOR variable = expr, statement, ., statement ENDThe columns of the expression are stored one at a

8、 time in the variable and then the following statements, up to the END, are executed. Some examples (assume N has already been assigned a value).FOR I = 1:N,FOR J = 1:N,A(I,J) = 1/(I+J-1);ENDEND同樣，如果想了解MATLAB中有關(guān)矩陣的操作運(yùn)算函數(shù)，可以鍵入：help matfunMatrix functions - numerical linear algebra.Matrix analysis.nor

9、m - Matrix or vector norm.normest - Estimate the matrix 2-norm.rank - Matrix rank.det - Determinant. trace - Sum of diagonal elements.null -Null space.orth - Orthogonalization.rref - Reduced row echelon form.subspace - Angle between two subspaces.Eigenvalues and singular values.eig - Eigenvalues and

10、 eigenvectors.svd - Singular value decomposition.gsvd - Generalized ingular value decomposition.eigs - A few eigenvalues.svds - A few singular values.poly - Characteristic polynomial.polyeig - Polynomial eigenvalue problem.condeig - Condition number with respect to eigenvalues.hess - Hessenberg form

11、.qz - QZ factorization for generalized eigenvalues.schur - Schur decomposition.Matrix functions.expm - Matrix exponential.logm - Matrix logarithm.sqrtm - Matrix square root.funm - Evaluate general matrix function.上面所列的都是有關(guān)矩陣的操作函數(shù)。如eig(A)可求出A的特征根及其特征向量，具體執(zhí)行方法為：A=0 1;-6 -5 輸入A矩陣A = 01-6-5E=eig(A) %求出方

12、陣A的特征根EE =-2-3V,D=eig(A) %求出方陣A的特征向量V及其A的對角型DV =0.4472-0.3162-0.89440.9487D =-200-33 MATLAB基本功能我們下面給出一些MATLAB的常用的功能，不過這只是MATLAB及其眾多TOOLBOX中的極少極少部分。用戶可以參閱有關(guān)MATLAB的手冊，或直接在MATLAB系統(tǒng)中用HELP命令查閱其它功能。31 MATLAB的主要線性代數(shù)運(yùn)算 如表A1所示為常用的矩陣和線性代數(shù)運(yùn)算函數(shù)，用戶可以用help matfun獲得更多內(nèi)容。表A1常用線性代數(shù)函數(shù)B=A矩陣轉(zhuǎn)置C=A+B矩陣相加C=A*B矩陣相乘C=Ak矩陣冪C

13、=A.*B矩陣點乘，即兩維數(shù)相同的矩陣各對應(yīng)元素相乘expm(A)指數(shù)矩陣，也就是inv(A)矩陣的逆矩陣det(A)矩陣的行列式的值rank(A)計算矩陣的秩eig(A)矩陣的特征值X,D=eig(A)矩陣的特征向量X和以特征值為元素的對角陣Dp=poly(A)矩陣的特征多項式r=roots(p)特征多項式方程的根conv(p1,p2)兩多項式相乘32常用的控制系統(tǒng)處理函數(shù)（1） TF2SS 將傳遞函數(shù)轉(zhuǎn)換到狀態(tài)空間表達(dá)式A,B,C,D = TF2SS(NUM,DEN) 將系統(tǒng)：轉(zhuǎn)換成：其中:NUM=bm,bm-1,b1,b0，DEN=1,an-1,an-2,a1,a0 （2）ZP2SS將零

14、極點型傳遞函數(shù)轉(zhuǎn)換到狀態(tài)空間表達(dá)式A,B,C,D = ZP2SS(Z,P,K) 除了以外，其它與TF2SS相同。（3） SS2TF 將狀態(tài)空間表達(dá)式轉(zhuǎn)換到傳遞函數(shù)NUM,DEN = SS2TF(A,B,C,D,iu) 即求第iu個輸入信號對輸出y(t)的傳遞函數(shù)，即：SS2TF的調(diào)用返回值為G(s)的分子多項式的系數(shù)矩陣NUM和分母多項式的系數(shù)向量DEN。（4） SS2ZP 將狀態(tài)空間表達(dá)式到零極點形式的傳遞函數(shù)的轉(zhuǎn)換Z,P,K = SS2ZP(A,B,C,D,iu) （5 ）TF2ZP 一般傳遞函數(shù)轉(zhuǎn)換到零極點型傳遞函數(shù)Z,P,K = TF2ZP（NUM,DEN） （6） ZP2TF 零極點

15、型傳遞函數(shù)轉(zhuǎn)換到一般傳遞函數(shù)NUM,DEN = ZP2TF(Z,P,K) （7） SS2SS 狀態(tài)空間表達(dá)式的線性變換A1,B1,C1,D1=ss2ss(A,B,C,D,T)其中T為變換矩陣。注意變換方程為：X1=TX，而不是常見的X=TX1。所以要與用戶習(xí)慣的變換方程一致，我們必須用T的逆代入上式，即：A1,B1,C1,D1=ss2ss(A,B,C,D,inv(T)（8）CANON 求狀態(tài)空間表達(dá)式的對角標(biāo)準(zhǔn)型As,Bs,Cs,Ds,Ts=canon(A,B,C,D,'mod')其中Ts為變換矩陣，注意變換方程為：Xs=TsX。 （9） CTRB 計算系統(tǒng)的可控判別矩陣MM=

16、ctrb(A,B)（10） OBSV 計算系統(tǒng)的可觀判別矩陣NN=obsv(A,C)（11） IMPULSE 求系統(tǒng)的單位脈沖響應(yīng)y,x=impulse(A,B,C,D,in,t)y,x=impulse(num,den,t)（12） STEP 求系統(tǒng)的單位階躍響應(yīng)y,x=step(A,B,C,D,in,t)y,x=step(num,den,t)（13） LSIM 求系統(tǒng)對任意輸入函數(shù)u(t)的響應(yīng)y,x=lsim(A,B,C,D,u,t)y,x=lsim(num,den,u,t)（14） C2D 連續(xù)系統(tǒng)狀態(tài)方程轉(zhuǎn)換為離散狀態(tài)方程,T為采樣周期G,H=c2d(A,B,T)相關(guān)的函數(shù)還有D2C,

17、D2D。（15） LYAP 求解如下形式的李雅普諾夫方程：所以求解用戶習(xí)慣的李雅普諾夫方程，我們必須用A的轉(zhuǎn)置A代入，即：P=Lyap（A，Q）。（16）PLACE 極點配置F=PLACE（A，B，P）（17） PLOT 繪圖函數(shù)plot(X,Y,str)可以用不同顏色、不同符號繪制曲線，其中str可以是下列三組選項的任意組合。y -yellow . -point - -solidm -magentao-circle: -dottedc -cyanx-x-mark-. -dashdot r -red + -plus- -dashed g -green* -starb -blues -squar

18、ew -white d -diamondk -black v -triangle (down) -triangle (up)< -triangle (left)> -triangle (right)p -pentagramh -hexagram4例題分析例1 給定某控制系統(tǒng)的狀態(tài)空間描述為：試求對角規(guī)范型和變換矩陣W，并根據(jù)其對角規(guī)范型繪制系統(tǒng)在初值，時的響應(yīng)曲線。解：close all;%狀態(tài)方程模型A=0 1 -1; -6 -11 6;-6 -11 5; B=0;0;1; C=1 0 0;%求取對角規(guī)范型W,lamda=eig(A); L=inv(W)*A*W; b=inv(W

19、)*B; c=C*W;%顯示結(jié)果disp('The Diagonal Canonical From of System is:'); L b c disp('Transformation Matrix is:');W%仿真數(shù)據(jù)初始化t=0:0.01:3; x0=5;10;15; xx0=inv(W)*x0; n=length(t);x=zeros(3,n); xx=zeros(3,n);%求解狀態(tài)變量for i=1:n xx(:,i)=expm(L*t(i)*xx0; x(:,i)=W*xx(:,i);end%繪圖plot(t,x)title('Syst

20、em Time Response with u=0');xlabel('Time-sec'); ylabel('Response-value');text(0.8,3.7,'x(1)'); text(0.8,1.5,'x(3)'); text(0.8,-0.4,'x(2)');>> The Diagonal Canonical From of System is:L = -1.0000 -0.0000 -0.0000 0.0000 -2.0000 -0.0000 -0.0000 0.0000 -

21、3.0000b = -2.8284 -13.7477 10.8628c = 0.7071 -0.2182 -0.0921Transformation Matrix is:W = 0.7071 -0.2182 -0.0921 0.0000 -0.4364 -0.5523 0.7071 -0.8729 -0.8285分析：通過對角規(guī)范型求借系統(tǒng)狀態(tài)方程的優(yōu)點是，不必計算矩陣指數(shù)，只需計算單個的標(biāo)量的指數(shù)即可。在本例具體的求解的過程中，需要注意的是求得對角規(guī)范型的狀態(tài)xx的時間響應(yīng)之后，還要用變換矩陣W將其轉(zhuǎn)化為原狀態(tài)變量x，然后畫圖。從系統(tǒng)的響應(yīng)曲線來看，各狀態(tài)變量最終都趨于零。這一方面是因為輸入

22、量u=0，另一方面是因為初值不為零，但矩陣A的特征值均為零，系統(tǒng)的零輸入響應(yīng)呈衰減趨勢。但狀態(tài)變量x（3）下降幅度和速度最大，x（1）下降幅度和速度最小。這是因為狀態(tài)變量x（3）對應(yīng)的特征值是-3，而x（1）對應(yīng)的特征值是-1，其衰減速率不同。 例2 給定線性定常系統(tǒng)的狀態(tài)空間描述為：試判斷系統(tǒng)的能控性。如果系統(tǒng)狀態(tài)完全能控，試求其能控標(biāo)準(zhǔn)型和變換矩陣。解：%系統(tǒng)狀態(tài)方程模型A=2 0 0 1;0 4 1 3;0 0 4 1;0 0 0 2;B=1;0;1;2;C=1 1 0 0;%系統(tǒng)階次n=length(A);%求解系統(tǒng)能控性矩陣Q=ctrb(A,B);%系統(tǒng)能控性矩陣的秩m=rank(Q

23、);%判斷系統(tǒng)是否狀態(tài)完全能控，并求解能控標(biāo)準(zhǔn)型if m=n Ac1=inv(Q)*A*Q; Bc1=inv(Q)*B; Cc1=C*Q; disp('System is Controllable.'); disp('System First Controllable Canonnical Form is:'); Ac1 Bc1 Cc1 disp('The Transformation Matrix is:'); Qelse disp('System State Variable cannot be totally Controlled.

24、'); disp('The rank of System Controllable Matrix is:'); mend運(yùn)行結(jié)果：System is Controllable.System First Controllable Canonnical Form is:Ac1 = 0 -0.0000 -0.0000 -64.0000 1.0000 -0.0000 -0.0000 96.0000 0 1.0000 -0.0000 -52.0000 0.0000 0.0000 1.0000 12.0000Bc1 =1 0 0 0Cc1 = 1 11 58 268The Tran

25、sformation Matrix is:Q = 1 4 12 32 0 7 46 236 1 6 28 120 2 4 8 16例3 對于上例的線性定常系統(tǒng)，試判斷系統(tǒng)的能觀性，如果系統(tǒng)狀態(tài)完全能觀的，試求其能觀標(biāo)準(zhǔn)型和變換矩陣。解：%系統(tǒng)狀態(tài)方程模型A=2 0 0 1;0 4 1 3;0 0 4 1;0 0 0 2;B=1;0;1;2;C=1 1 0 0;%系統(tǒng)階次n=length(A);P=obsv(A,C);%系統(tǒng)能觀性矩陣的秩m=rank(P);if m=n Ao1=P*A*inv(P); Bo1=P*B; Co1=C*inv(P); disp('System is Obse

26、rvable.'); disp('System First Observable Canonnical Form is:'); Ao1 Bo1 Co1 disp('The Transformation Matrix is:'); inv(P)else disp('System State Variable cannot be totally Observed.'); disp('The rank of System Observable Matrix is:'); mend運(yùn)行結(jié)果：ystem is Observable.

27、System First Observable Canonnical Form is:Ao1 = 0 1 0 0 0 0 1 0 0 0 0 1 -64 96 -52 12Bo1 =1 11 58 268Co1 = 1 0 0 0The Transformation Matrix is:ans = -14.0000 16.0000 -5.3750 0.5625 15.0000 -16.0000 5.3750 -0.5625 0 1.0000 -0.7500 0.1250 -8.0000 8.0000 -2.5000 0.2500例4某控制系統(tǒng)的狀態(tài)方程描述如下。試判斷其穩(wěn)定性并繪制其時間響應(yīng)來

28、驗證上述判斷。解：z =-1p =-1.4737 + 2.2638i -1.4737 - 2.2638i -0.0263 + 0.7399i -0.0263 - 0.7399ik =1System is stable分析：系統(tǒng)的極點均在復(fù)平面的左半平面，因而系統(tǒng)是穩(wěn)定的。從系統(tǒng)的時間響應(yīng)曲線上看也是如此，系統(tǒng)的輸出振蕩衰減，但可以看出，系統(tǒng)衰減響應(yīng)曲線的振蕩頻率很高，衰減速率很慢。這是因為系統(tǒng)有兩個極點的實部為-0.0263，非常靠近虛軸，在實際過程中，這種系統(tǒng)屬于性能非常差的一類。穩(wěn)定裕量很小，稍有干擾系統(tǒng)的極點就可能越過虛軸跑到復(fù)平面的右半面去。從而造成系統(tǒng)的不穩(wěn)定。對于這類系統(tǒng)，一般都要

29、加校正手段，使其極點遠(yuǎn)里離虛軸。例5含積分環(huán)節(jié)的伺服系統(tǒng)設(shè)計，設(shè)對象為設(shè)計控制器，使閉環(huán)系統(tǒng)具有極點。解：%Pole placementclose all;a=0 1 0;0 0 1;0 -2 -3;b=0;0;1;c=1 0 0;d=0;disp('The rank of controllability matrix')rc=rank(ctrb(a,b)p=-2+2*sqrt(3)*j -2-2*sqrt(3)*j -10;k=place(a,b,p)a1=a-b*k;b1=b*k(1);c1=c;d1=d;figure(1)step(a1,b1,c1,d1)title(&#

30、39;Step Responese of Designed Servo System')figure(2)y,x,t=step(a1,b1,c1,d1);plot(t,x)title('Step Responese Curves for x1,x2,x3')grid onThe rank of controllability matrixrc = 3k = 160.0000 54.0000 11.0000例6開環(huán)系統(tǒng)的狀態(tài)方程為設(shè)計全維狀態(tài)觀測器，使觀測器的閉環(huán)極點為。解：狀態(tài)觀測器的結(jié)構(gòu)見教材圖5-5，為求出狀態(tài)觀測器的反饋增益陣K，先為原系統(tǒng)構(gòu)造一對偶系統(tǒng)然后采用極點配置方法對對偶系統(tǒng)進(jìn)行閉環(huán)極點培植，得到反饋增益陣k，從而可由對偶原理得到原系統(tǒng)的狀態(tài)觀測器增益陣。程序如下：%Design of a full-order state observera=0 1 0;0 0 1;-6 -11 -6;b=0;0;1;c=1 0 0;%Check observabilitydisp('The rank of observability matrix')r0=rank(obsv(a,c)%Construct the daul s