控制系统仿真与CAD_实验报告

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

2、型与转换、线性控制系统的计 算机辅助分析1.考虑茗名的Mssol化学反应方程组= _y_ 金& = * I ny E = 6 1 （X C）3选定口 =6=02 c = 5.且工1（3=工式0）=工式0） =0,绘制仿真结果的三维相轨迹，并得 出其在K-y平面上的投影.学习参考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=od e45(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),g

3、rid学习参考2.4解二之,最优化:口:用1H1Dy=(x)x(1)A2-2*x(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,ABAeq,Beq,xm,xM,wzhfc1,ff)Warning: Options LargeScale = off and Algorithm =.学习参trust-region-reflective conflict.Igno

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

5、ss than twice the selected value of the step size tolerance and constraints are satisfied to within the selected value of the constraint tolerance.Active inequalities (to within options.TolCon = 1e-20): lower upper ineqlin ineqnonlin2 x =1.000001.0000 f =-1.0000 iterations: 5funcCount: 20Issteplengt

6、h: 1stepsize: 3.9638e-26algorithm: medium-scale: SQP , Quasi-Newton, line-searchfirstorderopt: 7.4506e-09constrviolation: 0message: 1x766 char3.请珞r面的传递由数模型输入乳MATLAB环境.3）的133十2言答：十53H=-十0四，丁1秒(a) s=tf(s);G=(sA3+4*s+2)/(sA3*(sA2+2)*(sA2+1)A3+2*s+5)sA3 + 4 s + 2sA11 + 5 sA9 + 9 sA7 + 2 sA6 + 12 sA5 + 4

7、 sA4 + 12 sA3Continuous-time transfer function.(b) z=tf(z,0.1);H=(zA2+0.568)/(z-1)*(zA2-0.2*z+0.99)H =zA2 + 0.568zA3 - 1.2 zA2 + 1.19 z - 0.99Sample time: 0.1 secondsDiscrete-time transfer function.4.假设描述系统的常微分方程为y十1时十的阴十的闵=2班小 请选择一纲状态变量， 并将此方程在MATLAE工作空间中表示出来,如果想得到系统的传递函数和零极点模型， 我门将如何求取“得出的结果又是怎样的？

8、由微分方程模型能否直接写出系统的传递函数模 型？A=010;001;-15-4-13;B=002;C=100;D=0;G=ss(A,B,C,D),Gs=tf(G),Gz=zpk(G)G =a =x1 x2 x3x1010x2001x3 -15-4 -13b =u1x1 0x2 0x3 2c =x1 x2 x3y1 100d =u1y1 0Continuous-time state-space model.Gs =2sA3 + 13 sA2 + 4 s + 15Continuous-time transfer function.Gz =2学习参考(s+12.78) (sA2 + 0.2212s

9、+ 1.174)Continuous-time zero/pole/gain model.5.已知某系统的差分方程模型为十？)十式/十1)十01办) = 十1)十2试即).试将其输入 到MATLAB工乍空间.设采样周期为0.01s z=tf(z,0.01);H=(z+2)/(zA2+z+0.16)H =z + 2zA2 + z + 0.16Sample time: 0.01 secondsDiscrete-time transfer function.6.假设某单位负反馈系统中，G=/十十区十5), &=+试用 MATLAB推导出闭环系统的传递函数模型.symsJKpKis;G=(s+1)/(

10、J*sA2+2*s+5);Gc=(Kp*s+Ki)/s;GG=feedback(G*Gc,1)GG =(Ki + Kp*s)*(s + 1)/(J*sA3 + (Kp + 2)*sA2 + (Ki + Kp + 5)*s + Ki) 7.从：面给出的典型反馈控制系统结构子模型中一求出总系统的状态方程与传递函数模型；并 得出各个模型的零极点模型表示.(Ai仃.1694 十III211.875 + 317.64+ 211)( + 94.34)(s + II.1684)35786.7Z-1 十 108444z-十 4j(二T 十 20)(-1 十 74.04)(a)s=tf(s);G=(211.87

11、*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 sA3 + 3.732e04 sA2 + 1.399e05 s + 1270560.01 sA6 + 2.185 sA5 + 142.1 sA4 + 2444 sA3 + 4.389e04 sA2 + 1.399e05s+127056Continuous-time transfer function.Gd =x1x2x3x4x5x6x1

12、-218.5-111.1-29.83-16.74-6.671-3.029x212800000x30640000x40032000x5000800x6000020b =u1x1 4x2 0x3 0x4 0x5 0x6 0yix1x2x3x4x5x6001.0973.5591.6680.7573u1y1 0Continuous-time state-space model.Gz =35933.152 (s+100) (s+2.358) (s+1.499)(sA2 + 3.667s + 3.501) (sA2 + 11.73s + 339.1) (sA2 + 203.1s + 1.07e04)Con

13、tinuous-time zero/pole/gain model.(b)设采样周期为0.1sz=tf(z,0.1);G=(35786.7*zA2+108444*zA3)/(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 zA5 + 1.844e04 zA4 + 1.789e04 *31.144e05 zA5 + 2.876e04 zA4 + 274.2 zA3 + 782.4 zA2 + 47.52 z + 0.5Sample time

14、: 0.1 secondsDiscrete-time transfer function.Gd =a =x1x2x3x4x5x1-0.2515-0.00959-0.1095-0.05318-0.01791x20.250000x300.25000x4000.12500x50000.031250b =u1x1 1x2 0x3 0x4 0x5 0x1 x2 x3 x4 x5y1 0.39960.63490.1038 0.05043 0.01698d =u1y1 -0.9482Sample time: 0.1 secondsDiscrete-time state-space model.Gz =-0.

15、94821 zA3 (z-0.5) (z+0.33)(z+0.3035) (z+0.04438) (z+0.01355) (zA2 - 0.11z + 0.02396)Sample time: 0.1 secondsDiscrete-time zero/pole/gain model.8.学习参考巳知系统的方枢图如图所示，试推导出从输入信号r(i)到输出信号y(t)的总系统模型.s=tf(s);g1=1/(s+1);g2=s/(sA2+2);g3=1/sA2;g4=(4*s+2)/(s+1)A2;g5=50;g6=(sA2+2)/(sA3+14);G1=feedback(g1*g2,g4);G

16、2=feedback(g3,g5);GG=3*feedback(G1*G2,g6)GG =3 sA6 + 6 sA5 + 3 sA4 + 42 sA3 + 84 sA2 + 42 ssA10 + 3 sA9 + 55 sA8 + 175 sA7 + 300 sA6 + 1323 sA5 + 2656 sA4 + 3715sA3 + 7732 sA2 + 5602 s + 1400Continuous-time transfer function.9.已知传递函数模型G(s)=(s + l)2(s2 + 2fi + 400)(s + 5)气岸-3s + 100) * (s2 十 3后十 2500

17、)对不同采样周期了=和T = 1秒对之进行离散化，比我原系统的阶跃喃应与各离散系统的阶跃嘀应曲线. 提示：后面将介绍，如果已知系统模型为G,则用smp(G)即可绘制出其阶跃响应曲援.s=tf(s);T0=0.01;T1=0.1;T2=1;G=(s+1)A2*(sA2+2*s+400)/(s+5)A2*(sA2+.学习参3*s+100)*(sA2+3*s+2500);Gd1=c2d(G,T0),Gd2=c2d(G,T1),Gd3=c2d(G,T2),st ep(G),figure,step(Gd1),figure,step(Gd2),figure,step(Gd3)Gd1 =4.716e-05

18、zA5 - 0.0001396 zA4 + 9.596e-05 zA3 + 8.18e-05 zA2 -0.0001289 z + 4.355e-05zA6 - 5.592 zA5 + 13.26 zA4 - 17.06 z3 + 12.58 zA2 - 5.032 z + 0.8521Sample time: 0.01 secondsDiscrete-time transfer function.Gd2 =0.0003982 zA5 - 0.0003919 zA4 - 0.000336 *3 + 0.0007842 zA2 -0.000766 z + 0.0003214zA6 - 2.644

19、 zA5 + 4.044 zA4 - 3.94 *3 + 2.549 zA2 - 1.056 z + 0.2019Sample time: 0.1 secondsDiscrete-time transfer function.Gd3 =8.625e-05 zA5 - 4.48e-05 zA4 + 6.545e-06 zA3 + 1.211e-05 zA2 - 3.299e-06 z + 1.011e-07zA6 - 0.0419 zA5 - 0.07092 zA4 - 0.0004549 *3 + 0.002495 zA2 -3.347e-05 z + 1.125e-07Sample time

20、: 1 secondsDiscrete-time transfer function.SIjbp Rea-DDn-xerune lS4I：Qrrd9110.判定下列差缴传递函数模型的稳定性.1 1 ,1“十以2一内十2 也 4十及：十2十。十J U）国4+匹一3M 营十2 (a) G=tf(1,1 2 1 2);eig(G),pzmap(G)ans =-2.0000-0.0000 + 1.0000i-0.0000 - 1.0000i系统为临界稳止o(b) G=tf(1,6 3 2 1 1);eig(G),pzmap(G) ans =-0.4949 + 0.4356i-0.4949 - 0.43

21、56i0.2449 + 0.5688i0.2449 - 0.5688i学习参考splkjubim) 0 4 r UAruj-Rea* Axe (MwneiB :1.学习参考有一对共腕复根在右半平面，所以系统不稳定(c) G=tf(1,1 1-3 -1 2);eig(G),pzmap(G) ans =-2.0000-1.00001.00001.0000PnhB-Zflroi Mur有两根在右半平面，故系统不稳定。11.3-3 - 0.39ff - 0.09(b)及=J 一 1,7 1 1.04;2 + 0-268;十 0,024判定下面采样系统的稳定柢-3二 十 2（力 =,二口,2二0.2二

22、| 0.田(1) H=tf(-3 2,1 -0.2 -0.25 0.05);pzmap(H),abs(eig(H) ans =0.50000.50000.2000-MMSUS. -sxy .ca_a_nE-Fleal Axn Iseconits1)系统稳定。(2) H=tf(3 -0.39 -0.09,1 -1.7 1.04 0.268 0.024);pzmap(H),abs(eig(H) ans =1.19391.19390.12980.1298系统不稳定12.治出连续系统的状态方程模型，请判定系统的稳定性-0.215000 0 10-0.51.600.T(i)=00-11.3S5.80x(

23、t)十0四）000-33.31000.0000-10-:打(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).学习参x1x2x3x4x5x1-0.20.5000x20-0.51.600x300-14.385.80x4000-33.3100x50000-10b =u1x1 0ans =-0.2000-0.5000-14.3000-33.3000-10.0000x2 0x3 0x4 0x5

24、30c =x1 x2 x3 x4 x5y1 00000d =u1y1 0Continuous-time state-space model.系统稳定。13.请求出下面自治系统状态方程的解析解-4 -1; -3 2 0 -4; A=sym(A);syms-5 出=_3 .-3 并和数值解得出的曲线比较一 A=-5 2 0 0; 0 -4 0 0; -3 2t;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*tA2*(exp(-4*t

25、)/(4*t)+exp(-4*t)/(2*tA2)+ 8*tA2*(exp(-4*t)/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(-520 0;0-400;-32-4-1;-320 -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(

26、tt,xx,tt,y,:)解析解和数值解的脉冲响应曲线如图所示，可以看出他们完全一致14.试绘制不列开环系统的根轨迹曲线，并大致确定使单位负反馈系统稳定的K值范围K(s 十 6)s(s+3)(A-F4-Lj)(s + 4-4j)学习参考(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学习参考Real Axa faeccnds 1),-3E一放大根轨迹图像，可以看到，根轨迹与虚轴交点处，K伯:为5.53,因此，

27、0K tau=2;n,d=paderm(tau,1,3);s=tf(s);G=tf(n,d)*(s-1)/(s+1)A5,rlocus(G)-1.5 sA2 + 4.5 s - 3sA8 + 8 sA7 + 29.5 sA6 + 65.5 sA5 + 95 sA4 + 91 sA3 + 55.5 sA2 + 19.5s + 3Continuous-time transfer function.Rx刈 LocusHedJ Axs7歪 Reo: Ldd筝Real Axis (leoondi-1)由图得0Ks=tf(s);G=8*(s+1)/(sA2*(s+15)*(sA2+6*s+10);bode

28、(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.2336Bode Oogr-am5*0sflfwMTO翁 I 1 工 21 - - MlSE省 nEnw 5 o 5 )-2-2-3用J B雪左Frquflz=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,nyqui

29、st(G),figure,nichols(G),Gm,y,wcg,wcp=margin(G),figure,step(feedback(G,1)Warning: The closed-loop system is unstable. In warning at 26In DynamicSystem.margin at 63Gm =0.9578 y =-1.7660wcg =1.0464wcp =1.0734FTtqMsney 什Mrt)270215135Opti-Laoe Fhaie (deojNgu&i DdigramReel Axd10001 如0200025W 3期3 如040M延帆I5

30、000Un* IBflGDnil系统不稳定17.假设受控对象模型为G=并假设由某种方法设计出串联控制器模 s( I 十 s/(L5)( I 十 s/50)型为G式身)= 绊坐共空加，试用频域响应的方法判定闭环系统的性能，并用时娥晌 (s 十 U.5)( 十 50)应检验得出的结论. 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+50);GG=G*Gc;nyquist(GG),grid,figure,bode(GG)figure,nichols(GG),grid,figure,step

31、(feedback(GG,1)由奈氏图可得，曲线不包围（-1,j0）点，而开环系统不含有不稳定极点，所以根据奈氏稳定判据闭环系统是稳定的。Sbep ResparrieL80.40.2Tm (sflcanda ji用阶跃响应来验证，可得系统是稳定的第二部分：Simulink在系统仿真中的应用、控制系统计算机辅助设计、控制工程中的仿真技术应用2.2考虑简单的践性微分手程峭如+ 5y+的十每十2y - b3t十唱4st疝in+开;，3)且子程的初 值为y(0) = 1,讽=式0) = l/2,y(0) = 02 试用Simulink搭建起系统的仿真模型，井绘 制出仿真结果曲线.由第2章介绢的知识，该

32、方程可以用微分方程数值解的形式进行分析， 试比较二者的分析结果， 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);endplot(tout,yout,tt,yy,:)3.建立起如图印7所示非线性系统网的Sinmlink框图，并观察在单位除跋信号输入下系统的输 出曲线和误差曲线.图齐了习题

33、5的系统方框图输出曲线及误差曲线4.已知某系统的Simulink仿真框密如图5-15所示，试由谟框图写出系统的数学横型公式.* 1/5 1/Slore grataeImegratcrl【口旧察皿口1加记富耳噂四p5in(u(1 )*exp(2.3*(-i(2)J)4FenPit-AoductClock图5-15习题7的Sim Mink仿直框图X1 = tsin(x2e2.3(f4)X2 = XiX3 = sin(x2e2.3( 4)5.假设已知直流电机拖动模型方框图如图二17所示，试利用Simulink提供的工具提取该系统的 总模型，并利用该工具绘制系统的阶跃响鹿、顽域响魔曲线” 图5-17习

34、题10的台流出机拖劫系统方椎.图 A,B,C,D=linmod(part2_5);G=ss(ABC,D)Warning: Using a default 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 the Diagnostics page

35、 of the configuration parameters dialog. In dlinmod at 172In linmod at 60a =x1x2x3x4x5.学习参考x6x1000000x20-1000000x31300-100000x40200-0.88-10000x50000-1000x6000294.1-29.41-149.3x70100-0.44000x8-27.5600001.045e+004x9000100-100x1000000x7x8x9x10x101.400x20000x30000x411.76000x501.400x60019.610x70000x80-6.

36、66700x90000x100000b =u1x10x21x30x40x50x60x70x80x9x100c =x1 x2 x3 x4 x5 x6 x7 x8 x9 x10y1 130000000000d =u1y1 0Continuous-time model.subplot(221),step(G),grid,subplot(222),bode(G),grid,subplot(223),nyquist(G),grid,subplot(224),nichols(G),grid阶跃响应和频率响应曲线6.Hyguifil Difigrflm口1口。 2DQ 300 4GQReaLAxitBode

37、 uagrm一岳B 强一 d (ej.eMiad 口rirluudo假设系统的对象模型为G(s)=-一 并已知我们可以设!t-个控 制器为&=请观察在该控制器式系统的动态特性-比较原系统和校正后 系统的幅值和相匝雷专堂给出进一步改进系统性能的建议.s=tf(s);G=210*(s+1.5)/(s+1.75)*(s+16)*(sA2+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(8595)口 tU4niM2月口与 Gm,garma,wcg,wcp=ma

38、rgin(G)Gm =4.8921garma =60.0634 wcg =7.9490wcp =3.9199 Gm,garma,wcg,wcp=margin(G*Gc)Warning: The closed-loop system is unstable. In warning at 26In DynamicSystem.margin at 63Gm =0.8090garma =-6.0615wcg =17.1659wcp =学习参考18.90297.若系统的状态方程模型为()3.选择加权矩阵Q =由电(1,2,3,4)及冗=心，则设计出这一线性二次型指标的最优控制擀及在 最优控制下的闭环系统极点位置，并绘制出闭环系统各个状态的曲线.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(ABQ,R),eig(A-B*K)-0.09781.21181.87670.7871-3.8819-0.46682.67131.0320P =5.44000.6152-2.31630.04520.61521.8354-0.0138-0.758