估计水塔的水流量

摘 要 :建立模型提供水流量及一天的总用水量的较为准确的估计值。本文采用多项式拟合的方法，利用数学软件MATLAB对水塔流量进行计算，两个供水时段水泵功率相当，可见计算结果与实际相吻合。

关键词：建模 流量 拟合

模型假设

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

(2)流量只取决于水位差。水塔中的水位不影响水流量的大小，气候条件、温度变化等也不影响水流量．因为物理学的Torricelli定律指出：水塔的最大水流量正比于水位高度

的平方根，由表中数据有，说明最高水位和最低水位的两个流量几乎相等，

可忽略水位对流速的影响。

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

(4)假设流量是对时间的连续函数。

(5)水塔的水流量与水泵是否工作的状态无关，并不因水泵工作而增加或减少水流量的大小。

(6)定义流量为单位时间流出的水的高度，即水位下降对时间变化率的绝对值，由于水塔截面积是常数s，最后给出结果时乘以s即可。任何时刻的流量：v(t)=-h’(t)*s

流量(水位变化率)估计方法

(1)拟合水位-时间函数，从表中的测量记录看，一天有两次供水时段和三次未供水时段，分别对第1，2未供水时段的测量数据直接作多项式函数拟合，可得到水位函数，根据多项式拟合的特点，为使拟合曲线比较光滑，此处拟合多项式的次数不宜过高，一般为3～6次。对第3未供水时段来说，3个测量记录的数据过少无法得到很好的拟合，故采用外推的方法解决。

T1=[0.00 0.92 1.84 2.95 3.87 4.98 5.90 7.01 7.93 8.97];

H1=[968 948 931 913 901 881 869 849 839 822];

p1=polyfit(T1,H1,5)%用三次多项式拟合第1时段水位，p1输出多项式的系数

a1=polyder(p1)%输出多项式(系数为p1)导数的系数

tp1=0:0.01:8.97;

h1=-polyval(a1,tp1);%输出多项式(系数为a1)在tp1点的函数值(取负后变为正值)，即tp1时刻的流量

hh1=-polyval(a1,T1);

plot(t1,hh1,'s',tp1,h1,'k')

执行命令后得到p1 = -0.0054 0.1380 -1.2945 5.5908 -26.7446 968.1889

a1 = -0.0271 0.5518 -3.8836 11.1817 -26.7446

第一未供水时间段流量图如下所示，

同时得到如下多项式：

H1=-0.0054*t^5+0.1380*t^4-1.2945*t^3+5.5908*t^2-26.7446*t+968.1889

h1=-0.0271*t^4+0.5518*t^3-3.8836*t^2+11.1817*t-26.7446

第二未供水时间段流量，类似可得，

(2)确定流量-时间函数，对于第1，2未供水时段只需将水位-时间函数求导数即可，而对于两个供水时段的流量，则用供水时段前后期(水泵不工作时段)的流量拟合得到，为使流量函数在t=9和t=11连续，只取4个点，用三次多项式拟合得到第一供水时段的时间-流量图。在第二供水时间段之前取t=19.96,20.84两点的流量，用第3未供水时段的3个记录做差分得到两个流量数据，然后用这4个数据点做三次多项式拟合得到第二供水时段与第3未供水时段的时间-流量图。

jt1=-polyval(a1,[7.93 8.97]); %取第1时段在t=7.93,8.97的流量

jt2=-polyval(a2,[10.95 12.03]);%取第2时段在t=10.95,12.03的流量

djt=[7.93 8.97 10.95 12.03];

djh=[jt1 jt2];

pjt=polyfit(djt,djh,5) % 4个记录拟合三次多项式

tpj=7.93:0.01:12.03;

hj=polyval(pjt,tpj);TT=[7.93 8.97 10.95 12.03];

plot(TT,djh,'o',tpj,hj,'b')

执行命令后得到：pjt =-0.0050 0.1530 -1.5324 5.2657 0 0

即得到一个五次多项式，djh=-0.0050*t^5+0.1530*t^4-1.5324*t^3+5.2657*t^2

第一供水时段的时间-流量图第二供水时段的时间-流量图

(3)一天总用水量的估计，总用水量等于两个供水时段和两个未供水时段的用水量之和，它们可由流量对时间的积分得到。

第一用水时段的用水量，其中h1通过梯

形公式计算，第二用水时段的用水量及第一，二供水时段的用水量同理可得。

HH1=0.01*trapz(h1)%第一用水时段的用水量(仍按高度计算)，0.01为步长

HH2=0.01*trapz(h2)%第二用水时段的用水量(仍按高度计算)，0.01为步长

HH3=0.01*trapz(hj)%第一供水时段的用水量

HH4=0.01*trapz(h3)%第二供水时段的用水量

s=pi*[((57/3.281)/2)^2]*10000;

V=(HH1+HH2+HH3+HH4)*s%一天的总用水量(立方厘米)

计算结果

五次多项拟合的时间-流量(以高度计)曲线图

从算法设计和分析可知,计算结果与各时段所用的拟合多项式的次数有关。下表给出的是对第1,2未供水时段分别用5次多项式拟合后得到的用水量结果。

时段

第1未

供水段

第2未

供水段

第1

供水段

第2

供水段

全天用水

用水高度

145.93

260.05

41.67

78.726

526.376

全天的用水量为V = 1.2478e+009立方厘米

流量及总用水量的检验

计算出各时刻的流量可用水位记录的数值微分来检验，各时段的用水高度可以用实际记录的水位下降高度来检验。例如，算得第1未供水段的用水量高度是145.93cm，而实际记录的水位下降高度为968-822=146cm，两者是吻合的；同样地，算得第2未供水段的用水量高度是260.05cm，而实际记录的水位下降高度为1082-822=260cm，两者也是吻合的。供水时段的用水量加上水位上升值260是该时段泵入的水量，除以时段长度得到水泵的功率(单位时间泵入的水量)，而两个供水时段水泵的功率应大致相等。第1、2时段水泵的功率可以计算如下：

P1=[V1+V11]/2 %第1供水时段水泵的功率

P2=(s*HH4+V22)/2.12 %第2供水时段水泵的功率

P1 =3.4600e+008，P2 = 3.7874e+008

模型分析

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

2.模型用到的只是简单易懂，模型容易完成；

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

缺点：

无法准确估计结果的误差。

参考文献：

(l)张德丰， 《MATLAB数值计算方法》 ，北京：机械工业出版社，2010年，p.162-168

附录：

表：某小镇某天的水塔水位(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

时间(小时)

水位(厘米)

时间(小时)

水位(厘米)

时间(小时)

水位(厘米)

0.00

968

10.95

1082

20.84

822

0.92

948

12.03

1050

23.88

1059

1.84

931

12.95

1021

24.99

1033

2.95

913

13.88

996

25.91

1018

3.87

901

14.98

962

4.98

881

15.90

939

5.90

869

16.83

918

7.01

849

17.93

892

7.93

839

19.04

866

8.97

822

19.96

843

详细程序：

T1=[0.00 0.92 1.84 2.95 3.87 4.98 5.90 7.01 7.93 8.97];

H1=[968 948 931 913 901 881 869 849 839 822];

p1=polyfit(T1,H1,5);

a1=polyder(p1); tp1=0:0.01:8.97;

h1=-polyval(a1,tp1); hh1=-polyval(a1,T1);

T2=[10.95 12.03 12.95 13.88 14.98 15.90 16.83 17.93 19.04 19.96 20.84];

H2=[1082 1050 1021 996 962 939 918 892 866 843 822];

p2=polyfit(T2,H2,5);

a2=polyder(p2);

tp2=10.95:0.01:20.84;

h2=-polyval(a2,tp2);

hh2=-polyval(a2,T2);

jt1=-polyval(a1,[7.93 8.97]);

jt2=-polyval(a2,[10.95 12.03]);

djt=[7.93 8.97 10.95 12.03];

djh=[jt1 jt2];

pjt=polyfit(djt,djh,5);

tpj=8.97:0.01:10.95;

hj=polyval(pjt,tpj);TT=[7.93 8.97 10.95 12.03];

TT3=[19.96 20.84 23.88 24.99];

T3=[23.88 24.99 25.91];

H3=[1059 1033 1018];

dt3=diff(T3); %最后三个时刻的两两之差

dh3=diff(H3); %最后三个水位的两两之差

dth3=-dh3/dt3; % T(22)和T(23)的流量

ht3=[-polyval(a2,TT3,dth3)]; % TT3各时刻的流量

p3=polyfit(TT3,ht3,5); %拟合三次多项式

tp3=19.96:0.01:24;

h3=polyval(p3,tp3); %输出第二供水时段外推至T=24各时刻的流量

hh3=polyval(p3,TT3); plot(T1,hh1,'s',tp1,h1,'k',TT,djh,'o',tpj,hj,'b',T2,hh2,'*',tp2,h2,'m',TT3,hh3,'p',tp3,h3,'k')

xlabel('时间-h')

ylabel('流速cm/h')

title('5次拟合的流量曲线')

legend('第一未供水时刻水流速','第一未供水时段','第一供水时刻水流速','第一供水时段','第二未供水时刻水流速','第二未供水时段','第二供水时刻水流速','第二供水时段')

HH1=0.01*trapz(h1)

HH2=0.01*trapz(h2)

HH3=0.01*trapz(hj)

HH4=0.01*trapz(h3)

s=pi*[((57/3.281)/2)^2]*10000;

>> V=(HH1+HH2+HH3+HH4)*s

V1=HH1*s

V2=HH2*s

V11=(968-822)*s

V22=(1082-822)*s

P1=[V1+V11]/2 %第1供水时段水泵的功率

P2=(s*HH4+V22)/2.12 %第2供水时段水泵的功率