連續(xù)系統(tǒng)數(shù)字仿真1_08

1、1第第3 3章章 連續(xù)系統(tǒng)的數(shù)字仿真連續(xù)系統(tǒng)的數(shù)字仿真 本章內(nèi)容本章內(nèi)容(1) 熟悉在數(shù)字計算機仿真技術中常用的幾種熟悉在數(shù)字計算機仿真技術中常用的幾種數(shù)值積分法數(shù)值積分法， 特別是四階龍格特別是四階龍格-庫塔法；庫塔法；(2) 典型環(huán)節(jié)典型環(huán)節(jié)及其系數(shù)矩陣的及其系數(shù)矩陣的確定確定；(3) 各各連接矩陣連接矩陣的確定；的確定；(4) 利用利用MATLAB在四階龍格在四階龍格-庫塔法的基礎上，對庫塔法的基礎上，對以狀態(tài)以狀態(tài) 空間表達式和方框圖描述的連續(xù)系統(tǒng)空間表達式和方框圖描述的連續(xù)系統(tǒng)進行進行仿真仿真；(5) 了解以了解以增廣矩陣法增廣矩陣法為基礎的連續(xù)系統(tǒng)的為基礎的連續(xù)系統(tǒng)的快速仿真快速仿

2、真方法方法。 2 用數(shù)字計算機來用數(shù)字計算機來仿真仿真或模擬一個連續(xù)控或模擬一個連續(xù)控制系統(tǒng)的目的就是制系統(tǒng)的目的就是求解系統(tǒng)的數(shù)學模型求解系統(tǒng)的數(shù)學模型。由控制理論知，一個由控制理論知，一個n n階連續(xù)系統(tǒng)可以被描階連續(xù)系統(tǒng)可以被描述成由述成由n n個積分器組成的模擬結構圖。因此個積分器組成的模擬結構圖。因此利用數(shù)字計算機來進行連續(xù)系統(tǒng)的仿真，利用數(shù)字計算機來進行連續(xù)系統(tǒng)的仿真，從本質(zhì)上講就是要在數(shù)字計算機上構造出從本質(zhì)上講就是要在數(shù)字計算機上構造出n n個數(shù)字積分器，也就是個數(shù)字積分器，也就是讓數(shù)字計算機進行讓數(shù)字計算機進行n n次數(shù)值積分運算次數(shù)值積分運算。可見，連續(xù)系統(tǒng)數(shù)字仿?？梢姡B

3、續(xù)系統(tǒng)數(shù)字仿真中的最基本的算法是數(shù)值積分算法。真中的最基本的算法是數(shù)值積分算法。343.1 3.1 數(shù)值積分法數(shù)值積分法 連續(xù)系統(tǒng)通常連續(xù)系統(tǒng)通常把數(shù)學模型化為狀態(tài)空間表達式把數(shù)學模型化為狀態(tài)空間表達式，為，為了對了對n階連續(xù)系統(tǒng)在數(shù)字計算機上仿真及求解，就要階連續(xù)系統(tǒng)在數(shù)字計算機上仿真及求解，就要采用數(shù)值積分法來求解系統(tǒng)數(shù)學模型中的采用數(shù)值積分法來求解系統(tǒng)數(shù)學模型中的n個一階微個一階微分方程。分方程。 設設n階連續(xù)系統(tǒng)由以下階連續(xù)系統(tǒng)由以下n個一階微分方程組成個一階微分方程組成 （3-13-1） 所謂數(shù)值積分法，就是要逐個求出區(qū)間所謂數(shù)值積分法，就是要逐個求出區(qū)間a ,b內(nèi)若內(nèi)若干個離散點干

4、個離散點a t0 t1 t0 0時，時，x( (t) )是未知的，因此式（是未知的，因此式（3-23-2）右端的積分是求不出的。為了解決這個問題，右端的積分是求不出的。為了解決這個問題，把積分間隔取得足夠小，使得在把積分間隔取得足夠小，使得在tk與與tk+1+1之間之間的的f( (t, ,x( (t)可以近似看作常數(shù)可以近似看作常數(shù)f( (tk,x( (tk)，這樣便得到用矩形公式積分的近似公式這樣便得到用矩形公式積分的近似公式或簡化為或簡化為 這就是這就是歐拉公式歐拉公式。htxtftxtxkkkk)(,()()(1hxtfxxkkkk),(17 以以x(t0)=x0作為初始值，應用歐拉公式

5、，作為初始值，應用歐拉公式，就可以一步步求出每一時刻就可以一步步求出每一時刻tk的的xk值，值， 即即 k =0,=0,x1 1x0 0+ +f( (t0 0, ,x0 0) )h k=1, =1, x2 2x1 1+ +f(t1 1, ,x1 1) )h k= =n-1-1, xn nxn-1n-1+ +f( (tn-1n-1, ,xn-1n-1) )h 這樣式這樣式(3-1)(3-1)的解的解x( (t) )就求出來了。歐拉法的計就求出來了。歐拉法的計算雖然比較簡單，但精度算雖然比較簡單，但精度較低。圖較低。圖3-13-1為歐拉法的為歐拉法的幾何解釋。幾何解釋。xx0 x2x1t0t2t1

6、(t1 ,x1)(t2 ,x2)t圖3 -1 歐 拉法的幾何解釋83.1.2 3.1.2 梯形法梯形法 由上可知由上可知歐拉公式歐拉公式中的積分是用中的積分是用矩形面積矩形面積f( (tk, ,xk) )h 來來近似近似的。的。 由圖由圖3-23-2知，用矩形面積知，用矩形面積tkabtk+1+1代替積分，其誤差就代替積分，其誤差就是圖中陰影部分。為了提高精是圖中陰影部分。為了提高精度現(xiàn)度現(xiàn)用梯形面積用梯形面積tkactk+1+1來代來代替積分替積分，即，即 于是可得梯形法的計算公式為于是可得梯形法的計算公式為 ),(),(2)(, (111kkkkttxtfxtfhdttxtfkk),(),

7、(2111kkkkkkxtfxtfhxxtdtdxf kt1ktabc圖3 -2 梯 形法的幾何意義9 由于上式由于上式右邊包含未知量右邊包含未知量xk+1+1, ,所以每一步都所以每一步都必須通過必須通過迭代求解迭代求解，每一步迭代的初值每一步迭代的初值xk+1+1(0)(0)通常采用歐拉公式來計算通常采用歐拉公式來計算，因此梯形法每一步，因此梯形法每一步迭代公式為迭代公式為 (3-3(3-3）式中式中 迭代次數(shù)迭代次數(shù)R=0,1,2,=0,1,2,),(),(2),()(11)1(1)0(1RkkkkkRkkkkkxtfxtfhxxxthfxx103.1.3 3.1.3 預估校正法預估校正

8、法雖然雖然梯形法梯形法比歐拉法精確，但是由于每比歐拉法精確，但是由于每一步都要進行多次疊代，計算量大，為了一步都要進行多次疊代，計算量大，為了簡化計算，有時只對式（簡化計算，有時只對式（3-33-3）進行）進行一次疊一次疊代代就可以了，因此可得就可以了，因此可得 通常稱這類方法為預估校正方法。它首通常稱這類方法為預估校正方法。它首先根據(jù)歐拉公式計算出先根據(jù)歐拉公式計算出xk+1的預估值的預估值xk+1(0)(0)，然后再對它進行校正，以得到更準確的，然后再對它進行校正，以得到更準確的近似值近似值xk+1(1)(1)。),(),(2),()0(11)1(1k)0(1kkkkkkkkkxtfxtf

9、hxxxthfxx113.1.4 3.1.4 龍格庫塔法龍格庫塔法 根據(jù)泰勒級數(shù)將式（根據(jù)泰勒級數(shù)將式（3-13-1）在）在tk+1=tk+h時刻的解時刻的解xk+1=x(tk+h) 在在tk附近展開，有附近展開，有 3-53-5） 可以看出，提高截斷誤差的階次，便可提高其精可以看出，提高截斷誤差的階次，便可提高其精度，但是由于計算各階導數(shù)相當麻煩，所以直接采度，但是由于計算各階導數(shù)相當麻煩，所以直接采用泰勒級數(shù)公式是不適用的，為了解決提高精度問用泰勒級數(shù)公式是不適用的，為了解決提高精度問題，龍格和庫塔兩人先后提出了間接使用泰勒級數(shù)題，龍格和庫塔兩人先后提出了間接使用泰勒級數(shù)公式的方法，即用函

10、數(shù)值公式的方法，即用函數(shù)值f (t,x)的線性組合來代替的線性組合來代替f (t,x)的導數(shù)，然后按泰勒公式確定其中的系數(shù)的導數(shù)，然后按泰勒公式確定其中的系數(shù), , 這這樣既能避免計算樣既能避免計算f (t,x)的導數(shù)，又可以提高數(shù)值計算的導數(shù)，又可以提高數(shù)值計算精度，其方法如下。精度，其方法如下。)( 0! 21 1)(21 ppkpkkkkhxphxhxhxx12因因故式（故式（3-53-5）可寫成）可寫成 （3-63-6） 為了避免計算式（為了避免計算式（3-63-6）中的各階導數(shù)項，可令）中的各階導數(shù)項，可令xk+1由由以下多項式表示。以下多項式表示。 （3-73-7）kxkxxttx

11、xttkkkkkfffdtdxxftfdtxtdfxfxtfxkkkkk )(),(),()(0)(! 2121pkxkkkkhfffhhfxxkvmmmkkkahxx1113 式中式中 am為待定因子，為待定因子，v為使用為使用f函數(shù)值的個數(shù)，函數(shù)值的個數(shù)，km滿足下列方程滿足下列方程 （3-83-8）即：即： 將式（將式（3-73-7）展開成）展開成h的冪級數(shù)并與微分方程式的冪級數(shù)并與微分方程式（3-13-1）精確解式（）精確解式（3-63-6）逐項比較）逐項比較, ,便可求得式便可求得式（3-73-7）和式（）和式（3-83-8）中的系數(shù)）中的系數(shù)am ,bmj和和cm等。等。vmchk

12、bxhctfkmjjmjkmkm, 2 , 1, 0),(111),(),(),(23213133121221hkbhkbxhctfkhkbxhctfkxtfkkkkkkk14現(xiàn)以現(xiàn)以v=2=2為例，來說明這些參數(shù)的確定方法。為例，來說明這些參數(shù)的確定方法。設設v=2=2，則有，則有 （3-93-9） 將將k1 1和和k k2 2在同一點（在同一點（tk ,xk）上用二元函數(shù)展開為）上用二元函數(shù)展開為),(),()(12122122111hkbxhctfkxtfkkakahxxkkkkkkkkkfxtfk),(1)( 0)( 0),(3212312122hff hbf hcfhxfhkbtfh

13、cxtfkkxkkxxttxxttkkkkkkk15將將k1 1和和k2 2代入式（代入式（3-93-9）整理后可得）整理后可得 （3-103-10）將上式與式（將上式與式（3-63-6）逐項進行比較，可得以下關系式）逐項進行比較，可得以下關系式若取若取 則則)( 0)()(3212222211hffbafcahfaahxxkxkkkkk21,21, 12122221bacaaa12c1,212121baa16于是可得于是可得 (3-11)(3-11) 由于式（由于式（3-113-11）只取到泰勒級數(shù)展開式的）只取到泰勒級數(shù)展開式的h2 2項，故稱這種方法為項，故稱這種方法為兩階龍格庫塔兩階龍

14、格庫塔法法，其截斷誤差為，其截斷誤差為0(0(h3 3) )。)(2211kkhxxkk),(),(121hkxhtfkxtfkkkkk17 同理當同理當v=4=4時，仿照上述方法可得如下時，仿照上述方法可得如下四階四階龍格龍格- -庫塔公式庫塔公式),()2,2()2,2(),()22(6342312143211hkxhtfkkhxhtfkkhxhtfkxtfkkkkkhxxkkkkkkkkkk18 通過上述龍格通過上述龍格- -庫塔法的介紹，可以把以上庫塔法的介紹，可以把以上介紹的幾種數(shù)值積分法統(tǒng)一起來，它們都是基介紹的幾種數(shù)值積分法統(tǒng)一起來，它們都是基于在初值附近展開成泰勒級數(shù)的原理，所

15、不同于在初值附近展開成泰勒級數(shù)的原理，所不同的是取泰勒級數(shù)多少項。的是取泰勒級數(shù)多少項。歐拉公式僅取到歐拉公式僅取到h項項，梯形法與二階龍格庫塔法相同均取到梯形法與二階龍格庫塔法相同均取到h2 2項，項，四四階龍格庫塔法取到階龍格庫塔法取到h4 4項項。從理論上講，取得的。從理論上講，取得的項數(shù)愈多，計算精度愈高，但計算量愈大，愈項數(shù)愈多，計算精度愈高，但計算量愈大，愈復雜，計算誤差也將增加，因此要適當?shù)倪x擇。復雜，計算誤差也將增加，因此要適當?shù)倪x擇。目前在數(shù)字仿真中，目前在數(shù)字仿真中，最常用的是四階龍格庫塔最常用的是四階龍格庫塔法，其截斷誤差為法，其截斷誤差為( (h5 5),), 已能滿足

16、仿真精度已能滿足仿真精度的要求。的要求。193.1.5 3.1.5 關于仿真數(shù)值積分法的幾點討關于仿真數(shù)值積分法的幾點討1 1單步法和多步法單步法和多步法 解解初值問題初值問題的數(shù)值解法的的數(shù)值解法的共同特點共同特點是是步進式步進式，即，即從最初一點或幾點出發(fā)，每一步根據(jù)從最初一點或幾點出發(fā)，每一步根據(jù)xk一點或前面一點或前面幾點幾點xk-1 , xk-2 ,來計算新的來計算新的xk+1的值，這樣逐步推的值，這樣逐步推進。進。 當從當從tk推進到推進到tk+1只需用只需用tk時刻的數(shù)據(jù)時，稱為時刻的數(shù)據(jù)時，稱為單步單步法，例如歐拉法和龍格庫塔法。法，例如歐拉法和龍格庫塔法。20 相反，需要用到

17、相反，需要用到tk以及以及過去時刻過去時刻tk-1 ,tk-2 ,的數(shù)據(jù)的數(shù)據(jù)時稱為時稱為多步法多步法。線性多步法的一般形式是。線性多步法的一般形式是 （3-133-13） 多步法不能從多步法不能從t=0=0自啟動，通常需要選用相同階次自啟動，通常需要選用相同階次精度的單步法來啟動，獲得所需前精度的單步法來啟動，獲得所需前k步數(shù)據(jù)后，方可步數(shù)據(jù)后，方可轉(zhuǎn)入相應多步法，因多步法利用信息量大，因而比單轉(zhuǎn)入相應多步法，因多步法利用信息量大，因而比單步法更精確。步法更精確。)(0111101nknkknknkkkfffhxxxx212 2顯式和隱式顯式和隱式 在計算在計算xk+1時公式時公式右端所用到

18、的數(shù)據(jù)均已知右端所用到的數(shù)據(jù)均已知時，稱時，稱為為顯式顯式算法。例如算法。例如歐拉法、龍格庫塔法歐拉法、龍格庫塔法和式（和式（3-133-13）中中-1-1=0=0的情況。相反，在算式右端中隱含有未知量的情況。相反，在算式右端中隱含有未知量xk+1時，稱為時，稱為隱式隱式算法。例如算法。例如梯形法、預估校正法梯形法、預估校正法和和式（式（3-133-13）中）中-1-100的情況。的情況。 顯式算法利用前幾步計算結果即可進行遞推求解下顯式算法利用前幾步計算結果即可進行遞推求解下步結果，因而易于計算。而隱式計算需要迭代法，先步結果，因而易于計算。而隱式計算需要迭代法，先用另一同階次顯式公式估計出

19、一個初值用另一同階次顯式公式估計出一個初值xk+1(0)(0), ,并求得并求得fk+1, ,然后再用隱式求得校正值然后再用隱式求得校正值xk+1(1)(1), , 若未達到所需若未達到所需精度要求，則再次迭代求解，直到兩次迭代值精度要求，則再次迭代求解，直到兩次迭代值xk+1 ( (i) )和和xk+1(i+1) 之間的誤差在要求的范圍內(nèi)為止，故之間的誤差在要求的范圍內(nèi)為止，故隱式隱式算法算法精度高精度高，對誤差有較強的抑制作用對誤差有較強的抑制作用。223.3.數(shù)值穩(wěn)定性與仿真誤差數(shù)值穩(wěn)定性與仿真誤差 仿真誤差與數(shù)值計算方法、計算機的精度以及仿真誤差與數(shù)值計算方法、計算機的精度以及計算步長

20、的選擇有關。當計算方程和計算機確定以計算步長的選擇有關。當計算方程和計算機確定以后，則僅與計算步長有關，所以在仿真中計算步長后，則僅與計算步長有關，所以在仿真中計算步長是一個重要的參數(shù)。仿真誤差一般有如下兩種：是一個重要的參數(shù)。仿真誤差一般有如下兩種：（1 1）截斷誤差截斷誤差由于仿真模由于仿真模型僅是原系統(tǒng)型僅是原系統(tǒng)模型模型的一種的一種逼近逼近，以及各種數(shù)值積分法的計算都以及各種數(shù)值積分法的計算都是近似的算法。通常計算步長是近似的算法。通常計算步長愈小，截斷誤差也愈小。愈小，截斷誤差也愈小。（2 2）舍入誤差舍入誤差由于由于計算機計算機的精度有限的精度有限（有限位數(shù)）所產(chǎn)（有限位數(shù)）所產(chǎn)生

21、。通常計算步長愈小，計算生。通常計算步長愈小，計算次數(shù)愈多次數(shù)愈多, , 舍入誤差愈大。舍入誤差愈大。圖3 -3 仿 真誤差曲線h 0h誤差總誤差截斷誤差舍入誤差023 對截斷誤差而言，計算步長愈小愈好，但太小不但對截斷誤差而言，計算步長愈小愈好，但太小不但會增加計算時間，而且由于舍入誤差的增加，不一定會增加計算時間，而且由于舍入誤差的增加，不一定能達到提高精度的目的，甚至可能出現(xiàn)數(shù)值不穩(wěn)情況。能達到提高精度的目的，甚至可能出現(xiàn)數(shù)值不穩(wěn)情況。顯然計算步長太大，不但精度不能滿足要求，而且計顯然計算步長太大，不但精度不能滿足要求，而且計算步長超過該算法的判穩(wěn)條件時，也會出現(xiàn)不穩(wěn)定情算步長超過該算法

22、的判穩(wěn)條件時，也會出現(xiàn)不穩(wěn)定情況。由此可見，計算步長只能在某一范圍內(nèi)選擇，圖況。由此可見，計算步長只能在某一范圍內(nèi)選擇，圖中的中的h0 0為最佳計算步長。為最佳計算步長。 一般一般控制系統(tǒng)控制系統(tǒng)的輸出的輸出動態(tài)響應動態(tài)響應在在開始段變化較快開始段變化較快，到到最后變化將會很緩慢最后變化將會很緩慢。這時，計算可以。這時，計算可以采用變步長采用變步長的方法的方法，即在開始階段步長取得小一些，在最后階段，即在開始階段步長取得小一些，在最后階段取得大一些，這樣即可以保證計算的取得大一些，這樣即可以保證計算的精度精度，也可以加，也可以加快計算的快計算的速度速度。24 對于一般工程計算，計算精度要求并不

23、對于一般工程計算，計算精度要求并不太高，故常用太高，故常用定步長定步長的方法。作為的方法。作為經(jīng)驗經(jīng)驗數(shù)據(jù)，當采用四階龍格庫塔法作數(shù)值數(shù)據(jù)，當采用四階龍格庫塔法作數(shù)值積分計算時，取計算步長積分計算時，取計算步長 h=tr r/10/10或或ts/40/40 式中式中 tr系統(tǒng)在階躍函數(shù)作用下的上升系統(tǒng)在階躍函數(shù)作用下的上升時間；時間；ts系統(tǒng)在階躍函數(shù)作用下的過渡系統(tǒng)在階躍函數(shù)作用下的過渡過程時間。若系統(tǒng)有多個回路，則應按過程時間。若系統(tǒng)有多個回路，則應按反應最快的回路考慮。反應最快的回路考慮。253.2 3.2 連續(xù)系統(tǒng)的數(shù)字仿真程序連續(xù)系統(tǒng)的數(shù)字仿真程序若系統(tǒng)的狀態(tài)空間表達式為若系統(tǒng)的狀態(tài)

24、空間表達式為 （3-143-14） （3-153-15）其中其中 A: nn; b: n1; C:1nCxybuAxx26 假設在仿真中，數(shù)值積分法采用四階龍格假設在仿真中，數(shù)值積分法采用四階龍格庫塔方法，因?qū)τ趲焖椒ǎ驅(qū)τ趎 n階系統(tǒng)，其狀態(tài)方程式階系統(tǒng)，其狀態(tài)方程式（3-143-14）可寫成以下）可寫成以下n個一階微分方程個一階微分方程 (3-16)(3-16)ubxaxaxaxininiii.2211), 2 , 1(), (321nixxxxtfni27 故根據(jù)式（故根據(jù)式（3-123-12）可得）可得求解以上一階微求解以上一階微分方程組的四階龍格庫塔公式分方程組的四階龍格庫塔公式

25、如下如下), 2 , 1()22(64321)()1(nikkkkhxxiiiikiki)()()()2()2()2()2()2()2()(33114221131111222111htubhkxahkxakhtubkhxakhxakhtubkhxakhxaktubxaxaxakkiikninikiikiikninikiikiikninikiikikninkikii式中式中 xik為為t=tk時刻的時刻的xi值，值，xik+1為為t=tk+h時時刻的刻的xi值。值。28令令TnTnTnTnTknkkkTknkkkkkkKkkkKkkkKkkkKxxxxxxxx424144323133222122

26、1211112111211129則式（則式（3-173-17）可寫成如下矩陣的形式）可寫成如下矩陣的形式根據(jù)式（根據(jù)式（3-153-15）可得）可得t=tk+1時刻的輸出時刻的輸出)()()2()2()2()2()()22(6342312143211htbuhKxAKhtbuKhxAKhtbuKhxAKtbuAxKKKKKhxxkkkkkkkkkk11kkCxy30給定系統(tǒng)參數(shù)并根據(jù)需要轉(zhuǎn)化為狀態(tài)空間表達式輸入仿真時間和計算步長根據(jù)龍格庫塔法求 K1,K2,K3,K4求 和 值1kx1kyt=Tf輸出結果Ny圖3 -4 連 續(xù)系統(tǒng)的數(shù)字仿真程序框圖開始給定輸入信號31例例3-13-1 假設單變

27、量系統(tǒng)如下圖所示。假設單變量系統(tǒng)如下圖所示。 試根據(jù)四階龍格庫塔法，試根據(jù)四階龍格庫塔法，求輸出量求輸出量y的動的動態(tài)響應。態(tài)響應。yr=20.40.002+ 圖3 -5-5 .638 .2132 .2045410562. 110875. 123466sssss32解解 仿真程序如下仿真程序如下ex3_1ex3_1 取仿真時間：取仿真時間：Tf f=5=5； 計算步長：計算步長：h=0.02 =0.02 在在MATLABMATLAB環(huán)境下執(zhí)行以上程序可得如圖環(huán)境下執(zhí)行以上程序可得如圖3-63-6所示仿真所示仿真曲線。曲線。33圖3-634 圖中帶數(shù)字的方框表示典型環(huán)節(jié)，圖中帶數(shù)字的方框表示典型

28、環(huán)節(jié)，, , ,表表示比例系數(shù)。示比例系數(shù)。V352r1r2r3+-y1y2x3 u4u3x2x1u2u1x4x5u5圖3 -9 多 變量系統(tǒng)結構圖 +353.3 3.3 面向系統(tǒng)結構圖的仿真面向系統(tǒng)結構圖的仿真 這種方法與上節(jié)介紹的方法相比，有以下這種方法與上節(jié)介紹的方法相比，有以下幾個主要優(yōu)點：幾個主要優(yōu)點：1）便于研究各環(huán)節(jié)參數(shù)對系統(tǒng)的影響；便于研究各環(huán)節(jié)參數(shù)對系統(tǒng)的影響；2）可以得到每個環(huán)節(jié)的動態(tài)響應；可以得到每個環(huán)節(jié)的動態(tài)響應；3）可對多變量系統(tǒng)進行仿真?？蓪Χ嘧兞肯到y(tǒng)進行仿真。下面具體介紹面向結構圖的仿真方法。下面具體介紹面向結構圖的仿真方法。363.3.1 3.3.1 典型環(huán)節(jié)的

29、確定典型環(huán)節(jié)的確定 一個控制系統(tǒng)可能由各種各樣的環(huán)節(jié)所一個控制系統(tǒng)可能由各種各樣的環(huán)節(jié)所組成，但比較常見環(huán)節(jié)有：組成，但比較常見環(huán)節(jié)有：（1）比例環(huán)節(jié)：比例環(huán)節(jié)： （2）積分環(huán)節(jié)：積分環(huán)節(jié)： （3）比例積分環(huán)節(jié)：比例積分環(huán)節(jié)：（4 4）慣性環(huán)節(jié)：）慣性環(huán)節(jié)：（5）超前滯后環(huán)節(jié)：超前滯后環(huán)節(jié)：（6 6）二階振蕩環(huán)節(jié)：）二階振蕩環(huán)節(jié)： 37 為了編制比較簡單而且通用的仿真程序必須恰當選為了編制比較簡單而且通用的仿真程序必須恰當選擇仿真環(huán)節(jié)。在這里選用圖擇仿真環(huán)節(jié)。在這里選用圖3-73-7所示的典型環(huán)節(jié)作為所示的典型環(huán)節(jié)作為仿真環(huán)節(jié)，即仿真環(huán)節(jié)，即式中式中 u為典型環(huán)節(jié)的輸入，為典型環(huán)節(jié)的輸入，x為

30、典型環(huán)節(jié)的輸出。為典型環(huán)節(jié)的輸出。 利用這個典型環(huán)節(jié)，只要利用這個典型環(huán)節(jié)，只要改變改變a ,b ,c 和和d 參數(shù)的參數(shù)的值值，便可分別表示以上所述的各一階環(huán)節(jié)，至于二階，便可分別表示以上所述的各一階環(huán)節(jié)，至于二階振蕩環(huán)節(jié)，則可用振蕩環(huán)節(jié)，則可用兩個一階環(huán)節(jié)等效連兩個一階環(huán)節(jié)等效連接得到，如圖接得到，如圖3-83-8所示。所示。bsadscsUsXsG)()()(xu圖3-7 典型環(huán)節(jié)bsadscx21TsTs1ku+-圖3-8 二階振蕩環(huán)節(jié)的等效結構圖383.3.2 3.3.2 連接矩陣連接矩陣一個控制系統(tǒng)用典型環(huán)節(jié)來描述時，必一個控制系統(tǒng)用典型環(huán)節(jié)來描述時，必須用連接矩陣把各個典型環(huán)節(jié)連

31、接起來。須用連接矩陣把各個典型環(huán)節(jié)連接起來。所謂連接矩陣，就是用矩陣的形式表示各所謂連接矩陣，就是用矩陣的形式表示各個典型環(huán)節(jié)之間的關系。下面介紹連接矩個典型環(huán)節(jié)之間的關系。下面介紹連接矩陣的建立方法陣的建立方法, , 假設多輸入多輸出系統(tǒng)的假設多輸入多輸出系統(tǒng)的結構圖如圖結構圖如圖3-93-9所示。所示。39 圖中帶數(shù)字的方框表示典型環(huán)節(jié)，圖中帶數(shù)字的方框表示典型環(huán)節(jié)，, , ,表表示比例系數(shù)。示比例系數(shù)。V352r1r2r3+-y1y2x3 u4u3x2x1u2u1x4x5u5圖3 -9 多 變量系統(tǒng)結構圖 +40 由圖可得由圖可得各環(huán)節(jié)的輸入各環(huán)節(jié)的輸入與與各環(huán)節(jié)的輸出各環(huán)節(jié)的輸出間的間

32、的關系關系以及系統(tǒng)的以及系統(tǒng)的輸出與環(huán)節(jié)的輸出輸出與環(huán)節(jié)的輸出間間關系關系分別為分別為和和 453433323222121551xuxurxxurxxurxu3241xyxy41寫成寫成矩陣矩陣形式形式和和32154321325543210000001000100010100000100001000010000rrrxxxxxuuuuu54321210010001000 xxxxxyy42 或?qū)懗苫驅(qū)懗?(3-20)(3-20) 定義式中的定義式中的W ，W0 和和 Wc陣為連接矩陣陣為連接矩陣，W反映了各典型環(huán)節(jié)輸入輸出間的連接關系，反映了各典型環(huán)節(jié)輸入輸出間的連接關系，W0反映了系統(tǒng)的參考

33、輸入與各環(huán)節(jié)輸入間的反映了系統(tǒng)的參考輸入與各環(huán)節(jié)輸入間的連接關系，連接關系， Wc c反映了系統(tǒng)的輸出與各環(huán)節(jié)輸反映了系統(tǒng)的輸出與各環(huán)節(jié)輸出間的關系。出間的關系。xWyrWWxuc043 一般也將系統(tǒng)中一般也將系統(tǒng)中各典型環(huán)節(jié)各典型環(huán)節(jié)的系數(shù)寫成如下的系數(shù)寫成如下矩陣的形式矩陣的形式（假設系統(tǒng)由（假設系統(tǒng)由n n個典型環(huán)節(jié)組成）個典型環(huán)節(jié)組成） (3-21)(3-21)nnnndcbadcbadcbaP22221111443.3.3 3.3.3 確定系統(tǒng)的狀態(tài)方程確定系統(tǒng)的狀態(tài)方程 典型環(huán)節(jié)和連接矩陣確定后典型環(huán)節(jié)和連接矩陣確定后，便可，便可求求得得系系統(tǒng)的統(tǒng)的狀態(tài)空間表達式狀態(tài)空間表達式，推

34、導過程推導過程如下。如下。 假設系統(tǒng)由假設系統(tǒng)由n n個典型環(huán)節(jié)組成，則根據(jù)典個典型環(huán)節(jié)組成，則根據(jù)典型環(huán)節(jié)的傳遞函數(shù)有型環(huán)節(jié)的傳遞函數(shù)有 （ i=1,2,）即即 sbasdcsUsXsGiiiiiii)()()()()()()(sUsdcsXsbaiiiiii45寫成矩陣形寫成矩陣形 （3-223-22）式中：式中：nnbbbBaaaA000000,0000002121nndddDcccC000000,0000002121)()()()(sUsDCsXsBA46 將式（將式（3-203-20）中的上式進行拉氏變換后代）中的上式進行拉氏變換后代入式（入式（3-223-22）中可得）中可得 對上

35、式兩邊取拉氏反變換得對上式兩邊取拉氏反變換得（3-233-23） 若參考輸入向量若參考輸入向量r=r1 1 r2 2 rm T T中中的的, ,r2 2,rm均為階躍函數(shù)，則上式可簡化為均為階躍函數(shù)，則上式可簡化為 （3-243-24）)()()()()()()()()()()(000ssRWDsRWCsXAWCssXWDBsRWsWXsDCsXsBArWDrWCxAWCxWDB00)()(rWCxAWCxWDB0)()(47令令則式（則式（3-243-24）可寫成）可寫成若若H的逆存在，則有的逆存在，則有再令再令可得可得 （3-253-25） 上式即為閉環(huán)系統(tǒng)的狀態(tài)方程，它是一個典型的上式即

36、為閉環(huán)系統(tǒng)的狀態(tài)方程，它是一個典型的狀態(tài)方程，利用前面介紹的求解方法可方便地求出各狀態(tài)方程，利用前面介紹的求解方法可方便地求出各典型環(huán)節(jié)的輸出響應，最后根據(jù)式（典型環(huán)節(jié)的輸出響應，最后根據(jù)式（3-203-20）中的第二）中的第二式便可求出系統(tǒng)的輸出響應。式便可求出系統(tǒng)的輸出響應。AWCQWDBH,rWCQxxH0rWCHQxHx011011,WCHBQHABrAxx483.3.4 3.3.4 面向結構圖的數(shù)字仿真程序面向結構圖的數(shù)字仿真程序 面向結構圖面向結構圖的數(shù)字仿真的數(shù)字仿真程序框圖如程序框圖如圖圖-10-10所所示，其程序示，其程序清單通過下清單通過下例給出。例給出。給定環(huán)節(jié)參數(shù)及連接

37、矩陣輸入仿真時間和計算步長求 和 Q 陣值1HH ，求A ,B 陣t =Tf輸出結果根據(jù)龍格庫塔法求狀態(tài)方程的解Ny圖3 -10 面 向結構圖的數(shù)字仿真程序框圖開始給定輸入信號49 例例3-23-2 假設某一系統(tǒng)由四個典型環(huán)節(jié)組成，假設某一系統(tǒng)由四個典型環(huán)節(jié)組成，如圖如圖3-113-11所示。求輸出量所示。求輸出量y的動態(tài)響應。的動態(tài)響應。r=10圖3-111 . 05 . 0sss122s1010sy+-x1 u2x2 u3x3 u4u1x450 解解 由圖可得各由圖可得各環(huán)節(jié)環(huán)節(jié)的的輸入與輸出輸入與輸出以及系統(tǒng)的以及系統(tǒng)的輸出與環(huán)節(jié)的輸出間關系為輸出與環(huán)節(jié)的輸出間關系為rxxxxuuuu0

38、00101000010000110004321432143211000 xxxxy51 根據(jù)以上兩式和各典型環(huán)節(jié)的系數(shù)值，根據(jù)以上兩式和各典型環(huán)節(jié)的系數(shù)值，可得如下可得如下連接矩陣和系數(shù)矩陣連接矩陣和系數(shù)矩陣1000,000101000010000110000cWWW0101100212011015 . 011 . 04444333322221111dcbadcbadcbadcbaP52取仿真時間取仿真時間: : Tf=10；計算步長；計算步長: : h=0.05 在在MATLABMATLAB環(huán)境下執(zhí)行以上程序可得如圖環(huán)境下執(zhí)行以上程序可得如圖3-123-12所示仿真曲線。所示仿真曲線。仿真程

39、序如下仿真程序如下ex3_2.mex3_2.m53圖3-12543.4 3.4 連續(xù)系統(tǒng)的快速仿真連續(xù)系統(tǒng)的快速仿真前面介紹過的兩種連續(xù)系統(tǒng)的數(shù)字仿真前面介紹過的兩種連續(xù)系統(tǒng)的數(shù)字仿真方法，當系統(tǒng)比較復雜并要求滿足較高的方法，當系統(tǒng)比較復雜并要求滿足較高的計算精度時，計算工作量較大，計算速度計算精度時，計算工作量較大，計算速度較慢，有時不能滿足實時仿真的要求，為較慢，有時不能滿足實時仿真的要求，為了解決這個問題，下面介紹一種連續(xù)系統(tǒng)了解決這個問題，下面介紹一種連續(xù)系統(tǒng)的快速數(shù)字仿真方法的快速數(shù)字仿真方法-增廣矩陣法增廣矩陣法。553.4.1 3.4.1 增廣矩陣法的基本原理增廣矩陣法的基本原理設連續(xù)系統(tǒng)的狀態(tài)方程為設連續(xù)系統(tǒng)的狀態(tài)方程為( (齊次齊次) ) （3-263-26）則它的解為則它的解為將將eAt展開成泰勒級數(shù)，即展開成泰勒級數(shù)，即則有則有)0()()()(0 xtxtAxtxt)0()(xetxAt3322! 31! 21tAtAAtIeAt)0(! 31! 21)(3322xtAtAAtItx56 可以證明，如果取可以證明，如果取eAt的泰勒級數(shù)的前五項，則式的泰勒級數(shù)的前五項，則式（3-273-27）的計算精度與四階龍格庫塔法相同，設計）的計算精度與四階龍格庫塔法相同，設計算步長為算步長為T, , 則式（則式（3-273-27）可寫成）可寫成 上式括號中只