发布时间 : 星期三 文章数学实验第二次作业 - 常微分方程数值求解更新完毕开始阅读
数学实验第二次作业——常微分方程数值求解
实验4常微分方程数值解
实验目的：
1. 练习数值积分的计算；
2. 掌握用MATLAB软件求微分方程初值问题数值解的方法； 3. 通过实例学习用微分方程模型解决简化的实际问题；
4. 了解欧拉方法和龙格——库塔方法的基本思想和计算公式，及稳定性等概念。 实验内容：
3.小型火箭初始质量为1400kg，其中包括1080kg燃料，火箭竖直向上发射是燃料燃烧率为18kg/s，由此产生32000N的推力，火箭引擎在燃料用尽时关闭。设火箭上升是空气阻力正比于速度的平方，比例系数为0.4kg/m，求引擎关闭瞬间火箭的高度，速度，加速度，及火箭到达最高点是的高度，速度和加速度，并画出高度，速度，加速度随时间变化的图形。 解答如下：
这是一个典型的牛顿第二定律问题，分析火箭受力情况； 先规定向上受力为正数
建立数学模型：
A燃料未燃尽前，在任意时刻（t<60s）
火箭受到向上的-F=32000N， 向下的重力G=mg，g=9.8，
向下的阻力f=kv^2, k=0.4, v表示此时火箭速度； 此时火箭收到的合力为F1=（F-mg-f）；
火箭的初始质量为1400kg，燃料燃烧率为-18kg/s; 此刻火箭质量为m=1400-18*t
根据牛顿第二定律知，加速度a=F1/m=(F-mg-f)/(m-r*t)= (32000-(0.4.*v.^2)-9.8.*(1400-18.*t)) 由此可利用龙格-库塔方法来实现，程序实现如下
Function [dx]=rocket[t,x] %建立名为rocket的方程 m=1400;k=0.4;r=-18;g=9.8; %给出题目提供的常数值 dx=[x(2);(32000-(k*x(2)^2)-g*(m+r*t))/(m+r*t)];
%以向量的形式建立方程
[a]=(32000-(k*x(2)^2)-g*(m+r*t))/(m+r*t); %给出a的表达式 End;
ts=0:60; %根据题目给定燃烧率计算出燃料燃尽的时间，确定终点 x0=[0,0]; %输入x的初始值 [t,x]=ode15s(@rocket,ts,x0); %调用ode15s计算
[t,x]; h=x(:,1); v=x(:,2);
plot(t,x(:,1)),grid; %绘出火箭高度与时间的关系曲线 title('h-t');
xlabel('t/s');ylabel('h/m'),pause;
plot(t,x(:,2)),grid ; %绘出火箭速度与时间的曲线关系 title('v-t');
xlabel('t/s');ylabel('v/m/s'),pause;
a=(32000-(0.4.*v.^2)-9.8.*(1400-18.*t))/(1400-18.*t); plot(t,a),grid; %绘出火箭加速度与时间的曲线关系 title('a-t');
xlabel('t/s'),ylabel('a/m^2/s'),pause 火箭高度随时间变化的曲线
火箭速度随时间变化的曲线
火箭加速度随时间变化的曲线
数据过多，故截取部分如下
第一列为时间，第二列为火箭高度，第三列为火箭速度