1、控制系统仿真及 CAD 试题（研 2012）1、4-2 某反馈控制系统的开环传递函数为 240KGsHss试绘制其根轨迹。 解：在 MATLAB 命令窗口中输入下列命令：ng=1.0;dg=poly(0,-4,-2+4j,-2-4j);rlocus(ng,dg)运行结果为：4-3 已知某系统传递函数为 2180()()1).3140sWs s试绘制其伯德图。解：分子分母同乘 100*200 得到 280(10)()2.5).32sWs在 Matlab 窗口中输入下列命令：k=80*200;num=1 100;den=conv(2.5 100, (1/200) 2*0.3 200);w=logs
2、pace(-1,1,100);m,p=bode(k*num,den,w);subplot(2,1,1);semilogx(w,20*log10(m);grid;xlabel(Frequency(rad/sec);ylabel(Gain(dB);subplot(2,1,2);semilogx(w,p);grid;xlabel(Frequency(rad/s);ylabel(Phase(deg);可绘制该系统的伯德图如下所示。4-4 设控制系统具有如下的开环传递函数 15KGsHs试求取当 K=10 时的相角裕度和幅值裕度，并画出其伯德图。解：在 MATLAB 命令窗口中输入下列命令：k=10;n
3、um=1;den=poly(0,-1,-5);m,p,w=bode(k*num,den);subplot(2,1,1);semilogx(w,20*log10(m);grid on;ylabel(Gain(dB);subplot(2,1,2);semilogx(w,p);grid on;xlabel(Frequency(rad/sec);ylabel(Phase(deg);gm,pm,wcg,wcp=margin(m,p,w)这里 gm,wcg 为幅值裕度值与相应的频率 pm,wcp 为相角裕度值与相应的频率，运行结果为：gm =3.0000，pm = 25.4489，wcg =2.2361，
4、wcp = 1.2241。因此，系统的幅值裕度和相角裕度分别为 3dB 和 25.4489。系统的伯德图如下所示。4-13 对于高阶系统的设计问题，往往要进行降阶近似处理，并要验证近似效果。已知某高阶系统模型为 234567819480264518768201380135()39417sssssGs经简化处理后，模型等效为 22.35()0sGs试比较两个模型在单位阶跃信号作用下的响应情况，并分析近似效果。解：在 Matlab 命令窗口中输入下列命令：num1=35,10861,13285,82402,278376,511812,422964,194480;den1=1,33,437,3017
5、,11870,27470,37492,28880,9600;sys1=tf(num1,den1);num2=35,284.98;den2=1,10.114,12.31;sys2=tf(num2,den2);step(sys1,-,sys2,:r,15)legend(原模型 , 简化后模型) ;grid on程序运行结果如下：从曲线中可以看出，降阶后系统响应无超调，调整时间缩短，系统响应变快，但是简化前后终值有差异。2-5 用四阶龙格-库塔法求解题 2-3 数值解，并与前两题结果相比较。解：四阶龙格-库塔法表达式 ，其截断误差为11234213243()6,)(),()kkkkhykftyhft
6、同阶无穷小，当 h 步距取得较小时，误差是很小的.5h(1) 程序如下：h=0.1;disp(四阶龙格- 库塔方法求解函数数值解为 );disp(y=);y=1;for t=0:h:1disp(y);k1=-y;k2=-(y+k1*h/2);k3=-(y+k2*h/2);k4=-(y+k3*h);y=y+(k1+2*k2+2*k3+k4)*h/6;end 得到结果 四阶龙格-库塔方法求解函数数值解为y= 1 0.9048 0.8187 0.7408 0.6703 0.6065 0.5488 0.4966 0.4493 0.4066 0.3679(2)比较这几种方法： 对于四阶龙格-库塔方法真值
7、1 0.9048 0.8187 0.7408 0.6703 0.6065 0.5488 0.4966 0.4493 0.4066 0.3679龙库1 0.9048 0.8187 0.7408 0.6703 0.6065 0.5488 0.4966 0.4493 0.4066 0.3679误差0 0 0 0 0 0 0 0 0 0 0显然四阶龙格-库塔法求解精度很高，基本接近真值。三种方法比较可以得到精度（四阶 ） 精度（二阶） 精度（欧拉）2-9 用题 2-8 仿真程序求解题 2-7 系统的闭环输出响应 y(t).解：题 2-7 单位反馈系统的开环传递函数已知如下2510()4.6)(3.6.
8、5)sGs已知开环传递函数，求得闭环传递函数为 2().)(.1.)0ss s在 matlab 命令行里键入 a=1 0; b=1 4.6; c=1 3.4 16.35; d=conv(a,b) ; e=conv(d,c)e = 1.0000 8.0000 31.9900 75.2100 0 f=0 0 0 5 100; g=e+fg = 1.0000 8.0000 31.9900 80.2100 100.0000%以上是计算闭环传递函数的特征多项式% m=5 100; z=roots(m)z = -20 %计算零点%那么 A= ,B= ,C= ,D=001018.23.980150m 文件为
9、：function y= hs(A,B,C,D,R,T,h) %T 为观测时间,h 为计算步长，R 为输入信号幅值%disp(数值解为);y=0;r=R;x=0;0;0;0;N=T/h;for t=1:N;k1=A*x+B*R;k2=A*(x+h*k1/3)+B*R;k3=A*(x+2*h*k2/3)+B*R;x=x+h*(k1+3*k3)/4;y(t)=C*x+D*R;end在命令行里键入 A=0 1 0 00 0 1 00 0 0 1-100 -80.21 -31.99 -8; B=0 0 0 1; C=-100 5 0 0; D=0; T=1; R=1; h=0.01; y= hs (A
10、,B,C,D,R,T,h)数值解为08.3333e-0075.8659e-0061.8115e-0053.9384e-0057.0346e-005。 %仅取一部分%3-5 图 3-61 中，若各环节的传递函数已知为：12345678910.71();();();0.80.2. .531.4();();();0.01sGssGsssss但 ；重新列写联接矩阵 ,W和非零元素矩阵 IJW，将程序 sp4_2.m10G=.2（ ）完善后，应用 sp4_2.m 求输出 7y的响应曲线。 G6(s)G8(s)G10(s)G4(s)G3(s) G7(s)G5(s)G1(s) G2(s)G9(s)- - -
11、y0 y7解： 根据图中 , 拓扑连结关系，可写出每个环节输入 受哪些环节输出iuy iu的影响，iy现列如入下:10293485677890.1uyuyuy把环节之间的关系和环节与参考输入的关系分别用矩阵表示出来， 0UWY（2） 00110001.2000 00W12934185670.2189IJW程序为：P=1,0.01,1,0;0,0.085,1,0.17;1,0.01,1,0;0,0.051,1,0.15;1,0.0067,70,0;1,0.15,0.21,0;0,1,130,0;1,0.01,0.1,0;1,0.01,0.0044,0;WIJ=1,0,1;2,1,1;2,9,-1
12、;3,2,1;4,3,1;4,8,-1;5,4,1;6,5,1;6,7,-0.212;7,6,1;8,6,1;9,7,1;n=9;Y0=10;Yt0=0,0,0,0,0,0,0,0,0;h=0.001;L1=10;T0=0;Tf=5;nout=7;%.形成闭环各系数阵.A=diag(P(:,1);B=diag(P(:,2);C=diag(P(:,3);D=diag(P(:,4);m=length(WIJ(:,1);W0=zeros(n,1);W=zeros(n,n);for k=1:mif (WIJ(k,2)=0);W0(WIJ(k,1)=WIJ(k,3);else W(WIJ(k,1),WI
13、J(k,2)=WIJ(k,3);end;end;Q=B-D*W;Qn=inv(Q);R=C*W-A;V1=C*W0;Ab=Qn*R;b1=Qn*V1;%.数值积分求解Y=Yt0;y=Y(nout);t=T0;N=round(Tf-T0)/(h*L1);R=Y0;for i=1:N;for j=1:L1;k1=Ab*Y+b1*R;k2=Ab*(Y+h*k1/2)+b1*R;k3=Ab*(Y +h*k2/2)+b1*R;k4=Ab*(Y+k3*h)+b1*R;Y=Y+(k1+2*k2+2*k3+k4)*h/6;end;y=y,Y(nout);t=t,t(i)+h*L1;end;t,yplot(t,y); xlabel(Time t/s)ylabel(幅度)
