清华大学数学实验报告2..docxVIP

数学实验报告2 常微分方程数值解 电13 苗键强 2011010645  一、实验目的 1.掌握MATLAB软件求微分方程初值问题数值解的方法； 2.通过实例学习用微分方程模型解决简化的实际问题； 3.了解欧拉方法和龙格-库塔方法的基本思想和计算公式，及稳定性等概念。 二、实验内容 题目1 【问题描述】 小型火箭初始重量为m(0)=1400kg，其中包括m1(0)=1080kg燃料。火箭竖直向上发射时燃料燃烧率为q=18kg/s，由此产生F=32000N的推力，火箭引擎在燃料用尽时关闭。设火箭上升时空气阻力正比于速度v的平方，比例系数为k=0.4kg/m，求引擎关闭瞬间火箭的高度h、速度v、加速度a，及达到最高点时的高度h和加速度a，并画出高度h、速度v、加速度a随时间变化的图形。 【分析与解】 将火箭的运动分为两个阶段：引擎开启阶段和引擎关闭阶段。 引擎开启阶段： 在这一阶段，火箭运动符合以下规律： 运动时间：t0=m1(0)/q=60 s ； 重力加速度：g=9.8 m/s2 即时重力：G=9.8*(1400-18t) ； 推力：F=32000N ； 空气阻力：f=0.4*v2 根据牛顿运动定律，火箭的运动方程为： F-f-G=ma ； 即： a=(32000- 且有： a=dv/ 及： h=v d 用龙格-库塔方法求t=60s时的速度v，在matlab中编程如下： function dv=rockets(t,v) F=32000; g=9.8; m=1400; q=18; dv=(F-0.4*v^2)/(m-q*t)-g; end 在实际计算时，取时间间隔为0.1s，运用上述的方法求解火箭在上升阶段结束时的运动学参数。写出运行函数的代码，在matlab中编程如下： tt=0:0.1:60; v0=0; h(1)=0; [t,v]=ode45(@rockets,tt,v0); [t,v] a=(32000-0.4*v.^2)./(1400-18*t)-9.8; for i=1:600 h(i+1)=h(i)+0.05*(v(i)+v(i+1)); end 通过运行以上代码，在matlab中求得火箭引擎关闭瞬间的速度v=267.2612 m/s加速度a=0.9143 m/s2,高度h=12190 引擎关闭阶段： 在此阶段，火箭的质量为m=1400-1080=320kg，推力F=0，初始速度v(0)=267.2612，h(0)=12190m。且根据牛顿运动定律，火箭的运动方程为： -f-G=ma ； 即： a=(- 且有： a=dv/ 及： h=h(0)+v 用龙格-库塔方法描述火箭在此阶段的运动方程，在matlab中编程如下： function dv=rockets2(t,v) dv=-0.4*v^2/320-9.8; end 接下来求火箭达到最高点时的高度h和加速度a，由于不知道结束时间的具体值，在这里先取一个比较大的时间t=100s，通过比较每段时间结束后的速度v若v0，则停止运行。因而，在matlab中编程如下： tt=60:0.1:100; v(1)=267.2612;h(1)= 12190; [t,v]=ode45(@rockets2,tt,267.2612); for i=1:400 if(v(i)=0v(i+1)=0) break end end; i 结果求得i=114，即t=60+114*0.1=71.4s时，火箭几乎已经到达最高点。所以通过设置时间区间为[60 , 71.4]，编程如下： tt=60:0.1:71.4; v(1)=267.2612;h(1)= 12190; [t,v]=ode45(@rockets2,ts,v(1)); [t,v] 求得在t=71.3s时，v=0.0206m/s，在t=71.4s时，v=-0.9595m/s。因此可认为在t=71.3s时，火箭达到最高点。根据t求解此阶段结束时的加速度a和高度h，在matlab中编程如下：