Matlab控制系統(tǒng)計(jì)算機(jī)輔助設(shè)計(jì)

1、實(shí) 驗(yàn) 目 錄實(shí)驗(yàn)一：Matlab環(huán)境熟悉與基本運(yùn)算（設(shè)計(jì)型）實(shí)驗(yàn)二：Matlab語言程序設(shè)計(jì)(設(shè)計(jì)型）實(shí)驗(yàn)三：控制系統(tǒng)模型的建立（設(shè)計(jì)型）實(shí)驗(yàn)四：Simulink仿真入門（驗(yàn)證型）實(shí)驗(yàn)五：控制系統(tǒng)時域仿真分析（設(shè)計(jì)型）實(shí)驗(yàn)六：Simulink環(huán)境下時域仿真實(shí)驗(yàn)七：控制系統(tǒng)根軌跡仿真分析實(shí)驗(yàn)八：控制系統(tǒng)頻域仿真分析（設(shè)計(jì)型）實(shí)驗(yàn)名稱：Matlab環(huán)境熟悉與基本運(yùn)算（設(shè)計(jì)型）1、矩陣運(yùn)算（1）矩陣的乘法A=1 2;3 4;B=5 5;7 8;y=A2*By = 105 115 229 251（2）矩陣除法A=1 2 3;4 5 6;7 8 9;B=1 0 0;0 2 0;0 0 3;y1=AB警

2、告: 矩陣接近奇異值，或者縮放錯誤。結(jié)果可能不準(zhǔn)確。RCOND = 1.541976e-18。 y1 = 1.0e+16 * -0.4504 1.8014 -1.3511 0.9007 -3.6029 2.7022 -0.4504 1.8014 -1.3511y2=A/By2 = 1.0000 1.0000 1.0000 4.0000 2.5000 2.0000 7.0000 4.0000 3.0000（3）矩陣的轉(zhuǎn)置及共軛轉(zhuǎn)置A=5+i,2-i,1;6*i,4,9-i;y1=A.'y1 = 5.0000 + 1.0000i 0.0000 + 6.0000i 2.0000 - 1.00

3、00i 4.0000 + 0.0000i 1.0000 + 0.0000i 9.0000 - 1.0000iy2=A'y2 = 5.0000 - 1.0000i 0.0000 - 6.0000i 2.0000 + 1.0000i 4.0000 + 0.0000i 1.0000 + 0.0000i 9.0000 + 1.0000i（4）使用冒號選出指定元素A=1 2 3;4 5 6;7 8 9;y1=A(1:2,3)y1 = 3 6y2=A(2:3,:)y2 = 4 5 6 7 8 9（5）復(fù)數(shù)矩陣的生成：syms a1 a2 a3 a4 b1 b2 b3 b4A=a1+b1*j a2+

4、b2*j;a3*exp(b3*j) a4*exp(b4*j) A = a1 + b1*1i, a2 + b2*1i a3*exp(b3*1i), a4*exp(b4*1i)a1=3;a2=-2;a3=9;a4=23;b1=5;b2=3;b3=6;b4=33;A=a1+b1*j a2+b2*j;a3*exp(b3*j) a4*exp(b4*j)A = 3.0000 + 5.0000i -2.0000 + 3.0000i 8.6415 - 2.5147i -0.3054 +22.9980i2、 多項(xiàng)式p = 1 0 2 4roots(p)%用roots函數(shù)求多項(xiàng)式的根ans = 0.5898 +

5、1.7445i 0.5898 - 1.7445i -1.1795 + 0.0000iA=1.2 3 5 0.9;5 1.7 5 6;3 9 0 1;1 2 3 4;P=poly(A);poly2sym(P) ans =x4 - (69*x3)/10 - (3863*x2)/50 - (8613*x)/100 + 12091/20polyval(P,20)ans = 7.2778e+04polyvalm(P,A)ans = 1.0e-11 * -0.4093 -0.4849 -0.3876 -0.4519 -0.5002 -0.8072 -0.6004 -0.5684 -0.4704 -0.61

6、96 -0.6480 -0.5230 -0.3297 -0.4455 -0.3595 -0.36383、基本繪圖命令（1）t=0:0.01:2*pi;y=cos(t);plot(t,y)運(yùn)行結(jié)果如右圖所示：（2）t=0:0.01:2*pi;y1=cos(t-0.25);y2=sin(t-0.5);plot(t,y1,t,y2)運(yùn)行結(jié)果如右圖所示： 4、基本繪圖控制t=0:0.01:4*pi;x1=10*sin(t);plot(t,x1,'-.b*')axis(-20 20 -20 20)grid ontitle('正弦曲線x1=10*sin(t)圖')xlabe

7、l('t=0:0.01:4*pi')ylabel('x1')text(pi/2,10,'最大值')運(yùn)行結(jié)果如右圖所示： 實(shí)驗(yàn)名稱：Matlab語言程序設(shè)計(jì)（設(shè)計(jì)型）1、編寫命令文件：計(jì)算 1+2+n<5000 時的最大 n 值；解： 建立M文件，文件名稱記為two1，編輯程序內(nèi)容如下所示：s=0;n=1;while(s<5000), s=s+n; n=n+1;end,n 保存，并在命令窗口中輸入two1，運(yùn)行結(jié)果如下所示：two1n = 1012、編寫函數(shù)文件：解：用for循環(huán)結(jié)構(gòu)編寫程序： 建立函數(shù)文件，文件名稱記為findsum1

8、，編輯程序內(nèi)容如下所示：function s1 =findsum1( )%findsum1 用for循環(huán)編寫的程序% 求2的0到30次冪的和,s表示和，i表示執(zhí)行次數(shù)s1=0;for i=0:30, s1=s1+2i;Endend保存，并在命令窗口中輸入s1=findsum1( )，運(yùn)行結(jié)果如下所示：s1=findsum1( )s1 = 2.1475e+09用while循環(huán)結(jié)構(gòu)編寫程序： 建立函數(shù)文件，文件名記為findsum2，編輯程序內(nèi)容如下所示：function s2 =findsum2( )%findsum2 用while循環(huán)編寫的程序% 求2的0到30次冪的和,s表示和，i表示執(zhí)行次

9、數(shù)s2=0;i=0;while(i<=30), s2=s2+2i; i=i+1;end end保存，并在命令窗口中輸入 s2 =findsum2( )，運(yùn)行結(jié)果如下所示： s2 =findsum2( )s2 = 2.1475e+093、解：建立M文件，文件名記為two3，編輯程序內(nèi)容如下：x=input('please input a character:n','s');if x='y'|x='Y', disp('x=1');else if x='n'|x='N', disp(

10、'x=0'); else disp('its erro'); endend保存，并在命令窗口中輸入two2，運(yùn)行結(jié)果如下所示：two2please input a character:nx=0two2please input a character:Sits errotwo2please input a character:Yx=11、在MATLAB環(huán)境中輸入下面的系統(tǒng)模型解：實(shí)驗(yàn)名稱：控制系統(tǒng)模型的建立和分析（設(shè)計(jì)型） 建立M文件，文件名稱記為thr1，編輯程序內(nèi)容如下所示：num=1,4,3,2; %系統(tǒng)模型的分子系數(shù)den=conv(1,0,0,conv(

11、1,1,1,8,20);%分母系數(shù)G1=tf(num,den);%得到開環(huán)系統(tǒng)模型G2=feedback(G1,1);%得到閉環(huán)系統(tǒng)模型p1,z1=pzmap(G1)%不加分號，在命令窗口可顯示出開環(huán)系統(tǒng)的零極點(diǎn)p2,z2=pzmap(G2)% 不加分號，在命令窗口可顯示出閉環(huán)系統(tǒng)的零極點(diǎn) 在命令窗口輸入thr1，則得出以下結(jié)果：thr1p1 = 0 0 -4.0000 + 2.0000i -4.0000 - 2.0000i -1.0000 z1 = -3.2695 -0.3652 + 0.6916i -0.3652 - 0.6916ip2 = -3.9174 + 2.1007i -3.917

12、4 - 2.1007i -1.1449 -0.0102 + 0.2972i -0.0102 - 0.2972iz2 = -3.2695 -0.3652 + 0.6916i -0.3652 - 0.6916i2、假設(shè)系統(tǒng)由兩個模塊和串聯(lián)連接而成，已知 且 解：建立M文件，文件名稱為thr2，編輯程序內(nèi)容如下所示：G3=tf(1,1,1,3,4);%實(shí)驗(yàn)三的第二題，題中的G1（s）模型G4=tf(1,3,5,1,4,3,2,1);% 題中的G2（s）模型G5=series(G3,G4)%用串聯(lián)函數(shù)，將G3和G4模型串聯(lián)，得到G5模型，在命令窗口顯示出模型Gss=ss(G5)%用ss得到G5模型的狀

13、態(tài)方程模型在命令窗口中輸入thr2，得到以下運(yùn)行結(jié)果：Transfer function: s3 + 4 s2 + 8 s + 5-s6 + 7 s5 + 19 s4 + 27 s3 + 19 s2 + 11 s + 4 a = x1 x2 x3 x4 x5 x6 x1 -7 -2.375 -1.688 -0.5938 -0.3438 -0.25 x2 8 0 0 0 0 0 x3 0 2 0 0 0 0 x4 0 0 2 0 0 0 x5 0 0 0 1 0 0 x6 0 0 0 0 0.5 0 b = u1 x1 1 x2 0 x3 0 x4 0 x5 0 x6 0 c = x1 x2 x

14、3 x4 x5 x6 y1 0 0 0.0625 0.125 0.25 0.3125 d = u1 y1 0 Continuous-time model.3、 假設(shè)系統(tǒng)的對象模型為解： 建立M文件，文件名稱為thr3，編輯程序內(nèi)容如下所示：s=tf('s');%拉氏算子G6=10/(s+1)3);%實(shí)驗(yàn)三第3題的對象模型Gpid=0.48*(1+1/(1.814*s)+0.4353*s/0.04353);%PID控制器G7=series(G6,Gpid);%控制器與對象模型進(jìn)行串聯(lián)得到的開環(huán)系統(tǒng)模型G8=feedback(G7,1)%閉環(huán)系統(tǒng)的傳遞函數(shù)模型z3,p3,k3=zp

15、kdata(G7,'v')%得到開環(huán)系統(tǒng)的零極點(diǎn)及增益z4,p4,k4=zpkdata(G8,'v') %得到閉環(huán)系統(tǒng)的零極點(diǎn)及增益在命令窗口中輸入thr3，得到以下運(yùn)行結(jié)果Transfer function: 3.79 s2 + 0.379 s + 0.2089-0.07896 s4 + 0.2369 s3 + 4.027 s2 + 0.458 s + 0.2089z3 = -0.0500 + 0.2294i -0.0500 - 0.2294ip3 = 0 -1.0000 -1.0000 + 0.0000i -1.0000 - 0.0000ik3 = 48z4

16、 = -0.0500 + 0.2294i -0.0500 - 0.2294ip4 = -1.4442 + 6.9670i -1.4442 - 6.9670i -0.0558 + 0.2217i -0.0558 - 0.2217ik4 = 484、解：三種典型模型建模函數(shù)： 傳遞函數(shù)模型：G=tf(num,den)； 零極點(diǎn)增益模型：G=zpk(z,p,k)； 狀態(tài)空間模型：G=ss(A,B,C,D)；各種模型之間的轉(zhuǎn)換語句：A,B,C,D=tf2ss(num,den);num,den=ss2tf(A,B,C,D,iu);z,p,k=tf2zp(num,den);num,den=zp2tf(z,

17、p,k);A,B,C,D=zp2ss(z,p,k);z,p,k=ss2zp(A,B,C,D,i)實(shí)驗(yàn)名稱：Simulink仿真環(huán)境入門（驗(yàn)證型）1、搭建下圖所示控制系統(tǒng)模型，自己設(shè)計(jì)輸入正弦曲線參數(shù)，觀察示波器及輸出數(shù)據(jù)。解：畫出的仿真圖如下，Sine Wave:O(t)=3*sin(2*t)，增益為5。實(shí)驗(yàn)名稱：Simulink仿真環(huán)境入門（驗(yàn)證型）出現(xiàn)以下波形圖：左圖為使用plot函數(shù)畫出的figure圖，右圖為示波器顯示的圖。實(shí)驗(yàn)名稱：控制系統(tǒng)時域仿真（設(shè)計(jì)型）1、時域分析（1）解：建立M文件，寫入以下程序：G1=tf(5,25,30,1,6,10,8);t=0:0.0005:25;y1

18、,t=step(G1,t);ymax,tp=max(y1);max_overshoot=ymax-1y2,t=impulse(G1,t);plot(t,y1,t,y2)figure圖如右圖所示：（2）解：建立M文件，寫入以下程序：wn=2:2:12;z=0.7;t=0:0.0005:25;hold onfor i=1:length(wn)G2=tf(wn(i)2,1,2*z*wn(i),wn(i)2);step(G2,t)endgrid on,hold offlegend('wn=2','wn=4','wn=6','wn=8',&

19、#39;wn=10','wn=12')figure圖如右圖所示：（3）解：建立M文件，寫入以下程序：wn=6;zetas=0.2 0.4 0.6 0.8 1.0 1.5 2.0;t=0:0.0005:25;hold onfor i=1:length(zetas)G3=tf(wn2,1,2*zetas(i)*wn,wn2);step(G3,t)endgrid on,hold offlegend('zetas=0.2','zetas=0.4','zetas=0.6','zetas=0.8','zetas

20、=1.0','zetas=1.5','zetas=2.0')figure圖如右上圖所示：2、系統(tǒng)穩(wěn)定性分析（1）解：零極點(diǎn)分布圖法程序：G4=zpk(-2,0;-1;-20,100)G5=feedback(G4,1)pzmap(G5)figure圖：由圖可知該系統(tǒng)的極點(diǎn)均位于左半平面，所以該系統(tǒng)穩(wěn)定。程序判定法程序：G4=zpk(-2,0;-1;-20,100)G5=feedback(G4,1)p=pole(G5);for i=1:length(p) if real(p(i)>0 disp('系統(tǒng)不穩(wěn)定') break endend

21、disp('系統(tǒng)穩(wěn)定')pzmap(G5)figure圖：實(shí)驗(yàn)名稱：Simulink環(huán)境下時域仿真1、 解：仿真圖： 階躍響應(yīng)： 方波響應(yīng)：2、解：仿真圖： 階躍響應(yīng)： 正弦信號輸入響應(yīng)： 斜坡輸入響應(yīng)：3、 解：仿真圖： 建立M文件，編寫以下程序： G1=tf(1 1,2 1); G2=tf(5,2 3 1); H1=tf(1,2 1); G3=series(G1,G2); G4=feedback(G3,H1) isstable(G4) 得到以下實(shí)驗(yàn)結(jié)果： g3 Transfer function: 10 s2 + 15 s + 5 - 8 s4 + 20 s3 + 18 s

22、2 + 12 s + 6 ans = 1 結(jié)論：編寫程序判斷該系統(tǒng)的穩(wěn)定性結(jié)果是該系統(tǒng)是穩(wěn)定的。 在Simulink環(huán)境下驗(yàn)證編程結(jié)果： 給出階躍信號，畫出以下仿真圖進(jìn)行仿真： 示波器顯示結(jié)果： 從上圖可以看出該系統(tǒng)是是穩(wěn)定的，所以編程結(jié)果是正確的。實(shí)驗(yàn)名稱：控制系統(tǒng)根軌跡分析1、（1）解：編寫程序如下所示：G=zpk(-0.5,0;-1;-2;-5,1)rlocus(G)k,p=rlocfind(G)Figure圖：系統(tǒng)的根軌跡圖：（2）解：編寫程序如下所示：G=zpk(,0;-2,1)rlocus(G)k,p=rlocfind(G)Figure圖：系統(tǒng)的根軌跡圖：2、根軌跡分析（1）解：編

23、寫程序如下所示：G=zpk(-3,0;-5;-6;-1+j;-1-j,1)rlocus(G)k,p=rlocfind(G)Figure圖：系統(tǒng)的根軌跡圖：由根軌跡圖可以看出點(diǎn)擊的這一點(diǎn)的增益k=35.2；閉環(huán)極點(diǎn)的位置p=-0.00587+1.35i，因?yàn)閜<0，所以該點(diǎn)系統(tǒng)閉環(huán)是穩(wěn)定的。（2）解：編寫以下程序：G=zpk(,0;-1;-2,1)rlocus(G)k=0:0.01:98;r,k=rlocus(G,k);m,n=size(r);for i=1:n if real(r(2,i)>0 break endendkg=k(i-1)程序的執(zhí)行結(jié)果：Zero/pole/gain: 1-s (s+1) (s+2) kg = 6由上述的執(zhí)行結(jié)果可以看出該系統(tǒng)穩(wěn)定的k值范圍是當(dāng)k<6時，系統(tǒng)是穩(wěn)定的。Figure圖如右圖所示：實(shí)驗(yàn)名稱：控制系統(tǒng)頻域仿真分析（設(shè)計(jì)型）1、（1）解：編寫以下程序：zeta=0.7;wn=2,4,6,8,10,12;for i=1:length(wn) num=wn(i)2; den=1,2*zeta*wn(i),wn(i)2; G1=tf(num,den) bode(G1)endlegend('wn=2','wn=4','wn=6