控制系統(tǒng)仿真與CAD-實驗報告(共33頁)

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

2、必備基礎知識、控制系統(tǒng)模型與轉換、線性控制系統(tǒng)的計算機輔助分析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.設采樣周期為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)設采樣周期為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有一對共軛復根在右半平面，所以系統(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,':')解析解和數值解的脈沖響

24、應曲線如圖所示，可以看出他們完全一致。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放大根軌跡圖像，可以看到，根軌跡與虛軸交點處，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)點，而開環(huán)系統(tǒng)不含有不穩(wěn)定極點，所以根據奈氏穩(wěn)定判據閉環(huán)系統(tǒng)是穩(wěn)定的。用階躍響應來驗證，可得系統(tǒng)是穩(wěn)定的。第二部分：Simulink 在系統(tǒng)仿真中的應用、控制系統(tǒ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階躍響應和頻率響應曲線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