常微分方程数值解(1)

《常微分方程数值解(1)》由会员上传分享，免费在线阅读，更多相关内容在应用文档-天天文库。

1、淮海工学院实验报告书课程名称：数学实验实验名称：常微分方程数值解班级数学091姓名：耿萍学号：090911107日期：2012.4.6地点数学实验室指导教师：曹卫平成绩：11数理科学系111.实验目的：1．用MATLAB软件掌握求微分方程数值解的方法。2．通过实例学习用微分方程模型解决简化的实际问题。2.实验内容：（1）容器盛满水后，底端直径为的小孔开启。根据水力学知识，当水面高度为h时，水从小孔中流出的速度为v=0.6（g为重力加速度，0.6为空口收缩系数）。1）若容器为倒圆锥形，现测得容器高和上底面直径均为1

2、.2m，小孔直径为3cm，问水从小孔中流完需要多少时间；2分钟时水面高度是多少。2）若容器为倒葫芦形，现测得容器高1.2m，小孔直径3cm，由底端向上每隔0.1m测出容器的直径如下表所示，问水从小孔中流完需要多少时间；2分钟水面高度是多少。x00.10.20.30.40.50.60.70.80.911.11.2D0.030.050.080.140.190.330.450.680.981.101.201.131.00（图略）（2）一只小船渡过宽为d的河流，目标是起点A正对着的另一岸B点。已知河水流速为v1与船在静水

3、中的速度v2之比为k。1）建立小船航线的方程，求其解析解。2）设d=100m，v1=1m/s，v2=2m11/s，用数值解法求渡河所需时间、任意时刻小船的位置及航行曲线，作图，并与解析解比较。1.实验步骤：（1）水面直径等于水深，设水深为h时，流量为0.6(π/4)d^2*(gh)^0.5，0.6*(g*h)^(0.5)*π*(d0/2)^2*dt=π/4*h^2*dh则水深下降dh所需时间：dt=-[(π/4)h^2*dh]/[0.6(π/4)d^2*(gh)^0.5]=-[h^1.5*dh]/[0.6d^2*

4、(g)^0.5]水深由1.2m至0定积分得水从小孔流完的时间：T（其中已知d=0.03m，g=9.8m*s(-2）设两分钟（120S）后水深为Xm，由dt=-[(π/4)h^2*dh]/[0.6*(π/4)*d^2*(gh)^0.5]=-[h^1.5*dh]/[0.6d^2*(g)^0.5]则263.93-120=X^2.5/[1.5*d^2*(g)^0.5]以d=0.03m，g=9.8m*s*-2代入上式得水深：X由知容器高1.2m，水深为h时，流量为0.6(π/4)d^2*(gh)^0.5，由于不同高度，倒葫

5、芦形半径不同，用欧拉方程和龙格—库塔方法则水深下降dh所需时间：dt=t(k+1)-t(k)=-[(π/4)h^2*dh]/[0.6(π/4)d^2*(gh)^0.5]=-[h^1.5*dh]/[0.6d^2*(g)^0.5]然后利用循环11fork=1：length(L)，t(k)=((h(k+1)-h(k))*(π/4)*d(k)^2)/(0.6*(π/4)*d^2*(g(1.2-h(k)))^0.5)，T=sum(t).可以求得水从小孔流完的总时间。④设两分钟（120S）后水深为Xm，由S=0,利用条件，当

6、120-s<0.0001时s=s+t(k)，x(k)把d=0.03m，g=9.8m*s(-2)代入上式得水深：X程序>>g=9.8;d=0.03;>>g=9.8;d=0.03;symshy=-h^(1.5)/(0.6*d^2*sqrt(2*g));T=int(y,h,1.2,0);t1=eval(T)x=((T-120)*(1.5*d^2*sqrt(2*g)))^(0.4);h1=eval(x)t1=263.9316h1=0.941611clearall,clc;g=9.8;d=0.03;k1=0;k2=0;h=

7、0.4;x(1)=1.2;forn=1:1000k1=0.15*sqrt(g*(x(n)))*d^2/(-43.6359*x(n)^8+213.0457*x(n)^7-414.873*x(n)^6+410.2075*x(n)^5-218.8936*x(n)^4+62.553*x(n)^3-8.3215*x(n)^2+0.49619*x(n)+.014892)^2;k2=0.15*sqrt(g*(x(n))-h*k1)*d^2/(-43.6359*(x(n)-h*k1)^8+213.0457*(x(n)-h*k1)

8、^7-414.873*(x(n)-h*k1)^6+410.2075*(x(n)-h*k1)^5-218.8936*(x(n)-h*k1)^4+62.553*(x(n)-h*k1)^3-8.3215*(x(n)-h*k1)^2+0.49619*(x(n)-h*k1)+.014892)^2;x(n+1)=x(n)-h*(k1+k2)/2;endh2=x(300)t=0:h:1