水塔流量估計

1、水塔流量估計水塔流量估計摘摘要要 .1關(guān)鍵字關(guān)鍵字 .2一一.問題重述問題重述 .2二二.問題分析問題分析 .2三三.模型假設(shè)模型假設(shè) .3四問題求解四問題求解 .3（1）數(shù)據(jù)差值法 .6（2）數(shù)值擬合法 .11五五.結(jié)果簡單分析結(jié)果簡單分析 .18六模型檢驗六模型檢驗 .18七模型改進七模型改進 .18八八.實驗小結(jié)實驗小結(jié) .18九參考文獻九參考文獻 .18附錄附錄 .181摘要針對水塔流量估計問題，由于水塔流量與水位變化密切相關(guān)，所以通過水位變化來確定其流量；用流速和時間的關(guān)系來確定水位的變化，即水位差；由物理學中 Torricelli 定律確定流速，繼而解決具體問題。我們對某一天水塔水

2、位統(tǒng)計數(shù)據(jù)進行分析，將其分為水泵供水時段和水泵不工作時段，然后通過插值與擬合模型分別進行計算。計算結(jié)果如下：圖 1.1 插值法算得各時段句全天的用水量（用高度計）第一未供水時段第二未供水時段第一供水時段混合時段全天拉格朗日插值法145.6231258.866454.0689127.8244586.3828分段線性插值法147.1430258.969749.6051115.0578570.7756三次樣條插值法145.6870258.654753.3334120.0178577.6929圖 1.2 數(shù)據(jù)擬合法算得各個時段及全天的用水量（用高度計）時段第一未供水時段第二未供水時段第一供水時段混合時

3、段全天用水高度146.5150257.760546.131776.3076526.7148關(guān)鍵字：流量 插值與擬合 torricelli 定律 一.問題重述:某一地區(qū)的用水管理機構(gòu)要求各社區(qū)提供以每小時多少加侖計的用水量以及每天所用水的總量。許多社區(qū)沒有測量流入或流出水塔的水量裝置，他們只能通過測量水塔每小時的水位高度來代替，其誤差不超過 0.05.但是，當水塔中的水位下降到最低水位 L 時水泵就自動啟動向水塔輸水直到最高水位 H，在此期間無法測量水泵的供水量。這樣，當水泵正在輸水時就不容易建立水塔中水位和用水量之間的關(guān)系，水泵每天輸水一次或兩次，每次約兩小時。當水位降至約 27.00ft（英

4、尺）時，水泵開始啟動向水塔供水，當水位升到35.50ft 時，水泵自動停止工作。已知某一小鎮(zhèn)的水塔高度為 40ft，直徑為 57ft 的正圓柱，下表記錄的是某一天的水塔水位的真實數(shù)據(jù)（注：1ft=0.3024m） 。表 1.1 某一天水塔水位統(tǒng)計表時間/s水位/0.01ft時間/s水位/0.01ft203175466363350331631104995332606635305453936316710619299457254308713937294760574301217921289264554292721240285068535284225223279571854276728543275275

5、021269732284269779254水泵開動35932水泵開動82649水泵開動39332水泵開動859683475394353550899533397433183445932703340試估計任何時刻（包括水泵正在輸水時間）從水塔流出的水流量和一天的用水總量。二.問題分析：流量是單位時間內(nèi)流出水的體積，由于水塔是正圓柱形，橫截面積時常數(shù)，所以在水泵不工作時段，流量很容易根據(jù)水位相對時間的變化率算出。問題的難點在于如何估計水泵供水時段的流量。水泵供水時段的流量只能依據(jù)供水時段前后的流量經(jīng)差值或擬合得到。作為用于差值或擬合的原始數(shù)據(jù)，我們希望水泵不工作時段的流量越準確越好。事實上，水泵不工

6、作時段的用水量可以由測量記錄直接得到嗎，由表中記錄的下降水位乘以水塔的橫截面積就該時段的用水量。這個數(shù)值可以用來檢驗數(shù)據(jù)差值或擬合的結(jié)果。為了表示方便，我們將表中的數(shù)據(jù)全部化為國際標準單位（如下表） ，時間用小時（h） ，高度用米（m）。表 2.1 一天內(nèi)水塔水位記錄時間（h）水位（m）時間（h）水位（m）09.6812.9510.210.929.4513.889.941.849.3114.989.652.959.1315.909.413.878.9816.839.184.988.8117.938.925.908.6919.048.667.008.5219.968.437.938.3920.8

7、48.228.978.2222.02水泵開動39.98水泵開動22.96水泵開動10.93水泵開動23.8810.5910.9510.8224.9910.3512.0310.5025.9110.18三.模型假設(shè)（1）假設(shè)流量只取決于水位差，與水位本身無關(guān)。由水塔的最低和最高水位分別是 8.1648m（27x0.3024）和 10.7352m（35.50 x0.3024） ，設(shè)出口的水位為零，從小孔流出的流速正比于水面高度的平方根，可忽略水位對流速的影響。（2）假設(shè)流量是時間的連續(xù)函數(shù)。為計算簡單，不妨將流量定義成單位時間流出水的高度，即水位對時間變化率的絕對值（水位時下降的） ，得到結(jié)果后在乘

8、以水塔橫截面積即可。（3）已知水塔橫截面積為 S=（57*0.3048）2*/4=237.8（平方米） 。四問題求解首先根據(jù)表 2.1 的數(shù)據(jù)，用 MATLAB 軟件作出水位時間散點圖 4.1。程序見附錄 1圖 4.1 水位時間散點圖下面計算流量與時間的關(guān)系。根據(jù)數(shù)據(jù)散點圖 1.0,一種簡單的處理方法是先將表 2.0 中的數(shù)據(jù)分為三段，然后對每一段數(shù)據(jù)做如下處理：設(shè)某段數(shù)句為（x0，y0） ， （x1,y1）,（xn，yn），相鄰數(shù)據(jù)中點的平均流速采用公式：流速=（左端點的水位-右端點的水位）/區(qū)間長度，即：v（）=(yi-y(i+1)/(x（i+1）-xi)2) 1(xiix4每段數(shù)據(jù)首尾點

9、的柳樹采用下面的公式計算：v(x0)=(3y0-4y1+y2)/(x2-x0),v(xn)=(-3yn+4y(n-1)-y(n-2).用以上公式算得流速與時間之間的數(shù)據(jù)表 4.0 如下：表 4.1 流速與時間關(guān)系數(shù)據(jù)表時間（h）流速（cm/h）時間（h）流速（cm/h）029.8913.4229.030.4621.7414.4326.361.3818.4815.4426.092.39516.2216.3724.733.4116.3017.3823.644.42515.3218.4923.425.4413.0419.5025.006.4515.4520.4023.867.46513.9820.8

10、422.178.4516.3522.02水泵開動8.9719.2922.96水泵開動9.98水泵開動23.8827.0910.93水泵開動24.4321.6210.9533.5025.4518.4811.4929.6325.9113.3012.4931.52由表 4.1 作出流速時間散點圖 4.2，程序見附錄 2：圖 4.2 流速時間關(guān)系散點圖5下面用數(shù)據(jù)差值和數(shù)據(jù)擬合兩種方法來估計水塔水流量：（1）數(shù)據(jù)差值法 由表 4.1，對水泵不工作時段 1 和水泵不工作時段 2 采用差值方法。可以得到任意時刻的流速，從而可以知道任意時刻的流量。對于水泵工作時段數(shù)據(jù)太少，我們將它與水泵工作時段 2 合并一

11、同進行差值。由于最后一段水泵不工作時段數(shù)據(jù)太少，我們將它與水泵工作時段 2 合并一同進行差值處理（該時段一下簡稱混合時段） 。這樣，總共需要對四段數(shù)據(jù)（第 1，2 未供水時段，第 1供水時段，混合時段）進行差值處理。A.第 1 未供水時段數(shù)據(jù)差值求解，分別用拉格朗日插值法，分段線性差值及三次樣條差值三種方法算出流量函數(shù)和用水量（用水高度） 。圖 4.3 中曲線 lglr，fdxx，scyt 分別表示用拉格朗日差值法，分段線性差值法及是三次樣條差值法對第 1 未供水時段數(shù)據(jù)差值得到的流速曲線，程序見附錄 3圖 4.3 第 1 未供水時段流速曲線B. 第 2 未供水時段數(shù)據(jù)差值求解，分別用拉格朗日

12、插值法，分段線性差值及三次樣條差值三種方法算出流量函數(shù)和用水量（用水高度） 。圖 4.4 中曲線 lglr，fdxx，scyt 分別表示用拉格朗日差值法，分段線性差值法及是那次樣條差值法對第 2 未供水時段數(shù)據(jù)差值得到的流速曲線， 程序代碼見附錄 4：6圖 4.4 第 2 未供水時段流速曲線C.第 1 供水時段數(shù)據(jù)差值求解，分別用拉格朗日插值法，分段線性差值及三次樣條差值三種方法算出流量函數(shù)和用水量（用水高度） 。圖 4.5 中曲線 lglr，fdxx，scyt 分別表示用拉格朗日差值法，分段線性差值法及是那次樣條差值法對第 2 未供水時段數(shù)據(jù)差值得到的流速曲線，程序代碼見附錄 5：圖 4.5

13、 第 2 未供水時段流速曲線7D.混合時段數(shù)據(jù)差值求解，分別用拉格朗日插值法，分段線性差值及三次樣條差值三種方法算出流量函數(shù)和用水量（用水高度） 。圖 4.2 中曲線 lglr，fdxx，scyt 分別表示用拉格朗日差值法，分段線性差值法及是那次樣條差值法對第 2 未供水時段數(shù)據(jù)差值得到的流速曲線:程序執(zhí)行后顯示結(jié)果，程序代碼見附錄 6：圖 4.6 混合時段流速曲線綜合得到：表 4.2 差值法算得個時段及全天的用水總量（以高度計）時段及方法第 1 未供水時段第 2 未供水時段第 1 供水時段混合時段全天拉格朗日插值法145.6231258.866454.0689127.8244586.3828

14、分段線性差值法147.1430258.969749.6051115.0578570.7756三次樣條差值法145.6870258.654753.3334120.0178577.6929圖 4.6 是用分段線性及三次樣條差值方法算得全天水塔水流速與時間的關(guān)系曲線，圖中曲線 scyt 表示三次樣條差值曲線，曲線 fdxx 表示分段線性差值的曲線，程序代碼見附錄 7：8圖 4.7 全天水塔流速與時間的關(guān)系曲線(2)數(shù)值擬合法1）擬合水位時間函數(shù)根據(jù)表 4.1 中的數(shù)據(jù)，分別對第 1、2 未供水時段的測量數(shù)據(jù)直接做三次多項式擬合，得到水位與實踐的關(guān)系函數(shù)。程序中，變量 x1 存放了以 0.1 為時間步

15、長第 1 未供水時段各個時刻的水位高度。圖 4.7 繪出的是第 1 未供水時段水位與時間的三次多項式擬合曲線，程序代碼見附錄 8：圖 4.8 第 1 未供水時段水位與時間的三次多項式擬合曲線9下圖4.9，Matlab程序?qū)Φ?未供水時段的數(shù)據(jù)進行擬合,得到的第2未供水時段水位與時間三次多項式擬合曲線，Matlab程序見附錄9.圖 4.9 第 2 未供水時段水位與時間的三次多項式擬合曲線2）確定流量時間函數(shù)對第 1,2 未供水時段的水位求導可得流量，用三次多項式擬合第 1,2 未供水時段的流速與時間關(guān)系曲線，程序代碼見附錄 10：圖 4.10 第 1,2 未供水時段的流速與時間關(guān)系曲線（三次多項

16、式擬合）10與圖 4.6 中相應部分段對比發(fā)現(xiàn)，用三次多項式擬合效果不是很好，若改用 5 次多項式擬合則得效果較好的第 1,2 未供水時段的流速與時間關(guān)系曲線，見圖 4.9(程序代碼見附錄 11):圖 4.11 第 1,2 未供水時段的流速與時間關(guān)系曲線（五次多項式擬合）第 1 供水時段的流速則用前后時刻的流速擬合得到。為使流速函數(shù)在 t=9和 t=11 處連續(xù)，只取 4 個點，用 3 次多項式擬合得第 1 供水時段流速與時間關(guān)系曲線圖 4.12，同圖 4.1 中相應部分段較為吻合，程序代碼見附錄 12：圖 4.12 第一供水時段流速與時間的三次多項式似合曲線11對于混合時段，在第 2 供水時

17、段之前取 t=19.96,20.84 兩點的流速，用第三末供水時段的 3 個記錄做差分得到兩個流速數(shù)據(jù) 21.62,18.48，然后用這 4 個數(shù)據(jù)做 3 次多項式似合得到混合時段的流速與時間的關(guān)系曲線如圖 4.13 所示，同數(shù)據(jù)插值法得到的曲線圖 4.1 中相應曲線段較為吻合，程序代碼見附錄 13：圖 4. 11 混合時段流速與時間的三次多項式擬合曲線3）估計一天的用水總量分別對供水的兩個時段和不供水的兩個時段積分（流量對時間積分）并求和得到一天的用水總量約為 526.8935（此數(shù)據(jù)是用水總高度 cm）.表 4.3.1 列出各段用水量，同數(shù)據(jù)插值法算得的表 4.1 數(shù)據(jù)相比，較為吻合。表

18、4.3 數(shù)據(jù)擬合算得各個時段的用水高度（一）時段第 1 未供水時段第 2 未供水時段第 1 供水時段混合時段全天用水高度145.67260.6646.6073.9635526.89354）流量及用水總量的檢驗計算出各時刻的流量可以用水位記錄數(shù)據(jù)的數(shù)值微分來檢驗，各時段的用水高度可以用實際記錄的水位下降高度來檢驗。例如，算得第 1 未供水時段的用水高度是 145.67，而實際記錄的水位下降高度為 968-822=146，兩者是吻合的同樣地，算得第 2 未供水時段的用水高度是 260.66，而實際記錄的水位下降高度為 1082-822=260，兩者也是吻合的。從算法設(shè)計和分析可知，計算結(jié)果與各時段

19、所用的擬合多項式的次數(shù)有關(guān)。表 4.4 給出的是對第 1，2 未供水時段分別用 5，6 次多項式擬合得到的用水高度。12表 4.4 數(shù)據(jù)擬合法算得各個時段的用水高度（二）時段第 1 未供水時段第 2 未供水時段第 1 供水時段混合時段全天用水高度146.5150257.760546.131776.3076526.7148五.結(jié)果簡單分析對于數(shù)據(jù)差值法，由表 4.1 可以看出，使用三次樣條差值法得到的第 1,2 未供水時段的用水高度結(jié)果 145.6870 和 258.6547 與表中記錄的下降高度 146 和260 相差不大，說明差值結(jié)果與原始數(shù)據(jù)比較吻合。用三次樣條差值法估計全天的用水量約為

20、539.4450 x237.8x10=1282800.21（升） 。對于數(shù)據(jù)擬合法，由表 4.3.2 可得全天的用水總量約為526.7148x237.8x10=1252527.7944（升） ，同數(shù)據(jù)差值法得到的結(jié)果也很接近。六模型檢驗七模型改進八.實驗小結(jié)本次實驗主要進行水塔水量的估算。第一種估計為數(shù)據(jù)差值方法，我們用了三種不同的差值方法進行估計。在求解的過程中，可以熟悉數(shù)據(jù)差值的理論和方法。第二種估計方法為數(shù)據(jù)擬合法， 用多項式進行擬合，得到水塔水流量的估計。九參考文獻附錄：下面程序是對上面問題的補充，均由 matlab 軟件編寫：附錄附錄 1：附錄附錄 2：x=xlsread(x=xls

21、read(stgj.xlsstgj.xls, ,sheet1sheet1, ,A2:A29A2:A29););y=xlsread(y=xlsread(stgj.xlsstgj.xls, ,sheet1sheet1, ,B2:B29B2:B29););plot(x,y,plot(x,y,b*b*) )xlabel(xlabel(xx軸軸 );ylabel();ylabel(yy軸軸 ););x=xlsread(x=xlsread(stgjls.xlsstgjls.xls, ,sheet1sheet1, ,A2:A32A2:A32););y=xlsread(y=xlsread(stgjls.xls

22、stgjls.xls, ,sheet1sheet1, ,B2:B32B2:B32););plot(x,y,plot(x,y,b*b*) )xlabel(xlabel(xx軸軸 );ylabel();ylabel(yy軸軸 ););13附錄補充：拉格朗日函數(shù)文件附錄附錄 3：functionunction y=lglrcz(x0,y0,x)y=lglrcz(x0,y0,x)n=length(x0);n=length(x0);m=length(x);m=length(x);forfor i=1:mi=1:m z=x(i);z=x(i); s=0.0;s=0.0; forfor k=1:nk=1:n

23、 p=1.0;p=1.0; forfor j=1:nj=1:n ifif j=kj=k p=p*(z-x0(j)/(x0(k)-x0(j);p=p*(z-x0(j)/(x0(k)-x0(j); endend endend s=p*y0(k)+s;s=p*y0(k)+s; endend y(i)=s;y(i)=s;endendt=0,0.46,1.38,2.395,3.41,4.425,5.44,6.45,7.465,8.45,8.97;v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.27;t0=0:0.1:8.9

24、7;lglr=lglrcz(t,v,t0);lglrjf=0.1*trapz(lglr)fdxx=interp1(t,v,t0);fdxxjf=0.1*trapz(fdxx)scyt=interp1(t,v,t0,spline);sancytjf=0.1*trapz(scyt)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)gtext(lglr)gtext(fdxx)gtext(scyt)14附錄附錄 4t=10.95,11.49,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.4t=10.95,11.49

25、,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.40,20.84;0,20.84;v=33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.8v=33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.86,22.17;6,22.17;t0=10.95:0.1:20.84;t0=10.95:0.1:20.84;lglr=lglrcz(t,v,t0);lglr=lglrcz(t,v,t0);lglrj

26、f=0.1*trapz(lglr)lglrjf=0.1*trapz(lglr)fdxx=interp1(t,v,t0);fdxx=interp1(t,v,t0);fdxxjf=0.1*trapz(fdxx)fdxxjf=0.1*trapz(fdxx)scyt=interp1(t,v,t0,spline);scyt=interp1(t,v,t0,spline);sancytjf=0.1*trapz(scyt)sancytjf=0.1*trapz(scyt)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)plot(t,v,*,t0,lglr,r,t0,fdxx,g

27、,t0,scyt,b)gtext(lglr)gtext(lglr)gtext(fdxx)gtext(fdxx)gtext(scyt)gtext(scyt)附錄附錄 5：t=0,0.46,1.38,2.395,3.41,4.425,5.44,6.45,7.465,8.45,8.97,10.95,11.t=0,0.46,1.38,2.395,3.41,4.425,5.44,6.45,7.465,8.45,8.97,10.95,11.49,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.40,20.84;49,12.49,13.42,14.43,

28、15.44,16.37,17.38,18.49,19.50,20.40,20.84;v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.2v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.29,33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.869,33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.

29、00,23.86,22.17;,22.17;t0=8.97:0.1:10.95;t0=8.97:0.1:10.95;lglr=lglrcz(t,v,t0);lglr=lglrcz(t,v,t0);lglrjf=0.1*trapz(lglr)lglrjf=0.1*trapz(lglr)fdxx=interp1(t,v,t0);fdxx=interp1(t,v,t0);fdxxjf=0.1*trapz(fdxx)fdxxjf=0.1*trapz(fdxx)scyt=interp1(t,v,t0,spline);scyt=interp1(t,v,t0,spline);sancytjf=0.1*tra

30、pz(scyt)sancytjf=0.1*trapz(scyt)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)gtext(lglr)gtext(lglr)gtext(fdxx)gtext(fdxx)gtext(scyt)gtext(scyt)15附錄附錄 6：t=10.95,11.49,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.40,20t=10.95,11.49,12.49,13.42,14.43,15.44,16.

31、37,17.38,18.49,19.50,20.40,20.84,23.88,24.43,25.45,25.91;.84,23.88,24.43,25.45,25.91;v=33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.86,22v=33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.86,22.17,27.09,21.62,18.48,13.30;.17,27.09,21.62,18.48,13.30;t0=20.84:0.1:22.96;

32、t0=20.84:0.1:22.96;lglr=lglrcz(t,v,t0);lglr=lglrcz(t,v,t0);lglrjf=0.1*trapz(lglr)lglrjf=0.1*trapz(lglr)fdxx=interp1(t,v,t0);fdxx=interp1(t,v,t0);fdxxjf=0.1*trapz(fdxx)fdxxjf=0.1*trapz(fdxx)scyt=interp1(t,v,t0,spline);scyt=interp1(t,v,t0,spline);sancytjf=0.1*trapz(scyt)sancytjf=0.1*trapz(scyt)plot(t,

33、v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)plot(t,v,*,t0,lglr,r,t0,fdxx,g,t0,scyt,b)gtext(lglr)gtext(lglr)gtext(fdxx)gtext(scyt)附錄附錄 7：t=0,0.46,1.38,2.395,3.41,4.425,5.44,6.45,7.465,8.45,8.97,10.95,11.49t=0,0.46,1.38,2.395,3.41,4.425,5.44,6.45,7.465,8.45,8.97,10.95,11.49,12.49,13.42,14.43,15.44,16.37,17.38,1

34、8.49,19.50,20.40,20.84,23.88,24,12.49,13.42,14.43,15.44,16.37,17.38,18.49,19.50,20.40,20.84,23.88,24.43,25.45,25.91;.43,25.45,25.91;v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.29,v=29.89,21.74,18.48,16.22,16.30,15.32,13.04,15.45,13.98,16.35,19.29,33.50,29.63,31.52,29.03,26.36,26

35、.09,24.73,23.64,23.42,25.00,23.86,22.33.50,29.63,31.52,29.03,26.36,26.09,24.73,23.64,23.42,25.00,23.86,22.17,27.09,21.62,18.48,13.30;17,27.09,21.62,18.48,13.30;t0=0:0.1:25.91;t0=0:0.1:25.91;fdxx=interp1(t,v,t0);fdxx=interp1(t,v,t0);fdxxjf=trapz(fdxx)fdxxjf=trapz(fdxx)scyt=interp1(t,v,t0,scyt=interp1

36、(t,v,t0,splinespline););sancytjf=trapz(scyt)sancytjf=trapz(scyt)plot(t,v,plot(t,v,*,t0,fdxx,t0,fdxx,gg,t0,scyt,t0,scyt,bb) )gtext(gtext(fdxxfdxx) )gtext(gtext(scytscyt) )附錄附錄 8：t=0,0.92,1.84,2.95,3.87,4.98,5.90,7.00,7.93,8.97,10.95,12.03,12.95,13.t=0,0.92,1.84,2.95,3.87,4.98,5.90,7.00,7.93,8.97,10.9

37、5,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.6688,14.98,15.90,16.83,17.93,19.04,19.96,20.84,23.88,24.99,25.66h=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,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.65,9.41,9.18

38、,8.92,8.66,8.43,8.22,10.59,10.35,10.18;9.94,9.65,9.41,9.18,8.92,8.66,8.43,8.22,10.59,10.35,10.18;c1=polyfit(t(1:10),h(1:10),3);c1=polyfit(t(1:10),h(1:10),3); % %輸入項輸入項 3 3 表示擬合多項式的次數(shù)表示擬合多項式的次數(shù)tp1=0:0.1:8.9;tp1=0:0.1:8.9;x1=polyval(c1,tp1);x1=polyval(c1,tp1);plot(tp1,x1)plot(tp1,x1)16附錄附錄 9：t=0,0.92,

39、1.84,2.95,3.87,4.98,5.90,7.00,7.93,8.97,10.95,12.03,12.95,t=0,0.92,1.84,2.95,3.87,4.98,5.90,7.00,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.6613.88,14.98,15.90,16.83,17.93,19.04,19.96,20.84,23.88,24.99,25.66h=9.68,9.48,9.31,9.13,8.98,8.81,8.69,8.52,8.3

40、9,8.22,10.82,10.50,10.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.65,9.41,9.18,8.92,8.66,8.43,8.22,10.59,10.35,10.18;21,9.94,9.65,9.41,9.18,8.92,8.66,8.43,8.22,10.59,10.35,10.18;c1=polyfit(t(1:10),h(1:10),5);c1=polyfit(t(1:10),h(1:10),5); % %輸入項輸入項 5 5 表示擬合多項式的次數(shù)表示擬合

41、多項式的次數(shù)tp1=0:0.1:8.9;tp1=0:0.1:8.9;x1=polyval(c1,tp1);x1=polyval(c1,tp1);plot(tp1,x1)plot(tp1,x1)附錄附錄 10：c1=polyfit(t(1:10),h(1:10),3);c1=polyfit(t(1:10),h(1:10),3); % %用三次多項式擬合用三次多項式擬合c2=polyfit(t(11:21),h(11:21),3);c2=polyfit(t(11:21),h(11:21),3);a1=polyder(c1);a1=polyder(c1);a2=polyder(c2);a2=poly

42、der(c2);tp1=0:0.01:8.97;tp1=0:0.01:8.97;tp2=10.95:0.01:20.84;tp2=10.95:0.01:20.84;x13=-polyval(a1,tp1);x13=-polyval(a1,tp1);x113=-polyval(a1,0:0.01:8.97);x113=-polyval(a1,0:0.01:8.97);wgsysl1=100*trapz(tp1,x113);wgsysl1=100*trapz(tp1,x113); % %計算第計算第 1 1 未供水時段的總用水量未供水時段的總用水量x14=-polyval(a1,7.93,8.97

43、);x14=-polyval(a1,7.93,8.97); % %為后面的程序準備數(shù)據(jù)為后面的程序準備數(shù)據(jù)x23=-polyval(a2,tp2);x23=-polyval(a2,tp2);x114=-polyval(a2,10.95:0.01:20.84);x114=-polyval(a2,10.95:0.01:20.84);wgsysl2=100*trapz(tp2,x114);wgsysl2=100*trapz(tp2,x114); % %計算第計算第 2 2 未供水時段的總用水量未供水時段的總用水量x24=-polyval(a2,10.95,12.03);x24=-polyval(a2

44、,10.95,12.03); % %為后面的程序準備數(shù)據(jù)為后面的程序準備數(shù)據(jù)x25=-polyval(a2,19.96,20.84);x25=-polyval(a2,19.96,20.84); % %為后面的程序準備數(shù)據(jù)為后面的程序準備數(shù)據(jù)subplot(1,2,1)subplot(1,2,1)plot(tp1,x13*100)plot(tp1,x13*100)subplot(1,2,2)subplot(1,2,2)plot(tp2,x23*100)plot(tp2,x23*100)17附錄附錄 11：c1=polyfit(t(1:10),h(1:10),3);c1=polyfit(t(1:1

45、0),h(1:10),3); % %用五次多項式擬合用五次多項式擬合c2=polyfit(t(11:21),h(11:21),3);c2=polyfit(t(11:21),h(11:21),3);a1=polyder(c1);a1=polyder(c1);a2=polyder(c2);a2=polyder(c2);tp1=0:0.01:8.97;tp1=0:0.01:8.97;tp2=10.95:0.01:20.84;tp2=10.95:0.01:20.84;x13=-polyval(a1,tp1);x13=-polyval(a1,tp1);x113=-polyval(a1,0:0.01:8.97);x113=-polyval(a1,0:0.01:8.97);wgsysl1=100*trapz(tp1,x113);wgsysl1=100*trapz(tp1,x113); % %計算第計算第 1 1 未供水時段的總用水量未供水時段的總用水量x14=