MATLAB--水塔流量的估计
水塔水流量的估计
摘要:数学建模方法是处理科学理论的一种经典方法，也是解决各类实际问题的常用方法。

本文采用曲线拟合的方法，并利用数学软件MATLAB对水塔流量进行计算，计算结果与实际记录基本吻合。

关键词:建模，流量，拟合，MATLAB
1.问题重述
美国某州的各用水管理机构要求各社区提供用水率（以每小时多少加仑计，英制单位下，1加仑=4.54596dm3，美制单位下，1加仑=3.78533dm3）以及每天所用的总用水量，但许多社区并没有测量流入或流出当地水塔的水量的设备，而只能以每小时测量水塔的水位代替，其精度在0.5%以内。

更为重要的是，无论什么时候，只要水塔中的水位下降到某一最低水位L时，水泵就启动向水塔重新充水直至某一最高水位H，但也无法得到水泵的供水量的测量数据。

因此，在水泵正在工作时，不容易建立水塔中水位与水泵工作时用水量之间的关系。

水泵每天向水塔充水一次或两次，每次大约2小时。

试估计在任何时候，甚至包括水泵正在工作的时间内从水塔流出的流量()
f t，并估计一天的总用水量。

水塔是一个垂直圆柱体，高为40英尺，直径为57英尺。

下表给出了某个小镇某一天的真实数据：
表1
某小镇某天的水塔水位（1m=3.281英尺）
时间（秒）水位
（英
尺）
时间
（秒）
水位
（英
尺）
时间
（秒）
水位
（英
尺）
0 31.75 35932 水泵工作68535 28.42 3316 31.10 39332 水泵工作71854 27.67 6635 30.54 39435 35.50 75021 26.97 10619 29.94 43318 34.45 79154 水泵工作13937 29.55 46636 33.50 82649 水泵工作17921 28.92 49953 32.67 85968 34.75 21240 28.50 53936 31.56 89953 33.89
25223 27.87 57254 30.81 93270 33.40 28543 27.52 60574 30.12
32284 26.97 64554 29.27
2.问题分析
数据的单位转换：
表2
时间（h）
水位
（m）
时间
（h）
水位
（m）
时间
（h）
水位
（m）
0 9.6769 9.98 水泵工作19.04 8.6620
0.92 9.4788 10.93 水泵工作19.96 8.4334
1.84 9.3081 10.95 10.8199 20.84 8.2201
2.95 9.1253 12.03 10.4998 22.01 水泵工作
3.87 8.9864 12.95 10.2103 22.96 水泵工作
4.98 8.8144 13.88 9.9573 23.88 10.5913
5.90 8.6864 14.98 9.6190 24.99 10.3292
7.01 8.5200 15.90 9.3904 25.91 10.1798
7.93 8.3877 16.83 9.1801
8.97 8.2201 17.94 8.9211
流量是单位时间流出的水的体积，由于水塔是正圆柱形，横截面积是常数，在水泵不工作的时段，流量很容易从水位对时间的变化率算出，问题是如何估计水泵供水时段的流量。

水泵供水时段的流量只能靠供水时段前后的流量拟合得到，作为用于拟合的原始数据，我们希望水泵不工作的时段流量越准确越好。

这些流量大体可由两种方法计算: 一是直接对表2中的水位用数值微分算出各时段的流量，用它们拟合其它时刻或连续时间的流量。

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

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

下面我们用第二种方法处理。

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

其实，水泵不工作时段的用水量可以由测量记录直接得到，如表2可知从t=0到t=8.97(h)水位下降了9.6769-8.2201=1.4568 (m)，乘以水塔的截面积就是这一时段的用水量。

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

3.模型假设
供水时段的假设
水泵第1次供水时段为t=9到t=11(h),第2次供水时段为t=20.8到t=23 (h)。

这是根据最低和最高水位分别是8.2201m 和10.8199m 及表2的水位测量记录作出的假设。

其中前3个时刻取自实测数据（精确到0.01h ）,最后1个时刻来自每次供水约两小时的已知条件（从记录看，第2次供水时段应在有记录的22.96h 之后不久结束）。

水泵工作时单位时间的供水量基本为常数，这个常数大于单位时间的平均流量。

流量是单位时间流出水的体积，这里假设流量是对时间的连续函数，即()t h h =。

为简化处理，不影响问题的解决，假设流量与水泵是工作无关。

由于水塔截面积是常数S ,为简单起见，计算中将流量定义为单位时间流出的水的高度，即水位对时间变化率的绝对值（水位是下降的），最后给出结果时再乘以S 即可。

即：
水位对时间的变化率(流量): 任何时刻的流量: ()()S t h t f ⋅'=
4.流量估计
4.1拟合水位-时间函数
从表2 测量记录看，一天有两个供水时段（以下称第1供水时段和第2供水时段）和3个水泵不工作时段(以下称第1用水时段t=0到t=8.97，第2用水时段t=10.95到t=20.48和第3用水时段t=23以后)。

对第1、2用水时段的测量数据分别作多项式拟合，得到水位函数()t h h 11=和()t h h 22=。

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

由于第3时段只有3个测量记录，无法对这一时段的水位作出比较好的拟合，可采用外推的办法解决。

4.2确定流量-时间函数
对于第1、2用水时段，只需将水位函数()2,1,==i t h h i i 求导数即可，对于两个供水时段的流量，则用供水时段前后（水泵不工作时段）的流量拟合得到，并且将拟合得到的第2供水时段流量外推，将第3用水时段流量包含在第2供水时段内，需要
dt
t dh h )
(=
'
拟合四个流量函数。

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

5.算法设计与计算结果
5.1拟合第1、2时段的水位，并得出流量 1第1时段的流速
设t 、h 为已输入的时刻和水位测量记录，实现如下：
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.94,19.04,19.96,20.84,23.88,24.99,25.91];
h=[ 9.6769,9.4788,9.3081,9.1253,8.9864,8.8144,8.6864,8.5200,8.3877,8.2201,10.8199,10.4998,10.2103,9.9573,9.6190,9.3904,9.1801,8.9211,8.6620,8.4334,8.2201,10.5913,10.3292,10.1798];
f1=polyfit(t(1:10),h(1:10),5); %用5次多项式拟合第1用水时段水位()t h h 11=，f1输出5次多项式的系数
b1=polyder(f1); % b1输出多项式（系数为f1）导数的系数，给出水位变化率 tm1=0:0.01:8.97; %将第一用水时段[0 , 8.97]细分
g1=-polyval(b1,tm1); %g1输出多项式b1在tm1点的函数值(取负后边为正值)，即tm1时刻的流量(水位下降的速率)。

2第2时段的流速 实现如下：
f2=polyfit(t(11:21),h(11:21),5); %用5次多项式拟合第2用水时段水位水位()t h h 22=，
⎰⎰'-='=t
t t t dt
h S dt V V 0
f2输出5次多项式的系数
b2=polyder(f2); %b2输出多项式（系数为f2）导数的系数,给出水位变化率
tm2=10.95:0.01:20.84; %将第二用水时段[10.95 , 20.84]细分
g2=-polyval(b2,tm2); %g2输出多项式（b2）在tm2点的函数值(取负后边为正值)，即tm2时刻的流量(水位下降的速率)
第1、2用水时段(水位变化率)曲线图:
5.2拟合供水时段的流量
1
在第1供水时段（t = 9~11）之前（即第1用水时段）和之后（第2用水时段）各取几点，其流量已经得到，用它们拟合第1供水时段的流量。

为使流量函数在t =9 和t =11连续，只取4个点，拟合5次多项式（即曲线必过这4个点）。

拟合5次多项式，实现如下：
q1=-polyval(b1,[7.93,8.97]); %取第1时段在t=7.93,8.97的流量
q2=-polyval(b2,[10.95,12.03]); %取第2时段在t=10.95,12.03的流量
dx=[7.93,8.97,10.95,12.03];。