控制系統(tǒng)仿真與CAD-實(shí)驗(yàn)報(bào)告(共33頁)

1、精選優(yōu)質(zhì)文檔-傾情為你奉上控制系統(tǒng)仿真與CAD實(shí)驗(yàn)課程報(bào)告 一、實(shí)驗(yàn)教學(xué)目標(biāo)與基本要求上機(jī)實(shí)驗(yàn)是本課程重要的實(shí)踐教學(xué)環(huán)節(jié)。實(shí)驗(yàn)的目的不僅僅是驗(yàn)證理論知識，更重要的是通過上機(jī)加強(qiáng)學(xué)生的實(shí)驗(yàn)手段與實(shí)踐技能，掌握應(yīng)用MATLAB/Simulink 求解控制問題的方法，培養(yǎng)學(xué)生分析問題、解決問題、應(yīng)用知識的能力和創(chuàng)新精神，全面提高學(xué)生的綜合素質(zhì)。通過對MATLAB/Simulink進(jìn)行求解，基本掌握常見控制問題的求解方法與命令調(diào)用，更深入地認(rèn)識和了解MATLAB語言的強(qiáng)大的計(jì)算功能與其在控制領(lǐng)域的應(yīng)用優(yōu)勢。上機(jī)實(shí)驗(yàn)最終以書面報(bào)告的形式提交，作為期末成績的考核內(nèi)容。二、題目及解答第一部分：MATLAB

2、必備基礎(chǔ)知識、控制系統(tǒng)模型與轉(zhuǎn)換、線性控制系統(tǒng)的計(jì)算機(jī)輔助分析1. >>f=inline('-x(2)-x(3);x(1)+a*x(2);b+(x(1)-c)*x(3)','t','x','flag','a','b','c');t,x=ode45(f,0,100,0;0;0,0.2,0.2,5.7);plot3(x(:,1),x(:,2),x(:,3),grid,figure,plot(x(:,1),x(:,2),grid2.>>y=(x)x(1)2-2*x(

3、1)+x(2);ff=optimset;ff.LargeScale='off'ff.TolFun=1e-30;ff.TolX=1e-15;ff.TolCon=1e-20;x0=1;1;1;xm=0;0;0;xM=;A=;B=;Aeq=;Beq=;x,f,c,d=fmincon(y,x0,A,B,Aeq,Beq,xm,xM,wzhfc1,ff)Warning: Options LargeScale = 'off' and Algorithm ='trust-region-reflective' conflict.Ignoring Algorithm

4、 and running active-set algorithm. To runtrust-region-reflective, setLargeScale = 'on'. To run active-set without this warning, useAlgorithm = 'active-set'. > In fmincon at 456 Local minimum possible. Constraints satisfied.fmincon stopped because the size of the current search dir

5、ection is less thantwice the selected value of the step size tolerance and constraints are satisfied to within the selected value of the constraint tolerance.<stopping criteria details>Active inequalities (to within options.TolCon = 1e-20): lower upper ineqlin ineqnonlin 2 x = 1.0000 0 1.0000f

6、 = -1.0000c = 4d = iterations: 5funcCount: 20lssteplength: 1stepsize: 3.9638e-26algorithm: 'medium-scale: SQP, Quasi-Newton, line-search'firstorderopt: 7.4506e-09constrviolation: 0message: 1x766 char3.(a) >> s=tf('s');G=(s3+4*s+2)/(s3*(s2+2)*(s2+1)3+2*s+5)G = s3 + 4 s + 2 - s11

7、 + 5 s9 + 9 s7 + 2 s6 + 12 s5 + 4 s4 + 12 s3 Continuous-time transfer function.(b) >> z=tf('z',0.1);H=(z2+0.568)/(z-1)*(z2-0.2*z+0.99)H = z2 + 0.568 - z3 - 1.2 z2 + 1.19 z - 0.99Sample time: 0.1 secondsDiscrete-time transfer function.4.>> A=0 1 0;0 0 1;-15 -4 -13;B=0 0 2'C=1

8、0 0;D=0;G=ss(A,B,C,D),Gs=tf(G),Gz=zpk(G)G = a = x1 x2 x3 x1 0 1 0 x2 0 0 1 x3 -15 -4 -13 b = u1 x1 0 x2 0 x3 2 c = x1 x2 x3 y1 1 0 0 d = u1 y1 0Continuous-time state-space model.Gs = 2 - s3 + 13 s2 + 4 s + 15 Continuous-time transfer function.Gz = 2 - (s+12.78) (s2 + 0.2212s + 1.174) Continuous-time

9、 zero/pole/gain model.5.設(shè)采樣周期為0.01s>> z=tf('z',0.01);H=(z+2)/(z2+z+0.16)H = z + 2 - z2 + z + 0.16 Sample time: 0.01 secondsDiscrete-time transfer function.6.>> syms J Kp Ki s;G=(s+1)/(J*s2+2*s+5);Gc=(Kp*s+Ki)/s;GG=feedback(G*Gc,1) GG = (Ki + Kp*s)*(s + 1)/(J*s3 + (Kp + 2)*s2 + (K

10、i + Kp + 5)*s + Ki)7.(a)>>s=tf('s');G=(211.87*s+317.64)/(s+20)*(s+94.34)*(s+0.1684);Gc=(169.6*s+400)/(s*(s+4);H=1/(0.01*s+1);GG=feedback(G*Gc,H),Gd=ss(GG),Gz=zpk(GG)GG = 359.3 s3 + 3.732e04 s2 + 1.399e05 s + - 0.01 s6 + 2.185 s5 + 142.1 s4 + 2444 s3 + 4.389e04 s2 + 1.399e05 s + Continu

11、ous-time transfer function.Gd = a = x1 x2 x3 x4 x5 x6 x1 -218.5 -111.1 -29.83 -16.74 -6.671 -3.029 x2 128 0 0 0 0 0 x3 0 64 0 0 0 0 x4 0 0 32 0 0 0 x5 0 0 0 8 0 0 x6 0 0 0 0 2 0 b = u1 x1 4 x2 0 x3 0 x4 0 x5 0 x6 0 c = x1 x2 x3 x4 x5 x6 y1 0 0 1.097 3.559 1.668 0.7573 d = u1 y1 0Continuous-time stat

12、e-space model.Gz = 35933.152 (s+100) (s+2.358) (s+1.499) - (s2 + 3.667s + 3.501) (s2 + 11.73s + 339.1) (s2 + 203.1s + 1.07e04)Continuous-time zero/pole/gain model.(b)設(shè)采樣周期為0.1s>>z=tf('z',0.1);G=(35786.7*z2+*z3)/(1+4*z)*(1+20*z)*(1+74.04*z);Gc=z/(1-z);H=z/(0.5-z);GG=feedback(G*Gc,H),Gd=

13、ss(GG),Gz=zpk(GG)GG = - z5 + 1.844e04 z4 + 1.789e04 z3 - 1.144e05 z5 + 2.876e04 z4 + 274.2 z3 + 782.4 z2 + 47.52 z + 0.5 Sample time: 0.1 secondsDiscrete-time transfer function.Gd = a = x1 x2 x3 x4 x5 x1 -0.2515 -0.00959 -0.1095 -0.05318 -0.01791 x2 0.25 0 0 0 0 x3 0 0.25 0 0 0 x4 0 0 0.125 0 0 x5 0

14、 0 0 0.03125 0 b = u1 x1 1 x2 0 x3 0 x4 0 x5 0 c = x1 x2 x3 x4 x5 y1 0.3996 0.6349 0.1038 0.05043 0.01698 d = u1 y1 -0.9482 Sample time: 0.1 secondsDiscrete-time state-space model.Gz = -0.94821 z3 (z-0.5) (z+0.33) - (z+0.3035) (z+0.04438) (z+0.01355) (z2 - 0.11z + 0.02396) Sample time: 0.1 secondsDi

15、screte-time zero/pole/gain model.8.>>s=tf('s');g1=1/(s+1);g2=s/(s2+2);g3=1/s2;g4=(4*s+2)/(s+1)2;g5=50;g6=(s2+2)/(s3+14);G1=feedback(g1*g2,g4);G2=feedback(g3,g5);GG=3*feedback(G1*G2,g6)GG = 3 s6 + 6 s5 + 3 s4 + 42 s3 + 84 s2 + 42 s - s10 + 3 s9 + 55 s8 + 175 s7 + 300 s6 + 1323 s5 + 2656

16、 s4 + 3715 s3 + 7732 s2 + 5602 s + 1400 Continuous-time transfer function.9.>>s=tf('s');T0=0.01;T1=0.1;T2=1;G=(s+1)2*(s2+2*s+400)/(s+5)2*(s2+3*s+100)*(s2+3*s+2500);Gd1=c2d(G,T0),Gd2=c2d(G,T1),Gd3=c2d(G,T2),step(G),figure,step(Gd1),figure,step(Gd2),figure,step(Gd3)Gd1 = 4.716e-05 z5 - 0

17、. z4 + 9.596e-05 z3 + 8.18e-05 z2 - 0. z + 4.355e-05 - z6 - 5.592 z5 + 13.26 z4 - 17.06 z3 + 12.58 z2 - 5.032 z + 0.8521 Sample time: 0.01 secondsDiscrete-time transfer function.Gd2 = 0. z5 - 0. z4 - 0. z3 + 0. z2 - 0. z + 0. - z6 - 2.644 z5 + 4.044 z4 - 3.94 z3 + 2.549 z2 - 1.056 z + 0.2019 Sample

18、time: 0.1 secondsDiscrete-time transfer function.Gd3 = 8.625e-05 z5 - 4.48e-05 z4 + 6.545e-06 z3 + 1.211e -05 z2 - 3.299e-06 z + 1.011e-07 - z6 - 0.0419 z5 - 0.07092 z4 - 0. z3 + 0. z2 - 3.347e-05 z + 1.125e-07 Sample time: 1 secondsDiscrete-time transfer function.10.(a)>> G=tf(1,1 2 1 2);eig(

19、G),pzmap(G)ans = -2.0000 -0.0000 + 1.0000i -0.0000 - 1.0000i系統(tǒng)為臨界穩(wěn)定。(b) >> G=tf(1,6 3 2 1 1);eig(G),pzmap(G)ans = -0.4949 + 0.4356i -0.4949 - 0.4356i 0.2449 + 0.5688i 0.2449 - 0.5688i有一對共軛復(fù)根在右半平面，所以系統(tǒng)不穩(wěn)定。(c) >> G=tf(1,1 1 -3 -1 2);eig(G),pzmap(G)ans = -2.0000 -1.0000 1.00001.0000有兩根在右半平面

20、，故系統(tǒng)不穩(wěn)定。11.(1) >> H=tf(-3 2,1 -0.2 -0.25 0.05);pzmap(H),abs(eig(H')ans = 0.5000 0.50000.2000系統(tǒng)穩(wěn)定。(2) >> H=tf(3 -0.39 -0.09,1 -1.7 1.04 0.268 0.024);pzmap(H),abs(eig(H')ans = 1.1939 1.1939 0.1298 0.1298系統(tǒng)不穩(wěn)定。12.(1)>> A=-0.2 0.5 0 0 0;0 -0.5 1.6 0 0;0 0 -14.3 85.8 0;0 0 0 -33

21、.3 100;0 0 0 0 -10;B=0 0 0 0 30'C=zeros(1,5);D=0;G=ss(A,B,C,D),eig(G)G = a = x1 x2 x3 x4 x5 x1 -0.2 0.5 0 0 0 x2 0 -0.5 1.6 0 0 x3 0 0 -14.3 85.8 0 x4 0 0 0 -33.3 100 x5 0 0 0 0 -10 b = u1 x1 0ans = -0.2000 -0.5000 -14.3000 -33.3000 -10.0000 x2 0 x3 0 x4 0 x5 30 c = x1 x2 x3 x4 x5 y1 0 0 0 0 0 d

22、 = u1 y1 0Continuous-time state-space model.系統(tǒng)穩(wěn)定。13.>> A=-5 2 0 0; 0 -4 0 0; -3 2 -4 -1; -3 2 0 -4; A=sym(A);syms t;x=expm(A*t)*1;2;0;1x =4*exp(-4*t) - 3*exp(-5*t)2*exp(-4*t)12*exp(-4*t) - 18*exp(-5*t) + 3*t*exp(-4*t) - 4*t2*(exp(-4*t)/(4*t) + exp(-4*t)/(2*t2) + 8*t2*(exp(-4*t)/2 - exp(-4*t)/(

23、2*t) - 16*t*(exp(-4*t) - exp(-4*t)/(2*t)6*exp(-4*t) - 9*exp(-5*t) - 8*t*(exp(-4*t) - exp(-4*t)/(2*t)>> G=ss(-5 2 0 0; 0 -4 0 0; -3 2 -4 -1; -3 2 0 -4,1;2;0;1,eye(4),zeros(4,1);tt=0:0.01:2; xx=;for i=1:length(tt) t=tt(i); xx=xx eval(x);endy=impulse(G,tt); plot(tt,xx,tt,y,':')解析解和數(shù)值解的脈沖響

24、應(yīng)曲線如圖所示，可以看出他們完全一致。14.(a) >> s=tf('s');G=(s+6)*(s-6)/(s*(s+3)*(s+4-4j)*(s+4+4j);rlocus(G),grid不存在K使得系統(tǒng)穩(wěn)定。(b) >> G=tf(1,2,2,1 1 14 8 0);rlocus(G),grid放大根軌跡圖像，可以看到，根軌跡與虛軸交點(diǎn)處，K值為5.53，因此，0<K<5.53時，系統(tǒng)穩(wěn)定。15.pade_app.mfunction Gr=pade_app(c,r,k)w=-c(r+2:r+k+1)'vv=c(r+1:-1:1)&#

25、39;zeros(k-1-r,1);W=rot90(hankel(c(r+k:-1:r+1),vv);V=rot90(hankel(c(r:-1:1);x=1 (Ww)'dred=x(k+1:-1:1)/x(k+1);y=c(1) x(2:r+1)*V'+c(2:r+1);nred=y(r+1:-1:1)/x(k+1);Gr=tf(nred,dred);paderm.mfunction n,d=paderm(tau,r,k)c(1)=1;for i=2:r+k+1,c(i)=-c(i-1)*tau/(i-1);endGr=pade_app(c,r,k);n=Gr.num1(k-

26、r+1:end);d=Gr.den1;>> tau=2;n,d=paderm(tau,1,3);s=tf('s');G=tf(n,d)*(s-1)/(s+1)5,rlocus(G)G = -1.5 s2 + 4.5 s - 3 - s8 + 8 s7 + 29.5 s6 + 65.5 s5 + 95 s4 + 91 s3 + 55.5 s2 + 19.5 s + 3Continuous-time transfer function.由圖得0<K<3.68能夠使得閉環(huán)系統(tǒng)穩(wěn)定。16.(a)>>s=tf('s');G=8*(s+1

27、)/(s2*(s+15)*(s2+6*s+10);bode(G),figure,nyquist(G),figure,nichols(G),Gm,y,wcg,wcp=margin(G),figure,step(feedback(G,1)Gm = 30.4686y = 4.2340wcg = 1.5811wcp =0.2336系統(tǒng)穩(wěn)定。(b)>>z=tf('z');G=0.45*(z+1.31)*(z+0.054)*(z-0.957)/(z*(z-1)*(z-0.368)*(z-0.99);bode(G),figure,nyquist(G),figure,nichols

28、(G),Gm,y,wcg,wcp=margin(G),figure,step(feedback(G,1)Warning: The closed-loop system is unstable. > In warning at 26 In DynamicSystem.margin at 63 Gm = 0.9578y = -1.7660wcg = 1.0464wcp =1.0734系統(tǒng)不穩(wěn)定。17.>>s=tf('s');G=100*(1+s/2.5)/(s*(1+s/0.5)*(1+s/50);Gc=1000*(s+1)*(s+2.5)/(s+0.5)*(s+

29、50);GG=G*Gc;nyquist(GG),grid,figure,bode(GG),figure,nichols(GG),grid,figure,step(feedback(GG,1)由奈氏圖可得，曲線不包圍(-1,j0)點(diǎn)，而開環(huán)系統(tǒng)不含有不穩(wěn)定極點(diǎn)，所以根據(jù)奈氏穩(wěn)定判據(jù)閉環(huán)系統(tǒng)是穩(wěn)定的。用階躍響應(yīng)來驗(yàn)證，可得系統(tǒng)是穩(wěn)定的。第二部分：Simulink 在系統(tǒng)仿真中的應(yīng)用、控制系統(tǒng)計(jì)算機(jī)輔助設(shè)計(jì)、控制工程中的仿真技術(shù)應(yīng)用2.>> syms y t;y=dsolve('D4y+5*D3y+6*D2y+4*Dy+2*y=exp(-3*t)+exp(-5*t)*sin(4*

30、t+pi/3)','y(0)=1','Dy(0)=1/2','D2y(0)=1/2','D3y(0)=1/5');tt=0:.05:10; yy=;for k=1:length(tt) ti=tt(k); yy=yy subs(y,'t',ti);endplot(tout,yout,tt,yy,':')3.輸出曲線及誤差曲線4.5.>> A,B,C,D=linmod('part2_5');G=ss(A,B,C,D)Warning: Using a default

31、value of 0.2 for maximum step size. The simulation stepsize will be equal to or less than this value. You can disable this diagnostic bysetting 'Automatic solver parameter selection' diagnostic to 'none' in theDiagnostics page of the configuration parameters dialog.> In dlinmod at

32、 172 In linmod at 60 a = x1 x2 x3 x4 x5 x6 x1 0 0 0 0 0 0 x2 0 -100 0 0 0 0 x3 130 0 -100 0 0 0 x4 0 200 -0.88 -100 0 0 x5 0 0 0 0 -100 0 x6 0 0 0 294.1 -29.41 -149.3 x7 0 100 -0.44 0 0 0 x8 -27.56 0 0 0 0 1.045e+004 x9 0 0 0 100 -10 0 x10 0 0 0 0 0 0 x7 x8 x9 x10 x1 0 1.4 0 0 x2 0 0 0 0 x3 0 0 0 0

33、x4 11.76 0 0 0 x5 0 1.4 0 0 x6 0 0 19.61 0 x7 0 0 0 0 x8 0 -6.667 0 0 x9 0 0 0 0 x10 0 0 0 0 b = u1 x1 0 x2 1 x3 0 x4 0 x5 0 x6 0 x7 0 x8 0 x9 0 x10 0 c = x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 y1 130 0 0 0 0 0 0 0 0 0 d = u1 y1 0 Continuous-time model.>>subplot(221),step(G),grid,subplot(222),bode(G),

34、grid,subplot(223),nyquist(G),grid,subplot(224),nichols(G),grid階躍響應(yīng)和頻率響應(yīng)曲線6.>>s=tf('s');G=210*(s+1.5)/(s+1.75)*(s+16)*(s2+3*s+11.25);Gc=52.5*(s+1.5)/(s+14.86);GG=feedback(G*Gc,1);step(feedback(G,1),figure,step(GG),xlim(85 95)>> Gm,garma,wcg,wcp=margin(G)Gm = 4.8921garma = 60.0634w

35、cg = 7.9490wcp = 3.9199>> Gm,garma,wcg,wcp=margin(G*Gc)Warning: The closed-loop system is unstable. > In warning at 26 In DynamicSystem.margin at 63 Gm = 0.8090garma = -6.0615wcg = 17.1659wcp = 18.90297.>>A=0 1 0 0;0 0 1 0;-3 1 2 3;2 1 0 0;B=1 0;2 1;3 2;4 3;Q=diag(1 2 3 4);R=eye(2);K,P=lqr(A,B,Q,R),eig(A-B*K)K = -0.0978 1.2118 1.8767 0.7871 -3.8819 -0.4668 2.6713 1.0320P = 5.4400 0.6152 -2.3163 0.0452 0.615