估计水塔水流量的求解模型要点.docx

估计水塔水流量的求解模型要点.docx

《估计水塔水流量的求解模型要点.docx》由会员分享，可在线阅读，更多相关《估计水塔水流量的求解模型要点.docx（13页珍藏版）》请在冰豆网上搜索。

估计水塔水流量的求解模型要点

估计水塔水流量的求解模型

摘要

由所给的题目可知，本问题是一个关于如何计算居民用水的问题，由题目给出的表格，可知不同时刻的水位，根据所要求的不同时刻水位的不同入手，此计算问题就可以转化为插值或拟合问题。

这里主要考虑采用插值的方法，可以利用MATLAB软件进行插值和曲线拟合计算并解决一些具体的实际问题。

根据题目建立模型并采用插值的方法进行求解，推算出任何时刻（包括水泵正供水时）从水塔流出的水流量，及一天的总用水量。

关键词：

用水规律与水泵的工作功率原始数据用水规律与水泵的工作功率

一、问题重述

1.1基本情况

某居民区有一供居民用水的圆柱形水塔，一般可以通过测量其水位来估计水的流量。

面临的困难是，当水塔水位下降到设定的最低水位时，水泵自动启动向水塔供水，到设定的最高水位的时候停止供水，这段时间无法测量水塔的水位和水泵的供水量。

通常水泵每天供水一两次，每次约3h.已知水塔是一个高为12.2m，直径为17.4m的正圆柱。

1.2所要解决的问题

现在需要了解该居民区用水规律与水泵的工作功率。

按照设计，当水塔的水位降至最低水位，约8.2米时，水泵自动启动加水；当水位升高到一个最高水位，约10.8米时，水泵停止工作。

可以考虑采用用水率（单位时间的用水量）来反映用水规律，并通过间隔一段时间测量水塔里的水位来估算用水率。

表1是某一天的测量记录数据，测量了28个时刻，但是由于其中有4个时刻遇到水泵正在向水塔供水，而无水位记录（表

中用符号//表示）。

所要解决的问题就是，要估计任何时刻（包括水泵正供水时）

从水塔流出的水流量，及一天的总用水量。

表1水位测量记录（符号//表示水泵启动）

二、问题背景

1991年的美国大学生数学建模竞赛A题（AMCM1991A）,由于它是水库调度、

自来水管理、公共场所的人流量估计等问题的代表,因此有许多文献对其进行了

研究,但一般都是采用差分与拟合的方法。

而由于居民何时用水是无法准确的预

报的,可能引起的水位的变化是随机事件,因此,可以以水容量作为随机变量,建

立一个随机数学模型,不仅可以给出了水塔流量函数,同时还可以讨论水容量函

数的数学期望。

现在越来越多的人开始研究居民区的用水问题，已知居民的用水是由一个圆柱形的水塔提供，要算出任何时刻（包括水泵正供水时）从水塔流出的水流量，及一天的总用水量，根据题目所给的信息和表格，需要建立合适的数学模型，推算任何时刻的用水率和水泵工作功率，从而计算出一天的总用水量。

三、问题分析

问题要求任意时刻的用水率，即求单位时间流出的水的体积，一般称为水流速度或流量。

由于水塔是一个圆柱体，体积V=

π

4Dh

2

可以很容易地通过水位高

度h计算出来，这样在水泵不工作的时间段，水流速度就可以从体积对时间的导数计算出来，由于没有水的体积关于时间的函数表达式，而只能利用问题中给定的原始数据表公式V=

π

4

Dh，计算出离散的在测量时刻的体积V，因此可以考

2

虑用差商代替微商，也就是用离散代替连续的思想。

为提高计算精度，采用二阶差商，即

f（ti）=-∇vi

2

由于所有数据被水泵两次供水分割成三组数据对每组数据的中间数据采用中心差商，前后两个数据不能采用中心差商，改用向前差商、向后差商或用中点公式进行差商。

1.中心差商公式：

∇2vi=

-vi+2+8vi+1-8vi-1+vi-2

12（ti+1-ti）

-vi+2+4vi+1-3vi

2（ti+1-ti）3vi-4vi-1+3vi-2

2（ti-ti-1）

2.向前差商公式：

∇2vi=

3.向后差商公式：

∇2vi=

4.中心公式：

∇vi=

v（i+1）-v（i-1）t（i+1）-t（i-1）

以上分析了水泵不工作的时段，用水率的计算。

对于水泵供水时段的用水率，计算难度较大，我们只好用供水时间段前后的用水率进行插值或拟合而得到。

有了任何时刻的用水率，可以采用数值积分计算一天的总用水量。

四、模型假设与约定

为建模的需要，给出如下假设：

（1）影响水从水塔中流出的流量的唯一因素是公众对水的传统要求．因为表１给出的数据没有提及任何其他的影响因素，我们假定所给数据反映了有代表性的一天，而不包括任何特殊情况，如自然灾害、火灾、水塔溢水、水塔漏水等对水的特殊要求．

（2）水塔中的水位不影响水流量的大小，气候条件、温度变化等也不影响水流量．因为物理学的Torricelli定律指出：

水塔的最大水流量与水位高度的平方根成正比，题目中给出水塔的最高和最低水位分别为10.82米和8.22米，所1

以对于这两种高度，最大水流速的比约为8.22

位的两个流量几乎相等.，说明最高水位和最低水

（3）水泵工作起止时间由它的水位决定，每次充水时间大约为,3个小时．水泵工作性能效率总是一定的，没有工作时需维修、使用次数多影响使用效率问题，水泵充水量远大于水塔水流量.

（4）水泵工作时单位时间的供水量大致为常数，这个常数大于单位时间内从水塔中流出的水流的最大流速，这是因为居民区内一直需要用水，不允许水塔中的水用光。

（5）水塔中水流量是时间的连续光滑函数，与水泵工作与否无关。

这是因为虽然就个别用户而言可能用水量有较大的变化，但由于个人的用水量与整个居民区用水量相比是非常小的，从统计意义上来讲，不太可能同时整个社区的用水量增长或减少。

（6）水塔的水流量曲线可以用一条光滑的曲线来逼近．这时，水流量曲线的两阶导数是连续的．

（7）水泵工作与否完全取决于水塔内水位的高度，且每次加水的工作时间为3小时，根据表1中的数据可知，水泵第一次供水时间段为[8.97,10.95]，第二次供水时间段为[20.84,23.88].

五、符号说明及名词定义

符号说明

t：

测量的时刻,单位为小时;

h：

水塔中水位的高度，是时间的函数，单位为米;

v：

水塔中水的体积,是时间的函数，单位为立方米;

f（t）：

水塔中水流速度，即水流量，是时间的函数，单位为米/小时；yi（i=1~3）：

第i时段用水量，即没有供水时间段用水量；

y12：

第1和第2供水时间段用水量之和；

y：

一天中总用水量；

p：

水泵工作时充水的水流量，是时间的函数，单位为米/小时。

六、模型建立与分析

通过以上对问题的分析，现在的问题已转化为根据某一天已测量的时刻水塔中水的流速，产生在整个区间（24小时）上的函数或函数值，一般来说插值和拟合是两种最常用的方法。

V=

（1）由体积公式：

π

4Dh2可以计算各时刻所对应的水塔中水的体积如下：

时刻

00.921.842.953.874.985.907.017.938.97

体积

2302225422142171213520952066202619951955

时刻

9.9810.9210.9512.0312.9513.8814.9815.9016.8317.93

体积

////25732497242823642295223821832121

时刻

19.0419.9620.8422.0122.9623.8824.9925.91

体积

205920051955////251824612421

整理成表格：

表2水塔中水的体积：

时刻（小时），体积（立方米）

（2）计算水流速度并画出散点图

MATLAB程序如下：

%估计水塔流量

%计算流速，并画流速散点图

%文件名：

ch1.m

t=[0,0.92,1.84,2.95,3.87,4.98,5.90,7.01,7.93,8.97,...

10.95,12.03,12.95,13.88,14.98,15.90,16.83,17.93,...

19.04,19.96,20.84,23.88,24.99,25.91]

h=[9.68,9.48,9.31,9.13,8.98,8.81,8.69,8.52,8.39,...

8.22,10.82,10.50,10.21,9.94,9.654,9.41,9.18,...

8.92,8.66,8.43,8.22,10.59,10.35,10.18]

%计算水塔中水的体积

v=（pi*17.4^2/4）*h

%对每组数据的中间数据计算中心差商

%对每组数据不能计算中心差商的计算向前或向后差商

fori=1:

2

f（i）=-（-v（i+2）+4*v（i+1）-3*v（i））/（2*（t（i+1）-t（i）））%计算向前差商

end

fori=3:

8

f（i）=-（-v（i+2）+8*v（i+1）-8*v（i-1）+v（i-2））/（12*（t（i+1）-t（i）））%计算

中心差商

end

fori=9:

10

f（i）=-（3*v（i）-4*v（i-1）+v（i-2））/（2*（t（i）-t（i-1）））%计算向后差商

end

fori=11:

12f（i）=-（-v（i+2）+4*v（i+1）-3*v（i））/（2*（t（i+1）-t（i）））%计

算向前差商

end

fori=13:

19

f（i）=-（-v（i+2）+8*v（i+1）-8*v（i-1）+v（i-2））/（12*（t（i+1）-t（i）））

%计算中心差商

end

fori=20:

21

f（i）=-（3*v（i）-4*v（i-1）+v（i-2））/（2*（t（i）-t（i-1）））;%计算向后差商

end

i=22;f（i）=-（-v（i+2）+4*v（i+1）-3*v（i））/（2*（t（i+1）-t（i）））%计算向前差商

i=23;f（i）=-（v（i+1）-v（i-1））/（（t（i+1）-t（i-1）））%用中点方法计算差商i=24;f（i）=-（3*v（i）-4*v（i-1）+v（i-2））/（2*（t（i）-t（i-1）））%计算向后差商disp（‘水塔流速’）

f

plot（t,f,‘b*’）

title（‘流速散点图’）;xlabel（‘时间（小时）’）;ylabel（‘流速（立方米/小时）’）

执行文件ch1.m后的结果如下：

t=

Columns1through10

00.92001.84002.95003.87004.98005.90007.0100

7.93008.9700

Columns11through20

10.950012.030012.950013.880014.980015.900016.830017.9300

19.040019.9600

Columns21through24

20.840023.880024.990025.9100

h=

Columns1through10

9.68009.48009.31009.13008.98008.81008.69008.5300

8.39008.9200

Columns11through16

8.66008.43008.220010.590010.350010.1800

v=

1.0e+003*

Columns1through10

2.30182.25422.2138

1.99502.1211

Columns11through16

2.05922.00451.9546

f=

55.5698

f=

55.569842.6466

f=

55.569842.646637.4890

f=

55.569842.646637.4890

2.17102.13532.51822.461142.43122.09492.06642.42072.0283

f=

55.569842.646637.489042.431234.6327

f=

55.569842.646637.489042.431234.632737.0465

f=

55.569842.6466

f=

55.569842.6466

f=

55.569842.6466

33.6004

f=

Columns1through9

55.569842.6466

33.6004

Column10

-197.7749

f=

Columns1through9

37.489042.431237.489042.431237.489042.431237.489042.431234.632737.046534.632737.046534.632737.046534.632737.046529.455629.455654.062129.455654.062129.455654.0621

55.569842.646637.489042.431234.632737.046529.455654.062133.6004

Columns10through11

-197.774952.8416

f=

Columns1through9

55.569842.646637.489042.431234.6327

33.6004

Columns10through12

-197.774952.8416387.6964

f=

Columns1through9

55.569842.646637.489042.431234.6327

33.6004

Columns10through13

-197.774952.8416387.6964-332.1776

f=

Columns1through9

55.569842.646637.489042.431234.6327

33.6004

Columns10through14

-197.774952.8416387.6964-332.1776-275.4368

37.046529.455637.046529.455637.046529.455654.062154.062154.0621

a）由上述实验结果可得水塔中水的流速，经整理成表格如下：

表3水塔中水的流速（水流量）：

时刻（小时），流速（立方米/小时）

b）散点图如下所示：

七、模型的进一步讨论

水泵供水时段的流量只能靠供水时段前后的流量拟合得到。

作为用于拟合的

原始数据，我们希望水泵不工作时段的流量越准确越好。

这些流量大体上可由两

种方法计算，一是直接对表1中的水位用数值微分算出各时段的流量，用他们拟

合其他时刻或连续时间的流量；二是先用表中的数据拟合水位~时间函数，求导数即可得到连续时间的流量。

一般来说数值微分的精度不高，何况测量记录还是不等距的，数值微分的计算尤为麻烦。

因此我们用第二种方法处理。

有了任何时刻的流量，就不难计算一天的总用水量。

其实，水泵不工作时段的用水量可以由测量记录直接得到，如由表4-2-1可知从0到8.97（小时）水位下降了968-822=146（厘米），乘以水塔的截面积就是这一时段的用水量。

这个数值可以用来检验拟合的结果。

1.水流量的估计

（1）拟合水位~时间函数从表4-2-1测量记录看，一天有两个水泵供水时段（以下称第1供水时段和第2供水时段），和3个水泵不工作时段（以下称第1时段0到8.97，第2时段10.95到20.84，和第3时段23.88以后）。

对第1、2时段的测量数据直接分别做多项式函数拟合，得到水位函数。

为使拟合曲线比较光滑，多项式次数不要太高，一般在3~6。

由于第3时段只有3个测量记录，无法对这一时段的水位做出较好的拟合。

（2）确定流量~时间函数对于第1、2时段只需将水位函数求导数即可，对于两个供水时段的流量，则用供水时段前后（水泵不工作时段）的流量拟合得到，并且将拟合得到的第2供水时段流量外推，将第3时段流量包含在第2供水时段内。

（3）一天总用水量的估计总用水量等于两个水泵不工作时段和两个供水时段用水量之和，他们都可以由流量对时间的积分得到。

2.算法设计与程序

（1）拟合第1、2时段的水位，并导出流量设为已输入的时刻和水位测量记录（水泵启动的4个时刻不输入），第1时段各时刻的流量可如下得到：

c1=ployfit（t（1：

10），h（1：

10），3）；%用3次多项式拟合第1时段水位，c1输出3次多项式的系数

a1=ployder（c1）；%a1输出多项式（系数为c1）导数的系数

tp1=0：

0.1：

9；

x1=-ployval（a1，tp1）；%x1输出多项式（系数为a1）在tp1点的函数值（取负后变为正值），即tp1时刻的流量类似地计算第2时段各时刻的流量。

（2）拟合水泵工作时段的流量在第1供水时段（）之前（即第1时段）和之后（即第2时段）各取几点，其流量已经得到，用他们拟合第1供水时段的流量。

为使流量函数在和连续，我们简单地只取4个点，拟合3次多项式（即曲线必过这4个点），实现如下：

xx1=-polyval（a1，[89]）；%取第1时段在t=8，9的流量

xx2=-polyval（a2，[1112]）；%取第2时段在t=11，12的流量

xx12=[xx1xx2]；

c12=polyfit（[891112]，xx12，3）；%拟合3次多项式

tp12=9：

0.1：

11；

x12=polyval（c12，tp12）；%x12输出第1供水时段各时刻的流量在第2供水时段之前取两点的流量，在该时段之后（第3时段）仅有3个水位记录，我们用差分得到流量，然后用这4个数值拟合第2供水时段的流量如下：

dt3=diff（t（26：

28））；%用最后3个时刻的两两之差

dh3=diff（h（26：

28））；%用最后3个时刻的两两之差

dht3=[-dh3./dt3；%t（26）和t（27）的流量

t3=[2020.8t（26）t（27）]；

xx3=[-polyval（a2，t3（1：

2）），dht3]；%取t3各时刻的流量

c3=polyfit（t3，xx3，3）；%拟合3次多项式

tp3=20.8：

0.1：

24；

x3=polyval（c3，tp3）；%x3输出第2供水时段（外推至t=24）各时刻的流量

（3）一天总用水量的估计第1、2时段和第1、2供水时段流量的积分之和，就是一天总用水量。

虽然诸时段的流量已表为多项式函数，积分可以解析的算出，这里仍用数值积分计算如下：

y1=0.1*trapz（x1）；%第1时段用水量（仍按高度计），0.1为积分步长y2=0.1*trapz（x2）；%第2时段用水量

y12=0.1*trapz（x12）；%第1供水时段用水量

y3=0.1*trapz（x3）；%第2供水时段用水量

y=（y1+y2+y12+y3）*237.8*0.01%一天总用水量（）

（4）流量及总用水量的检验计算出的各时刻的流量可用水位记录的数值微分来检验。

用水量y1可用第1时段的水位测量记录中的下降高度968-822=146来检验，类似地，y2用1082-822=260检验。

供水时段流量的一种检验办法如下：

供水时段的用水量加上水位升高值260是该时段水泵注入的水量，除以时段长度得到水泵的功率（单位时间泵入的水量），而两个供水时段水泵的功率应大致相等。

第1、2时段水泵的功率可计算如下：

p1=（y12+260）/2；%第1供水时段水泵的功率（水量仍以高度计）

tp2=20.8：

1：

23；

xp2=ployval（c3，tp2）；%xp2输出第2供水时段各时刻的流量

p2=（0.1*trapz（xp2）+260）/2.2；第2供水时段水泵的功率（水量仍以高度计）

§2.6计算结果从上面的分析和算法设计过程看，计算结果与各时段所用拟合的多项式的次数有关，下面只给出拟合第1、2时段水位函数时，采用不同次数的多项式所得流量及总用水量的结果。

表4-1-2是各时段的用水量，一天总用水量及两个供水时段水泵的功率。

表4-1-22各时段和一天总用水量及两个供水时段水泵的功率（符号和单位见上）

（n1,n2）

y1y2y12y3

y

p1p2

（3,4）

146.18258.1048.5078.50

1263.4

154.25143.59

（5,6）

146.52257.7646.1376.30

1252.5

153.07142.67

§2.7分析及改进由表4-1-2可以看出，第1时段用水量与水位测量记录中的下降高度146相差无几，第2时段用水量与记录中的下降高度260相差无几，所以数据拟合、数值积分的精度是足够的。

对不同次数的拟合多项式第1、2供水时段用水量相差稍大，两供水时段水泵的功率亦有差别，都说明供水时段用3次曲线通过4点的做法不够好，应该多取几点拟合，但要注意让流量曲线在不同时段相接处保持连续。

由图可以看出，流量曲线与原始记录基本上相吻合，零点到十点钟流量很低，十点到下午三点是用水高峰，全天流量平均在22（cm）/h）左右。

若按这个平均流量计算，一天总用水量应为，与表2中的结果非常接近。

八、模型的分析与改进

模型评价

1、模型的优缺点

优点：

1.模型可以用于有一个标准水箱的小镇使用，容易推广；

2.模型用到的知识简单易懂,模型容易完成；

3.模型不仅提供了水流量及一天用水量的较为准确的估计,还可以估计任何时刻的水流量,包括水泵工作时的水流量.

缺点:

无法准确估计结果的误差

九、参考文献

[1]飞思科技产品研发中心编著MATLAB基础与提高电子工业出版社2005年4月

[2]：

姜启源、谢金星、叶俊编，《数学模型—3版》北京：

高等教育出版社，2003

年8月

[3]:

施吉林、刘淑珍、陈桂芝编，《计算机数值方法-第三版》北京：

高等教育出版社，2009年4月