控制系統(tǒng)計(jì)算機(jī)輔助設(shè)計(jì)實(shí)驗(yàn)報(bào)告

1、控制系統(tǒng)計(jì)算機(jī)輔助設(shè)計(jì)實(shí)驗(yàn)報(bào)告姓名：學(xué)號(hào)：學(xué)院：自動(dòng)化學(xué)院專業(yè)：自動(dòng)化2013-11實(shí)驗(yàn)一、實(shí)驗(yàn)要求:1、用matlab語(yǔ)言求下列系統(tǒng)的狀態(tài)方程、傳遞函數(shù)、零極點(diǎn)增益、和部分分 式形式的模型參數(shù)，并分別寫(xiě)出其相應(yīng)的數(shù)學(xué)模型表達(dá)式：s3 + 752 + 24s + 24 兒町=卩+麗+ 35辭+莎+可'f 2.25-5-1.25-0.54225-4.25-1.250,252i(i)=力+叩)0.25-0,51.25-121.25-1.750衛(wèi)5一 0.750 i y(Z) =02 (1 2x(/)2、 用歐拉法求下面系統(tǒng)的輸出響應(yīng) y(t)在0< t < 1 上, h=0.1

2、時(shí)的數(shù)值。 y ' = -y, y(0) =1要求保留4位小數(shù)，并將結(jié)果與真解y(t) = e -t比較。3、用二階龍格庫(kù)塔法求解2的數(shù)值解，并于歐拉法求得的結(jié)果比較。.、實(shí)驗(yàn)步驟：1、求(1)的M文件如下:clear;n um=1 7 24 24;den=1 10 35 50 24;sys=tf( num,de n)A,B,C,D=tf2ss( num,de n)Z,P,K=tf2zp( num,de n)R,P,H=residue( num,de n)1.1系統(tǒng)系數(shù)矩陣A,系統(tǒng)輸入矩陣B,系統(tǒng)輸出矩陣C,直接傳輸矩陣D分別為:-1035-50-241000Q1QQ00I0B =L0

3、CQD =所以系統(tǒng)的狀態(tài)方程為：x(t)=A x(t) +B u (t) ; y (t) =C x (t)1.2P =l 二-4.QOOOK 二-2_ 7306 4- 2.S5311-3.0000-2.7306 - 2.3531i-2.OOCO1-1_ 5388'1.0000零極點(diǎn)增益模型：(s)=【( s+2.7306-2.8531i )(s+2.7306+2.8531i ) (s+1.5388 )】 /【(s+4)( s+3)( s+2)( s+1)l1.3系統(tǒng)零點(diǎn)向量乙 極點(diǎn)向量P,系數(shù)H分別為：F =R 二4.0000-3-0000TT-6_ OCOOH =2-OCOO-2.Q

4、QQQ1.0000-L0000口部分分式形式：G(s)=4/ (s+4) -6/ (s+3) +2/ (s+2) +1/ (s+1)2.求(2)的M文件如下:clear;a=2.25,-5,-1.25,-0.5;2.25,-4.25,-1.25,-0.25;0.25,-0.5,-1.25,-1;1.25,-1.75,-0.25,-0.75; b=4;2;2;0;c=0,2,0,2;d=0;sys=ss(a,b,c,d)nu m,de n=ss2tf(a,b,c,d)Z,RK=ss2zp(a,b,c,d)R,PH=residue (nu m,de n)2.1nixm.=04.000314*000

5、022. OO'OO15.OQOOden.=LOCOO1 00006.2B005. 25002. 2500傳遞函數(shù)模型參數(shù)：G(S)=(4 sA3 + 14 s A2+ 22 s + 15)/(sA4 + 4 sA3 + 6.25 s A2+ 5.25 s + 2.25)2.2系統(tǒng)零點(diǎn)向量乙 極點(diǎn)向量P,系數(shù)K分別為:-0.5000 + 0. 86601-L 0000 十 h 224 71-1.0000 - 1. 22471-1. 5000-O.SOOO - 0.secoi-1.50004.COCO)(s+1+1.2247i1.5000零極點(diǎn)增益模型參數(shù)：G(s)=【4 (s+1-1.

6、2247i(s+0.5-0.866i)( s+0.5+0.866i s+1.5)】2.3R =p 二4. 0090'.5000-0.oooo-1.5000H =0.0000 - Z 3O94i-0. 5000 + 0. 8660i10.0000 十 2.30941-0,5000 - 0.86601部分分式形式的模型參數(shù)：G ( s)=4/( s+1.5)-2.3094i/(s+0.5-0.866i)+2.3094i/(s+0.5+0.866i)3原理：把f(t,y)M文件如下：在t k, yk區(qū)間內(nèi)的曲邊面積用矩形面積近似代替cleary=1;h=0.1;j=0;for i=1:11j

7、=j+1; a(j)=y y=y+h*(-y); end j=0； for i=0:0.1:1 f=exp(-i); j=j+1； b(j)=f；end figure(1) x=0:0.1:1; a b plot(x,a,'y-*') hold on plot(x,b,'-ro') 得到圖形：0,0.0.0.0.0.00. 20. 40. 60飛使用歐拉法得到的結(jié)果和真值對(duì)比:歐 拉10.90000.81000.72900.65610.59050.53140.47830.43050.38740.3487真10.90480.81870.74080.67030.60

8、650.54880.49660.44930.40660.3679值誤0-0.004-0.000-0.011-0.014-0.016-0.017-0.018-0.018-0.019-0.019差8782043822顯然誤差與h2為同階無(wú)窮小，歐拉法具有一階計(jì)算精度，精度較低，但算法簡(jiǎn)單。4. 原理：把f(t,y) 在t k，yQ區(qū)間內(nèi)的曲邊面積用上下底為f k和fk+i、高為h的 梯形面積近似代替。M文件如下：clear;y=i;h=0.1;j=0；for i=1:11j=j+1；a(j)=yk仁 _y;k2=-(y+0.5*h*k1);y=y+h*k2;endj=0;for i=0:0.1:1

9、f=exp(-i);j=j+1;b(j)=f;endfigure(2)x=0:0.1:1;abplot(x,a,'y-*')hold onplot(x,b,'-ro')得到圖形：0,90.0.a0,0.0. 60.8比較歐拉法與二階龍格-庫(kù)塔法求解.直/、值10.90480.81870.74080.67030.60650.54880.49660.44930.40660.3679龍 庫(kù)10.90500.81900.74120.67080.60710.54940.49720.45000.40720.3685誤 差00.00020.00030.00040.00050.

10、00060.00060.00060.00070.00060.0006明顯誤差為h3得同階無(wú)窮小，具有二階計(jì)算精度，而歐拉法具有以階計(jì)算精度,二階龍格-庫(kù)塔法比歐拉法計(jì)算精度咼。三、實(shí)驗(yàn)總結(jié)：此次實(shí)驗(yàn)只要平時(shí)上課認(rèn)真聽(tīng)過(guò)課，參考課件和書(shū)本便能順利完成實(shí)驗(yàn)。由 此實(shí)驗(yàn)也可以總結(jié)出很多問(wèn)題都會(huì)有多種解法，我們要通過(guò)實(shí)踐總結(jié)出最佳解 法。實(shí)驗(yàn)二實(shí)驗(yàn)內(nèi)容:1、用四階龍格-庫(kù)塔法求解題2-3數(shù)值解，并與前兩題結(jié)果相比較2、已知二階系統(tǒng)狀態(tài)方程為(1)(2)ii(f)帀(0)吃(0)61 012(721° 2220寫(xiě)出取計(jì)算步長(zhǎng)為h時(shí)，該系統(tǒng)狀態(tài)變量忑=xj (£), X2(t)Hl(t

11、)2(t)咱);的四階龍格-庫(kù)塔法遞推關(guān)系式令上式中u （t）=0，用試探法選取參數(shù)L一S'一a12石1021 2220 帶入（a）所得公式，給出仿真圖形。要求選取兩組參數(shù)，一組使系統(tǒng)穩(wěn) 定，一組使系統(tǒng)發(fā)散。（注：系統(tǒng)穩(wěn)定從仿真圖形上看，可視為系統(tǒng)的狀 態(tài)曲線x（ t）趨于一定的值，發(fā)散可視為系統(tǒng)的狀態(tài)曲線 x（ t）趨于無(wú) 窮，當(dāng)時(shí)間t趨于無(wú)窮時(shí)。）實(shí)驗(yàn)步驟:1. 求四階龍格-庫(kù)塔方法求解函數(shù)數(shù)值解:M文件：clear;y=1；h=0.1;j=0；for i=1:11j=j+1;a(j)=yk仁 _y;k2=-(y+0.5*h*k1);k3=-(y+0.5*h*k2);k4=-(y+

12、h*k3);y=y+h/6*(k1+2*k2+2*k3+k4);endj=0；for i=0:0.1:1 f=exp(-i); j=j+1; b(j)=f；end figure(3) x=0:0.1:1; a b plot(x,a,'y-*') hold on plot(x,b,'-ro')10.90.30,70.60.5a. 4得到圖形:對(duì)于四階龍格-庫(kù)塔方法:真 值10.90480.81870.74080.67030.60650.54880.49660.44930.40660.3679龍庫(kù)10.90480.81870.74080.67030.60650.54

13、880.49660.44930.40660.3679誤差00000000000四階龍格-庫(kù)塔法得到的結(jié)果與真值完全重合，所以四階龍格庫(kù)塔法求解精度高 于二階龍格庫(kù)塔法，二階龍格庫(kù)塔法求解精度高于歐拉法。2 當(dāng) u（ t）=0 時(shí):M源程序：clear;h=0.1;i=1;j=1；x=2;1;A=-1,0；2,-2for t=O:h:1Odisp(x);k1=A*x;k2=A*(x+k1*h/2);k3=A*(x+k2*h/2);k4=A*(x+k3*h);M(i)=x(1,:);T(i)=x(2,:);i=i+1j=j+1m=xX=m+h/6*(k1+2*k2+3*k3+k4);endeig(

14、A)x=0:h:10;plot(x,M)hold onplot(x,T)0.5616得到結(jié)果：特征根 ans =-3.5616圖像：將A改為卜1,0；2,-2得到：特征根ans = -2, -1圖形為：三、實(shí)驗(yàn)總結(jié):此次實(shí)驗(yàn)需要耐心調(diào)整矩陣 A的值，并且h需要設(shè)置合適的大小，才能保證圖形的圓滑。實(shí)驗(yàn)三實(shí)驗(yàn)內(nèi)容：1、針對(duì)2-6中問(wèn)題（b）,對(duì)所選取的使系統(tǒng)發(fā)散的一組參數(shù)，設(shè)置控制 u（t）=Kx（t）使系統(tǒng)穩(wěn)定，其中K可以設(shè)計(jì)為一個(gè)常數(shù)（一般而言是個(gè) 負(fù)數(shù)）或者為一個(gè)2*2的矩陣（一般而言其特征值均為負(fù)）。2、將上述控制系統(tǒng)在 Matlab/Simulink 平臺(tái)上進(jìn)行仿真，并選取不同的 仿真算

15、法，比較所得的結(jié)果。（注：這里的不同仿真算法是指，在Simulink 仿真參數(shù)配置對(duì)話框中分別選?。憾ú介L(zhǎng)和變步長(zhǎng)進(jìn)行仿真，在定步長(zhǎng)中 又可以分為歐拉法，或其他，變步長(zhǎng)中也可以選擇其他算法，并比較不同 的仿真算法對(duì)仿真結(jié)果的影響。）二、實(shí)驗(yàn)步驟：。在Simulink下建立系統(tǒng)框圖如下：X2; 1；A-1, 2； 2, -2；B1; 1；K=-1, -1在Simulink仿真參數(shù)配置對(duì)話框中分別選取不同算法：定步長(zhǎng)的Euler 法、Runge-Kutta 法;變步長(zhǎng)的 Adams 法、Bogacki-Shampine 法、 Dormand-Prinee法。其中定步長(zhǎng)時(shí)步長(zhǎng)為 0.2。變步長(zhǎng)模式可以在仿真的過(guò)程中改變步長(zhǎng)，提供誤差控制和過(guò)零檢 測(cè)。固定步長(zhǎng)模式在仿真過(guò)程中提供固定的步長(zhǎng)，不提供誤差控制和過(guò)零 檢測(cè)。1. 1定步 長(zhǎng)Euler法如圖：1. 2定步 長(zhǎng) Runge-Kutta 法:對(duì)于定步長(zhǎng)分析可知，定步長(zhǎng)Runge -K