機(jī)電系統(tǒng)動態(tài)性能的計算機(jī)仿真

1、4機(jī)電系統(tǒng)動態(tài)性能的計算機(jī)仿真41 概述 機(jī)電系統(tǒng)計算機(jī)仿真是目前對復(fù)雜機(jī)電系統(tǒng)進(jìn)行分析的重要手段與方法。在進(jìn)行機(jī)電系統(tǒng)分析綜合與設(shè)計工作過程中，除了需要進(jìn)行理論分析外，還要對系統(tǒng)的特性進(jìn)行實驗研究。系統(tǒng)性能指標(biāo)與參數(shù)是否達(dá)到預(yù)期的要求？它的經(jīng)濟(jì)性能如何？這些都需要在系統(tǒng)設(shè)計中給出明確的結(jié)論。對于那些在實際調(diào)試過程中存在很大風(fēng)險或?qū)嶒炠M用昂貴的系統(tǒng)，一般不允許對設(shè)計好的系統(tǒng)直接進(jìn)行實驗，然而沒有經(jīng)過實驗研究是不能將設(shè)計好的系統(tǒng)直接放到生產(chǎn)實際中去的，因此就必須對其進(jìn)行模擬實驗研究。當(dāng)然在有些情況下可以構(gòu)造一套物理模擬裝置來進(jìn)行實驗，但這種方法十分費時而且費用又高，而在有的情況下物理模擬幾乎是不

2、可能的。近年來隨著計算機(jī)的迅速發(fā)展，采用計算機(jī)對機(jī)電系統(tǒng)進(jìn)行數(shù)學(xué)仿真的方法已被人們采納。所謂機(jī)電系統(tǒng)計算機(jī)仿真就是以機(jī)電系統(tǒng)的數(shù)學(xué)模型為基礎(chǔ)，借助計算機(jī)對機(jī)電系統(tǒng)的動靜態(tài)過程進(jìn)行實驗研究。這里講的機(jī)電系統(tǒng)計算機(jī)仿真是指借助數(shù)字計算機(jī)實現(xiàn)對機(jī)電系統(tǒng)的仿真分析。這種實驗研究的特點是：將實際系統(tǒng)的運(yùn)動規(guī)律用數(shù)學(xué)表達(dá)式加以描述，它通常是一組常微分方程或差分方程，然后利用計算機(jī)來求解這一數(shù)學(xué)模型，以達(dá)到對系統(tǒng)進(jìn)行分析研究的目的。 對機(jī)電系統(tǒng)進(jìn)行計算機(jī)仿真的基本過程包括：首先建立系統(tǒng)的數(shù)學(xué)模型，因為數(shù)學(xué)模型是系統(tǒng)仿真的基本依據(jù)，所以數(shù)學(xué)模型極為重要。然后根據(jù)系統(tǒng)的數(shù)學(xué)模型建立相應(yīng)的仿真模型，一般需要通過一

3、定的算法或數(shù)值積分方法對原系統(tǒng)的數(shù)學(xué)模型進(jìn)行離散化處理，從而建立起相應(yīng)的仿真模型，這是進(jìn)行機(jī)電系統(tǒng)仿真分析的關(guān)鍵步驟；最后根據(jù)系統(tǒng)的仿真模型編制相應(yīng)的仿真程序，在計算機(jī)上進(jìn)行仿真實驗研究并對仿真結(jié)果加以分析。 機(jī)電系統(tǒng)計算機(jī)仿真的應(yīng)用與發(fā)展已經(jīng)過了近40年的歷程，進(jìn)入20世紀(jì)80年代以來，隨著微型計算機(jī)技術(shù)以及軟件技術(shù)的飛速發(fā)展與廣泛應(yīng)用，使得機(jī)電系統(tǒng)計算機(jī)仿真獲得了實質(zhì)性的發(fā)展，并使其走進(jìn)廣大的機(jī)電系統(tǒng)生產(chǎn)、設(shè)計、研究的第一線。目前在該領(lǐng)域應(yīng)用最為廣泛的軟件包就是美國Mathworks公司開發(fā)的MATLAB語言軟件。雖然它最初并非是為機(jī)電系統(tǒng)仿真與設(shè)計開發(fā)的，但是它的強(qiáng)大的矩陣運(yùn)算功能與圖形

4、處理及繪制能力，以及在MATLAB平臺上開發(fā)出來的各種工具箱，和面向結(jié)構(gòu)圖的SIMULINK 系統(tǒng)分析環(huán)境，為機(jī)電系統(tǒng)計算機(jī)仿真提供了強(qiáng)有力的軟件工具，從而為機(jī)電系統(tǒng)計算機(jī)仿真與輔助設(shè)計開辟新的方法與途徑提供了充分的條件。 考慮到部分讀者不具備MATLAB基礎(chǔ)知識，因此，在本書的附錄一中對MATLAB語言基礎(chǔ)知識作了簡單的介紹。42面向微分方程的連續(xù)系統(tǒng)仿真與MATLAB實現(xiàn) 對于控制系統(tǒng)的最基本數(shù)學(xué)描述就是微分方程，因此有必要在介紹機(jī)電系統(tǒng)計算機(jī)仿真時對這種方法加以介紹。實際上在前面章節(jié)已經(jīng)介紹了幾種典型數(shù)值積分方法以及它們的MATLAB實現(xiàn)。數(shù)值積分的數(shù)學(xué)基礎(chǔ)就是微分方程（狀態(tài)方程），如果

5、一個系統(tǒng)以微分方程的形式加以描述，就可以利用前面的知識對該系統(tǒng)進(jìn)行仿真研究，有關(guān)這部分內(nèi)容本章就不再加以討論了。這里將介紹如何利用MATLAB中為用戶所提供的常微分方程解函數(shù)ode和對微分方程進(jìn)行求解計算，以及在面向微分方程的系統(tǒng)仿真研究中如何應(yīng)用M函數(shù)。421 基于ode函數(shù)的面向微分方程的系統(tǒng)仿真1. 常微分方程解函數(shù)ode的基本格式 （l）ode函數(shù)的基本分類。在 MATLAB 5.x中 ode函數(shù)可以分為兩類五種形式。一類是應(yīng)用于非剛性微分方程求解，如ode45；另一類屬于求解剛性微分方程的ode解函數(shù)，如ode23s.他們的具體形式如下：1ode45函數(shù)用于求解非剛性微分方程它屬于中

6、階龍格一庫塔法；2ode23函數(shù)用于求解非剛性微分方程，它屬于底階龍格-庫塔法；3ode113函數(shù)用于求解非剛性微分方程的變階法；4ode15s函數(shù)用于求解剛性微分方程的變階法；5ode23s函數(shù)用于求解剛性微分方程的低階法。 （2）ode函數(shù)的基本格式。上還這五種 ode函數(shù)的格式基本上一樣這里主要以最常用的ode45函數(shù)為例加以介紹，其結(jié)果多數(shù)也可用于其他ode函數(shù)。ode45的基本格式主要有以下幾種。 t，y=ode45('F'，tspan，y0) 其中：'F'是ode文件名字符串；tspan= T0 TFINAL為一個向量，它指定了積分的啟始時刻和結(jié)束時

7、刻； y0為積分運(yùn)算的初始條件。 該函數(shù)對微分方程系統(tǒng)從T0到TFINAL時間段進(jìn)行積分。函數(shù)返回一個列向量。在解向量Y中的每一行與列向量T的返回時間相對應(yīng)。為了在指定的時刻點T0，T1，TFINAL得到解，則使用tspan=T0，T1，TFINAL。 t , yode45（'F'， tspan， y0，options）該函數(shù)的解如同前一格式，只是在這里用選項options中的參數(shù)來代替前一格式中的缺省值，options中的參數(shù)變量由函數(shù)odeset予以建立。最常使用的選項參數(shù)是標(biāo)量相對容差'RelTol'（缺省值為1e-3）與絕對容差向量'AbsTol

8、'（缺省的全部向量元素均為1e-6）。 （3）odeset函數(shù)的基本格式。odeset函數(shù)用來建立或改變options選項結(jié)構(gòu)，其基本格式為 options=odeset（'namel'，valuel，'name2'，value2，）通過該odeset函數(shù)建立積分器選項options的結(jié)構(gòu)，在options中所命名的屬性具有指定值，對于任何未指定的屬性均為缺省值。 optionsodeset（oldopts，'namel'， valuel，）通過該odeset函數(shù)改變現(xiàn)存選項結(jié)構(gòu)oldopts。 options=odeset（oldopt

9、， newopts）通過該odeset函數(shù)將現(xiàn)有選項結(jié)構(gòu)oldopts與新選項結(jié)構(gòu)newopts結(jié)合在一起，新選項中的屬性將覆蓋對應(yīng)原有屬性。2機(jī)電系統(tǒng)仿真的ode函數(shù)實現(xiàn)由ode函數(shù)可知，當(dāng)已知系統(tǒng)由微分方程加以描述時，并建立以'F'作為ODE文件名，則可應(yīng)用 t，y=ode45('F'，tspan，y0)實現(xiàn)對系統(tǒng)的求解分析。這里所講的ODE文件實際就是M函數(shù)文件，因此對于以高階微分方程加以描述的系統(tǒng)，首先應(yīng)將其寫成一組一階微分方程的形式，當(dāng)然如果系統(tǒng)是以狀態(tài)方程描述的，即可直接對其編寫M函數(shù)文件。函數(shù)文件的特征就是文件第一行一定是以function開始，而

10、且在函數(shù)文件中的變量均為局部變量，同時由函數(shù)文件構(gòu)造的功能函數(shù)可以像MAT-LAB的其他函數(shù)一樣進(jìn)行調(diào)用。實際上在ode函數(shù)t，yode45（'F'，tspan，y0）中，對以'F'作為文件名的ODE文件就是這類函數(shù)文件。下面通過實例說明對系統(tǒng)仿真過程中的具體應(yīng)用。 例4-1 之已知一三階系統(tǒng)，其狀態(tài)方程表達(dá)式為 其中 ， 輸入為u= 100，試應(yīng)用ode函數(shù)求解該系統(tǒng)。 解： 由于該系統(tǒng)是以狀態(tài)方程形式給出的，因此即可直接對其編寫M函數(shù)文件。 function xd= rr（t，x） u100； xd=-14*x(1)+9*x(2)+10*x(3)+2.8*u

11、；12*x(1)-9*x(2)+10*x(3)+4*u；24*x(1)-24*x(2)-18*x(3)+12*u；該函數(shù)文件對所研究系統(tǒng)的狀態(tài)方程加以描述并形成rr.m文件以便在ode解函數(shù)中使用后即可利用ode解函數(shù)對系統(tǒng)進(jìn)行求解分析?？紤]狀態(tài)的零初值條件，即 圖4-1 系統(tǒng)狀態(tài)曲線 clear x0=0；0；0； t，x=ode45('rr'，0 3，x0)； plot(t，x)對其同樣可存一M文件，如li32.m文件，這樣即可在命令窗口執(zhí)行該文件。在運(yùn)行過程中它會自動調(diào)用相關(guān) rr函數(shù)文件，實現(xiàn)預(yù)定的工作，圖 4-1給出了相應(yīng)的運(yùn)行結(jié)果。根據(jù)ode函數(shù)用途可知，該函數(shù)應(yīng)用

12、于求解常微分方程，因此對于含有輸出方程的系統(tǒng)狀態(tài)空間表達(dá)式就不能直接利用ode函數(shù)求解系統(tǒng)輸出響應(yīng)?？紤]到ode函數(shù)返回變量的向量格式，可以根據(jù)輸出方程的具體要求，由狀態(tài)變量的線性組合直接得到所要求的系統(tǒng)輸出。利用ode函數(shù)可以方便地對以微分方程形式描述的系統(tǒng)進(jìn)行仿真分析，然而ode解函數(shù)并非是萬能的，在應(yīng)用ode函數(shù)時對其進(jìn)行干預(yù)一般是不方便的，因此在很多情況下，對以微分方程形式給出的系統(tǒng)，還必須采用前面所介紹的各種數(shù)值積分方法進(jìn)行仿真分析。422 基于 M函數(shù)的面向微分方程的系統(tǒng)仿真 M函數(shù)是MATLAB中所特有的一類函數(shù)，前面介紹ode函數(shù)時已經(jīng)使用了M函數(shù)。這里將介紹采用M函數(shù)建立描述

13、系統(tǒng)的數(shù)學(xué)模型并在此基礎(chǔ)上編制相應(yīng)的程序。根據(jù)需要調(diào)用該M文件，實現(xiàn)所要求的運(yùn)算。下面通過實例加以說明。例4-2已知系統(tǒng)狀態(tài)方程為 應(yīng)用數(shù)值積分方法（預(yù)估校正梯形法）實現(xiàn)對系統(tǒng)狀態(tài)的仿真研究。 解：首先建立描述該系統(tǒng)數(shù)學(xué)模型的M函數(shù)文件，相應(yīng)的程序如下。 function df=rrl（X，u） A= -0.5572 -0.7814； 0.7814 0； B=1；0 U= 20； df=A*X+B*u 將該函數(shù)文件存盤取作rrlm文件。然后根據(jù)預(yù)估校正梯形公式編制相應(yīng)的程序。 h=0.0l；t0=0；tf=15； X=0 0'； X1=0 0； T=(tf-t0)/h for i=l：

14、T K1=rrl（X，1）； X= X+h*K1 K2=rrl（X，1）； X= X+h*（K1+K2）/2 X1=X1；X'；endplot（X1）在該程序段中兩次調(diào)用了rr1.m函數(shù)文件，在第一次調(diào)用rrl.m后，對各狀態(tài)的計算將產(chǎn)生新值，因此K2就是經(jīng)過預(yù)估計算所得到的相關(guān)值,從而即可應(yīng)用梯形法實現(xiàn)對系統(tǒng)狀態(tài)的仿真分析。圖4-2給出了系統(tǒng)狀態(tài)響應(yīng)仿真曲線。 圖4-2 系統(tǒng)狀態(tài)響應(yīng)43 面向傳遞函數(shù)的連續(xù)系統(tǒng)仿真一個控制系統(tǒng)的動態(tài)模型經(jīng)常是一傳遞函數(shù)的形式給出，特別是當(dāng)給出開環(huán)傳遞函數(shù)， 在進(jìn)行系統(tǒng)仿真時就要將開環(huán)傳遞函數(shù)轉(zhuǎn)換為開環(huán)狀態(tài)方程，再由開環(huán)狀態(tài)方程求出對應(yīng)的閉環(huán)狀態(tài)方程，

15、然后，再應(yīng)用數(shù)值積分方法進(jìn)行分解和分析仿真。下面將就這類問題加以介紹并給出MATLAB的實現(xiàn)。431 面向傳遞函數(shù)的系統(tǒng)仿真考慮一個連續(xù)系統(tǒng)，它的開環(huán)傳遞函數(shù)的一般描述為 其相應(yīng)的狀態(tài)空間表達(dá)式為 其中A、B、C為相應(yīng)的系統(tǒng)系數(shù)陣，這些系數(shù)陣可以寫作能控標(biāo)準(zhǔn)型形式或能觀標(biāo)準(zhǔn)型形式，當(dāng)然也可以寫成其他的形式，這些都不會影響系統(tǒng)分析的結(jié)果。現(xiàn)如果假定系統(tǒng)結(jié)構(gòu)圖如圖4-3所示，對于SISO系統(tǒng)，其v、u、y 、f 均為標(biāo)量。由結(jié)構(gòu)圖可知u = v-fy，將其代入狀態(tài)空間表達(dá)式得到 圖 4-3 系統(tǒng)結(jié)構(gòu)圖式中 （4-1）輸出方程為： （4-2）（4-1）式即為系統(tǒng)閉環(huán)狀態(tài)方程。因此可以方便地利用各種數(shù)

16、值積分方法對其進(jìn)行求解，實現(xiàn)對以傳遞函數(shù)形式給出的系統(tǒng)進(jìn)行仿真分析研究。顯然，如果系統(tǒng)的數(shù)學(xué)模型是以微分方程或狀態(tài)方程形式給出的，一般就可以直接利用前述數(shù)值積分方法或ode函數(shù)方法對系統(tǒng)進(jìn)行仿真研究。432 面向傳遞函數(shù)的系統(tǒng)仿真的MATLAB實現(xiàn)下面通過一個例子說明如何利用MATLAB語言進(jìn)行以傳遞函數(shù)形式給出的系統(tǒng)的仿真分析。例4-3 設(shè)單位反饋系統(tǒng)的開環(huán)傳遞函數(shù)為 試給出系統(tǒng)閉環(huán)的階躍響應(yīng)。 解 根據(jù)開環(huán)傳遞函數(shù)求出系統(tǒng)開環(huán)狀態(tài)方程，應(yīng)用MATLAB控制工具箱中的函數(shù)，可以方便地求出系統(tǒng)開環(huán)狀態(tài)方程。 num =2*1,2; den = conv(conv(1,0,1,3),1,2,2)

17、 G1 = tf(num,den); G2 = ss(G1); a,b,c,d = ssdata (G);這樣即可以得到系統(tǒng)系數(shù)陣分別為 a = -5.0000 -2.0000 -0.7500 0 4.0000 0 0 0 0 2.0000 0 00 0 1.0000 0 b = 1 0 0 0 c = 0 0 0.2500 0.5000 d = 0進(jìn)而就可以方便地求出系統(tǒng)閉環(huán)狀態(tài)方程 Ab = a-b*c; Bb = b; Cb = c Db = 0Ab = -5.0000 -2.0000 -1.0000 -0.5000 4.0000 0 0 0 0 2.0000 0 00 0 1.0000

18、 0 Bb = 1 0 0 0 Cb = 0 0 0.2500 0.5000 Db = 0從而可以應(yīng)用數(shù)值積分方法求取系統(tǒng)的階躍響應(yīng)。下面給出了應(yīng)用MATLAB語言編寫的四階龍格-庫塔數(shù)值積分方法的仿真程序，相應(yīng)的仿真曲線由圖4-4 給出。X = 0;0;0;0; Y = 0; t = 0; U = 25; Tf = 15; h=0.01; N = Tf/h; Ab = -5 -2 -1 -0.5; 4 0 0 0; 0 2 0 0; 0 0 1 0;圖4-4 四階龍格-庫塔數(shù)值積分方法的仿真曲線Bb = 1; 0; 0;0;Cb = 0 0 0.25 0.5;Db = 0;for i = 1

19、:nK1 = Ab * X + Bb*U;K2 = Ab * (X + h*K1/2)+Bb*U;K3 = Ab * (X + h*K2/2)+Bb*U;K4 = Ab * (X + h*K3)+Bb*U;X = X+ h*(K1 + 2*K2+ 2*K3 +K4)/6;Y = Y, Cb*X;t = t, t(i)+h;endplot (t,Y) Grid當(dāng)由已知系統(tǒng)開環(huán)傳遞函數(shù)求取系統(tǒng)開環(huán)狀態(tài)方程時，除了可以利用ss函數(shù)直接求取外，還可以使用狀態(tài)方程相似變換函數(shù)ss2ss進(jìn)行求取，下面給出上例中應(yīng)用相似變換函數(shù)ss2ss求取系統(tǒng)開環(huán)狀態(tài)方程的基本用法。這里主要利用矩陣水平翻轉(zhuǎn)函數(shù)flipl

20、r構(gòu)造出反對角陣作為變換陣，從而得到相應(yīng)的狀態(tài)方程。num = 2 *1,2; den = conv(conv(1,0,1,3),1,2,2);G1= tf (num, den);G2 = ss(G1);T = fliplr (eye(4);G = ss2ss(G2,T)運(yùn)行結(jié)果： a = x1 x2 x3 x4 x10 1.000000 0x200 2.000000x30004.00000x40 -0.75000-2.00000-5.00000b =u1x1 0x20x30x4 1.00000 c = x1 x2 x3 x4 y10.500000.2500000d = u1y10從計算結(jié)果可

21、以看出，該狀態(tài)方程十分類似通常定義的系統(tǒng)可控標(biāo)準(zhǔn)型形式。 當(dāng)然如果在研究系統(tǒng)時，需要狀態(tài)方程以可控標(biāo)準(zhǔn)型形式給出，則可應(yīng)用MATLAB編程得到。但不論哪種形式的狀態(tài)方程，最終所要研究的系統(tǒng)階躍響應(yīng)都是一樣的。44面向結(jié)構(gòu)圖的連續(xù)系統(tǒng)仿真441 基于典型環(huán)節(jié)的系統(tǒng)仿真 由于習(xí)慣上往往采用結(jié)構(gòu)圖形式對系統(tǒng)進(jìn)行描述，因此在對系統(tǒng)進(jìn)行仿真時，有必要研究如何以結(jié)構(gòu)圖為基礎(chǔ)進(jìn)行。在此將介紹以典型環(huán)節(jié)為基本環(huán)節(jié)的面向結(jié)構(gòu)圖的系統(tǒng)仿真方法，這種方法的主要優(yōu)點是： 便于研究各環(huán)節(jié)參數(shù)對系統(tǒng)的影響 可以得到每個環(huán)節(jié)的動態(tài)響應(yīng)； 對于存在系統(tǒng)非線性環(huán)節(jié)較容易處理； 不但適合單人單出系統(tǒng)，而且對于多入多出系統(tǒng)也適用。

22、 該方法只要建立起各個典型環(huán)節(jié)的以及各環(huán)節(jié)之間的連接關(guān)系，就能方便地求取任意環(huán)節(jié)和系統(tǒng)的動態(tài)響應(yīng)。由于直接輸入各環(huán)節(jié)的參數(shù)，因此研究參數(shù)變化對系統(tǒng)的影響也十分方便。1典型環(huán)節(jié)的確定 控制系統(tǒng)常見的環(huán)節(jié)有：比例環(huán)節(jié) K、積分環(huán)節(jié)、慣性環(huán)節(jié)、一階超前（滯后）環(huán)節(jié)、二階振蕩環(huán)節(jié)等。為了仿真的方便，必須選擇一種能夠代表這些環(huán)節(jié)的典型環(huán)節(jié)。經(jīng)過分析可以看出，只要采用下述形式的環(huán)節(jié)作為典型環(huán)節(jié)，就可以滿足對系統(tǒng)研究的要求，即典型模型可表示為 （4-3）顯然，當(dāng)?shù)湫铜h(huán)節(jié)中時即為慣性環(huán)節(jié)；時即是積分環(huán)節(jié)；至于振蕩環(huán)節(jié)則只要用兩個典型環(huán)節(jié)串聯(lián)并加上一個負(fù)反饋即可得到。其中典型環(huán)節(jié)中的、均為對角陣，其對角線元素值

23、由各個環(huán)節(jié)中的相關(guān)參數(shù)確定。根據(jù)典型環(huán)節(jié)可以得到第i個環(huán)節(jié)的輸出為 因此得到 （i=1，2，n） 將其寫為矩陣方程形式，得到 （4-4） 其中 ， B= ，2系統(tǒng)連接矩陣 前面給出了從微觀角度對系統(tǒng)環(huán)節(jié)的描述，然而為了研究一個系統(tǒng)，還必須從宏觀角度對系統(tǒng)加以描述，這就需要將各個環(huán)節(jié)之間的連接關(guān)系用一種數(shù)學(xué)形式加以表示。考察圖4-5所示系統(tǒng)，該系統(tǒng)由4個典型環(huán)節(jié)構(gòu)成，根據(jù)系統(tǒng)結(jié)構(gòu)圖可以得到各環(huán)節(jié)的輸入輸出關(guān)系方程組。 圖4-5系統(tǒng)結(jié)構(gòu)圖 將其寫作矩陣形式 （4-5）其中 ， '， '， 一般稱與為系統(tǒng)連接矩陣。為一方陣，元素表示結(jié)構(gòu)圖中第j 個環(huán)節(jié)對第i個環(huán)節(jié)的作用；若第j個環(huán)節(jié)

24、對第i個環(huán)節(jié)沒有作用，則。表示外部輸入信號對系統(tǒng)的連接情況。3周環(huán)系統(tǒng)的狀態(tài)方程 將(4-5)式代入(4-4)式，消去U得 對上式展開，整理并取拉氏反變換得 考慮在外輸入為階躍信號時，避免對應(yīng)的初始導(dǎo)數(shù)無窮大，則要求。如果令 （4-6） （4-7）可以得到 如果的逆存在，則得到 若令 ， （4-8）則上式可寫作 （4-7）顯然上式即為一典型的狀態(tài)方程，從而可利用前面介紹的數(shù)值積分方法對其進(jìn)行求解，實現(xiàn)對系統(tǒng)的仿真研究。下面以一個結(jié)構(gòu)圖形式的系統(tǒng)為例，介紹基于MATLAB的面向結(jié)構(gòu)圖的系統(tǒng)仿真實現(xiàn)。 例4-4 已知系統(tǒng)結(jié)構(gòu)圖如圖4-5所示，采用MATLAB實現(xiàn)面向結(jié)構(gòu)圖的系統(tǒng)仿真 。 解：根據(jù)系

25、統(tǒng)結(jié)構(gòu)圖首先利用MATLAB建立各個環(huán)節(jié)參數(shù)構(gòu)成的系數(shù)陣（此處的系數(shù)陣均為對角陣）與系統(tǒng)連接陣，在此基礎(chǔ)上再依（4-6）式、（4-7）式、（4-8）式建立閉環(huán)系統(tǒng)的相應(yīng)系數(shù)陣。下面給出在MATLAB環(huán)境下的相應(yīng)程序。 AA=diag（1 6.5 1 1）； BB=diag（4 1 3 2 ）; CC=diag（4 1 2.4 5）； DD=diag（0 0 0 1）； W=0 0 -4 0;1 0 0 -1;0 1 0 0 ;0 0 1 0; W0=1 0 0 0 ； Q1BB-DD*W； Qinv(Q1）； RCC*WAA； AQ*R； BQ*CC*W0； 一旦得到了閉環(huán)系統(tǒng)的系數(shù)陣，即可應(yīng)

26、用數(shù)值積分方法實現(xiàn)對系統(tǒng)的分析仿真。應(yīng)注意，在該狀態(tài)方程中各個狀態(tài)變量的實際含義。這里應(yīng)用四階龍格一庫塔法對其進(jìn)行仿真分析，下面給出了相應(yīng)的仿真程序。圖4-6給出了采用面向結(jié)構(gòu)圖系統(tǒng)階躍響應(yīng)曲線。圖4-6面向結(jié)圖仿真（階躍響應(yīng)） 圖4-7環(huán)節(jié)參數(shù)改變后面向結(jié)構(gòu)圖的仿真 X0；0；0；0；Y0；t=0；U25； Tf=15；h001； n=Tfh； Cb=0 01 0； for i 1:n K1A*X+B*U； K2A*（Xh*K1/2） B*U； K3=A*（X h*K2/2） B*U； K4=A*（X h*K3） B*U； XX+h*(K12*K22*K3K4）/6； Y=Y，Cb*X； t

27、=t,t(i)+h; end plot(t,y) grid在這里由于采用的是面向結(jié)構(gòu)圖形式的系統(tǒng)仿真，因此研究因各環(huán)節(jié)參數(shù)發(fā)生變化對系統(tǒng)的影響就十分方便，僅對環(huán)節(jié)系數(shù)陣加以修改即可。如果將第三個環(huán)節(jié)修改為2.4/( s +1)時，只需修改B陣的相關(guān)元素即可，其他部分無需改動，圖4-7給出了改變環(huán)節(jié)三參數(shù)后的仿真結(jié)果。由響應(yīng)曲線看出，系統(tǒng)的響應(yīng)快速性明顯提高了，同時超調(diào)量也有所降低。另外采用面向結(jié)構(gòu)圖仿真，還可十分方便地研究各個環(huán)節(jié)輸出的響應(yīng)情況，只要根據(jù)需要修改相應(yīng)系數(shù)陣即可。圖4-8（a）、（b）和（c） 圖4-8 環(huán)節(jié) 1（a）、環(huán)節(jié)2（b）、環(huán)節(jié)4（c）的響應(yīng)分別給出了環(huán)節(jié)1、環(huán)節(jié)2和環(huán)

28、節(jié)4的輸出響應(yīng)。442 基于Connect連接函數(shù)的系統(tǒng)仿真 除了前面介紹的面向結(jié)構(gòu)圖的系統(tǒng)仿真方法外，還可利用MATLAB與控制工具箱所提供的相關(guān)函數(shù)與命令，實現(xiàn)以結(jié)構(gòu)圖形式給出的系統(tǒng)仿真。下面通過實例說明如何以這種方法進(jìn)行仿真工作。 例4-5 已知系統(tǒng)如圖4-9所示，研究該系統(tǒng)的階躍響應(yīng)。 解： 顯然該系統(tǒng)是以結(jié)構(gòu)圖形式給出的較復(fù)雜系統(tǒng)；利用MATLAB和控制工具箱可以方便地實現(xiàn)系統(tǒng)仿真。就基于Connect連接函數(shù)的仿真而言，一般可分為四步進(jìn)行。1建立系統(tǒng)中間狀態(tài)模型 首先應(yīng)對系統(tǒng)各環(huán)節(jié)進(jìn)行編號（圖4-9），然后即可據(jù)此與給定系統(tǒng)參數(shù)建立有關(guān)中間模型，以便運(yùn)行時在MATLAB的工作空間建

29、立相應(yīng)的系統(tǒng)環(huán)節(jié)模型。其基本語句和格式為 nblocks10； nl=1；dl= 001 1； n2=017 1；d2=0085 0; n31；d3=001 1； n4=015 1；d4=0051 0； n5=70；d5=00067 1； 圖4-9系統(tǒng)結(jié)構(gòu)圖 n6= 021；d6= 015 1； n7= 130；d7=1 0; n8=-0212；d81； n9= - 01 ；d90011 ； n10=- 00044；d10= 0 01 1； blkbuild； 其中ni、di表示第i個環(huán)節(jié)傳遞函數(shù)的分子、分母多項式的系數(shù)向量，如果所研究的系統(tǒng)是以狀態(tài)方程形式給出的各個環(huán)節(jié)，則應(yīng)采用ai、bi、

30、ci、di來表示環(huán)節(jié)狀態(tài)方程的系數(shù)陣。對于每個環(huán)節(jié)只能取二者之一，MATLAB會自動根據(jù)所輸入的陣列元素數(shù)量加以判定并進(jìn)行相應(yīng)的運(yùn)算。nblocks需根據(jù)系統(tǒng)環(huán)節(jié)個數(shù)對其進(jìn)行賦值， blkbuild語句將對nblocks語句后的各環(huán)節(jié)輸入情況進(jìn)行校驗，并自動生成系統(tǒng)的隱含中間狀態(tài)方程模型（a、b、c、d系數(shù)陣）。2建立系統(tǒng)連接矩陣 在MATLAB工作空間建立了系統(tǒng)各環(huán)節(jié)的模型后，還應(yīng)該建立反映各環(huán)節(jié)之間相互連接關(guān)系的Q陣，以及輸入、輸出向量（inputs、outputs)，下面結(jié)合實例將Q陣與inputs、outputs列寫規(guī)則加以介紹。根據(jù)已知系統(tǒng)可以得到系統(tǒng)連接陣Q和輸入、輸出向量： Q1

31、 0 0； 2 1 10； 3 2 0； 4 3 9； 5 4 0； 6 5 8； 7 6 0； 8 7 0; 9 6 0; 10 7 0; inputs=1; outputs=7;Q陣的行數(shù)與系統(tǒng)環(huán)節(jié)數(shù)相同，Q陣的第一列應(yīng)為相應(yīng)環(huán)節(jié)編號，第i行的第二列及以后的所有元素為進(jìn)入第i個環(huán)節(jié)的所有環(huán)節(jié)編號。輸入變量inputs為輸入信號進(jìn)入的環(huán)節(jié)編號；輸出變量ouputs為系統(tǒng)信號輸出所對應(yīng)的環(huán)節(jié)。如果是多輸入多輸出系統(tǒng)，輸入變量和輸出變量應(yīng)為一向量。3建立整體系統(tǒng)模型 完成了上述工作之后，就可以應(yīng)用connect函數(shù)建立整體系統(tǒng)的數(shù)學(xué)模型。其基本格式為 A， B， C，D=connect（a，b，

32、c，d，Q，inputs，outputs）或 G= connect（a， b， c， d， Q， inputs， outputs）在connect函數(shù)中的a、b 、c、d系數(shù)陣是由blkbuild函數(shù)自動生成的，其返回變量既可是系數(shù)陣，也可是系統(tǒng)陣G，不論哪種情況都將生成完整的系統(tǒng)模型，從而就可利用它們對該復(fù)雜系統(tǒng)進(jìn)行分析研究。4系統(tǒng)仿真 由于已經(jīng)求出了系統(tǒng)的系數(shù)陣，為了獲得系統(tǒng)階躍響應(yīng)，可以利用Step（A，B，C，D）求取系統(tǒng)的階躍相應(yīng)，圖4-10即為系統(tǒng)的階躍響應(yīng)曲線。當(dāng)然也可應(yīng)用數(shù)值積分方法對其進(jìn)行研究，可參考前面內(nèi)容。 圖4-10系統(tǒng)階躍響應(yīng)45 Matlab中的系統(tǒng)仿真工具SIMU

33、LINK應(yīng)用451 Simulink應(yīng)用簡介 Simulink是MATLAB下的面向結(jié)構(gòu)圖方式的仿真環(huán)境，對于結(jié)構(gòu)復(fù)雜的控制系統(tǒng)，要快速地建立系統(tǒng)模型是較為困難的，Mathworks公司所提供的這種采用圖形化方式的對系統(tǒng)建模和仿真的工具為用戶帶來了極大的方便。Simulink與用戶交互接口是基于Windows的圖形編程方法，因此非常易于為用戶所接受，使用十分靈活方便。特別是它的SFunction使其更增強(qiáng)了功能擴(kuò)展能力。目前Simulink作為MATLAB的重要組成部分，已經(jīng)成為仿真研究的重要工具。下面將簡要介紹有關(guān)Simulink的基本應(yīng)用。 在進(jìn)入MATLAB后，鍵人Simulink命令，

34、就可以打開一個名為Simulink的模型庫窗口，如圖4-11所示。圖4-11 Simulin 模型庫窗口 圖4-12 信號源Source 庫 圖4-13 輸出方式(Sinks)庫 Simulink模型庫包含功能齊全的子模型庫：Sources (信號源庫）、Sinks輸出方式庫）、Discrete(離散模型庫）、Linear（線性環(huán)節(jié)庫）、Nonlinear（非線性環(huán)節(jié)庫）、Connections（連接及接口庫）、Blocksets and toolboxes（模塊建立和工具箱庫）和 Demos（示例庫）。每個子模型庫包含很多的功能模塊，雙擊這些模型庫即可打開它們，圖4-12、圖4-13、圖4-

35、14、圖4-15分別給出了雙擊打開的信號源庫、輸出方式庫、線性環(huán)節(jié)庫以及連接模塊庫等。每個子模塊都可以從相應(yīng)的庫被鼠標(biāo)拖（即復(fù)制）到用戶建立的模型窗口中。在用戶窗口內(nèi)，可以用鼠標(biāo)對各個模塊進(jìn)行連接，并相應(yīng)設(shè)置子模塊參數(shù)、仿真算法參數(shù)，以便進(jìn)行仿真。設(shè)計好的系統(tǒng)模型應(yīng)加以存盤，在 MATLAB 5.x中 Simulink的文件名以mdl作為后綴，這一點是與 MATLAB 4.x不同的。為了說明 Simulink的使用，下面以簡單示例介紹在 Simulink環(huán)境中建立系統(tǒng)模型以及對其進(jìn)行仿真分析的方法。 例4-6已知系統(tǒng)開環(huán)傳遞函數(shù)為試分析閉環(huán)系統(tǒng)的階躍響應(yīng)。 解 首先在Simulink的Libr

36、ary窗口中的File中選擇New，打開一個新工作平臺，然后分別在信號源庫、輸出方式庫以及線性環(huán)節(jié)庫中找到 Step、Scope、Transfer Fun以及 Sum各功能塊，將其拖至工作平臺并按要求加以連接。在對相關(guān)參數(shù)進(jìn)行正確設(shè)置后，即可啟動仿真命令，如圖4-16所示，從而自動得到由圖4-17給出的系統(tǒng)仿真結(jié)果。 有關(guān)Simulink的詳細(xì)應(yīng)用請參閱Mathworks的simulink使用手冊。 圖4-14 線性環(huán)節(jié)模塊(Linear)庫 圖4-15 連接及接口模塊(Connection)庫 圖4-16 在Simulink 中建立仿真系統(tǒng) 圖4-17 Simulink仿真結(jié)果452 S-F

37、unction （S函數(shù)）應(yīng)用 前面已經(jīng)介紹過Simulink是一類面向結(jié)構(gòu)圖的仿真環(huán)境，然而由于Simulink提供了SFunction（S函數(shù)）的功能，從而極大地提高和豐富了Simulink對系統(tǒng)進(jìn)行分析與研究的能力。SFunction通過MATLAB語言進(jìn)行編程與Simulink中的SFunction功能模塊有機(jī)地結(jié)合在一起，實現(xiàn)了用戶所需的各種運(yùn)算和分析功能，因此說 S-Function是Simulink中極為重要的一部分內(nèi)容。SFunction的使用必須嚴(yán)格按照所規(guī)定的格式加以進(jìn)行，SFunction的運(yùn)行走向是由Flag標(biāo)志加以控制的，其實現(xiàn)的功能主要由sys參數(shù)及其相關(guān)的功能語句

38、塊所決定。在高版本（5.x）的MATLAN中，用戶可以利用模板程序進(jìn)行編程，保留所需功能，刪掉多余的內(nèi)容，從而簡化了編程工作量，保證了程序編制質(zhì)量，SFunction的模版程序可以在Blocks目錄中找到。 在S-Function中根據(jù)Flag提供的標(biāo)志，執(zhí)行6種不同的任務(wù)（case0，case1，case2，case3，case4與case9)，其中最為重要的是下述三種： Case 0， sys，x0，str， tsmdllnitializeSizes； ease l， sysmdlDerivatives（t，x，u）； Case 3， sysmdloutputs（t，x，u）； 它們分別為

39、SFunction提供進(jìn)行結(jié)構(gòu)初始化設(shè)置、狀態(tài)微分和狀態(tài)輸出功能。因此，欲靈活地使用SFunction，就必須首先掌握它們的應(yīng)用。在mdllnitializeSizes的返回變量中sys是根據(jù)Flag的標(biāo)志運(yùn)算得到的解，x0為初始狀態(tài)，str對于M類文件設(shè)置為空陣（S-Function也可由其他類文件構(gòu)成），ts為列向量，定義了取樣時間與偏移量。下面通過一個示例對SFunction的具體應(yīng)用加以說明 例46 已知系統(tǒng)狀態(tài)方程的系數(shù)陣為 試應(yīng)用SFunction求取系統(tǒng)的階躍響應(yīng)。 解：應(yīng)用S-Function一般可分為三步進(jìn)行。在Simulink窗口構(gòu)造所需結(jié)構(gòu)圖，其中必須包括SFunctio

40、n功能塊。根據(jù)用戶的需要編寫SFunction程序，這是最為重要的一步，當(dāng)然可以借助模版程序加以進(jìn)行。編寫相關(guān)的主運(yùn)行程序。 圖4-18 系統(tǒng)結(jié)構(gòu)圖首先建立結(jié)構(gòu)圖，如圖4-18所示。在結(jié)構(gòu)圖中S-Function模塊是關(guān)鍵，該模塊名應(yīng)與所編制的S-function程序名相同；To Workspace模塊將計算數(shù)據(jù)寫入工作空間；Fcn模塊為自定義函數(shù)模塊，用戶可根據(jù)需要自行定義。然后可利用模版根據(jù)要求編制如下相應(yīng)程序。 S-function程序:ssq.m functionsys,x0,str,tsssq（t， x，u，flag） A-18 19 -24；19 -21 20；44 -24 -36

41、； B=0 l l； C=1 0 2； switch flag， case 0， sys, x0, str, ts= mdlInitializeSizes（A,B,C）; case 1， sys= mdlDerivatives(t ，x、u, A, B, C); case 2， sys= ; case 3， sys=mdlOutputs(t，x，u，A，B，C); sys=simsizes(sizes); x0=0 0 0; str ； ts= 0 0 ； function sys=mdlDerivatives（t,x,u, A，B,C） ys A* x + B*u； function sys

42、= mdlUpdate（t，x，u） sys ； function case 4， sys= ; case 9， sys= ; otherwise error（'Unhandled flag，num2str（flag）； end functionsys，x0,str, ts =mdlInitializeSizes（A，B，C） sizes=simsizes; sizes.NumContstates =3； sizes.NumDiscStates =0； sizes.NumOutputs =1; sizes.NumInputs =1； sizes.DirFeedthrough =1; s

43、izes.NumSampleTimes =1； at least one sample time is needed sys=mdlOuputs（t，x，u，A，B，C) sys=C* x； function sys=mdlGetTimeOfNextVarHit(t,x,u) sampleTime 1； sys t sampleTime； function sys=mdlTerminate（t， x， u） sys= ; 無關(guān)函數(shù)一律賦空陣，函數(shù)是根據(jù)用戶的具體情況加以設(shè)置的，因為本系統(tǒng)有三個狀態(tài)變量，所以對 NumContStates設(shè)置 3 ；又因為是單輸入單輸出系統(tǒng)，所以 NumOutP

44、uts=1, NumInputs1。最后編制主程序，在這里主程序的功能就是運(yùn)行該SFunction并將結(jié)果繪制出來。主程序為 sim（sspaceq） plot（t，simout）運(yùn)行主程序后，圖417示出了系統(tǒng)運(yùn)行結(jié)果。這里僅給出了一個簡單的示例說明SFunction的應(yīng)用，顯然這種功能可以推廣到各種復(fù)雜的系統(tǒng)，請讀者自行加以嘗試。表41至表47分別給出了Simulink的基本模塊庫及其相應(yīng)的功能說明。 無關(guān)函數(shù)一律賦空陣，函數(shù)是根據(jù)用戶的具體情況加以設(shè)置的，因為本系統(tǒng)有三個狀態(tài)變量，所以對 NumContStates設(shè)置 3 ；又因為是單輸入單輸出系統(tǒng)，所以 NumOutPuts=1, NumInputs1。最后編制主程序，在這里主程序的功能就是運(yùn)行該SFu