控制系統(tǒng)數(shù)字仿真與CAD習(xí)題

1、控制系統(tǒng)計(jì)算機(jī)仿真及輔助設(shè)計(jì)課程上機(jī)作業(yè)1月5日上機(jī)作業(yè)一：2-1 思考題：（1）數(shù)學(xué)模型的微分方程，狀態(tài)方程，傳遞函數(shù)，零極點(diǎn)增益和部分分式五種形式，各有什么特點(diǎn)？（2）數(shù)學(xué)模型各種形式之間為什么要互相轉(zhuǎn)換？（3）控制系統(tǒng)建模的基本方法有哪些？他們的區(qū)別和特點(diǎn)是什么？（4）控制系統(tǒng)計(jì)算機(jī)仿真中的“實(shí)現(xiàn)問題”是什么含意？（5）數(shù)值積分法的選用應(yīng)遵循哪幾條原則？（1）微分方程是直接描述系統(tǒng)輸入和輸出量之間的制約關(guān)系，是連續(xù)控制系統(tǒng)其他數(shù)學(xué)模型表達(dá)式的基礎(chǔ)。狀態(tài)方程能夠反映系統(tǒng)內(nèi)部各狀態(tài)之間的相互關(guān)系，適用于多輸入多輸出系統(tǒng)。傳遞函數(shù)是零極點(diǎn)形式和部分分式形式的基礎(chǔ)。零極點(diǎn)增益形式可用于分析系統(tǒng)的

2、穩(wěn)定性和快速性。利用部分分式形式可直接分析系統(tǒng)的動(dòng)態(tài)過(guò)程。（2）不同的控制系統(tǒng)的分析和設(shè)計(jì)方法，只適用于特定的數(shù)學(xué)模型形式。（3）控制系統(tǒng)的建模方法大體有三種：機(jī)理模型法，統(tǒng)計(jì)模型法和混合模型法。機(jī)理模型法就是對(duì)已知結(jié)構(gòu)，參數(shù)的物理系統(tǒng)運(yùn)用相應(yīng)的物理定律或定理，經(jīng)過(guò)合理的分析簡(jiǎn)化建立起來(lái)的各物理量間的關(guān)系。該方法需要對(duì)系統(tǒng)的內(nèi)部結(jié)構(gòu)和特性完全的了解，精度高。統(tǒng)計(jì)模型法是采用歸納的方法，根據(jù)系統(tǒng)實(shí)測(cè)的數(shù)據(jù)，運(yùn)用統(tǒng)計(jì)規(guī)律和系統(tǒng)辨識(shí)等理論建立的系統(tǒng)模型。該方法建立的數(shù)學(xué)模型受數(shù)據(jù)量不充分，數(shù)據(jù)精度不一致，數(shù)據(jù)處理方法的不完善，很難在精度上達(dá)到更高的要求。混合法是上述兩種方法的結(jié)合。（4）“實(shí)現(xiàn)問題”

3、就是根據(jù)建立的數(shù)學(xué)模型和精度，采用某種數(shù)值計(jì)算方法，將模型方程轉(zhuǎn)換為適合在計(jì)算機(jī)上運(yùn)行的公式和方程，通過(guò)計(jì)算來(lái)使之正確的反映系統(tǒng)各變量動(dòng)態(tài)性能，得到可靠的仿真結(jié)果。（5）數(shù)值積分法應(yīng)該遵循的原則是在滿足系統(tǒng)精度的前提下，提高數(shù)值運(yùn)算的速度和并保證計(jì)算結(jié)果的穩(wěn)定。2-2.用matlab語(yǔ)言求下列系統(tǒng)的狀態(tài)方程、傳遞函數(shù)、零極點(diǎn)增益、和部分分式形式的模型參數(shù)，并分別寫出其相應(yīng)的數(shù)學(xué)模型表達(dá)式：(1) G（s）= (2) =y=0 2 0 2 X解：（1）程序代碼如下：clear;clc;num=1 7 24 24;den=1 10 35 50 24;Z,P,K=tf2zp(num,den);A,B

4、,C,D=zp2ss(Z,P,K);R,P1,H=residue(num,den);G1=tf(num,den)G2=zpk(Z,P,K)G3=ss(A,B,C,D)程序的輸出結(jié)果如下： Transfer function: s3 + 7 s2 + 24 s + 24-s4 + 10 s3 + 35 s2 + 50 s + 24 Zero/pole/gain:(s+1.539) (s2 + 5.461s + 15.6)-(s+4) (s+3) (s+2) (s+1)a = x1 x2 x3 x4 x1 -3 -1.414 0 0 x2 1.414 0 0 0 x3 1 1.088 -7 -3.

5、464 x4 0 0 3.464 0 b = u1 x1 1 x2 0 x3 0 x4 0 c = x1 x2 x3 x4 y1 1 1.088 -1.539 1.038 d = u1 y1 0 Continuous-time model.（2）程序代碼如下：clear;clc;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;num,den=ss2tf(A,B,C,D);Z,P,K=ss2zp(A,B,C,D);R,P1,H

6、=residue(num,den);G1=tf(num,den)G2=zpk(Z,P,K)G3=ss(A,B,C,D)程序輸出結(jié)果如下：Transfer function: 4 s3 + 14 s2 + 22 s + 15-s4 + 4 s3 + 6.25 s2 + 5.25 s + 2.25 Zero/pole/gain:4 (s+1.5) (s2 + 2s + 2.5)- (s+1.5)2 (s2 + s + 1) a = x1 x2 x3 x4 x1 2.25 -5 -1.25 -0.5 x2 2.25 -4.25 -1.25 -0.25 x3 0.25 -0.5 -1.25 -1 x4

7、 1.25 -1.75 -0.25 -0.75 b = u1 x1 4 x2 2 x3 2 x4 0 c = x1 x2 x3 x4 y1 0 2 0 2 d = u1 y1 0 Continuous-time model.2-3.用歐拉法求下面系統(tǒng)的輸出響應(yīng)y(t)在0t1上，h=0.1時(shí)的數(shù)值。 要求保留4位小數(shù)，并將結(jié)果與真解比較。解：程序代碼如下：clear;clc;h=0.1;y=1;disp('歐拉法求得值為：');for i=0:h:1 m=y; disp(y); y=m-m*h;enddisp('真解為：');for i=0:h:1 y=exp(

8、-i); disp(y);end輸出結(jié)果如下：歐拉法值10.90000.81000.72900.65610.59050.53140.47830.43050.38740.3487真值10.90480.81870.74080.67030.60650.54880.49660.44930.40660.36792-4用二階龍格庫(kù)塔法求解2-3的數(shù)值解，并于歐拉法求得的結(jié)果比較。解：程序代碼如下：clear;clc;h=0.1;y=1;disp('二階龍格庫(kù)塔法求得值為：');for i=0:h:1 disp(y); k1=-y; k2=-(y+k1*h); y=y+(k1+k2)*h/2

9、;enddisp('歐拉法求得值為：'); y=1;for i=0:h:1 m=y; disp(y); y=m-m*h;end程序輸出結(jié)果如下：歐拉法值10.90000.81000.72900.65610.59050.53140.47830.43050.38740.3487龍格庫(kù)塔法值10.90500.81900.74120.67080.60710.54940.49720.45000.40720.36852-5用四階龍格-庫(kù)塔法求解題2-3數(shù)值解，并與前兩題結(jié)果相比較。解：程序代碼如下：clear;clc;h=0.1;y=1;disp('四階龍格庫(kù)塔法求得值為：'

10、;);for i=0:h:1 disp(y); k1=-y; k2=-(y+k1*h/2); k3=-(y+k2*h/2); k4=-(y+k3*h); y=y+(k1+2*k2+2*k3+k4)*h/6;end程序結(jié)果如下：歐拉法值10.90000.81000.72900.65610.59050.53140.47830.43050.38740.3487二階龍格庫(kù)塔法值10.90500.81900.74120.67080.60710.54940.49720.45000.40720.3685四階龍格庫(kù)塔法值10.90480.81870.74080.67030.60650.54880.49660.

11、44930.40660.36792-6已知二階系統(tǒng)狀態(tài)方程為寫出取計(jì)算步長(zhǎng)為h時(shí)，該系統(tǒng)狀態(tài)變量X=的四階龍格-庫(kù)塔法遞推關(guān)系式。2-7單位反饋系統(tǒng)的開環(huán)傳遞函數(shù)已知如下用matlab語(yǔ)句 、函數(shù)求取系統(tǒng)閉環(huán)零極點(diǎn)，并求取系統(tǒng)閉環(huán)狀態(tài)方程的可控標(biāo)準(zhǔn)型實(shí)現(xiàn)。解：程序代碼如下：clear;clc;num=5 100;den1=1 0;den2=1 4.6;den3=1 3.4 16.35;den4=0 0 0 5 100;den=conv(den1,den2);den=conv(den,den3);den=den+den4Z,P,K=tf2zp(num,den);P1=roots(den)Z1=r

12、oots(num)ZP輸出結(jié)果如下：den = 1.0000 8.0000 31.9900 80.2100 100.0000P1 = -0.9987 + 3.0091i -0.9987 - 3.0091i -3.0013 + 0.9697i -3.0013 - 0.9697iZ1 = -20Z = -20P = -0.9987 + 3.0091i -0.9987 - 3.0091i -3.0013 + 0.9697i -3.0013 - 0.9697i2-8用matlab語(yǔ)言編制單變量系統(tǒng)三階龍格-庫(kù)塔法求解程序，程序入口要求能接收狀態(tài)方程各系數(shù)陣（A,B,C,D）,和輸入階躍函數(shù)r(t)=R

13、*1(t);程序出口應(yīng)給出輸出量y（t）的動(dòng)態(tài)響應(yīng)數(shù)值解序列。解：函數(shù)程序代碼如下：function y = hs( A,B,C,D,R,T,h ) y=0; r=R; x=0;0;0;0; N=T/h; for t=1:1:N k1=A*x+B*R; k2=A*(x+h*k1/3)+B*R; k3=A*(x+2*h*k2/3)+B*R; x=x+h*(k1+3*k3)/4; y=C*x+D*R;end 將程序保存為名為hs的文件，用以求解2-1求解出的A,B,C,D進(jìn)行驗(yàn)證如下：>> hs(A,B,C,D,1,1,1)ans = 2.66672-10.用式（2-34）梯形法求解試

14、驗(yàn)方程，分析對(duì)計(jì)算步長(zhǎng)h有何限制，說(shuō)明h對(duì)數(shù)值穩(wěn)定性的影響。1月6日上機(jī)作業(yè)二：3-1 求解下列線性方程，并進(jìn)行解得驗(yàn)證：(1),(2) 由A*X=B得：X=AB解：（1）程序代碼如下：clear;clc;A=7 2 1 -2 9 15 3 -2 -2 -2 11 5 1 3 2 13;B=4 7 -1 0;x=inv(A)*B輸出結(jié)果如下：x = 0.4979 0.1445 0.0629 -0.0813（2）程序代碼如下：clear;clc;A= 5 7 6 5 1 7 10 8 7 2 6 8 10 9 3 5 7 9 10 4 1 2 3 4 5;B=24 96 34 136 36 14

15、4 35 140 15 60;x=inv(A)*B輸出結(jié)果如下：x = 1.0000 4.0000 1.0000 4.0000 1.0000 4.0000 1.0000 4.00001.0000 4.00003-2進(jìn)行下列計(jì)算，給出不使用for和while等循環(huán)語(yǔ)句的計(jì)算方法。（1） (2)求出y=x*sin(x) 在0<x<100條件下的每個(gè)峰值解：（1）程序代碼如下：clear;clc;a1=1;q=2;n=64;sum=a1*(1-qn)/(1-q)輸出結(jié)果如下：sum = 1.8447e+019（2）先畫出0<x<100條件下y的圖形代碼如下：x=0:0.1:1

16、00;y=x.*sin(x);plot(x,y)xlabel('x');ylabel('y');圖形如下：可以看出sin(x)在峰值時(shí)y也是在峰值位置，即每個(gè)點(diǎn)的峰值為3-3繪制下面的圖形。（1）sin(1/t),-1<t<1 (2) -1<t<1解：（1）代碼如下：t=-1:0.01:1;y=sin(1./t);plot(t,y);xlabel('t');ylabel('y');輸出圖形如下：（2）代碼如下：clear;clc;t=-1:0.01:1;y=1-(cos(7*t).3;plot(t,y);x

17、label('t');ylabel('y');輸出圖形如下：3-4已知元件的實(shí)驗(yàn)數(shù)據(jù)如下，擬合這一數(shù)據(jù)，并嘗試給出其特性方程。 X 0.0100 1.0100 2.0100 3.0100 4.0100 Y 2.5437 7.8884 9.6242 11.6071 11.9727 X 5.0100 6.0100 7.0100 8.0100 9.0100 y 13.2189 14.2679 14.6134 15.4045 15.0805 解：程序代碼如下：clear;clc;x= 0.0100 1.0100 2.0100 3.0100 4.0100 5.0100 6

18、.0100 7.0100 8.0100 9.0100;y= 2.5437 7.8884 9.6242 11.6071 11.9727 13.2189 14.2679 14.6134 15.4045 15.0805 ;P=polyfit(x,y,3);xi=0:0.01:9.01;yi=polyval(P,xi);plot(x,y,'o',xi,yi,'*');plot(x,y,xi,yi);xlabel('x');ylabel('y');輸出結(jié)果如下：3-5分別使用解微分方程方法、控制工具箱、simulink求解具有如下閉環(huán)傳遞函

19、數(shù)的系統(tǒng)的階躍響應(yīng)。解：Simulink中建立如下系統(tǒng)階躍響應(yīng)圖像如下：輸入以下代碼：num=6 26 6 20;den=1 3 4 2 2;step(num,den,10);輸出波形如下3-6已知系統(tǒng)的閉環(huán)傳遞函數(shù)，試分析該系統(tǒng)的穩(wěn)定性。解：畫出系統(tǒng)零極點(diǎn)圖代碼如下：num=6 26 6 20;den=1 8 36 40 10;pzmap(num,den)輸出結(jié)果為：可以看出系統(tǒng)的幾點(diǎn)都在左半平面中，所以該系統(tǒng)為穩(wěn)定系統(tǒng)。3-8某小功率隨動(dòng)系統(tǒng)動(dòng)態(tài)結(jié)構(gòu)如下圖所示，已知： 若系統(tǒng)輸入分別為，適用simulink分析系統(tǒng)的輸出分別如何？（1）輸入為1（t）:輸出為：（2）輸入為 t 時(shí)： 輸出：

20、（3） 輸入為l(t)-l(t-1.5)時(shí)：1月7日上機(jī)作業(yè)三：4-2設(shè)典型閉環(huán)結(jié)構(gòu)控制系統(tǒng)如圖4-47所示，當(dāng)階躍輸入幅值 時(shí)，用sp4_1.m求取輸出的響應(yīng)。解：代碼為：a=0.016 0.864 3.27 3.42 1;b=30 25;X0=0 0 0 0; %系統(tǒng)狀態(tài)向量初值為零V=2; %反饋系數(shù) n=4;T0=0;Tf=10;h=0.01;R=20 ; sp4_1;plot(t,y)輸出結(jié)果為：4-8求圖4-49非線性系統(tǒng)的輸出響應(yīng)y(t),并與無(wú)非線性環(huán)節(jié)情況進(jìn)行比較。解：程序代碼為：P=0.1 1 0.5 1;0 1 20 0;2 1 1 0;10 1 1 0;WIJ=1 0 1;1 4 -1;2 1 1;3 2 1;4 3 1 ;Z=0 0 0 0;S=0 0 0 0;h=0.01;L1=25;n=4;T0=0Tf=20;nout=4;Y0=10;sp4_4;plot(t,y,'r') hold on