1、实验一 连续系统的数字仿真一、实验目的1 熟悉 Matlab 中 m 文件的编写；2 掌握龙格库塔法的基本原理。二、实验设备计算机、MATLAB 软件三、实验内容假设单变量系统如图所示。试根据四阶龙格库塔法，求系统输出 y 的动态响应。1首先把原系统转化为状态空间表达式： ，根据四阶龙格库塔公式，CXybuA可得到： 114321 )(6kkk KKh（1）其中： )()(22)()(343121htbuhKXAtbukkkk（2）根据（1） 、 （2）式编写仿真程序。2在 Simulink 环境下重新对上述系统进行仿真，并和 1 中结果进行比较。四、实验结果及分析要求给出系统输出响应曲线，并
2、分析计算步长对龙格库塔法的影响。计算步长对龙格库塔法的影响：单从每一步看，步长越小，截断误差就越小，但随着步长的缩小，在一定求解范围内所要完成的步数就增加，不但引起计算量的增大，而且可能导致舍入误差严重积累，因此同积分的数值计算一样，微分方程的解法也有选择步长的问题。源程序：r=5;numo=1;deno=1 4 8 5;numh=1;denh=1;num,den=feedback(numo,deno,numh,denh);A,b,C,d=tf2ss(num,den);_+ 1s5312sr=5 yTf=input(仿真时间 Tf= );h=input(计算步长 h=);x=zeros(len
3、gth(A),1);y=0;t=0;for i=1:Tf/h;K1=A*x+b*r;K2=A*(x+h*K1/2)+b*r;K3=A*(x+h*K2/2)+b*r;K4=A*(x+h*K3)+b*r;x=x+h*(K1+2*K2+2*K3+K4)/6;y=y;C*x;t=t;t(i)+h;endplot(t,y)Tf=5 h=0.02五、思考题1 试说明四阶龙格库塔法与计算步长关系，它与欧拉法有何区别。计算步长对龙格库塔法的影响：单从每一步看，步长越小，截断误差就越小，但随着步长的缩小，在一定求解范围内所要完成的步数就增加，不但引起计算量的增大，而且可能导致舍入误差严重积累，因此同积分的数值计
4、算一样，微分方程的解法也有选择步长的问题。区别：四阶龙格库塔法与欧拉法都是基于在初值附近展开成泰勒级数的原理，所不同的是取泰勒级数的项数，欧拉公式仅取到 h 项，四阶龙格库塔法取到 h4 项。实验二 面向结构图的仿真一、实验目的1 掌握连接矩阵及系统状态方程的确定方法；2 掌握面向结构图的仿真方法。二、实验设备计算机、MATLAB 软件三、实验内容假设某系统由三个典型环节组成，如下图所示，求输出量 y 的动态响应。仿真基本步骤：1 给定输入信号，确定典型环节及环节参数；2 确定连接矩阵；3 输入仿真时间和计算步长；4 求 H,H-1 和 Q 阵，确定 A、B 阵；5 根据龙格库塔法求状态方程的
5、解；6 根据 15 编写仿真程序。四、实验结果及分析_+ s12sr=10 y10s0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 500.10.20.30.40.50.60.70.80.9源程序：r=10;P=0 1 1 0;2 1 2 0;10 1 10 0;W=0 0 -1;1 0 0; 0 1 0;W0=1;0;0;Wc=0 0 1;Tf=input(仿真时间 Tf =);h=input( 计算步长 h=);A1=diag(P(:,1);B1=diag(P(:,2);C1=diag(P(:,3);D1=diag(P(:,4);H=B1-D1*W;Q=C1*W-A1;A=in
6、v(H)*Q;B=inv(H)*C1*W0;x=zeros(length(A),1);y=zeros(length(Wc(:,1),1);t=0;for i=1:Tf/h;K1=A*x+B*r;K2=A*(x+h*K1/2)+B*r;K3=A*(x+h*K2/2)+B*r;K4=A*(x+h*K3)+B*r;x=x+h*(K1+2*K2+2*K3+K4)/6;y=y;Wc*x;t=t;t(i)+h;endplot(t,y) 仿真时间Tf=10 计算步长h=0.050 2 4 6 8 10 12024681012图一 仿真曲线五、思考题1 典型环节的确定必须满足什么条件？答： G(S)= 式中
7、u 为典型环节的输入，x 为典型环节的输出。bsadcU(S)X1.为了保证 H 的逆 存在，应严格按照 0 的原则确定每个典型环节，既避免以纯比1- i例、纯微分环节作为典型环节。2.在输入向量不全为阶跃函数的情况下，只要在确定典型环节时，注意使含有微分项系数即（ 0）的环节不直接与参考输入连接。id实验三 连续系统的快速仿真一、实验目的1 熟悉增广矩阵的构建方法；2 掌握连续系统的快速仿真基本原理。二、实验设备计算机、MATLAB 软件三、实验内容假设某系统结构图如下，要求采用快速仿真方法求系统输出响应。仿真基本步骤：1 给定输入信号，确定典型环节及环节参数；2 确定连接矩阵；3 输入仿真
8、时间和计算步长；4 求 H,H-1 和 Q 阵，确定 A、B 阵；5 构建增广矩阵；6 采用增广矩阵法求齐次方程的解和系统输出响应。根据 16 编写仿真程序。4、实验结果和分析源程序：r=10;num,den=series(1,1,0,2,1,2);num,den=series(num,den,10,1,10);num,den=cloop(num,den);A1,B1,C1,D1=tf2ss(num,den);A=A1,B1;0,0,0,0;C=C1,0;Tf=10;h=0.05;k1=eye(size(A);_+ s12sr=10 y10sk2=A*h;k3=k2*k2/2;k4=k3*k2
9、/3;k=k1+k2+k3+k4;x=zeros(size(A1,1),1);10;y=0;t=0;for i=1:Tf/hx=k*x;y=y,C*x;t=t,i*h;endplot(t,y)五、思考题1增广矩阵法和四阶龙格库塔法有何区别？答： .tA!31t2tIeA！t如果取 t的泰勒级数前五项，则增广矩阵的计算与四阶龙格库塔法相同。四阶龙格库塔法具有很高的精度，但运行速度很慢，而增广矩阵法既满足了精度又满足了速度。实验四 离散相似法仿真一、实验目的1 掌握离散化的基本原理；2 掌握非线性系统的离散化仿真方法。二、实验设备计算机、MATLAB 软件三、实验内容已知非线性系统结构图如下，求系
10、统输出响应。1 给定参考输入，输入仿真时间，采样周期（T=0.1s） ；2 将被控对象离散化；y_+r=10 10ss12滞环非线性特性 .05s y3 给定参数初始值；4 计算误差和非线性环节输出；5 计算系统输出；6 参数值更新；7 若仿真时间到，则结束，否则转 1。编写仿真程序，求解系统输出响应。4、实验结果及分析源程序：%backlash.mfunction x,u1=backlash(u1,u,x1,s)if (uu1)if (u-s)=x1) x=u-s;else x=x1;endelse if (u50，幅值裕量 kg10db，试确定串联校正装置。设计基本步骤：1 根据性能指标确
11、定开环增益 k；2 利用确定的开环增益 k 画出未校正系统的 Bode 图，并求其相位裕量 r0 和幅值裕量 kg；3 确定所需增加的超前相位角 ， 一般取 515；0rc4 令超前校正装置最大超前角 ，则 ；cmcsin15 将校正装置的最大超前相位角处的频率 作为校正后系统的剪切频率 ，则：wcw；1|)(|0cjwG6 根据 ，可求得： ；cmcwT17 画出校正后系统 Bode 图，检验指标是否满足，如不满足，增大 重新设计。四、实验结果及分析源程序见课本 P220num0=40;den0=conv(1,0,1,2);Gm1,Pm1,Wcg1,Wcp1=margin(num0,den0
12、);r=50;r0=Pm1;w=logspace(-1,3);mag1,phasel=bode(num0,den0,w);for epsilon=5:15phic=(r-r0+epsilon)*pi/180;alpha=(1+sin(phic)/(1-sin(phic);i1,ii=min(abs(mag1-1/sqrt(alpha);wc=w(ii);T=1/(wc*sqrt(alpha);numc=alpha*T,1;denc=T,1;num,den=series(num0,den0,numc,denc);Gm,Pm,Wcg,Wcp=margin(num,den);if(Pm=r);bre
13、ak;endendprintsys(numc,denc)printsys(num,den)mag2,phase2=bode(numc,denc,w);mag,phase=bode(num,den,w);subplot(2,1,1);semilogx(w,20*log10(mag),w,20*log10(mag1),-,w,20*log10(mag2),-.);grid;ylabel(dB);title(-Go,-.Gc,-GoGc);subplot(2,1,2);semilogx(w,phase,w,phasel,-,w,phase2,-.,w,(w-180-w),:);grid;ylabel
14、(相位（度）);xlabel(频率（radc）)title(校正后：幅值裕量=,num2str(20*log10(Gm),dB ，,相位裕量=,num2str(Pm),度);disp(校正前：幅值裕量=,num2str(20*log10(Gm1),dB,相位裕量=,num2str(Pm1),度);disp(校正后：幅值裕量=,num2str(20*log10(Gm),dB, 相位裕量=,num2str(Pm),度);五、思考题1简述串联超前校正的基本思想。答：串联校正装置的主要作用是通过其相位超前效应来改变频率响应曲线的形状，产生足够大的相位超前脚，以补偿原来系统中元件造成的过大相位滞后。因此
15、校正时应使校正装置的最大超前相位出现在校正后系统的开环剪切频率。实验七 控制系统的计算机辅助分析一、实验目的1掌握控制系统的稳定性分析；2掌握控制系统的时域分析；3学会利用 MATLAB 软件求解控制系统性能指标。二、实验设备计算机、MATLAB 软件三、实验内容1已知闭环系统的开环传递函数为 125342)(34sssG判断系统的稳定性，并给出不稳定极点。num=3 2 1 4 2;den=3 5 1 2 2 1;z,p=tf2zp(num,den);ii=find(real(p)0);n1=length(ii);if(n10)disp(The Unstable Poles are:);di
16、sp(p(ii);else disp(System is Stable);endThe Unstable Poles are:0.4103 + 0.6801i0.4103 - 0.6801i2对于典型二阶系统 22)(nssG试绘制出自然振荡频率 n=6，阻尼比 0.2 时系统的单位阶跃响应，并求出稳态误差、上升时间，最大超调量及调节时间。wn=6;zeta=0.2;num=wn2;den=1,2*zeta*wn,wn2;t=0:0.01:10;y,x,t=step(num,den,t);plot(t,y)%求稳态误 差num0=wn2;den0=1,2*zeta*wn,0;sys=tf(nu
17、m0,den0); sys1=1+sys;sys2=tf(sys1.den,sys1.num); num1=1,0;den1=1;sys3=tf(num1,den1);num2=1;den2=1 0;r=tf(num2,den2); dcg=dcgain(sys2*sys3*r);disp(稳态误 差 dcg= num2str(dcg);%求超调 量M,lab=max(y);M=(M-1)/1)*100;disp(最大超 调量 M= num2str(M) %)%求上升时间t1=t(lab);Val,lab1=min(abs(y(1:lab)-1); Tr=t(lab1);disp(上升 时间
18、Tr= num2str(Tr);%求调节时间s=1;L=length(t);while slab2=find(abs(y(L:end)-1)0.05);if length(lab2)0;s=0;else L=L-1;endendL=L+1;Ts=t(L);disp(调节时间 Ts= num2str(Ts);四、实验结果及分析实验八 积分分离 PID 控制算法及仿真一、实验目的1掌握积分分离 PID 控制算法原理；2掌握积分分离 PID 控制器设计。二、实验设备计算机、MATLAB 软件三、实验内容已知某被控对象传递函数为 ，采样时间为 20s，延迟时间为 4 个采样160)(8seG时间，试设
19、计积分分离 PID 控制器，并求出系统输出响应。积分分离算法： kjdip TkeTeeku0 )1()()()(式中， 为积分项开关系数， )(1ke设计基本步骤：1、给定参考输入，输入仿真时间，采样周期（T=0.1s） ，积分分离 PID 控制器参数；2、被控对象离散化：y(k)=-den(2)y(k-1)+num(2)u(k-5)；3、给定输入信号 r(k)=40；4、设定阈值 0，计算当前误差 e(k)，若 时，采用 PD 控制，避免超调过大，)(ke同时使系统响应较快；当 时，采用 PID 控制，保证精度。)(ke5、计算系统响应 y(k)。四、实验结果及分析%Integration
20、 Separation PID Controllerclear all;close all;ts=20;%Delay plantsys=tf(1,60,1,inputdelay,80);%构造exp(80)，inputdelay不能少dsys=c2d(sys,ts,zoh);num,den=tfdata(dsys,v);u_1=0;u_2=0;u_3=0;u_4=0;u_5=0;y_1=0;y_2=0;y_3=0;error_1=0;ei=0;for k=1:1:200time(k)=k*ts;rin(k)=40;kp=0.80;ki=0.005;kd=3.0;%Delay plantyout
21、(k)=-den(2)*y_1+num(2)*u_5;%I separationerror(k)=rin(k)-yout(k);M=2;if M=1 %Using integration separationif abs(error(k)=30endif u(k)=-110u(k)=-110;endu_5=u_4;u_4=u_3;u_3=u_2;u_2=u_1;u_1=u(k); y_3=y_2;y_2=y_1;y_1=yout(k); error_2=error_1;error_1=error(k);endfigure(1);plot(time,rin,b,time,yout,r);xlabel(time(s);ylabel(rin,yout);figure(2);plot(time,u,r);xlabel(time(s);ylabel(u);五、思考题1、试比较 PID 与积分分离 PID 控制算法。
