黄河小浪底调水调沙问题
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
黄河小浪底调水调沙问题
内容摘要:为了确定排沙量与时间、排沙量与水流量的函数关系,我们可以用SAS软件做
线性回归得到排沙量与时间的函数关系式,再利用所求函数在区间[0,24]上进行积分得到总排沙量1.93962亿吨。对于排沙量与水流量之间的关系,按时间分为两段进行拟合,最终用MATHLAB软件来画出图像,确定排沙量与排水量之间的函数关系式。
关键词:调水调沙实验,sas,排沙量,排水量,matlab,拟合,线性回归
问题的提出:在小浪底水库蓄水后,黄河水利委员会进行了多次试验,特别是2004年
6月到7月进行的黄河第三次调水调沙试验具有典型的意义。这次试验首次由小浪底、三门峡和万家寨三大水库联合调度,进行接力式防洪预泄放水,形成人造洪峰进行调水调沙试验成功。这次试验的一个重要目的就是由小浪底上游的三门峡和万家寨水库泄洪,在小浪底形成人造洪峰,冲刷小浪底的库区的沉积泥沙。在小浪底开闸泄洪以后,从6月27日开始三门峡水库和万家寨水库陆续开闸放水,人造洪峰于29日先后到达小浪底,7月3日达到最大流量2720m^3/s,使小浪底水库的排沙量也不断的增加。表一是由小浪底观测站从6月29日到7月10日检测到的试验数据。
问题分析:1、对于问题一,所给数据中水流量x 和含沙量h 的乘积即为该时刻的排沙量
y 即:y=hx 。
2、对于问题二,研究排沙量与排水量的关系,从实验数据中可以看出,开始排沙量随水量增加而增加,而后随水流量的增加而减少,显然变化关系并非线性的关系,为此,把问题分为两部分,从水流量增加到最大值为第一阶段,从水流量最大值到结束为第二阶段,分别来研究水流量与排沙量之间的函数关系。
模型假设:1、水流量和排沙量都是连续的,不考虑上游泄洪所带来的含沙量和外界带来
的含沙量。
2、时间是连续变化的,所取时间点依次为1,2,3,…,24,单位时间为12h.
模型的建立与求解:<一>对于问题一,因为排沙量与时间的散点图基本符合正态曲
线,如图二所示。所以,排沙量的对数与时间的函数关系就应该符合二次函数关系,因而排沙量取对数后,再与时间t 进行二次回归,排沙量取自然后的数据见表2. 假设排沙量与时间函数关系的数学模型是
两边取对数得 Lny=at^2+bt+c
先由表二做出排沙量的自然对数lny 与时间t 的散点图见图一,并利用SAS 软件进行拟合,得到排沙量的自然对数与时间的回归方程为: Lny=-0.0209t^2+0.4298t+10.6321
由回归拟合参数表可知回归方程是显著的,因为相关系数人R^2=0.9629,误差均方S^2=0.0543,说明回归曲线拟合效果很好。 所以排沙量与时间之间的函数关系式为
e
c
bt at y ++=
2^e
t t y 6312
.104289.02^0209.0++-=
图二:排沙量对时间的曲线图
最后对所求出的函数关系在区间[0,24]之间进行积分
结果为总排沙量1.93962亿吨,此与媒体报道的排沙量几乎一样。 <二>对于第二个问题,两个阶段的数据如表三表四所示 表三:第一阶段试验数据
序号 1 2 3 4 5 6 7 8 水流量x
1800 1900 2100 2200 2300 2400 2500 2600 含沙量h 32 60 75 85 90 98 100 102 表四:第二阶段的试验观测数据
序号 1 2 3 4 5 6 7 8 水流量x
2650 2600 2500 2300 2200 2000 1850 1820 含沙量h
116
118
120
118
105
80
60
50
对于第一阶段,有表四用MATLAB 作图(如图三)可以看出其变化趋势,我们用多项式做最小二乘拟合。
设三次拟合函数关系h=a0+a1x+a2x^2+a3x^3 其中a0,a1,a2,a3,为待定系数。
四次拟合函数关系h= a0+a1x+a2x^2+a3x^3+a4x^4 其中a0,a1,a2,a3,a4为待定系数。
dt
e
t t ⎰
++-24
6312
.104289.02^0209.0*60*60*12
图三:第一阶段水流量与排沙量之间的关系图
三次多项式拟合由MATLAB拟合函数求解出a0=a1=0,a2=0.0032,a3=-2.4929.则拟合函数h=0.0032x^2-2.4929x^3,拟合效果如图四所示
图四:三次多项式拟合效果,红线为拟合曲线
类似的四次多项式拟合由MATLAB拟合函数求解出a0=a1=a2=0,a3=0.0121,a4= -7.4347则拟合函数
h=0.0121x^3-7.4347x^4,拟合效果如图五所示
图五:四次多项式拟合效果,蓝线线为拟合曲线
对于第二阶段,有表五用MATLAB作图可以看出其变化趋势,我们用多项式做最小二乘拟合。设三次拟合函数关系h=a0+a1x+a2x^2+a3x^3
其中a0,a1,a2,a3,为待定系数。
四次拟合函数关系h= a0+a1x+a2x^2+a3x^3+a4x^4
其中a0,a1,a2,a3,a4为待定系数。
三次多项式拟合由MATLAB拟合函数求解出a0=a1=0,a2=-0.9475,a3= 464.9601.则拟合函数
h=-0.9475x^2+464.9601x^3,拟合效果如图图六所示
图六:三次拟合函数拟合效果
类似的四次多项式拟合由MATLAB拟合函数求解出a0=a1=0,a2=-0.0013,a3= 1.1219 a4=-354.5952则拟合函数
h=-0.0013x^2+1.1219x^3-354.5952x^4,拟合效果如图七所示
图七:四次拟合函数拟合效果