资源描述：
《实验一 经典的连续系统仿真建模方法.doc》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。
1、实验一经典的连续系统仿真建模方法一实验目的1．了解和掌握利用仿真技术对控制系统进行分析的原理和步骤。2．掌握机理分析建模方法。3．深入理解一阶常微分方程组数值积分解法的原理和程序结构，学习用Matlab编写数值积分法仿真程序。4．掌握和理解四阶Runge-Kutta法，加深理解仿真步长与算法稳定性的关系。二实验内容1.编写四阶Runge_Kutta公式的计算程序，对非线性模型（3）式进行仿真。（1）将阀位u增大10％和减小10％，观察响应曲线的形状;（2）研究仿真步长对稳定性的影响，仿真步长取多大时RK4算
2、法变得不稳定？（3）利用MATLAB中的ode45()函数进行求解，比较与（1）中的仿真结果有何区别。2.编写四阶Runge_Kutta公式的计算程序，对线性状态方程（18）式进行仿真（1）将阀位增大10％和减小10％，观察响应曲线的形状;（2）研究仿真步长对稳定性的影响，仿真步长取多大时RK4算法变得不稳定？（4）阀位增大10％和减小10％，利用MATLAB中的ode45()函数进行求解阶跃响应，比较与（1）中的仿真结果有何区别。三实验代码及结果编写四阶Runge_Kutta公式的计算程序，对非线性模型（
3、3）式进行仿真。functiondY=f(Y,u)%f函数k=0.2;Qd=0.15;A=2;a1=0.20412;a2=0.21129;dY=zeros(2,1);dY(1)=1/A*(k*u+Qd-a1*sqrt(Y(1)))dY(2)=1/A*(a1*sqrt(Y(1))-a2*sqrt(Y(2)))%RK4文件clccloseY=[1.2,1.4]';u=0.45;Y=0.5;TT=[];XX=[];fori=1:Y:100k1=f(Y,u);k2=f(Y+Y*k1/2,u);k3=f(Y+Y*k2
4、/2,u);k4=f(Y+Y*k3,u);Y=Y+Y*(k1+2*k2+2*k3+k4)/6;TT=[TTi];XX=[XXY];end;plot(TT,XX,'-.')xlabel('x')ylabel('Y')仿真曲线如下u=0.5稳态值Y1=1.2，Y2=1.4。u增大10％和减小10％之仿真曲线如下u=0.45u=0.55步长为60，u=0.5仿真曲线可见步长越大越不稳定采用ode45算法程序如下functiondY=f(Y,u)k=0.2;u=0.5;Qd=0.15;A=2;a1=0.20412;
5、a2=0.21129;dY=zeros(2,1);dY(1)=1/A*(k*u+Qd-a1*sqrt(Y(1)))dY(2)=1/A*(a1*sqrt(Y(1))-a2*sqrt(Y(2)))[T,Y]=ode45('f',[1,200],[1.2,1.4]);%在命令窗口运行以下程序plot(T,Y(:,1),'-',T,Y(:,2),'--')仿真曲线与四阶Runge_Kutta公式的计算一致编写四阶Runge_Kutta公式的计算程序，对线性状态方程（18）式进行仿真functiondY=f1(Y,u
6、)%f1函数k=0.2;Qd=0.00001;A=2;a1=0.20412;a2=0.21129;R1=2*sqrt(1.5)/a1;R2=2*sqrt(1.4)/a2;dY=zeros(2,1);dY(1)=k/A*u+1/A*Qd-1/(A*R1)*Y(1)dY(2)=1/(A*R1)*Y(1)-1/(A*R2)*Y(2)%RK4clccloseY=[0.001,0.001]';u=0.00001;Y=0.1;TT=[];XX=[];fori=1:Y:200k1=f1(Y,u);k2=f1(Y+Y*k1
7、/2,u);k3=f1(Y+Y*k2/2,u);k4=f1(Y+Y*k3,u);Y=Y+Y*(k1+2*k2+2*k3+k4)/6;TT=[TTi];XX=[XXY];end;Yoldonplot(TT,XX(1,:),'--',TT,XX(2,:));xlabel('time')ylabel('Y')gtext('Y1')gtext('Y2')Yoldon仿真曲线如下u=0.00001u增大10％和减小10％之仿真曲线如下u=0.000011u=0.000009当步长40时，曲线如下，可见步长越大越不稳定
8、采用ode45函数求解程序如下functiondY=f1(Y,u)k=0.2;Qd=0.00001;A=2;a1=0.20412;a2=0.21129;R1=2*sqrt(1.5)/a1;R2=2*sqrt(1.4)/a2;dY=zeros(2,1);dY(1)=k/A*u+1/A*Qd-1/(A*R1)*Y(1)dY(2)=1/(A*R1)*Y(1)-1/(A*R2)*Y(2)[T,Y]=ode45('f',[
