控制系统仿真与CAD-实验报告

1、控制系统仿真与 CAD-实验报告控制系统仿真与 CAD实验课程报告一、实验教学目标与基本要求 上机实验是本课程重要的实践教学环节。实验的目的不仅仅是验证理论知 识，更重要的是通过上机加强学生的实验手段与实践技能，掌握应用 MATLAB/Simulink 求解控制问题的方法，培养学生分析问题、解决问题、应用 知识的能力和创新精神，全面提高学生的综合素质。通过对 MATLAB/Simulink进行求解，基本掌握常见控制问题的求解方法与命 令调用，更深入地认识和了解 MATLA语B 言的强大的计算功能与其在控制领域的 应用优势。上机实验最终以书面报告的形式提交，作为期末成绩的考核内容。二、题目及解答

2、第一部分： MATLAB必 备基础知识、控制系统模型与转换、线性控制系统的计 算机辅助分析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=o de45(f,0,100,0;0;0,0.2,0.2,5.7);plot3(x(:,1),x(:,2),x(:,3),grid,figure,plot(x(:,1)2.y=(x)x(1)2-2*x(1)+x(2);ff=optimset;ff.LargeScale=off;ff.TolFun=1e-30;ff. TolX=1e-15;ff.TolCon=1e-20

3、;x0=1;1;1;xm=0;0;0;xM=;A=;B=;Aeq=;Be q=;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 and running active-set algorithm. To run trust-region-reflective, setLargeScale = on. To run active-set with

4、out this warning, use Algorithm = active-set. In fmincon at 456Local minimum possible. Constraints satisfied.fmincon stopped because the size of the current search direction is less than twice the selected value of the step size tolerance and constraints are satisfied to within the selected value of

5、 the constraint tolerance.Active inequalities (to within options.TolCon = 1e-20): lower upper ineqlin ineqnonlin2x =1.000001.0000-1.00004d =iterations: 5 funcCount: 20 lssteplength: 1 stepsize: 3.9638e-26algorithm: medium-scale: SQP, Quasi-Newton, line-search firstorderopt: 7.4506e-09 constrviolatio

6、n: 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 + 2s11 + 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.568z3 - 1.2 z2 + 1.19 z - 0.99Sample time: 0.1 secondsDiscrete-time

7、 transfer function.4. A=0 1 0;0 0 1;-15 -4 -13;B=0 0 2;C=1 0 0;D=0;G=ss(A,B,C,D),Gs=tf(G),Gz=zpk(G)G =a =x1x2 x3x101 0x200 1x3-15-4 -13b =u1x10x20x32c =x1x2 x3y110 0d =u1y1 0Continuous-time state-space model.Gs =s3 + 13 s2 + 4 s + 15Continuous-time transfer function. Gz =(s+12.78) (s2 + 0.2212s + 1.

8、174) Continuous-time zero/pole/gain model.5.设采样周期为 0.01s z=tf(z,0.01);H=(z+2)/(z2+z+0.16) H =z2 + z + 0.16Sample time: 0.01 seconds Discrete-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 + (Ki + Kp + 5)*s

9、+ 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 + 1270560.01 s6 + 2.185 s5 + 142.1 s4 + 2444 s3 + 4.389e04 s2 + 1.399e05 s + 127056 Continuous-time transfer fun

10、ction.Gd =x1x2x3x4x5x6x1-218.5-111.1 -29.83-16.74-6.671-3.029x212800000x30640000x40032000x5000800x6000020b =u1x1 4x2 0x3 0x4 0x5 0x6 0c =x1 x2 x3 x4 x5 x6y1 0 0 1.097 3.559 1.668 0.7573d =u1y1 0 Continuous-time state-space model. Gz =35933.152 (s+100) (s+2.358) (s+1.499)(s2 + 3.667s + 3.501) (s2 + 1

11、1.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+108444*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=ss(GG),Gz=zpk(GG) GG =-108444 z5 + 1.844e04 z4 + 1.789e04 z31.144e05 z5 + 2.876e04 z4 + 274.2 z3

12、 + 782.4 z2 + 47.52 z + 0.5Sample time: 0.1 secondsDiscrete-time transfer function.Gd =a =x1x2x3x4x5x1-0.2515-0.00959-0.1095 -0.05318-0.01791x20.250000x300.25000x4000.12500x50000.031250b =u1x11x2 0x3 0x4 0x5 0c =x1 x2 x3 x4 x5 y1 0.3996 0.6349 0.1038 0.05043 0.01698d =u1y1 -0.9482Sample time: 0.1 se

13、condsDiscrete-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 secondsDiscrete-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=feedbac

14、k(g3,g5);GG=3*feedback(G1* G2,g6)GG =3 s6 + 6 s5 + 3 s4 + 42 s3 + 84 s2 + 42 ss10 + 3 s9 + 55 s8 + 175 s7 + 300 s6 + 1323 s5 + 2656 s4 + 3715 s3 + 7732 s2 + 5602 s + 1400Continuous-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

15、,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.0001396 z4 + 9.596e-05 z3 + 8.18e-05 z2 - 0.0001289 z +4.355e-05z6 - 5.592 z5 + 13.26 z4 - 17.06 z3 + 12.58 z2 - 5.032 z + 0.8521Sample time: 0.01 seconds Discrete-time transfer function.

16、Gd2 =0.0003982 z5 - 0.0003919 z4 - 0.000336 z3 + 0.0007842 z2 - 0.000766 z + 0.0003214z6 - 2.644 z5 + 4.044 z4 - 3.94 z3 + 2.549 z2 - 1.056 z + 0.2019Sample time: 0.1 seconds Discrete-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-07z6 -

17、0.0419 z5 - 0.07092 z4 - 0.0004549 z3 + 0.002495 z2 - 3.347e-05 z +1.125e-07Sample time: 1 seconds Discrete-time transfer function.10.(a) G=tf(1,1 2 1 2);eig(G),pzmap(G) ans =-2.0000-0.0000 + 1.0000i-0.0000 - 1.0000i系统为临界稳定。(b) G=tf(1,6 3 2 1 1);eig(G),pzmap(G) ans =-0.4949 + 0.4356i-0.4949 - 0.4356

18、i0.2449 + 0.5688i0.2449 - 0.5688i有一对共轭复根在右半平面，所以系统不稳定(c) G=tf(1,1 1 -3 -1 2);eig(G),pzmap(G) ans =-2.0000-1.00001.00001.0000有两根在右半平面，故系统不稳定11.(1) H=tf(-3 2,1 -0.2 -0.25 0.05);pzmap(H),abs(eig(H) ans =0.50000.50000.2000系统稳定。(2) H=tf(3 -0.39 -0.09,1 -1.7 1.04 0.268 0.024);pzmap(H),abs(eig(H) ans =1.19

19、391.19390.12980.1298系统不稳定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.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)a =x1x2x3x4x5x1-0.20.5000x20-0.51.600x300-14.385.80x4000-33.3100x50000-10u1x1 0 ans =-0.2000-0.5000-14.3000-33.3000-10.0000x2 0x3 0x4 0x5 30c

20、 =x1 x2 x3 x4 x5y1 0 0 0 0 0d =u1y1 0Continuous-time state-space model. 系统稳定。13. A=-5 2 0 0; 0 -40 0; -3 2-4 -1; -3 20 -4;A=sym(A);symst;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)/

21、2 - exp(-4*t)/(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); end y=impulse(G,tt); plot(tt,xx,tt,y,:)解析解和数值解的脉冲响应曲

22、线如图所示，可以看出他们完全一致。14.(a) s=tf(s);G=(s+6)*(s-6)/(s*(s+3)*(s+4-4j)*(s+4+4j);rlocus(G),grid不存在 K 使得系统稳定(b) G=tf(1,2,2,1 1 14 8 0);rlocus(G),grid放大根轨迹图像， 可以看到，根轨迹与虚轴交点处， K 值为 5.53，因此，0K tau=2;n,d=paderm(tau,1,3);s=tf(s);G=tf(n,d)*(s-1)/(s+1)5,rlocus(G)-1.5 s2 + 4.5 s - 3s8 + 8 s7 + 29.5 s6 + 65.5 s5 + 95

23、 s4 + 91 s3 + 55.5 s2 + 19.5 s + 3Continuous-time transfer function.由图得 0Ks=tf(s);G=8*(s+1)/(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.46864.2340 wcg =1.5811 wcp =0.2336系统稳定。(b)z=tf(z);G=0.45*(z+1.31)*(z+0.054)*(z-0.957)

24、/(z*(z-1)*(z-0.368)*(z-0.99);b ode(G),figure,nyquist(G),figure,nichols(G),Gm,y,wcg,wcp=margin(G),figure,ste p(feedback(G,1)Warning: The closed-loop system is unstable. In warning at 26In DynamicSystem.margin at 63Gm =0.9578y =-1.7660 wcg =1.0464 wcp =Ki-oei-m-ilC-io!ilg-Lnop 畫 n 一一*rri*itj系统不稳定17.s=

25、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+50);GG=G*Gc;nyquist(GG),grid,figure,bode(GG),figure,nichols(GG),grid, figure,step(feedback(GG ,1)由奈氏图可得，曲线不包围 (-1,j0) 点，而开环系统不含有不稳定极点，所以根据 奈氏稳定判据闭环系统是稳定的。可得系统是稳定的用阶跃响应来验证，第二部分： Simulink 在系统仿真中的应用、控制系统计算机辅助设计、控制工程中的 仿真技术应用2

26、. syms y t; y=dsolve(D4y+5*D3y+6*D2y+4*Dy+2*y=exp(-3*t)+exp(-5*t)*sin(4*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);end plot(tout,yout,tt,yy,:)3.输出曲线及误差曲线4.x1x2x3x4t sin( x2e2.3( x4) ) x1 sin(x2e2.3( x4) x35. A,B,C,D=linmod(part2_5)

27、;G=ss(A,B,C,D)Warning: Using a default value of 0.2 for maximum step size.The simulationstepsize will be equal to or less than this value. You can disable this diagnostic by setting Automatic solver parameter selection diagnostic to none in the Diagnostics page of the configuration parameters dialog

28、. In dlinmod at 172In linmod at 60a =x1x2x3x4x5x6x1000000x20-1000000x31300-100000x40200-0.88-10000x50000-1000x6000294.1-29.41x701000x8-27.5601.045e+004x9000x10000x7x8x101.4x200x300x411.760x501.4x600x700x80-6.667x900x1000b=u1x10x21x30x40x50x60x70x80x90x100c =x1x2x3x4 x5x6y1130000 00d=u1y10-149.3x7x8

29、x9 x10000 0Continuous-time model.-0.44000000100-10000x9x10000000000019.61000000000subplot(221),step(G),grid,subplot(222),bode(G),grid,subplot(223),nyquist(G),g rid,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

30、,1);step(feedback(G,1),figure,step(GG),xlim(85 95) Gm,garma,wcg,wcp=margin(G) Gm =4.8921 garma =60.0634 wcg =7.9490 wcp =3.9199 Gm,garma,wcg,wcp=margin(G*Gc) Warning: The closed-loop system is unstable. In warning at 26In DynamicSystem.margin at 63 Gm =0.8090 garma =-6.0615 wcg =17.1659 wcp =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)-0.09781.21181.87670.7871-3.8819-0.46682.67131.03205.44000.6152-2.3163