资源预览内容
第1页 / 共15页
第2页 / 共15页
第3页 / 共15页
第4页 / 共15页
第5页 / 共15页
第6页 / 共15页
第7页 / 共15页
第8页 / 共15页
第9页 / 共15页
第10页 / 共15页
亲,该文档总共15页,到这儿已超出免费预览范围,如果喜欢就下载吧!
资源描述
实验四 常微分方程的数值解 化21华晓【实验目的】1. 掌握用MATLAB 软件求微分方程初值问题数值解的方法;2. 通过实例学习用微分方程模型解决简化的实际问题;3. 了解欧拉方法和龙格-库塔方法的基本思想和计算公式,及稳定性等概念。【实验内容】五、概要:官司争论。将装满放射性废物的圆桶扔到水深300ft的海底,圆桶体积55gal,装满废料的桶重为527.436lbf,在海中浮力为470.327lbf。此外,下沉时受到的阻力与速度成正比,比例系数为0.08lbf/s。实验发现当圆桶速度超过40ft/s时,就会因与海底冲撞而破裂。(1)建立解决上述问题的微分方程模型(2)用数值和解析两种方法求解微分方程,并回答谁赢得了官司。建立模型:下落的圆桶在水中受到三个力,分别是重力,浮力和摩擦力。取向下为正方向,由牛顿第二定律可知:mg-F-f=ma 对于摩擦力,有:f=kv 同时,有:v=dxdta=d2xdt2得到微分方程:md2xdt2+kdxdt+F-mg=0对于该微分方程,有初始条件:x0=0x0=0其中F:圆桶所受浮力F=470.327*0.4536*9.82090.7Nv:圆桶速度,单位m/sK:阻力与速度的比例系数k=0.08*0.4536*9.8/0.30481.1667Ns/mm:圆桶质量 m=527.436*0.4536239.24kg可得,当圆桶速度超过12.192m/s 时就会破裂;水深为91.44米。模型求解1. 数值求解首先编写微分方程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,xsubplot(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变换,设 Lxt=Fs 有ms2Fs+ksFs+F-mgs=0于是Fs=mg-Fs(ms2+ks)=g-Fms+kms2作Laplace逆变换,得-F - mgmexpktm- m + ktk2代入数据解方程得v=13.6366m/s这个速度大于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两个方向的分解,可列出常微分方程如下:dxdt=v1-v2xx2+y2dydt=-v2yx2+y2其初值条件为x(0)=0,y(0)=-d.2问题分析:1解析解分析两式相除,得dxdy=-v1x2+y2v2y+xy=-kxy2+12+xy变量代换,令u=xy, 则有dxdy=ydudy+u代入有ydudy=-ku2+12分离变量积分可得y-k=c(u2+12+u)代入初值,可解出c=-d-k整理得解析解最终表达式为x=-d2(y-d)1-k-(y-d)1+k2数值求解首先建立微分方程:function dx=gh(t,x,v1,v2)d=100;v1=1;v2=2;s=(x(1)2+x(2)2)0.5;dx=v1-x(1)/s*v2;-x(2)/s*v2;Matlab程序如下:ts=0:0.1:100; x0=0,-100; t,x=ode23s(gh,ts,x0);t xsubplot(1,3,1),plot(t,x(:,1),grid%画出x的坐标关于时间的关系图hold on ,z=0*t;plot(t,z,r),hold off%在此关系图上增加一条显示最终x的坐标的线subplot(1,3,2),plot(t,x(:,2),grid%画出小船y坐标关于时间的关系图hold on ,z=0*t;plot(t,z,r),hold off%在此图的基础上增加一条显示终点y坐标的线subplot(1,3,3),plot(x(:,1),x(:,2),g),grid;%画出x y坐标的对应关系线得到如下图形:部分重要数据为:0000.01000.00010.01060.02000.00020.02120.03000.00050.03180.04000.00080.04250.05000.00130.0531.66.30000.3902-0.005966.40000.2903-0.003366.50000.1903-0.001466.60000.0903-0.000366.70000.0000066.80000.0000066.90000.0000067.0000-0.00000由数据可得,小船在66.7s时到达彼岸与解析解相比较:x与y的如下:x+12ckyk+1-12c-ky1-k=0其中c为未知常数,而k=v1v2另外由题目中已知可以解得k=0.5,再由初始条件x=0,y=-100,解得c=-0.01利用matlab做出相应图形,程序如下:y=-100:0.1:0;%给定初值k=0.5 ;x=-y.*(-0.5*(-0.01*y).k)+0.5*(-0.01*y).(-k);%给出相应的函数关系式plot(x,y),gridxlabel(x);ylabel(y);得到图形为:发现数值解与解析解图形重合度极好,故可得66.7s到达彼岸。(一)当V1=0m/s时,function dx=gh(t,x,v1,v2)d=100;v1=0;v2=2;s=(x(1)2+x(2)2)0.5;dx=v1-x(1)/s*v2;-x(2)/s*v2;ts=0:0.1:100;x0=0,-100;t,x=ode15s(gh,ts,x0);plot(t,x(:,1),r),gridhold on,plot(t,x(:,2),g),hold offgtext(x(t);gtext(y(t);pauseplot(x(:,1),x(:,2),g),grid;xlabel(x);ylabel(y);得到图形为:重要数据为: 0 0100.0000 0.1000 099.8000 0.2000 099.6000 0.3000 099.4000. 49.6000 00.8000 49.7000 00.6000 49.8000 00.4000 49.9000 00.2000 50.0000 00.0000故知,小船会在50s时到达彼岸。与解析解相比较:由数学计算可知:此时t=d/v2=100/2=50s结果一致。(二)当V1=0.5m/s时,主程序不变,只把gh.m脚本文件中v1数值改变,即得如下数据:00-1000.10.049947-99.80.20.099796-99.60.30.149545-99.40.40.199194-99.20.50.248742-99 . . .52.90.604925-0.30225530.474462-0.2161853.10.340737-0.1372453.20.203201-0.0678453.30.060869-0.01339得到图形为:水速为0.5m/s时,时间为53.3s时船的距河岸0.060869m,可以认为船到了河对岸
网站客服QQ:2055934822
金锄头文库版权所有
经营许可证:蜀ICP备13022795号 | 川公网安备 51140202000112号