自动控制原理4控制系统数字仿真

自控原理实验四：控制系统数字仿真一、实验目的 通过本实验掌握利用四阶龙格-库塔法进行控制系统数字仿真的方法，并分析系统参数改变对系统性能的影响。二、实验方法1、四阶龙格库塔法四、实验设备硬件： PC机一台软件：MATLAB软件，Microsoft Windows XP。五、实验报告1、由公式计算出分别为5%、25%、50%时的分别为0.690、0.404、0.216。画出在K为0到时的闭环根轨迹，如下图再画出分别为0.690、0.404、0.216的阻尼线，求出阻尼线与根轨迹的3个交点。则可求出K分别为30.88、59.25和103.55。K也可以这样算：若系统有超调量，则由主导极点法可知原系统可简化为二阶系统，两个闭环极点共轭靠近虚轴，另一个闭环极点远离虚轴，分别设为，则，故，即可算出K。2、根据仿真结果，绘制阶跃响应曲线并求出(the settling time)和%(the over shoot)当K30.88时，以矩阵形式输入a:10,25,30.88以矩阵形式输入c:0,0,30.88请输入步长h:0.025请输入打印步长mh之m:8请输入迭代次数N*m之N:60the over shoot =4.425886 %the settling time = 3.000000 s.当K=59.25时，以矩阵形式输入a:10,25,59.25以矩阵形式输入c:0,0,59.25请输入步长h:0.025请输入打印步长mh之m:8请输入迭代次数N*m之N:80the over shoot =23.117615 %.the settling time =3.150000 s. 当K=103.55时，以矩阵形式输入a:10,25, 103.55以矩阵形式输入c:0,0,103.55请输入步长h:0.025请输入打印步长mh之m:8请输入迭代次数N*m之N:80the over shoot is =45.722310 %the settling time =4.950000 s六、实验结论可以看出，当系统其它参数不变，根轨迹增益K的值增加时，一对主导极点起作用；调节时间增大，响应速度变慢，快速性降低；超调量增加，振荡加剧，稳定性变坏。附录程序：function y=runge_kutta(f)%输入ai,ci,h,m,N.a=input(以矩阵形式输入a:n);%闭环传递函数分母的系数，除最高项系数1外。cc=input(以矩阵形式输入c:n);%闭环传递函数分子的系数，元素数与a的相同。h=input(请输入步长h:);m=input(请输入打印步长mh之m:);N=input(请输入迭代次数N*m之N:);%计算A,b,cA=0 1 0;0 0 1;-a(3) -a(2) -a(1);b=0 0 1;c=cc(3) cc(2) cc(1);u=1; %时域形式的u(t)，单位阶跃输入时u(t)=1.x=zeros(3,N*m);%给X初值t0=0; %t初值t=t0+0:N*m*h;y0=0; %y初值y=zeros(N*m,1);y(1,:)=y0;fprintf(t的值为%fn,t0); %输出t,y的初值。fprintf(y的值为%fn,y0);f=fun;%函数句柄for i=1:N for j=1:mk1=h*feval(f,t(j+m*(i-1),x(:,j+m*(i-1),u,A,b); k2=h*feval(f,t(j+m*(i-1)+h/2,x(:,j+m*(i-1)+k1/2,u,A,b); k3=h*feval(f,t(j+m*(i-1)+h/2,x(:,j+m*(i-1)+k2/2,u,A,b); k4=h*feval(f,t(j+m*(i-1)+h,x(:,j+m*(i-1)+k3,u,A,b); x(:,j+1+m*(i-1)=x(:,j+m*(i-1)+(k1+2*k2+2*k3+k4)/6; y(j+1+m*(i-1),:)=c*x(:,j+1+m*(i-1); end %四阶龙格库塔算法 fprintf(t的值为%fn,t(m*i+1); fprintf(y的值为%fn,y(m*i+1); %打印(m*i+1)*h步长时的t,yendY,k=max(y); %将一系列y中的最大值赋给Ypercentovershoot=100*(Y-1)/1; %计算超调量（对稳定系统的单位阶跃响应，终值为1）fprintf(the over shoot is %f %.n,percentovershoot); %打印超调量值i=length(t);while (y(i,:)0.98)&(y(i,:)1.02) i=i-1;endsettlingtime=t(i); %找出2%误差带时的调节时间fprintf(the settling time is %f s.n,settlingtime); %打印调节时间值plot(t,y); %利用数值解绘制单位阶跃响应曲线grid;title(step response);xlabel(t/s);ylabel(y/v);function f=fun(t,x,u,A