现代控制工程及测试技术

《现代控制工程及测试技术》作业班级：硕911姓名：学号：用MATLAB求解微分方程的不同命令求解如下微分方程。性+2性+3业+2y=0.5，y(0)=0,y(0)=0.4,y(0)=-0.2dt3 dt2 dt至少选用两种求解微分方程的命令；在同一幅图上，用不同属性、颜色的曲线表示y(t)和y(t)；解：编写m文件程序如下：% 第一种方法采用函数ode23或ode45解----clearall;closeall;t0=0;tf=15;y0=[0,0.4,-0.2]';[t,y]=ode23('vdpl',t0,tf,y0);figure(1)plot(t,y(:,1),'g-',t,y(:,2),'r--')title('用ode23函数实现微分方程的数值解')xlabel('time/sec')ylabel('value')legend('y','y''')grid% 第二种方法采用dsolve函数求解——t1=0:0.05:15y=dsolve('D3y+2*D2y+3*Dy+2*y=0.5','y(0)=0,Dy(0)=0.4,D2y(0)=-0.2')s=subs(y,t1);dy=diff(y);s1=subs(dy,t1);figure(2)plot(t1,s,'g-',t1,s1,'r--')title('用dsolve函数实现微分方程的符号解')xlabel('time/sec')ylabel('value')legend('y','y''')

grid在方法一中将高阶微分方程等效表达成一阶微分方程组的程序如下：functionxdot=vdpl(t,x)xdot=zeros(3,1);xdot(1)=x(2);xdot(2)=x(3);xdot(3)=-2.*x(3)-3.*x(2)-2.*x(1)+0.5;程序运行的结果及输出图形如图1.1，图1.2所示：y=1/4-3/20*exp(-t)+2/35*7A(1/2)*exp(-1/2*t)*sin(1/2*7N1/2)*t)-1/10*exp(-1/2*t)*cos(1/2*7A(1/2)*t)dy=3/20*exp(-t)+3/140*7人(1/2)*exp(-1/2*t)*sin(1/2*7N1/2)*t)+1/4*exp(-1/2*t)*cos(1/2*7A(1/2)*t)0.40.15用□快23函数实现微分方程的数值解0.350.30.250.20.10.050-0.05-0.11015time/sec图1.10.40.15用□快23函数实现微分方程的数值解0.350.30.250.20.10.050-0.05-0.11015time/sec图1.1用ode23函数求微分方程的解及解的一阶导数用,晶函数实I见微君方程的笋耳解0.35nn0.25U.Jn侦o.IU.U50-nnF用,晶函数实I见微君方程的笋耳解0.35nn0.25U.Jn侦o.IU.U50-nnF61U16liliih/^iii.u.iu图1.2用dsolve函数及diff函数求微分方程的解及解的一阶导数负反馈系统的前向通道和反馈通道传递函数分别为G(s)= 4s+16 ；丑(s)=-^±1^s3+s2+5s+20 (s+3)(s+5)求闭环系统的标准传递函数模型，零极点增益模型，状态空间模型；并将状态空间表达模型转换成可控标准型和可观测标准型。用传递函数模型求系统的单位阶跃响应；用零极点增益模型求单位斜坡响应；用状态空间表达模型求单位脉冲响应。解：编写m文件程序如下：% 第二题 clearall;closeall;sys1=tf([4,16],[1,1,5,20]); %前向通道传递函数sys2=zpk([-1],[-3,-5],2); %反馈通道传递函数%——生成标准传递函数模型、零极点增益模型、状态空间模型disp('闭环系统的零极点增益模型为:')zpksys=feedback(sys1,sys2)[num,den]=tfdata(zpksys,'v');disp('闭环系统的标准传递函数模型为:')tfsys=tf(num,den)disp('闭环系统的状态空间模型为:')[A,B,C,D]=ssdata(zpksys)abcdsys=ss(zpksys)%——下面将状态控制模型转换成可控标准型和可观测标准型%——判断系统是否可控M=ctrb(A,B);r1=rank(M);l1=length(A);ifr1<l1disp('系统是状态不完全可控的！’)；elsedisp('系统是状态完全可控的！’)；disp('将状态空间模型转换为可控标准型：')JA=poly(A);a4=JA(2);a3=JA(3);a2=JA(4);a1=JA(5);a0=JA(6);W=[a1a2a3a41;a2a3a410;a3a4100;a41000;10000];%计算变换矩阵TT=M*W;Ac=inv(T)*A*TBc=inv(T)*BCc=C*TDc=Dend% 判断系统是否可观V=[C'A'*C'A'*A'*C'(A')A3*C'(A')A4*C'];r2=rank(V);l2=size(A,1);ifr2<l2disp('系统是不完全可观的');elsedisp('系统是状态完全可观的');disp('将状态空间模型转换为可观测标准型：')%计算变换矩阵QQ=inv(W*V');Ag=inv(Q)*A*QBg=inv(Q)*BCg=C*QDg=Dend% 求系统的单位阶跃响应，单位斜坡响应，单位脉冲响应t1=0:0.2:5;figure(1)step(tfsys,t1) %传递函数模型求系统的单位阶跃响应title('传递函数模型求系统的单位阶跃响应');grid%——零极点增益模型求单位斜坡响应% 转换为求zpksys与1/s乘积的单位阶跃响应zpk2sys=zpk([],[0],1); %zpk2sys=1/sG=series(zpksys,zpk2sys);t2=0:0.2:5;figure(2)step(G,t2);title('零极点增益模型求单位斜坡响应');grid%------用状态空间模型求单位脉冲响应t3=0:0.2:5;figure(3)impulse(A,B,C,D,1,t3)title('状态空间模型求单位脉冲响应');grid程序运行的结果如下：1)闭环系统的零极点增益模型为：Zero/pole/gain:4(s+3)(s+4)(s+5)(s+3.332)(s+5.153)(s+2.029)(sA2-1.514s+9.532)闭环系统的标准传递函数模型为：Transferfunction:4sA3+48sA2+188s+240sA5+9sA4+28sA3+83sA2+275s+332闭环系统的状态空间模型为:A=0.75682.9931-0.79161.47410-2.99310.7568-0.35220.6558000-5.15332.14790000-3.33162.00000000-2.0288B=00004C=3.2269 0 -0.2685 0.5000 0

系统是状态完全可控的！将状态空间模型转换为可控标准型:Ac=-0.00001.0000-0.0000-0.0000-0.0000-0.00000.00001.00000.00000.0000-0.00000.00000.00001.00000.00000.00000.00000.00000.00001.0000-332.0000Bc=0.0000-0.0000-0.00000.00001.0000Cc=-275.0000-83.0000-28.0000-9.0000240.0000Dc=188.000048.00004.000000系统是状态完全可观的将状态空间模型转换为可观测标准型：Ag=0.0000-0.00000.0000-0.0000-332.00001.0000-0.00000.0000-0.0000-275.00000.00001.00000.00000.0000-83.00000.0000-0.00001.00000.0000-28.00000.0000-0.00000.00001.0000 -9.0000Bg=240.0000188.000048.00004.0000-0.0000Cg=0.0000Dg=0-0.00000.0000-0.0000 1.00002）运行结果如图2.1、图2.2、图2.3所示:骂SIL-t件厚函数模型求系梳的单位所跃响应0.52.5Timefsecj图2.1传递函数模型求系统的单位阶跃响应lllwur咨翅点噌益枳型审中应斜匠□问应0.5 1 1.5 2 2.5 3 3.5 4 4.5 5Time图2.2零极点增益模型求单位斜坡响应状志空间担型末华在脉冲响应0.5706050告'喜-iTimereeclon-图2.3状体空间模型求单位脉冲响应设反馈系统的开环传递函数为G(G(s)H(s)=K(s+1)s(s一1)(s2+4s+16)1）绘制系统的根轨迹。2）确定系统稳定时的K值范围。解：编写m文件程序如下：% 第二题 clearallcloseallnum=[11];den=conv(conv([1,0],[1,-1]),[1,4,16]);G=tf(num,den)rlocus(G) %求系统的根轨迹gridxlabel('实轴');ylabel('虚轴');title('根轨迹图');%-----下面求根轨迹与虚轴交点的增益及频率[K,Wcg]=plzy(G)程序运行结果如图3所示：根轨迹囹642-20-6-4 -2实轴0 2 4-6图3系统的根轨迹图K=23.3163 35.6837Wcg=1.5616 2.5615从图3中可以看出，使系统稳定的K值的取值范围为：23.3163<K<35.6873此时对应的频率为1.5616rad/s和2.5615rad/s。单位反馈系统具有如下的开环传递函数。G(s)H(s)= 10(S+1) (0.1s+1)(5s+1)(8s+1)绘制系统的Nyquist图，判别系统的稳定性。绘制系统的Bode图。确定系统的幅值交界频率，相位交界频率，相位裕量和幅值裕量。解：编写m文件程序如下：% 第四题 clearallcloseallnum=[10,10];den=conv(conv([0.1,1],[5,1]),[8,1]);sys=tf(num,den);% 绘制Nyquist图figure(1)nyquist(sys,'-b');title('系统Nyquist图')；v=[-510-88];axis(v);xlabel('实轴');ylabel('虚轴');% 画系统的Bode图figure(2)bode(sys)grid;title('系统Bode图')%——Wcm为幅值交界频率，Wcg为相位交界频率，Pm为相位裕量，Gm为幅值裕量[Gm,Pm,Wcg,Wcm]=margin(sys)运行结果如图4.1，图4.2所示：

8系统Nyquist图4.1系统的Nyquist图50-15Q□-90系统日况己匡］□8系统Nyquist图4.1系统的Nyquist图50-15Q□-90系统日况己匡］□-50-100-1B0-3 -2 -110 10 10a 1 210 W 10103图4.2系统的Bode图从图4.1可以看出，极坐标图不包围（-1,j0）点，而系统开环传递函数没有右半s平面的极点，因此闭环系统是稳定的。运行程序可知，系统的幅值交界频率Wcm=0.5018，相位交界频率Wcg=Inf，相位裕量Pm=59.4929和幅值裕量Gm=Inf。如题图所示的系统。题5图在rltool环境设计一个PID控制器Gc（s），使闭环系统一对共轭复数极点位于s=_1土顶石。分析系统给定值阶跃扰动时系统的瞬态响应过程和性能指标；并说明极点s=—1±顶另能否作为闭环系统的主导极点，为什么？解：编写m程序如下：% 第五题程序 clearallcloseallnum=1;den=[101];sys=tf（num,den）产生的受控对象模块sys为：Transferfunction:1sA2+1在MATLAB的CommandWindow窗口键入rltool打开根轨迹设计的GUI窗口，导入受控对象模块sys。首先观察C（s）=K时的根轨迹图，可知无论如何调整增益K，根轨迹都不可能经过s=—1土应的希望极点，可见单纯的P控制器不能满足控制要求，因此还需要添加零点或极点。在，平面上添加一实数零点，用鼠标拖动零点观察此时根轨迹的变化趋势，当零点移至-1.49后，根轨迹就经过希望的闭环极点s=—1土志，如图5-1所示，此时采用理想微分作用的PD控制器，控制器的传递函数为c（5）=3（0.67s+1），闭环系统的零点为-1.49，故系统的闭环传递函数为：—、 2(s+1.49) (s)一(s+1-Fs+1+E3 °打开LTIViewerforSISODesignTool窗口，显示闭环阶跃响应，如图5.2所示。性能参数：超调量为41.2%，峰值时间为1.06S，调整时间为3.72S。由于此时的零点为-1.49，极点为s=—1±必，在极点附近有闭环零点，故s=—1±八7不能作为闭环主极点。

4E暨ME_RoodLacjfiEdHor(C)Dragihsclosed-loopmleaicnathelocustoadjustthelocpsan.Ci.rreHtlocstiDh：-14E暨ME_RoodLacjfiEdHor(C)Dragihsclosed-loopmleaicnathelocustoadjustthelocpsan.Ci.rreHtlocstiDh：-1i-1.73iOampinq:0.502NMuraiFrequenDy:2rad/sea图5.1c(.)=k时的根轨迹图Time(sec；图5.2采用PD控制器的阶跃响应控制系统具有如下控制对象传递函数Y(s) 10 U(s)=(s+1)(s+2)(s+3)定义状态变量尤=V，尤=x，x=x1J2 1 3 2利用状态反馈U=-Kx，把闭环极点配置到

，1=—2+顶2点，%=-2-点，%=-1。试用MATLAB求所需的状态反馈增益矩阵k。解：根据系统的传递函数可以得到系统的线型微分方程为：y⑶+6y(2)+11》+6y=10ur.气=%整理得：\X2=x3 写成矩阵向量形式即为：X=Ax+Bu。式中:X=-6y-11y-6y(2)+10uxvxv1X=X2x21-3」010_「0一001，B=0-6-11-610编写m程序如下:% 第六题 %----依题意可知系统的状态方程模型——A=[010;001;-6-11-6];B=[0;0;10];P=[-2+j*2*1.732-2-j*2*1.732-10];K=acker(A,B,P)程序运行结果如下：K=15.3993 4.4999 0.8000D为蒸汽流量信号;题图所示为电站锅炉三冲量给水控制系统典型示意图。图中D为蒸汽流量信号;md为蒸汽流量测量变送装置放大系数，=0.2=0.2：互）。W为给水流量信号;mw为给水流量测量变送装置放大系数，mwmw=0.2（cn、:未）。H为汽包水位信号;mh为汽包水位测量变送装置放大系数，m^=0.5gm^=0.5g1。Imm)Gw(s)=心0.3W(s)=s(11s+1)为给水流量对汽包水位传递函数。G,（s）鱼-0J0035为蒸汽流量对汽包水位传递函数。dD（s）8s1sHzf1.8$为执行器构传递函数。Gc1（s）,Gc2（s）分别为主、副控制器。nd0.6,nw 0.4分别为蒸汽、给水流量分流系数。mhHmhH题7图锅炉三冲量给水控制系统试用SIMULINK及LTIViewer对控制器参数进行设计和整定。基本要求如下：1） 系统的超调量不大于20%；2） 主控制器选择PI或PID控制器；副控制器选择P,PD或PI控制器；3） 选择不同的主、副控制器类型组合，进行仿真和参数阐述整定研究，画出阶跃响应曲线，综合评价比较性能指标。解：（1）搭建SIMULINK下控制系统的仿真模型，如下图7.1所示：图7.1图7.1锅炉三冲量给水控制的Simulink模型(2)按照两步法来整定串级控制系统的参数：将主控制器设置成K=1的工作状态，用衰减曲线法整定副回路，当副回路75%衰减率时的比例增益"=2.5，震荡周期r=30S。然后置副控制器比例增益K2=2.5，求得当主回路75%衰减率时主控制器的比例增益K=1，瞬态响应周期r=；0s。即：pls 1sK2=2.5；气=30s；K]=1；二=30s当副控制器采用P作用，主控制器采用PI作用时K=2.5；K=0.83x1=0.83，T=0.5T=15s，K=Kp1/=0.055P2 P1 1 1 11 i1此时系统单位阶跃响应曲线图如图7.2所示：，%=0.055，%=0.055时的阶跃响应图7.2K2=2.5，K1=0.83此时不满足超调量在20%范围内的要求，对控制器参数进行调整如下:K2=2.5；K1=0.1；K.=0.00001此时系统的单位阶跃响应如图7.3所示:图7.3K2=2.5，K1=0.1，K]=0.00001时的阶跃响应

超调量为6%，响应时间为117s，满足要求。当副控制器采用P作用，主控制器采用PID作用时K=2.5；K=1.25x1=1.25；K=KpgT=°-139；K=KpixTd1=3・75。P2 P1 " i1 d1此时系统的阶跃响应如图7.4所示：图7.4K-2.5，K.=1.25，K=°.139，K=3.75时的阶跃响应由图可知，超调量为57%，不符合要求。经调整后的PID参数如下：K2=2.5，K=1.1，K=°.°39，K11=9此时系统的阶跃响应如图7.5所示：图7.5K2=2.5，K1=1.1，K1=°.°39，K11=9时的阶跃响应超调量为17.5%，调整时间为100s，满足要求。

当副控制器采用PI作用，主控制器采用PI作用时K2=2.5x