第三章連續(xù)系統(tǒng)仿真方法學(xué)

1、第三章第三章 連續(xù)系統(tǒng)仿真方法學(xué)連續(xù)系統(tǒng)仿真方法學(xué)本章主要內(nèi)容本章主要內(nèi)容n連續(xù)系統(tǒng)建模方法連續(xù)系統(tǒng)建模方法n模型變換模型變換n連續(xù)系統(tǒng)仿真算法連續(xù)系統(tǒng)仿真算法n采樣控制系統(tǒng)仿真采樣控制系統(tǒng)仿真*n分布參數(shù)系統(tǒng)仿真分布參數(shù)系統(tǒng)仿真*第一節(jié)第一節(jié) 連續(xù)系統(tǒng)建模方法連續(xù)系統(tǒng)建模方法n先驗知識建模先驗知識建模n機(jī)理建模方法機(jī)理建模方法知識模型知識模型 n常見形式：各種學(xué)科的公理、定理、定律等常見形式：各種學(xué)科的公理、定理、定律等 n專家系統(tǒng)方法專家系統(tǒng)方法邏輯關(guān)系模型邏輯關(guān)系模型 n符號、關(guān)系式、專家知識庫、推理規(guī)則等符號、關(guān)系式、專家知識庫、推理規(guī)則等 n模糊系統(tǒng)方法模糊系統(tǒng)方法模糊模型模糊模型

2、n高矮、大小等模糊語言量化成定量的表示形式，高矮、大小等模糊語言量化成定量的表示形式，按某種算法得到定量的結(jié)果后再轉(zhuǎn)換為模糊語言按某種算法得到定量的結(jié)果后再轉(zhuǎn)換為模糊語言 以上方法有時也用于離散事件系統(tǒng)建模以上方法有時也用于離散事件系統(tǒng)建模 n系統(tǒng)辨識建模系統(tǒng)辨識建模 n經(jīng)驗方法經(jīng)驗方法 n直接觀察數(shù)據(jù)曲線得出模型方程，如線性系統(tǒng)，直接觀察數(shù)據(jù)曲線得出模型方程，如線性系統(tǒng)，一階對象等一階對象等 n表格插值，一種靜態(tài)建模技術(shù)，主要用于計算動表格插值，一種靜態(tài)建模技術(shù)，主要用于計算動態(tài)方程中的系數(shù)態(tài)方程中的系數(shù) n統(tǒng)計建模（數(shù)理統(tǒng)計的方法）統(tǒng)計建模（數(shù)理統(tǒng)計的方法） n最小二乘法及其改進(jìn)形式、極大

3、似然估計法等最小二乘法及其改進(jìn)形式、極大似然估計法等 n神經(jīng)網(wǎng)絡(luò)神經(jīng)網(wǎng)絡(luò) n混合建模方法混合建模方法 n若干種模型形式（輸出）互相補(bǔ)充若干種模型形式（輸出）互相補(bǔ)充n 給定輸入后，從機(jī)理模型中產(chǎn)生輸出，與辨給定輸入后，從機(jī)理模型中產(chǎn)生輸出，與辨識模型的輸出按某種方式得到系統(tǒng)輸出，反識模型的輸出按某種方式得到系統(tǒng)輸出，反過來可以用輸出誤差繼續(xù)修正辨識模型過來可以用輸出誤差繼續(xù)修正辨識模型 第二節(jié)第二節(jié) 模型變換模型變換n連續(xù)系統(tǒng)常用的模型表示形式連續(xù)系統(tǒng)常用的模型表示形式n連續(xù)時間模型連續(xù)時間模型n系統(tǒng)的輸入量系統(tǒng)的輸入量u(t)，輸出量，輸出量y(t)及內(nèi)部狀態(tài)變量及內(nèi)部狀態(tài)變量x(t)均為時

4、間的連續(xù)函數(shù)均為時間的連續(xù)函數(shù)n微分方程微分方程n傳遞函數(shù)傳遞函數(shù)n權(quán)函數(shù)權(quán)函數(shù)n狀態(tài)空間表達(dá)式狀態(tài)空間表達(dá)式微分方程微分方程其中其中mnn用古典方法求解時非常復(fù)雜，用古典方法求解時非常復(fù)雜，n高階系統(tǒng)通常沒有封閉解或解析解高階系統(tǒng)通常沒有封閉解或解析解1111011011nnmmnmmnnmmd ydydyd uduaaa ybbbub udtdtdtdtdt1111011011nnmmnmmnnmmd ydydyd uduaaa ybbbub udtdtdtdtdt傳遞函數(shù)傳遞函數(shù) n當(dāng)初始條件為零時，對上述微分方程式作拉氏變換，當(dāng)初始條件為零時，對上述微分方程式作拉氏變換，可得傳遞函數(shù)形

5、式可得傳遞函數(shù)形式 n求解時可先用部分分式展開求解時可先用部分分式展開 n再進(jìn)行反變換即得時間解再進(jìn)行反變換即得時間解 nS域（復(fù)頻域）內(nèi)求解較為簡便但對多變量、時變或高域（復(fù)頻域）內(nèi)求解較為簡便但對多變量、時變或高階系統(tǒng)仍求解困難階系統(tǒng)仍求解困難 11101110( )( )( )mmmmnnnb sbsb sbY sW sU ssasa sa12( )( ) ( )( )( )( )lY sW s U sW sW sW s112( ) ( )( )( )( )ly tLY sy ty ty t權(quán)函數(shù)權(quán)函數(shù) 權(quán)函數(shù)權(quán)函數(shù) g(t)指初始條件為指初始條件為0時系統(tǒng)在理想脈沖函時系統(tǒng)在理想脈沖函

6、數(shù)數(shù)(t)作用下的響應(yīng)，又稱脈沖過渡函數(shù)作用下的響應(yīng)，又稱脈沖過渡函數(shù)系統(tǒng)對任意輸入的響應(yīng)可由卷積積分公式求出系統(tǒng)對任意輸入的響應(yīng)可由卷積積分公式求出0( )( ) ()ty tug td權(quán)函數(shù)與傳遞函數(shù)有如下關(guān)系：權(quán)函數(shù)與傳遞函數(shù)有如下關(guān)系： ( )( )L g tG s 為為n維狀態(tài)向量，維狀態(tài)向量，u為為r維輸入向量，維輸入向量，y為為m維輸出向量維輸出向量狀態(tài)空間表達(dá)式狀態(tài)空間表達(dá)式 動態(tài)系統(tǒng)的狀態(tài)是指能完全描述系統(tǒng)行為的最小一組變量動態(tài)系統(tǒng)的狀態(tài)是指能完全描述系統(tǒng)行為的最小一組變量若知道若知道t=t0時刻的初始狀態(tài)向量時刻的初始狀態(tài)向量x0及及tt0時的輸入時的輸入u就能完全就能完全

7、確定系統(tǒng)在確定系統(tǒng)在tt0時刻的行為時刻的行為狀態(tài)空間表達(dá)式由狀態(tài)方程和輸出方程組成狀態(tài)空間表達(dá)式由狀態(tài)方程和輸出方程組成xAxBuyCxDu1,Tnxxxn nAA為參數(shù)矩陣（或稱動態(tài)矩陣），為參數(shù)矩陣（或稱動態(tài)矩陣），n rBB為輸入矩陣，為輸入矩陣，m nCC為輸出矩陣為輸出矩陣m rDD運(yùn)用矩陣計算方法且借助于計算機(jī)很容易運(yùn)用矩陣計算方法且借助于計算機(jī)很容易對狀態(tài)空間方程求解對狀態(tài)空間方程求解 為關(guān)聯(lián)矩陣（輸入和輸出直接關(guān)聯(lián)）為關(guān)聯(lián)矩陣（輸入和輸出直接關(guān)聯(lián)）n離散時間模型離散時間模型n系統(tǒng)的輸入量，輸出量及內(nèi)部狀態(tài)量均為時間的系統(tǒng)的輸入量，輸出量及內(nèi)部狀態(tài)量均為時間的離散函數(shù)，即時間序

8、列離散函數(shù)，即時間序列u(kT),y(kT),x(kT)n差分方程差分方程nZ傳遞函數(shù)傳遞函數(shù)n權(quán)序列權(quán)序列n離散狀態(tài)空間模型離散狀態(tài)空間模型差分方程差分方程nT為采樣周期為采樣周期 11010() (1) (1) ()() (1) ()()nmmy kn Tay knTa y kTa y kTc u km Tcu kmTc u kTmnZ傳遞函數(shù)傳遞函數(shù)n對差分方程作對差分方程作Z變換，設(shè)所有初值為變換，設(shè)所有初值為0，則有，則有11101110( )( )( )mmmmnnnb zbzb zbY zW zU zzczc zc權(quán)序列權(quán)序列n權(quán)序列權(quán)序列h(k)為對初始條件為為對初始條件為0的

9、系統(tǒng)施加單位脈沖的系統(tǒng)施加單位脈沖序列序列(k)所得到的響應(yīng)所得到的響應(yīng) n系統(tǒng)對于任意輸入系統(tǒng)對于任意輸入u(k)的響應(yīng)為一卷積的響應(yīng)為一卷積 n與與Z傳遞函數(shù)間關(guān)系傳遞函數(shù)間關(guān)系 0( )( ) ()kiy ku i h ki ( )( )Z h kW z離散狀態(tài)空間模型離散狀態(tài)空間模型連續(xù)連續(xù)-離散混合模型離散混合模型n如計算機(jī)控制系統(tǒng)，對連續(xù)對象進(jìn)行控制時，狀態(tài)量如計算機(jī)控制系統(tǒng)，對連續(xù)對象進(jìn)行控制時，狀態(tài)量中既有連續(xù)的也有離散的中既有連續(xù)的也有離散的 (1)( )( )( )( )x kFx kGu ky kx k n連續(xù)系統(tǒng)模型之間的變換連續(xù)系統(tǒng)模型之間的變換 n微分方程、傳遞函數(shù)

10、、權(quán)函數(shù)模型描述系統(tǒng)微分方程、傳遞函數(shù)、權(quán)函數(shù)模型描述系統(tǒng)的輸入與輸出關(guān)系，稱為系統(tǒng)的外部模型的輸入與輸出關(guān)系，稱為系統(tǒng)的外部模型 n狀態(tài)方程則稱為系統(tǒng)的內(nèi)部模型狀態(tài)方程則稱為系統(tǒng)的內(nèi)部模型 n通常在仿真時，需要將系統(tǒng)的各種描述形式通常在仿真時，需要將系統(tǒng)的各種描述形式轉(zhuǎn)換成內(nèi)部模型，稱為模型結(jié)構(gòu)變換轉(zhuǎn)換成內(nèi)部模型，稱為模型結(jié)構(gòu)變換n化微分方程為狀態(tài)方程化微分方程為狀態(tài)方程 n化連續(xù)狀態(tài)方程為離散狀態(tài)方程化連續(xù)狀態(tài)方程為離散狀態(tài)方程化微分方程為狀態(tài)方程化微分方程為狀態(tài)方程設(shè)有微分方程設(shè)有微分方程 先考慮右邊僅含先考慮右邊僅含u的形式，令的形式，令 = 則有則有 1111011011nnmmnm

11、mnnmmd ydydyd uduaaa ybbbub udtdtdtdtdt12,Tnxx xx(1) , ,nTy yy 1223101121nnnnnxxxxxxxa xa xaxu 寫成矩陣形式：寫成矩陣形式： 輸出矩陣寫為輸出矩陣寫為 其中其中 xAxBu101100nnIAaaa001B 100yCxC當(dāng)右式包含導(dǎo)數(shù)項時，狀態(tài)方程形式為當(dāng)右式包含導(dǎo)數(shù)項時，狀態(tài)方程形式為 A,C與前相同，與前相同， 其中其中 xAxBuyCxDu10,TnBD1,(0,1,)iin injijjbain現(xiàn)代控制理論中還介紹了其它形式的轉(zhuǎn)換方程，現(xiàn)代控制理論中還介紹了其它形式的轉(zhuǎn)換方程，如能控標(biāo)準(zhǔn)型、

12、能觀標(biāo)準(zhǔn)型等如能控標(biāo)準(zhǔn)型、能觀標(biāo)準(zhǔn)型等 化連續(xù)狀態(tài)方程為離散狀態(tài)方程化連續(xù)狀態(tài)方程為離散狀態(tài)方程 連續(xù)狀態(tài)方程對應(yīng)的離散狀態(tài)表達(dá)式為連續(xù)狀態(tài)方程對應(yīng)的離散狀態(tài)表達(dá)式為 T為采樣周期或者計算步長為采樣周期或者計算步長 為確定為確定(T)和和H(T)，可利用連續(xù)狀態(tài)方程解，可利用連續(xù)狀態(tài)方程解 其中其中Ate為系統(tǒng)的矩陣指數(shù)或狀態(tài)轉(zhuǎn)移矩陣，為系統(tǒng)的矩陣指數(shù)或狀態(tài)轉(zhuǎn)移矩陣，x(0)為初始狀態(tài)向量為初始狀態(tài)向量 當(dāng)采用零階保持器時，當(dāng)采用零階保持器時， 即認(rèn)為即認(rèn)為u(t)在每個采樣周期內(nèi)保持常值在每個采樣周期內(nèi)保持常值 u(t)=u(kT)，( kTt(k+1)T ) (1) ( ) ()( ) (

13、)(0,1,)x kTT x kTH T u kTk 0( )(0)( )tAtAtAx te xeeBud(1)(1)(1)00(1)(1)0(1) (0)( )(1)(0)( )(2)(1)(2)(1) ()( )()()( ()()kTA kTA kTAkTAkTAkTAATkTATA kTAkTTATATAtATAx kTexeeBudx kTexeeBudex kTex kTeeBudex kTeeBu kT dtu kTex kTe與得：不變0()( ) ()( ) ()TBu kT dtT x kTH T u kT則有則有 其中其中和和H與與T有關(guān)，當(dāng)有關(guān)，當(dāng)T確定后，確定后，和

14、和H為常值矩陣為常值矩陣 離散化公式的核心在于計算矩陣指數(shù)及其積分，離散化公式的核心在于計算矩陣指數(shù)及其積分，常用級數(shù)展開的算法，即常用級數(shù)展開的算法，即 其余離散化的表示形式與連續(xù)形式之間的轉(zhuǎn)換其余離散化的表示形式與連續(xù)形式之間的轉(zhuǎn)換在計算機(jī)控制中介紹在計算機(jī)控制中介紹 22112!ATkkeIATA TA Tk第三節(jié)第三節(jié) 連續(xù)系統(tǒng)的仿真算法連續(xù)系統(tǒng)的仿真算法算法的基本概念算法的基本概念n系統(tǒng)模型系統(tǒng)模型計算機(jī)模型：二次建模，算法計算機(jī)模型：二次建模，算法是核心問題是核心問題n算法：解題方案的準(zhǔn)確而完整的描述，一算法：解題方案的準(zhǔn)確而完整的描述，一般采用文字、算式以及框圖的形式般采用文字、

15、算式以及框圖的形式n需要關(guān)注：需要關(guān)注：n算法性能分析：誤差、收斂性、計算效率等算法性能分析：誤差、收斂性、計算效率等n算法的比較與選擇算法的比較與選擇 n浮點數(shù)運(yùn)算浮點數(shù)運(yùn)算 n計算機(jī)上進(jìn)行數(shù)值計算時，實數(shù)計算機(jī)上進(jìn)行數(shù)值計算時，實數(shù)x用用t位十進(jìn)位十進(jìn)制浮點數(shù)表示：制浮點數(shù)表示：n其中其中m為為t 位十進(jìn)制小數(shù)，且位十進(jìn)制小數(shù)，且-1m1，c為為十進(jìn)制整數(shù)，若十進(jìn)制整數(shù)，若0.1m1，則稱此浮點數(shù)，則稱此浮點數(shù)系統(tǒng)為規(guī)格化的，系統(tǒng)為規(guī)格化的，t 稱該數(shù)的精度，特定的稱該數(shù)的精度，特定的計算機(jī)有固定的浮點數(shù)精度計算機(jī)有固定的浮點數(shù)精度 10cxm采用浮點數(shù)運(yùn)算存在的常見問題采用浮點數(shù)運(yùn)算存在

16、的常見問題n舍入誤差舍入誤差n計算機(jī)有一組操作浮點數(shù)的指令，用以模擬加、減、計算機(jī)有一組操作浮點數(shù)的指令，用以模擬加、減、乘、除運(yùn)算，但不可能精確。如乘法運(yùn)算時，乘積乘、除運(yùn)算，但不可能精確。如乘法運(yùn)算時，乘積應(yīng)有應(yīng)有2t位精度，但實際僅能保留位精度，但實際僅能保留t位，即存在舍入位，即存在舍入誤差。復(fù)雜計算（迭代等）中舍入誤差的累積可能誤差。復(fù)雜計算（迭代等）中舍入誤差的累積可能會影響結(jié)果，應(yīng)在算法分析中考慮會影響結(jié)果，應(yīng)在算法分析中考慮 n溢出溢出 n計算機(jī)對指數(shù)計算機(jī)對指數(shù)c范圍有限制，乘、除時可能會上溢、范圍有限制，乘、除時可能會上溢、下溢，也應(yīng)進(jìn)行處理下溢，也應(yīng)進(jìn)行處理 數(shù)值穩(wěn)定性問

17、題數(shù)值穩(wěn)定性問題 若運(yùn)算過程中計算誤差不斷增長，稱算法為數(shù)值不穩(wěn)定的若運(yùn)算過程中計算誤差不斷增長，稱算法為數(shù)值不穩(wěn)定的反之則為穩(wěn)定的反之則為穩(wěn)定的 例：計算例：計算 由分部積分得遞推公式：由分部積分得遞推公式： 用用Taylor展開計算：展開計算： 若取若取k=7，并保留，并保留4位小數(shù)，可得位小數(shù)，可得 截斷誤差：截斷誤差： 110(0,1,)nxnIex e dxn1111001,1xnnInIIee dxe 21( 1)( 1)1( 1)2!kek 10.3679e147110.3679108!4Re只考慮初值誤差，對只考慮初值誤差，對I00.6321，遞推計算，結(jié)果如表中第一行，遞推計

18、算，結(jié)果如表中第一行該積分不可能為負(fù)值，顯然算法有問題該積分不可能為負(fù)值，顯然算法有問題 分析計算誤差：分析計算誤差： 滿足關(guān)系：滿足關(guān)系： 誤差增長迅速誤差增長迅速 n01234567890.63210.36790.26420.20740.17040.14800.11200.2160-0.7287.5520.63210.36790.26430.20730.17080.14550.12680.11210.10350.0684nnnEII1nnEnE 0( 1) ( !)nnEn E 如果換一種算法，可以減小誤差，考慮積分估計值：如果換一種算法，可以減小誤差，考慮積分估計值： 逆向計算：逆向計算

19、： 取取n=9時，時， ，結(jié)果如表中第二行，結(jié)果如表中第二行 分析計算誤差，滿足關(guān)系：分析計算誤差，滿足關(guān)系： 顯然誤差一直減小顯然誤差一直減小 111110001011(min)(max)11xnxnnxxeeex dxIeex dxnn *11(1)nnIIn1911()0.06842 1010eI11nnEEn 病態(tài)問題病態(tài)問題 如線性方程組如線性方程組 精確解：精確解： 12123.0004.12715.41(1)1.0001.3745.147(2)xxxx1213.66586.2xx 如果用如果用4位有效數(shù)字進(jìn)行運(yùn)算：位有效數(shù)字進(jìn)行運(yùn)算：(2)-(1)/3.000，逐步計算后，可得逐

20、步計算后，可得x2=-5.000 ，誤差很大，誤差很大當(dāng)方程特征根相差太大時出現(xiàn)病態(tài)問題，也稱剛性當(dāng)方程特征根相差太大時出現(xiàn)病態(tài)問題，也稱剛性(Stiff)問題問題需要設(shè)計有效的算法需要設(shè)計有效的算法 二、數(shù)值積分法二、數(shù)值積分法n實際系統(tǒng)模型多為低階微分方程形式，求解實際系統(tǒng)模型多為低階微分方程形式，求解時本質(zhì)上應(yīng)用積分運(yùn)算時本質(zhì)上應(yīng)用積分運(yùn)算n對高階方程，可先轉(zhuǎn)換為多個一階方程，因?qū)Ω唠A方程，可先轉(zhuǎn)換為多個一階方程，因此，最終問題轉(zhuǎn)化為求解一階微分方程此，最終問題轉(zhuǎn)化為求解一階微分方程n常見算法：常見算法：nEuler法法 nRunge-Kutta法法(R-K法法)nAdams法法(多步法

21、多步法)Euler法法設(shè)有模型方程設(shè)有模型方程 ，初始條件，初始條件 歐拉法用歐拉法用tk點切線近似該點附近的曲線點切線近似該點附近的曲線f(t,y)，則有則有 其中其中是曲線上的點，是曲線上的點，是切線上的點是切線上的點 稱為第稱為第k 步的計算步長步的計算步長 此類方法稱為此類方法稱為“微分方程初值問題的數(shù)值計算法微分方程初值問題的數(shù)值計算法”，也稱也稱“數(shù)值積分法數(shù)值積分法” 優(yōu)點：簡單易行優(yōu)點：簡單易行 缺點：缺點：h取得大時，計算速度快，單步誤差大；取得大時，計算速度快，單步誤差大； h取得小，計算速度慢，且累計誤差大取得小，計算速度慢，且累計誤差大1()ky t1ky( , )dy

22、f t ydt00( )y ty111()( ,)()kkkkkkky tyyf tytt1kkkhttRunge-Kutta法法(R-K法法) 二階形式：二階形式： 其中其中 迭代公式由迭代公式由Taylor展開并保留展開并保留h2項獲得項獲得 注意：注意：R-K方法實質(zhì)是用均差代替導(dǎo)數(shù)，方法實質(zhì)是用均差代替導(dǎo)數(shù)， 其中其中k項的加權(quán)系數(shù)可任選項的加權(quán)系數(shù)可任選 1112()()2kkkhy tyykk121( ,)(,)kkkkkf tykf th yk h四階形式（固定步長）：四階形式（固定步長）：其中其中 四階形式在精度和復(fù)雜度方面都有較好的表現(xiàn)，最常用四階形式在精度和復(fù)雜度方面都有較

23、好的表現(xiàn)，最常用 Euler法與法與R-K法計算時僅用到前一步的結(jié)果，法計算時僅用到前一步的結(jié)果，稱單步法，已知初值后可自啟動稱單步法，已知初值后可自啟動 111234()(22)6kkkhy tyykkkk1213243( ,)(/2,/2)(/2,/2)(,)kkkkkkkkkf tykf thyk hkf thyk hkf th yk hAdams法法(多步法多步法) Euler法是用矩形公式（面積）近似定積分，在曲線與法是用矩形公式（面積）近似定積分，在曲線與矩形的邊之間有較大誤差；矩形的邊之間有較大誤差；Adams法考慮用梯形公式（面積）代替矩形公式法考慮用梯形公式（面積）代替矩形公

24、式 稱二階隱式稱二階隱式Adams法公式法公式 因公式右端包含未知項，不能直接求解，可以用迭代因公式右端包含未知項，不能直接求解，可以用迭代法求解，設(shè)有初值法求解，設(shè)有初值 迭代公式：迭代公式： 直至規(guī)定的精度直至規(guī)定的精度 1111( , ) (,)( ,)()22kktkkkkkkthhf t y dtf tyf tyff11()2kkkkhyyff(0)1ky( )(1)111 (,)( ,)2nnkkkkkkhyyf tyf ty隱式迭代計算步數(shù)太多，為此可降低精度，設(shè)計顯式隱式迭代計算步數(shù)太多，為此可降低精度，設(shè)計顯式Adams法，公式：法，公式： Adams法的統(tǒng)一形式為：法的統(tǒng)一

25、形式為： 算法特點：多步法，不能自啟動，隱式方法還需迭代求解算法特點：多步法，不能自啟動，隱式方法還需迭代求解實際應(yīng)用中，先用顯式法計算初值，再用隱式法校正一次，實際應(yīng)用中，先用顯式法計算初值，再用隱式法校正一次，稱預(yù)報稱預(yù)報-校正法校正法與與R-K法比較，同樣階次和精度下法比較，同樣階次和精度下Adams法計算次數(shù)較少法計算次數(shù)較少 1111()3 ( ,)(,)2kkkkkkkhy tyyf tyf ty111011kkkkNkNyyh B fB fBf算法分析算法分析n穩(wěn)定性分析穩(wěn)定性分析試驗方程：試驗方程： 若數(shù)值積分公式為若數(shù)值積分公式為 其中其中 是一個高階多項式函數(shù)是一個高階多項

26、式函數(shù) 僅當(dāng)僅當(dāng) 時算法才穩(wěn)定時算法才穩(wěn)定 如如Euler法穩(wěn)定條件為法穩(wěn)定條件為 隱式一階、二階隱式一階、二階Adams法恒穩(wěn)定，更高階條件穩(wěn)定法恒穩(wěn)定，更高階條件穩(wěn)定R-K法正好相反，階次越高穩(wěn)定域越大法正好相反，階次越高穩(wěn)定域越大 ,Re0dyyidt1()kkyp hy()p h()1p h11hn積分步長的選擇與控制積分步長的選擇與控制 兩個原則：保證穩(wěn)定性，要求一定的計算精度兩個原則：保證穩(wěn)定性，要求一定的計算精度受穩(wěn)定性限制，受穩(wěn)定性限制，h應(yīng)在系統(tǒng)中最小時間常數(shù)量級應(yīng)在系統(tǒng)中最小時間常數(shù)量級如如R-K4，要求步長小于系統(tǒng)中最小時間常數(shù)的，要求步長小于系統(tǒng)中最小時間常數(shù)的2.78

27、倍倍實際應(yīng)用中對大量的仿真計算可以考慮采用變步長法實際應(yīng)用中對大量的仿真計算可以考慮采用變步長法自動改變步長自動改變步長Matlab中中Simulink仿真時缺省的數(shù)值積分法為變步長法仿真時缺省的數(shù)值積分法為變步長法如如RKM3-4法法 首先進(jìn)行誤差估計：首先進(jìn)行誤差估計： 分別找一個三階和一個四階分別找一個三階和一個四階R-K公式：公式： 113451145(39126)30(4)6kkkkhyykkkkhyykkk其中其中1213124135134( ,)(/3,/3)(/3,/6/6)(/2,/83/8)(,/23/22)kkkkkkkkkkkf tykf thyk hkf thyk h

28、k hkf thyk hk hkf th yk hk hk h則誤差為則誤差為 1111345(298)30kkkhEyykkkk變步長策略：變步長策略： 設(shè)定一個最小誤差限設(shè)定一個最小誤差限 ，一個最大誤差限，一個最大誤差限 每一步的局部誤差取為每一步的局部誤差取為 ，第，第k+1步有效，下一步用步有效，下一步用2h積分；積分； ，保持，保持h不變；不變； ，第，第k+1步無效，步長變?yōu)椴綗o效，步長變?yōu)閔/2 minmax111/(| 1)kkkeEy1minkemin1maxke1maxke三、離散相似法三、離散相似法原理：將連續(xù)模型離散化后再仿真計算原理：將連續(xù)模型離散化后再仿真計算 形

29、式：形式：傳遞函數(shù)傳遞函數(shù)Z傳遞函數(shù)：傳遞函數(shù)：Z域離散相似模型，域離散相似模型， 狀態(tài)空間模型狀態(tài)空間模型離散狀態(tài)方程：時域離散相似模型，離散狀態(tài)方程：時域離散相似模型，右邊第三項表示使用一階保持器增加的項右邊第三項表示使用一階保持器增加的項( )( ) ( )hG zZ G s G s(1) ( ) ()( ) ()( ) ()mmx kTT x kTT u kTT u kT xAxBu 由采樣定理，為使離散相似模型中重構(gòu)的信號能精由采樣定理，為使離散相似模型中重構(gòu)的信號能精確表示原信號，應(yīng)有采樣時間小于系統(tǒng)最小時間常確表示原信號，應(yīng)有采樣時間小于系統(tǒng)最小時間常數(shù)的一半，或者采樣頻率是最大

30、信號頻率的兩倍數(shù)的一半，或者采樣頻率是最大信號頻率的兩倍離散相似法的優(yōu)點：不易受模型方程特性的影響，離散相似法的優(yōu)點：不易受模型方程特性的影響，尤其對特征根相差較大的系統(tǒng)十分有效，計算速度尤其對特征根相差較大的系統(tǒng)十分有效，計算速度也更快；但有的模型不易離散化也更快；但有的模型不易離散化 的計算誤差；的計算誤差；u(t)在采樣間隔中的近似處理在采樣間隔中的近似處理對對后者后者，當(dāng)輸入是典型函數(shù)時，可通過增廣矩陣，當(dāng)輸入是典型函數(shù)時，可通過增廣矩陣法消除誤差法消除誤差 誤差處理誤差處理考慮時域離散相似法，誤差來源：考慮時域離散相似法，誤差來源：ATe狀態(tài)方程狀態(tài)方程 的解為：的解為： 對典型函數(shù)

31、，考慮將輸入對典型函數(shù)，考慮將輸入u(t)增廣到狀態(tài)向量增廣到狀態(tài)向量x(t)中，中，得得 齊次解，避免積分近似出現(xiàn)的誤差齊次解，避免積分近似出現(xiàn)的誤差0( )(0)( )tAtAtAx te xeeBudxAxBu( )(0)Atx te x增廣矩陣法示例增廣矩陣法示例假設(shè)系統(tǒng)為假設(shè)系統(tǒng)為n階，模型階，模型xAxBuyCx0(0)xx對階躍輸入對階躍輸入0( )1( )u tUt定義定義1( )( )nxtu t，則，則1( )0nxt增廣后的狀態(tài)方程及輸出方程為增廣后的狀態(tài)方程及輸出方程為1100nnxxABxx10nxyCx初始條件初始條件010(0)(0)nxxxU增廣矩陣法示例增廣矩

32、陣法示例系統(tǒng)模型同上系統(tǒng)模型同上，對斜坡輸入，對斜坡輸入0( )u tU t定義：定義：1( )( )nxtu t210( )( )( )nnxtxtu tU則則2( )0nxt，增廣后的狀態(tài)方程及輸出方程為，增廣后的狀態(tài)方程及輸出方程為11220001000nnnnxxABxxxx1200nnxyCxx初始條件初始條件0102(0)(0)0(0)nnxxxUx課堂作業(yè)第四節(jié)第四節(jié) 采樣控制系統(tǒng)仿真采樣控制系統(tǒng)仿真n典型采樣控制系統(tǒng)結(jié)構(gòu)圖典型采樣控制系統(tǒng)結(jié)構(gòu)圖 其中信號比較環(huán)節(jié)可在控制器內(nèi)、外進(jìn)行其中信號比較環(huán)節(jié)可在控制器內(nèi)、外進(jìn)行 n本結(jié)構(gòu)與離散相似法得到的系統(tǒng)相似本結(jié)構(gòu)與離散相似法得到的系

33、統(tǒng)相似（被控對象連續(xù)，有采樣器、保持器）（被控對象連續(xù)，有采樣器、保持器）n但但n前者采樣周期、采樣器位置、保持器類型均前者采樣周期、采樣器位置、保持器類型均為實際存在，而后者均為虛擬的為實際存在，而后者均為虛擬的n前者在仿真時需要考慮仿真步距與實際采樣前者在仿真時需要考慮仿真步距與實際采樣周期的關(guān)系，還需要處理離散和連續(xù)部分所周期的關(guān)系，還需要處理離散和連續(xù)部分所得差分模型之間的聯(lián)系，后者直接離散化即得差分模型之間的聯(lián)系，后者直接離散化即可可一、采樣周期與仿真步距一、采樣周期與仿真步距采樣控制系統(tǒng)方塊圖采樣控制系統(tǒng)方塊圖其中其中G(s)為被控對象傳遞函數(shù)，為被控對象傳遞函數(shù)，H(s)為保持器

34、為保持器傳遞函數(shù)，傳遞函數(shù)，D(z)為數(shù)字控制器的為數(shù)字控制器的z傳遞函數(shù)，傳遞函數(shù)，Ts為實際采樣周期為實際采樣周期仿真步距仿真步距T的選擇必須根據(jù)被控對象結(jié)構(gòu)、采的選擇必須根據(jù)被控對象結(jié)構(gòu)、采樣周期大小、保持器類型及仿真精度和仿真速樣周期大小、保持器類型及仿真精度和仿真速度的要求綜合考慮度的要求綜合考慮 仿真步距的選擇（仿真步距的選擇（1）n仿真步距仿真步距T采樣周期采樣周期Tsn要求：要求：Ts較小，系統(tǒng)階次較低，仿真要求不較小，系統(tǒng)階次較低，仿真要求不高高n此時連續(xù)部分此時連續(xù)部分H(s)G(s)部分不再增加虛擬部分不再增加虛擬采樣保持器采樣保持器n為此必須計算為此必須計算G(z)=Z

35、H(s)G(s)（常用狀態(tài)空間表達(dá)式表示，易于計算）（常用狀態(tài)空間表達(dá)式表示，易于計算）仿真步距的選擇（仿真步距的選擇（2）n仿真步距仿真步距T采樣周期采樣周期Ts 更常見，因為：更常見，因為：nTs往往受軟硬件限制必須取較大往往受軟硬件限制必須取較大n連續(xù)部分有硬非線性時，系統(tǒng)往往分割處連續(xù)部分有硬非線性時，系統(tǒng)往往分割處理，需要引入更多采樣保持器，產(chǎn)生較大理，需要引入更多采樣保持器，產(chǎn)生較大的幅值和相位誤差，為保證精度，必須使的幅值和相位誤差，為保證精度，必須使TTs n一般為計算方便取一般為計算方便取Ts=NT，N為正整數(shù)，則仿為正整數(shù)，則仿真時每一次大循環(huán)中連續(xù)部分應(yīng)計算真時每一次大循

36、環(huán)中連續(xù)部分應(yīng)計算N次次n系統(tǒng)中還可能有不同頻率的實際采樣開關(guān)系統(tǒng)中還可能有不同頻率的實際采樣開關(guān) 如智能汽車自動駕駛系統(tǒng)，其中駕駛盤執(zhí)行轉(zhuǎn)如智能汽車自動駕駛系統(tǒng)，其中駕駛盤執(zhí)行轉(zhuǎn)角的局部反饋為內(nèi)反饋，汽車位置的偏差反饋角的局部反饋為內(nèi)反饋，汽車位置的偏差反饋為外反饋，前者因執(zhí)行機(jī)構(gòu)固定頻率高，采樣為外反饋，前者因執(zhí)行機(jī)構(gòu)固定頻率高，采樣周期短，后者則較長周期短，后者則較長 二、改變數(shù)字控制器的采樣間隔二、改變數(shù)字控制器的采樣間隔 n有時實際有時實際Ts較小，為分析系統(tǒng)取較大較小，為分析系統(tǒng)取較大T，須重，須重新求差分模型新求差分模型n原理：原理：S平面上，有相同零極點和穩(wěn)態(tài)值的系平面上，有相

37、同零極點和穩(wěn)態(tài)值的系統(tǒng)等價，故可由統(tǒng)等價，故可由Z域脈沖傳遞函數(shù)映射到域脈沖傳遞函數(shù)映射到S平面平面上，再按新的上，再按新的Ts*映射到映射到z平面平面 例：有數(shù)字控制器例：有數(shù)字控制器0.98( )2.620.64zD zzTs=0.04s，取Ts*=0.1s仿真，求D(z)00000.64ln()/11.16ln()/5.5050.98pppszzszzszTszTz 0.32770.9508s ps zT spT szzeze 0.9508( )0.3277zzD zkz再根據(jù)穩(wěn)態(tài)值相等原則確定再根據(jù)穩(wěn)態(tài)值相等原則確定kz，注意與輸入信號有關(guān)，注意與輸入信號有關(guān)本例中若輸入單位階躍信號，

38、由終值定理本例中若輸入單位階躍信號，由終值定理11( )lim( )0.145561zzzuD zzz 同樣同樣10.9508lim0.145561.9890.3277zzzzkkz第五節(jié)第五節(jié) 分布參數(shù)系統(tǒng)仿真分布參數(shù)系統(tǒng)仿真n物理系統(tǒng)宏觀均有空間分布特性，需要物理系統(tǒng)宏觀均有空間分布特性，需要用偏微分方程描述，不宜用解析法求解，用偏微分方程描述，不宜用解析法求解，通常用離散化模型描述，采用數(shù)值計算通常用離散化模型描述，采用數(shù)值計算方法。方法。n分布參數(shù)系統(tǒng)的求解方法較少，僅限于分布參數(shù)系統(tǒng)的求解方法較少，僅限于有限差分法和有限元法，另外還可以用有限差分法和有限元法，另外還可以用物理集總參數(shù)

39、法分隔物理空間并建立一物理集總參數(shù)法分隔物理空間并建立一組聯(lián)立的常微分方程組組聯(lián)立的常微分方程組 n對于常見的連續(xù)系統(tǒng)來說，通常在使用對于常見的連續(xù)系統(tǒng)來說，通常在使用數(shù)學(xué)模型時人們更關(guān)注的是模型對對象數(shù)學(xué)模型時人們更關(guān)注的是模型對對象的外部特性的描述，因此，雖然復(fù)雜的的外部特性的描述，因此，雖然復(fù)雜的對象往往是分布參數(shù)系統(tǒng)，但可以按照對象往往是分布參數(shù)系統(tǒng)，但可以按照對象各測點的位置、結(jié)構(gòu)特點及其物理對象各測點的位置、結(jié)構(gòu)特點及其物理化學(xué)等特性進(jìn)行適當(dāng)分區(qū)，分區(qū)內(nèi)采用化學(xué)等特性進(jìn)行適當(dāng)分區(qū)，分區(qū)內(nèi)采用集中參數(shù)系統(tǒng)模型描述，從而可以避開集中參數(shù)系統(tǒng)模型描述，從而可以避開非線性偏微分方程的迭代求

40、解過程。非線性偏微分方程的迭代求解過程。 例：扭振桿系統(tǒng)例：扭振桿系統(tǒng)n當(dāng)在自由端施加當(dāng)在自由端施加輸入扭矩輸入扭矩T(L,t)后，后，y點產(chǎn)生的輸出轉(zhuǎn)點產(chǎn)生的輸出轉(zhuǎn)角角(y,t)，求輸，求輸入扭矩與輸出轉(zhuǎn)入扭矩與輸出轉(zhuǎn)角間的表達(dá)式角間的表達(dá)式 =0(y,t)ydyLy=0T(L,t)考慮厚度為考慮厚度為dy的一段扭桿，其力矩平衡方程為：的一段扭桿，其力矩平衡方程為：222222222PPPPR dyRI GI GI GdyIdyyyytt 上段力矩下段力矩其中其中42PIR為圓截面的極慣性矩為圓截面的極慣性矩2222tGy其中其中是桿的線密度，是桿的線密度，G G是材料的剪切彈性模量。是材料

41、的剪切彈性模量。這是一個物理學(xué)經(jīng)典方程，直接求解非常繁瑣。這是一個物理學(xué)經(jīng)典方程，直接求解非常繁瑣。 得一維波動方程：得一維波動方程：1解析法解析法 采用頻域分析，先對采用頻域分析，先對t進(jìn)行拉氏變換，有：進(jìn)行拉氏變換，有：2222()()ttLLyGt222( , )( , )( , )( , )y ssy ssy sy syG取初始條件為取初始條件為0，則有，則有2220dsdyG解的形式為解的形式為12( , )expexp(1)y scsGycsGy 由于桿的固定端有由于桿的固定端有(0,t)=0，則，則(0,s)=0， 可得可得 c1= - c2 (2)在在y=L處，處，( , )(

42、 , )( , ), ( , )(3)PPL tL sT L tI GT L sI Gyy(1)式對式對y求偏導(dǎo)：求偏導(dǎo)：12expexp(4)scsyscsyyGGGG令令y=L，由，由(2),(3),(4)可求出可求出c1，c2（略），最終得（略），最終得expexp( , )( , )(expexp)sinh()1cosh()PPsGysGyy sT L sIGssGysGysGyIGssGy令令s=jw可得可得y=L處角運(yùn)動的頻率響應(yīng)處角運(yùn)動的頻率響應(yīng)sin()( ,)1(5)( ,)cos()PLGL jT L jIGLG顯然其中有無限多個顯然其中有無限多個固有頻率：固有頻率：,(1

43、,3,5,)2kkGkL若結(jié)構(gòu)以其中某一固有頻率振動，則此時的動態(tài)若結(jié)構(gòu)以其中某一固有頻率振動，則此時的動態(tài)扭轉(zhuǎn)曲線稱為振型扭轉(zhuǎn)曲線稱為振型(5)式可用來求振型式可用來求振型( ,)( ,)( ,)sin()( ,)( ,)( ,)2kkkkkky jy jT L jkyL jT L jL jL 2物理集總參數(shù)法物理集總參數(shù)法 42122/iKKKr G L扭振系統(tǒng)可用物理集總參數(shù)法離散形式（集總塊）扭振系統(tǒng)可用物理集總參數(shù)法離散形式（集總塊）近似，集總塊數(shù)目可通過經(jīng)驗或試驗的方法確定。近似，集總塊數(shù)目可通過經(jīng)驗或試驗的方法確定。設(shè)為設(shè)為2個集總塊個集總塊選擇每個集總塊質(zhì)心作為集總選擇每個集總

44、塊質(zhì)心作為集總慣量慣量J 所在位置，并以此確定所在位置，并以此確定K（質(zhì)心之間桿長的彈性），有（質(zhì)心之間桿長的彈性），有412/4JJLr 原系統(tǒng)可以表示為一原系統(tǒng)可以表示為一個兩段集總模型：個兩段集總模型：24242141421442)(4)( LrLGrLGrLrLGrTi可求得兩個固有頻率可求得兩個固有頻率121.533.70GLGL而準(zhǔn)確值為而準(zhǔn)確值為1.57和和4.71當(dāng)集總數(shù)量增大時，可預(yù)測出更多的固有頻率，數(shù)當(dāng)集總數(shù)量增大時，可預(yù)測出更多的固有頻率，數(shù)值也更精確值也更精確 3. 有限差分法有限差分法 對應(yīng)物理離散法有數(shù)學(xué)離散法，如有限差分法，對應(yīng)物理離散法有數(shù)學(xué)離散法，如有限差分

45、法，其中常用的為中心差分法其中常用的為中心差分法設(shè)設(shè)y=f(x,t)，當(dāng)，當(dāng)t為常值時，中心差分法即用切點為常值時，中心差分法即用切點Pn上的中心差分近似函數(shù)曲線的斜率上的中心差分近似函數(shù)曲線的斜率112nnnyyyxh二階導(dǎo)數(shù)類似二階導(dǎo)數(shù)類似111122,nnnnnnyyyyyyxhxh即一階導(dǎo)數(shù)為即一階導(dǎo)數(shù)為2111122222/nnnnnnyyyyyyhxxxh對于具有兩個位置變量的函數(shù)對于具有兩個位置變量的函數(shù)f(x,y)，偏導(dǎo)數(shù)為，偏導(dǎo)數(shù)為(,)(,)2iiiif xh yf xh yfxh簡記簡記f(xi+h,yi)=fi+1,j則：則：21,1,1,1,222,1,1,1,12221,11,11,11,122,22,24ijijiji jiji ji ji ji ji jijijijijfffffffxhxhfffffffyhyhfffffx yh 利用中心差分公式可將空間變量離散化，利用中心差分公式可將空間變量離散化，結(jié)果是把偏微分方程轉(zhuǎn)化為常微分方程結(jié)果是把偏微分方程轉(zhuǎn)化為常微分方程 扭振桿分析扭振桿分析