实验四 常微分方程的数值解.docxVIP

实验四 常微分方程的数值解化21华晓【实验目的】1. 掌握用MATLAB 软件求微分方程初值问题数值解的方法；2. 通过实例学习用微分方程模型解决简化的实际问题；3. 了解欧拉方法和龙格-库塔方法的基本思想和计算公式，及稳定性等概念。【实验内容】五、概要：官司争论。将装满放射性废物的圆桶扔到水深300ft的海底，圆桶体积55gal，装满废料的桶重为527.436lbf，在海中浮力为470.327lbf。此外，下沉时受到的阻力与速度成正比，比例系数为0.08lbf/s。实验发现当圆桶速度超过40ft/s时，就会因与海底冲撞而破裂。（1）建立解决上述问题的微分方程模型（2）用数值和解析两种方法求解微分方程，并回答谁赢得了官司。建立模型：下落的圆桶在水中受到三个力，分别是重力，浮力和摩擦力。取向下为正方向，由牛顿第二定律可知：对于摩擦力，有：同时，有：得到微分方程：对于该微分方程，有初始条件：其中F：圆桶所受浮力F=470.327*0.4536*9.8＝2090.7Nv：圆桶速度，单位m/sK：阻力与速度的比例系数k==0.08*0.4536*9.8/0.3048＝1.1667N·s/mm：圆桶质量m=527.436*0.4536＝239.24kg可得，当圆桶速度超过12.192m/s 时就会破裂；水深为91.44米。模型求解数值求解首先编写微分方程function dx=yt(t,x)F=2092.11;m=239.41;k=1.1674;g=9.8;dx=[x(2); -(k*x(2)+F-m*g)/m];我们推测在20S内桶就已经沉到海底，现在我们来看看圆桶的速度与深度和时间的关系Matlab程序如下：ts=0:0.01:20;x0=[0,0];[t,x]=ode45(@yt,ts,x0);[t,x]subplot(1,3,1),plot(ts,x(:,1))%画出深度关于时间的坐标图hold on ,z=91.44+0*ts;plot(ts,z ,r),hold off%在此图的基础上增加一条显示深度为91.44m的线subplot(1,3,2),plot(ts,x(:,2))%画出速度关于时间的坐标图hold on ,z=12.192+0*ts;plot(ts,z,r),hold off%在此图的基础上增加一条显示极限速度的的线subplot(1,3,3),plot(x(:,1),x(:,2)),grid; %画出速度与深度的关系图得到如下三个图形：取其中重要的数据如下：0000.10000.00530.10610.20000.02120.21220.30000.04770.31820.40000.08490.42410.50000.13260.5300 . . .13.000087.822213.369913.100089.164213.469513.200090.516113.569013.300091.878013.668513.400093.249813.768013.500094.631613.867413.600096.023313.966713.700097.424914.066013.800098.836514.1653由数据可知，当圆桶到达海底时，速度已经接近13.6m/s,故会破裂。2.求解析解首先，我们直接求解这个微分方程初值问题对x作Laplace变换，设有于是作Laplace逆变换，得代入数据解方程得这个速度大于12.19m/s的临界速度，故圆桶会破裂。Matlab程序设计如下：syms v x %设定字母变量dsolve(D2x=(2344.6-2090.7-1.1667*(Dx))/239.24,x(0)=0,Dx(0)=0)得到结果为：x=217.6244*t+ 44625/exp(0.0048767*t)-44625因此可解得，当x=91.44时，t为13.2835s.而再由x的导数可以求得这个时刻速度为13.636m/s.已经超过极限速度，故圆桶会撞击破裂。综上，数值解与解析解很好的符合，可得出圆桶会撞击海底而破裂，工程师们会赢得这场官司。六、一只小船度过宽为d的河流，目标是起点A正对着的B点，已知河水流速与船在静止的水中的速度之比为k。（1）建立描述小船航线的微分方程模型。（2）设d=100m,=1m/s,=2m/s，用数值方法求渡河所需时间，任意时刻小船位置及航行曲线，作图并与解析解比较；（3）若流速=0, 0.5, 1.5, 2（m/s），结果将如何？1模型建立：此种模型的前提是船并不事先知道水速，否则只要调整一个合适的角度，直接沿直线通过即可。现小船不知道水速，则它的策略应为始终使船头瞄准B点。对速度进行xy两个方向的分解，可列出常微分方程如