常微分方程数值解54633.docx
《常微分方程数值解54633.docx》由会员分享，可在线阅读，更多相关《常微分方程数值解54633.docx（7页珍藏版）》请在冰点文库上搜索。
常微分方程数值解54633
淮海工学院
实验报告书
课程名称：
数学实验
实验名称：
常微分方程数值解
班级数学091
姓名：
耿萍学号：
090911107
日期：
2012.4.6地点数学实验室
指导教师：
曹卫平成绩：
数理科学系
1.实验目的：
1．用MATLAB软件掌握求微分方程数值解的方法。
2．通过实例学习用微分方程模型解决简化的实际问题。
2.实验内容：
（1）容器盛满水后，底端直径为
的小孔开启。
根据水力学知识，当水面高度为h时，水从小孔中流出的速度为v=0.6
（g为重力加速度，0.6为空口收缩系数）。
1）若容器为倒圆锥形，现测得容器高和上底面直径均为1.2m，小孔直径为3cm，问水从小孔中流完需要多少时间；2分钟时水面高度是多少。
2）若容器为倒葫芦形，现测得容器高1.2m，小孔直径3cm，由底端向上每隔0.1m测出容器的直径如下表所示，问水从小孔中流完需要多少时间；2分钟水面高度是多少。
x
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
1.1
1.2
D
0.03
0.05
0.08
0.14
0.19
0.33
0.45
0.68
0.98
1.10
1.20
1.13
1.00
（图略）
（2）一只小船渡过宽为d的河流，目标是起点A正对着的另一岸B点。
已知河水流速为v1与船在静水中的速度v2之比为k。
1）建立小船航线的方程，求其解析解。
2）设d=100m，v1=1m/s，v2=2m/s，用数值解法求渡河所需时间、任意时刻小船的位置及航行曲线，作图，并与解析解比较。
3.实验步骤：
（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*（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，由于不同高度，倒葫芦形半径不同，用欧拉方程和龙格—库塔方法则水深下降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]
然后利用循环
fork=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,利用条件，当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;
symsh
y=-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.9316
h1=
0.9416
clearall,clc;
g=9.8;
d=0.03;
k1=0;
k2=0;
h=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）^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;
end
h2=x（300）
t=0:
h:
1000*h;
plot（t,x）;
axis（[0,400,0,1.21]）;
h2=
1.0278
（2）functiondx=key（t,x）
d=100;v1=1;v2=2;
s=sqrt（（d-x
（2））^2+（（v1*t）^2））;
dx=[v2*x
（1）/s;v2*（d-x
（2））/s];
ts=1:
1:
50;
x0=[0,0];
[t,x]=ode45（@key,ts,x0）;
[t,x]
plot（t,x）,grid,
gtext（'x（t）'）,gtext（'y（t）'）,pause
plot（x（:
1）,x（:
2））,grid,
gtext（'x'）,gtext（'y'）
ans=
2.000000
4.000003.9980
6.000007.9923
8.0000011.9801
10.0000015.9582
12.0000019.9225
14.0000023.8680
16.0000027.7883
18.0000031.6757
20.0000035.5209
22.0000039.3125
24.0000043.0372
26.0000046.6794
28.0000050.2218
30.0000053.6456
32.0000056.9316
34.0000060.0612
36.0000063.0182
38.0000065.7902
40.0000068.3690
42.0000070.7524
44.0000072.9424
46.0000074.9458
48.0000076.7726
50.0000078.4348
52.0000079.9454
54.0000081.3178
56.0000082.5652
58.0000083.6999
60.0000084.7333
62.0000085.6759
64.0000086.5371
66.0000087.3254
68.0000088.0482
70.0000088.7123
4.实验数据记录及分析（或程序及运行结果）：
通过这次试验，我学会了用MATLAB软件求微分方程数值解的方法，并且通过实例学习了微分方程模型解决简化的实际问题，我觉得学的知识只有运用到实际生活中才发挥了它的作用。
升级会员 