控制系统仿真实验报告

控制系统仿真实验实验报告实验一经典的连续系统仿真建模方法一实验目的:1了解和掌握利用仿真技术对控制系统进行分析的原理和步骤。2掌握机理分析建模方法。3深入理解阶常微分方程组数值积分解法的原理和程序结构，学习用Matlab编写数值积分法仿真程序。4掌握和理解四阶Runge-Kutta法，加深理解仿真步长与算法稳定性的关系。二实验原理:1非线性模型仿真221122112111HHAdtdHHQukAdtdHduduQuAAkHHARARARH00111012121212三实验内容:1.编写四阶Runge_Kutta公式的计算程序，对非线性模型（3）式进行仿真。（1）将阀位u增大10％和减小10％，观察响应曲线的形状;2.编写四阶Runge_Kutta公式的计算程序，对线性状态方程（18）式进行仿真（1）将阀位增大10％和减小10％，观察响应曲线的形状;（2）研究仿真步长对稳定性的影响，仿真步长取多大时RK4算法变得不稳定？四程序代码:龙格库塔:%Rk4文件function[x,y]=Rk4(f,y0,h,a,b,u)%四阶龙格库塔公式%微分方程组的函数名称，初始值，步长，时间起点，时间终点,阀门开度n=floor((b-a)/h);%求步数x(1)=a;%时间起点y(1,:)=y0;%赋初值fori=1:nx(i+1)=x(i)+h;k1=f(x(i),y(i,:),u);k2=f(x(i)+h/2,y(i,:)+h*k1/2,u);k3=f(x(i)+h/2,y(i,:)+h*k2/2,u);k4=f(x(i)+h,y(i,:)+h*k3,u);y(i+1,:)=y(i,:)+h*(k1+2*k2+2*k3+k4)/6;end状态方程：functiondy=f(t,y,u)dy=zeros(2,1);dy(1)=0.5*(0.2*u+0.15-0.20412*sqrt(y(1)));dy(2)=0.5*(0.20412*sqrt(y(1))-0.21129*sqrt(y(2)));dy=dy'水箱模型:clearu=0.6;%阀位大小y0=[1.51.4];%初始值a=0;b=600;%仿真时间h=1;%仿真步长[T,Y]=Rk4(@f,y0,h,a,b,u);%调用PK4plot(T,Y(:,1),'b-',T,Y(:,2),'m-.')title('非线性模型曲线图');legend('第一个水箱液位H1','第二个水箱液位H2');xlabel('时间t');ylabel('液位高度h');1、阀值u对仿真结果的影响U=0.4;h=1;U=0.6;h=1;2编写四阶Runge_Kutta公式的计算程序，对线性状态方程（18）式进行仿真：%RK_4文件function[x,y]=RK_4(f2,y0,h,a,b,u)%参数表顺序依次是微分方程组的函数名称，初始值向量，步长，时间起点，时间终点,阀门开度n=floor((b-a)/h);%求步数x(1)=a;%时间起点y(1,:)=y0;%赋初值，可以是向量，但是要注意维数fori=1:nx(i+1)=x(i)+h;k1=f2(x(i),y(i,:),u);k2=f2(x(i)+h/2,y(i,:)+h*k1/2,u);k3=f2(x(i)+h/2,y(i,:)+h*k2/2,u);k4=f2(x(i)+h,y(i,:)+h*k3,u);y(i+1,:)=y(i,:)+h*(k1+2*k2+2*k3+k4)/6;%按照龙格库塔方法进行数值求解end状态方程：functiondy=f2(t,y,u)A=2;Ku=0.1/0.5;dy=zeros(2,1);%acolumnvectora12=0.25/sqrt(1.5);a2=0.25/sqrt(1.4);R12=2*sqrt(1.5)/a12;R2=2*sqrt(1.4)/a2;AX=[-1/(A*R12)0;1/(A*R12)-1/(A*R2)];BX=[Ku/A1/A;00];dy=AX*y'+BX*u';dy=dy'end线性函数：clcclearallu=[0.50];y0=[00];a=0;b=500;h=1;[T,Y]=RK_4(@f2,y0,h,a,b,u);Y(:,1)=Y(:,1)+1.5;Y(:,2)=Y(:,2)+1.4;plot(T,Y(:,1),'g-',T,Y(:,2),'r-.')title('线性模型曲线图');legend('第一个水箱的液位H1','第二个水箱的液位H2');xlabel('时间t');ylabel('液位高度h');1阀值u对仿真结果的影响:U=0.4；h=1；0501001502002503003504004505001.41.61.822.22.42.62.8线性模型响应曲线图时间t液位高度hH1H2U=0.6；h=1；0501001502002503003504004505001.41.61.822.22.42.62.83线性模型响应曲线图时间t液位高度hH1H22步长h对仿真结果的影响:U=0.5；h=5；0501001502002503003504004505001.41.61.822.22.42.62.8线性模型响应曲线图时间t液位高度hH1H2U=0.5；h=20；0501001502002503003504004505001.41.61.822.22.42.62.8线性模型响应曲线图时间t液位高度hH1H2U=0.5；h=35；0501001502002503003504004505001.41.61.822.22.42.62.8线性模型响应曲线图时间t液位高度hH1H2U=0.5；h=45；0501001502002503003504004505001.41.61.822.22.42.62.8线性模型响应曲线图时间t液位高度hH1H2当步长为45时仿真结果开始不稳定，仿真步长越大,仿真结果越不稳定。实验二面向结构图的仿真一、实验目的1．掌握理解控制系统闭环仿真技术。2．掌握理解面向结构图的离散相似法的原理和程序结构。3．掌握MATLAB中C2D函数的用法，掌握双线性变换的原理二、实验内容根据上面的各式，编写仿真程序，实现无扰动时给定值阶跃仿真实验1.取Kp=1.78，Ti=85s，T=10s，H2s=80,Qd=0,tend=700，进行仿真实验，绘制响应曲线。三、实验程序及代码响应曲线的绘制：1）线性程序代码：clccleartend=700;T=10;Qd=0;Kp=1.78;Ti=85;A=2;%参数Ku=0.1/0.5;H2set_percent=80;alpha12=0.25/sqrt(1.5);alpha2=0.25/sqrt(1.4);R12=2*sqrt(1.5)/alpha12;R2=2*sqrt(1.4)/alpha2;H1spanLo=0;H1spanHi=2.52;H2spanLo=0;H2spanHi=2.52;H1=1.5;H2=1.4;u=0.5;k=2;Kc=Kp/Ti;bc=Ti;Kd=1/A;ad=1/(A*R12);K1=Ku/A;a1=1/(A*R12);K2=1/(A*R12);a2=1/(A*R2);uc(1)=0;ud(1)=0;u1(1)=0;u2(1)=0;xc(1)=0;xd(1)=0;x1(1)=0;x2(2)=0;yc(1)=0;yd(1)=0;y1(1)=0;y2(2)=0;fort=10:T:tendud(k)=Qd;xc(k)=xc(k-1)+Kc*T*uc(k-1);xd(k)=exp(-ad*T)*xd(k-1)+Kd*(1-exp(-ad*T))/ad*ud(k-1);x1(k)=exp(-a1*T)*x1(k-1)+K1*(1-exp(-a1*T))/a1*u1(k-1);x2(k)=exp(-a2*T)*x2(k-1)+K2*(1-exp(-a2*T))/a2*u2(k-1);y2(k)=x2(k);uc(k)=H2set_percent/100-(y2(k-1)+H2-H2spanLo)/(H2spanHi-H2spanLo);yc(k)=xc(k)+bc*Kc*uc(k);yd(k)=xd(k);y1(k)=x1(k);u1(k)=yc(k);u2(k)=y1(k);k=k+1;endH1_present=(y1+H1-H1spanLo)/(H1spanHi-H1spanLo)*100;%h1的高度H2_present=(y2+H2-H2spanLo)/(H2spanHi-H2spanLo)*100;%h2的高度yc=(yc+u)*100;%控制器的输出plot(0:T:tend,[H1_present',H2_present',yc'])title('线性模型阶跃响应曲线图');legend('H1','H2','yc');xlabel('时间t');ylabel('液位高度h/控制器输出（100%）');曲线图01002003004005006007005060708090100110线性模型阶跃响应曲线图时间t液位高度h/控制器输出（100%）H1H2yc2）非线性程序代码：clccleartend=700;T=10;Qd=0;Kp=1.78;Ti=85;A=2;c=0.5;%参数Ku=0.1/0.5;H2set_percent=80;alpha12=0.25/sqrt(1.5);alpha2=0.25/sqrt(1.4);R12=2*sqrt(1.5)/alpha12;R2=2*sqrt(1.4)/alpha2;H1spanLo=0;H1spanHi=2.52;H2spanLo=0;H2spanHi=2.52;H1=1.5;H2=1.4;u=0.5;k=2;Kc=Kp/Ti;bc=Ti;Kd=1/A;ad=1/(A*R12);K1=Ku/A;a1=1/(A*R12);K2=1/(A*R12);a2=1/(A*R2);uc(1)=0;ud(1)=0;u1(1)=0;u2(1)=0;xc(1)=0;xd(1)=0;x1(1)=0;x2(2)=0;yc(1)=0;yd(1)=0;y1(1)=0;y2(2)=0;fort=10:T:tendud(k)=Qd;xc(k)=xc(k-1)+Kc*T*uc(k-1);xd(k)=exp(-ad*T)*xd(k-1)+Kd*(1-exp(-ad*T))/ad*ud(k-1);x1(k)=exp(-a1*T)*x1(k-1)+K1*(1-exp(-a1*T))/a1*u1(k-1);x2(k)=exp(-a2*T)*x2(k-1)+K2*(1-exp(-a2*T))/a2*u2(k-1);y2(k)=x2(k);uc(k)=H2set_percent/100-(y2(k-1)+H2-H2spanLo)/(H2spanHi-H2spanLo);yc(k)=xc(k)+bc*Kc*uc(k);ifyc(k)cyc(k)=c;elseifyc(k)-cyc(k)=-c;endyd(k)=xd(k);y1(k)=x1(k);u1(k)=yc(k);u2(k)=y1(k);k=k+1;endH1_present=(y1+H1-H1spanLo)/(H1spanHi-H1spanLo)*100;H2_present=(y2+H2-H2spanLo)/(H2spanHi-H2spanLo)*100;yc=(yc+u)*100;plot(0:T:tend,[H1_present',H2_present',yc']);title('含有非线性环节的模型阶跃响应曲线图');legend('H1','H2','yc');xlabel('时间t');ylabel('液位高度h/控制器输出（100%）');曲线图：0