任剛紅 杜坤 和麗蓉 徐冰峰 杜雨

摘要：供水管網(wǎng)阻力系數(shù)識別是指通過調(diào)整管網(wǎng)水力模型中管道阻力系數(shù)使模型計(jì)算值與監(jiān)測值相符的過程。由于實(shí)際中監(jiān)測點(diǎn)數(shù)量有限，管網(wǎng)阻力系數(shù)識別為欠定的優(yōu)化問題。現(xiàn)行方法通常采用管道分組這一參數(shù)化方法將欠定問題轉(zhuǎn)換為超定，應(yīng)用遺傳算法或其它隨機(jī)搜索算法求解。提出了基于先驗(yàn)信息的供水管網(wǎng)阻力系數(shù)識別算法，所提出算法根據(jù)管道管材、管齡等先驗(yàn)信息對管道阻力系數(shù)進(jìn)行估計(jì)，并將估計(jì)值作為偽觀測值引入目標(biāo)函數(shù)將欠定優(yōu)化問題轉(zhuǎn)換為超定，采用高斯-牛頓算法進(jìn)行求解。與現(xiàn)有方法相比，所提出算法避免了管道分組不唯一的問題；再者，推導(dǎo)了供水管網(wǎng)阻力系數(shù)雅克比矩陣解析式用于搜索向量構(gòu)造，提高了參數(shù)識別計(jì)算效率。采用小型管網(wǎng)闡明了雅克比矩陣計(jì)算及搜索向量構(gòu)造，利用大型管網(wǎng)驗(yàn)證了算法的實(shí)用性。

關(guān)鍵詞：供水管網(wǎng)；阻力系數(shù)識別；先驗(yàn)信息；雅克比矩陣解析式

中圖分類號：TU99.3 文獻(xiàn)標(biāo)志碼：A文章編號：16744764（2018）02004607

收稿日期：20170314

基金項(xiàng)目：國家自然科學(xué)基金（51608242）；云南省應(yīng)用基礎(chǔ)研究青年項(xiàng)目（2017FD094）；云南省人才培養(yǎng)計(jì)劃項(xiàng)目（14118943）；云南省教育廳基金項(xiàng)目（2015Y077）

作者簡介：任剛紅（1992），女，主要從事市政工程研究，Email：554769994@qq.com。

杜坤（通信作者），男，博士，Email ：250977426@qq.com。

Received：20170314

Foundation item：National Natural Science Foundation of China（No.51608242）； Applied Basic Research Youth Project of Yunnan Province（No.2017FD94）； Personnel training program of Yunnan Province （No.14118943）； Education Department Fund Project of Yunnan Province（No. 2015Y077）

Author brief：Ren Ganghong（1992），main research interest：Municipal engineering， Email：554769994@qq.com.

Du Kun（corresponding author）， doctor， lecturer，Email：250977426@qq.com.Pipe resistance coefficient identification of water distribution

system based on prior information

Ren Ganghong1， Du Kun1， He Lirong1， Xu Bingfeng1， Du Yu2

（1. Faculty of Civil Engineering and Mechanics， Kunming University of Science and Technology， Kunming， 650500， China；

2. The third construction engineer company LTD， of China construction second engineer bureau， Hubei， Wuhan， 430022， China）

Abstract：Pipe resistance coefficients （PRCs） identification of water distribution systems （WDSs） is a process of adjusting the PRCs in hydraulic model of WDSs to make its predictions consisting with measurements. Because the number of monitoring sensors is limited in practice， the identification of PRCs of WDSs is an underdetermined optimization problem. Existing methods trend to use a parametric method of pipe grouping to convert the underdetermined problem to overdetermined， and then solve it using GA or other stochastic searching algorithms. This paper presents a prior information based algorithm for PRCs identification of WDSs. In the proposed method， the PRCs are estimated previously according to prior information of pipe material and pipeage， and then used as pseudo observations introduced into objective function to convert the underdetermined optimization problem to overdetermined one， and the Gauss Newton algorithm is utilized to solve it. Compared to existing method， the proposed algorithm avoids the nonuniqueness problem of pipe grouping； in addition， the analytic formula of Jacobian matrix of PRC is deduced for searching vector construction， which improves the calculation efficiency of parameter identification. A simple network was used to illustrate the calculation of Jacobian matrix and the construction of search vector， and a larger network was utilized to validate the practicability of the method.

Keywords：Water distribution system； Resistance coefficient identification； Prior information； Jacobian matrix formula

管網(wǎng)水力模型被越來越多的水廠用于優(yōu)化供水調(diào)度、指導(dǎo)運(yùn)營管理，如何使水力模型比較準(zhǔn)確的反映管網(wǎng)實(shí)際運(yùn)行狀態(tài)，保證決策結(jié)果的可靠性，是目前許多水廠面臨的難題。在管網(wǎng)水力模型中，相對于管道長度、管徑等參數(shù)，管道阻力系數(shù)具有較大不確定性，需要根據(jù)實(shí)測的節(jié)點(diǎn)水壓及管道流量進(jìn)行識別，以保證管網(wǎng)水力模型精度。

與傳統(tǒng)水力平差計(jì)算正好相反，管網(wǎng)阻力系數(shù)識別以監(jiān)測的節(jié)點(diǎn)水壓及管道流量作為已知量反算模型中管道阻力系數(shù)，國內(nèi)外學(xué)者通常將該反問題轉(zhuǎn)換為優(yōu)化問題進(jìn)行求解。袁一星等[1]提出了CGADFP混合優(yōu)化算法進(jìn)行管網(wǎng)阻力系數(shù)識別。王卓然[2]將SCEMUA算法與EPANET水力計(jì)算模塊相結(jié)合識別管網(wǎng)阻力系數(shù)。詹書俊等[3]通過建立以模型計(jì)算值與監(jiān)測值差的多目標(biāo)優(yōu)化問題，采用NSGA算法求解優(yōu)化問題實(shí)現(xiàn)管網(wǎng)阻力系數(shù)識別。Dini等[4]提出了基于蟻群算法的管網(wǎng)參數(shù)識別方法。劉永鑫等[5]利用遺傳算法求解管網(wǎng)連續(xù)性及能量方程識別管網(wǎng)阻力系數(shù)。信昆侖等[6]運(yùn)用全局靈敏度法進(jìn)行管道摩阻靈敏度分析，采用NSGAII算法對靈敏度較大的管道進(jìn)行參數(shù)識別。

筆者注意到，由于管網(wǎng)中監(jiān)測點(diǎn)數(shù)量遠(yuǎn)少于管道數(shù)（即已知量個(gè)數(shù)少于未知量個(gè)數(shù)），管網(wǎng)阻力系數(shù)識別是欠定的優(yōu)化問題。針對該問題的處理，Kang等[7]根據(jù)管道的管材及管齡對管道分組，并假設(shè)同組管道阻力系數(shù)相等，將欠定問題轉(zhuǎn)換為超定進(jìn)行求解。Wu等[8]根據(jù)管道在管網(wǎng)中位置對管道分組減少未知量個(gè)數(shù)，使管道阻力系數(shù)識別結(jié)果唯一、可靠。Mallick等[9]從識別結(jié)果穩(wěn)健性的角度出發(fā)探討了管道分組問題，其研究表明，管道分組這一參數(shù)化方法能有效降低識別結(jié)果方差，但由于實(shí)際管網(wǎng)拓?fù)浣Y(jié)構(gòu)不同，及監(jiān)測點(diǎn)數(shù)量、布置差異，不存在唯一準(zhǔn)則適用于所有管網(wǎng)。

本文提出了基于先驗(yàn)信息的供水管網(wǎng)阻力系數(shù)識別算法，其基本思路是根據(jù)管道的管材、管齡等先驗(yàn)信息，對管道阻力系數(shù)進(jìn)行估計(jì)，并將估計(jì)值作為偽觀測值引入目標(biāo)函數(shù)，采用加權(quán)最小二乘法求解優(yōu)化問題識別管道阻力系數(shù)。與現(xiàn)有方法相比，所提出算法無需對管道分組，通過利用先驗(yàn)信息將欠定問題轉(zhuǎn)化為超定，克服了現(xiàn)有方法中管道分組不唯一的缺點(diǎn)；再者，通過權(quán)重系數(shù)權(quán)衡先驗(yàn)信息與測量信息避免了參數(shù)過擬合問題。此外，推導(dǎo)供水管網(wǎng)雅克比矩陣解析式用于構(gòu)造搜索向量，提高了管道阻力系數(shù)識別計(jì)算效率。

第2期 任剛紅，等：基于先驗(yàn)信息的供水管網(wǎng)阻力系數(shù)識別1基于先驗(yàn)信息的供水管網(wǎng)阻力系數(shù)

識別框架實(shí)際中管道阻力系數(shù)除了與管道的管材及管齡相關(guān)外，還與管道內(nèi)壁涂料厚度、腐蝕程度及管網(wǎng)水力狀態(tài)等隨機(jī)因素相關(guān)。所提出識別算法一方面承認(rèn)基于先驗(yàn)信息的管道阻力系數(shù)包含一定有用信息，另一方面要求參數(shù)識別結(jié)果應(yīng)盡量減少模型計(jì)算值與實(shí)測值間差異?；诖?，優(yōu)化問題的目標(biāo)函數(shù)可構(gòu)建為f（C）=∑nHi=1wHHoi-Hi（C）2+

∑mqj=1wqqoj-qj（C）2+∑mk=1wCCok-Ck2（1）式中：nH為水壓監(jiān)測點(diǎn)數(shù)；n為管網(wǎng)中節(jié)點(diǎn)數(shù)；mq為管道流量監(jiān)測點(diǎn)數(shù)，m為管道數(shù)；Co為管道阻力系數(shù)經(jīng)驗(yàn)值；wH、wq、wC分別為節(jié)點(diǎn)水壓、管道流量及管道阻力系數(shù)權(quán)重系數(shù)（分別為水壓與流量監(jiān)測值誤差方差、管道阻力系數(shù)估計(jì)值方差的倒數(shù)）；Ho為水壓監(jiān)測值、H（C）為對應(yīng)的模型計(jì)算值；qo為管道流量監(jiān)測值、q（C）為對應(yīng)的模型計(jì)算值。優(yōu)化問題的約束條件為供水管網(wǎng)質(zhì)量與能量守恒方程。為便于推導(dǎo)，將式（1）化為矩陣形式：f（C）=Ho-H（C）

qo-q（C）

Co-CTWHo-H（C）

qo-q（C）

Co-C（2）由于管網(wǎng)能量方程為非線性，故需采用迭代法求解優(yōu)化問題。本文采用高斯牛頓算法求解上述優(yōu)化問題，相對于廣泛使用的遺傳算法，其具有計(jì)算效率高的優(yōu)點(diǎn)，且不需要進(jìn)行額外的參數(shù)設(shè)置。若第k次迭代的解為：f（Ck+ΔCk）=

Ho-H（Ck+ΔCk）

qo-q（Ck+ΔCk）

Co-（Ck+ΔCk）TWHo-H（Ck+ΔCk）

qo-q（Ck+ΔCk）

Co-（Ck+ΔCk）（3）式（3）的線性展開式為：f（Ck+ΔCk）≈

ΔHk-JH（Ck）ΔCk

Δqk-Jq（Ck）ΔCk

ΔCok-ΔCkTWΔHk-JH（Ck）ΔCk

Δqk-Jq（Ck）ΔCk

ΔCok-ΔCk（4）式中：JH（C）、Jq（C）為梯度向量，為與監(jiān)測值對應(yīng)的雅克比矩陣的行向量；ΔHk=HoΔH（Ck）；Δqk=qoΔq（Ck）、ΔCok=Co-Ck。當(dāng)目標(biāo)函數(shù)取得極小解時(shí)，有：f（Ck+ΔCk）Ck=

-2JH（Ck）

Jq（Ck）

ITWΔHk-JH（Ck）ΔCk

Δqk-Jq（Ck）ΔCk

ΔCok-ΔCk=0（5）式中：JH（Ck）為nH×m矩陣，Jq（Ck）為mq×m矩陣，I為m×m的單位矩陣，[JH（Ck）Jq（Ck） I]T為（nH+mq+m）×m矩陣。因?yàn)閚H+mq+m>m，即矩陣[JH（Ck） Jq（Ck） I]T的行向量個(gè)數(shù)大于列向量個(gè)數(shù)，管道阻力系數(shù)修正值為：ΔCk=JH（Ck）

Jq（Ck）

ITWJH（Ck）

Jq（Ck）

I-1·

JH（Ck）

Jq（Ck）

IWΔHk

Δqk

ΔCok （6）圖1給出了供水管網(wǎng)阻力系數(shù)識別流程圖，其中識別模塊與正計(jì)算模塊將對方的輸出作為輸入進(jìn)行反復(fù)運(yùn)算直到ΔC達(dá)到規(guī)定精度（本文取001）。理論上ΔC可取任意小值，由于實(shí)際中C的取值范圍為90～150，則當(dāng)ΔC=0.01時(shí)，ΔCC<0.001，即最終計(jì)算誤差小于千分之一，滿足實(shí)際工程需要。

圖1供水管網(wǎng)阻力系數(shù)識別流程圖

Fig.1Water supply pipe network resistance

coefficient identification flow chart2供水管網(wǎng)雅克比矩陣計(jì)算

如圖1所示，在管網(wǎng)阻力系數(shù)識別過程中需要反復(fù)計(jì)算雅克比矩陣用于構(gòu)造搜索向量，但筆者注意到，目前大多數(shù)學(xué)者[10～13]采用有限差分法估算管網(wǎng)雅克比矩陣，其需要逐個(gè)擾動(dòng)參數(shù)反復(fù)進(jìn)行管網(wǎng)水力平差計(jì)算，故會導(dǎo)致巨大計(jì)算量，這不利于大型管網(wǎng)阻力系數(shù)識別。鑒于此，本文推導(dǎo)了供水管網(wǎng)雅克比矩陣解析式，以提高參數(shù)識別計(jì)算效率。管網(wǎng)質(zhì)量與能量守恒方程為Aq-Q=0

ATH+h=0 （7）式中：A為管網(wǎng)銜接矩陣，q為管道流量向量；Q為節(jié)點(diǎn)流量向量；H為節(jié)點(diǎn)水壓向量；h為管道水頭損失向量。式（7）的微分式為AΔq-ΔQ=0

ATΔH+Δh=0 （8）配水管網(wǎng)水頭損失多采用海澄威廉公式計(jì)算：h=KqC1.852Ld4.871（9）式中：K為單位換算系數(shù)；q、C為管道流量及海澄威廉系數(shù)；L、d為管長及管徑。管道水頭損失對管道流量的偏微分式為hq=KuqC1.8521.852Ld4.871q=h1.852q（10）根據(jù)式（10）還可得：qC=hd4.871KuL11.852=qC（11）

hC=-1.852hC（12）根據(jù)式（10）、（11）及（12），管道水頭損失對管道阻力系數(shù)的向量微分方程可寫為Δh=ΔB-1SΔC（13）其中B=q11.852h10…0

0q11.852h1…0

00…q11.852h1

S=q1C10…0

0q2C2…0

00…qmCm當(dāng)管網(wǎng)中存在水泵時(shí)，矩陣B中對應(yīng)元素為（cb）-1|q|1-c，設(shè)水泵方程為hpump=a-bqc，a、b及c為水泵性能參數(shù)。根據(jù)式（7），可得Δh=-ATΔH（14）將式（14）帶入式（13），可得BATΔH=SΔC（15）根據(jù)式（15），可得ABATΔH=ASΔC（16）根據(jù)式（16），可得ΔH=（ABAT）-1ASΔC（17）同樣地，管道流量的向量微分方程為Δq=SΔC+BΔh（18）將式（14）帶入式（18），可得Δq=SΔCΔBATΔH（19）將式（17）帶入式（19），可得Δq=SΔCΔBAT（ABAT）-1ASΔC（20）根據(jù)式（17）、（20），節(jié)點(diǎn)水壓及管道流量對管道阻力系數(shù)的雅克比矩陣的解析式為HC=（ABAT）-1AS

qC=S-BAT（ABAT）-1AS（21）3基于數(shù)值仿真的監(jiān)測值生成

供水管網(wǎng)參數(shù)識別存在補(bǔ)償誤差問題，例如：當(dāng)調(diào)整管道阻力系數(shù)或節(jié)點(diǎn)流量都能使模型計(jì)算值與監(jiān)測值相符時(shí)，則無法分辨模型誤差源于節(jié)點(diǎn)流量或管道阻力系數(shù)，即二者間存在補(bǔ)償誤差。如何獲得有用監(jiān)測值是在利用優(yōu)化算法進(jìn)行參數(shù)識別前首先應(yīng)回答的問題。針對管網(wǎng)阻力系數(shù)識別，Ostfeld等[14]表明，應(yīng)通過消火栓放水并記錄放水量以減小節(jié)點(diǎn)流量補(bǔ)償誤差、獲得有用監(jiān)測值。Ormsbee等[15]指出消火栓放水至少應(yīng)保證管網(wǎng)供水壓力下降大于3.5m，使管網(wǎng)處于高負(fù)荷水力“緊繃”狀態(tài)，以加大監(jiān)測值對管道阻力系數(shù)敏感度，否則收集的監(jiān)測值是無用的。

從算法驗(yàn)證角度來說，通過開展實(shí)地消火栓放水獲得監(jiān)測值，需要投入大量人力、財(cái)力，且影響管網(wǎng)正常運(yùn)行，顯得代價(jià)過高。再者，由于實(shí)際管網(wǎng)阻力系數(shù)未知，而監(jiān)測值又存在誤差，這導(dǎo)致參數(shù)識別結(jié)果及模型準(zhǔn)確性都失去參照，不利于算法驗(yàn)證。鑒于上述兩方面原因，為便于算法驗(yàn)證，本文參考文獻(xiàn)[11]所采用的數(shù)值仿真法產(chǎn)生監(jiān)測值。

1）管網(wǎng)水力模型構(gòu)建；EPANET是目前使用最廣泛的管網(wǎng)水力計(jì)算引擎，故在EPANET中構(gòu)建管網(wǎng)水力模型開展相關(guān)研究；

2）管道阻力系數(shù)“真值”生成；考慮到實(shí)際管道阻力系數(shù)不僅與管齡及管材相關(guān)，還與管道內(nèi)壁涂料厚度、腐蝕程度，及管網(wǎng)水力狀態(tài)等隨機(jī)因素相關(guān)。為準(zhǔn)確反映實(shí)際情況，采用隨機(jī)抽樣法生成管道阻力系數(shù)；不失一般性，假定管道阻力系數(shù)真值服從N（Co，σ2）的正態(tài)分布，其中Co為根據(jù)先驗(yàn)信息估計(jì)的管道阻力系數(shù)；

3）監(jiān)測值生成；應(yīng)用EPANET進(jìn)行管網(wǎng)水力計(jì)算，采用隨機(jī)抽樣法產(chǎn)生隨機(jī)誤差添加到計(jì)算的管道流量及節(jié)點(diǎn)水壓中作為“真實(shí)”監(jiān)測值。這里添加的隨機(jī)誤差可包括水壓、流量本身的監(jiān)測誤差，同時(shí)還可包括節(jié)點(diǎn)流量的補(bǔ)償誤差；

4）管道阻力系數(shù)識別及結(jié)果評判；根據(jù)產(chǎn)生的“真實(shí)”監(jiān)測值，應(yīng)用所提出算法識別管網(wǎng)中各管道阻力系數(shù)。在評判識別結(jié)果時(shí)，一方面可將識別結(jié)果與管道阻力系數(shù)“真值”進(jìn)行比較，另一方面可觀察模型計(jì)算精度的改善情況。

4案例分析

4.1案例1

案例1的主要目的是闡明雅克比矩陣計(jì)算及搜索向量構(gòu)造。為便于闡明，選取圖2的小型管網(wǎng)作為例子，其中各管道管長均為500 m，管徑均為200 mm。假設(shè)根據(jù)管齡、管材等先驗(yàn)信息估計(jì)的管道阻力系數(shù)為90（海澄威廉系數(shù)），真實(shí)管道阻力系數(shù)服從N（90，102）的正態(tài)分布，隨機(jī)抽樣所得管道阻力系數(shù)真值如圖2所示。此外，假定節(jié)點(diǎn)1的水壓及水泵供水量被監(jiān)測，通過添加隨機(jī)誤差產(chǎn)生監(jiān)測值，如圖2所示。

圖2舉例管網(wǎng)1

Fig.2Example pipe network1表1給出了節(jié)點(diǎn)水壓對管道阻力系數(shù)的雅克比矩陣（ABAT）-1（AS），表2給出了管道流量對管道阻力系數(shù)的雅克比矩陣SBAT（ABAT）-1（AS）。表1雅克比矩陣（ABAT）-1（AS）

Table 1Jacobian matrix（ABAT）-1（AS）管道1管道2管道3管道4管道5節(jié)點(diǎn)10.072 1 -0.001 3 0.003 2 0.003 1 -0.000 6 節(jié)點(diǎn)20.066 6 0.004 7 0.008 0 0.000 2 0.000 0 節(jié)點(diǎn)30.067 6 -0.000 3 0.000 7 0.005 4 0.001 9 節(jié)點(diǎn)40.058 4 0.001 8 -0.004 4 -0.004 2 0.000 8 表2雅克比矩陣SBAT（ABAT）-1（AS）

Table 2Jacobian matrixSBAT（ABAT）-1（AS）管道1管道2管道3管道4管道5管道10.132 30.004 0-0.009 9-0.009 50.001 9管道2-0.057 8-0.020 70.050 9-0.030 50.006 2管道3-0.057 8-0.020 70.050 9-0.030 50.006 2管道40.074 4-0.016 70.041 0-0.040 00.008 1管道50.074 4-0.016 70.041 0-0.040 00.008 1 水泵-0.132 3-0.004 00.009 90.009 5-0.001 9

根據(jù)表1及表2給出的雅克比矩陣，所提出算法第一次迭代時(shí)的搜索向量能構(gòu)建為：表3第一次迭代時(shí)搜索向量[JH（C1） Jq（C1） I]T

Table 3Search for vectors for the first

iteration[JH（C1） Jq（C1） I]TJH（C1）Jq（C1）I0.072 1-0.001 30.003 20.003 1-0.000 6-0.132 3-0.0040.009 90.009 5-0.001 91000001000001000001000001其中JH（C1）為矩陣（ABAT）-1（AS）的第一行（詳表1），Jq（C1）為矩陣SBAT（ABAT）-1（AS）的最后一行（詳表2）。如前述，權(quán)重矩陣W中元素wH、wq、wC分別為監(jiān)測誤差方差及管道阻力系數(shù)估計(jì)值方差的倒數(shù)。不失一般性，認(rèn)為監(jiān)測值誤差服從正態(tài)分布，其中水壓監(jiān)測值均方差σH=0.3 m，流量監(jiān)測值均方差σq=2 L/s，則權(quán)重矩陣W為W=0.3-2000000

02-200000

0010-20000

00010-2000

000010-200

0000010-20

00000010-2第一次迭代修正計(jì)算值為：ΔC1=JH（C1）

Jq（C1）

ITWJH（C1）

Jq（C1）

I-1·

JH（C1）

Jq（C1）

IWΔH1

Δq1

ΔCo1=5.91

-0.25

0.62

0.59

-0.12表4給出了迭代過程中ΔC值及管道阻力系數(shù)。

最終識別結(jié)果，表5給出了參數(shù)識別前后模型誤差。 表4代過程中ΔC值及最終識別結(jié)果

Table 4ΔC value in the process and final identification results初值C1C2C3C4C5||ΔC||2ΔCk=15.91 -0.25 0.62 0.59 -0.12 5.98 ΔCk=20.16 -0.03 -0.02 -0.03 -0.02 0.16 ΔCk=30.00 0.00 0.00 0.00 0.00 0.00 終值96.06 89.72 90.60 90.56 89.86 真值981058075102

表5參數(shù)識別前后模型計(jì)算值與真實(shí)值間的差

Table 5The difference between the value of the model and

the real value before and after the parameter identification添加的

監(jiān)測誤差識別前

模型誤差識別后

模型誤差節(jié)點(diǎn)1水壓0.11 m0.43 m 0.01 m節(jié)點(diǎn)2水壓0.48 m 0.09 m節(jié)點(diǎn)3水壓0.43 m 0.04 m節(jié)點(diǎn)4水壓0.59 m0.27 m平均絕對誤差0.48 m0.10 m管道1流量1.36 L/s 0.61 L/s 續(xù)表5添加的

監(jiān)測誤差識別前

模型誤差識別后

模型誤差管道2流量-0.67 L/s -0.36 L/s 管道3流量-0.67 L/s -0.36 L/s 管道4流量0.69 L/s 0.25 L/s 管道5流量0.69 L/s 0.25 L/s 水泵流量1.5 L/s-1.36 L/s -0.61 L/s平均絕對誤差0.91 L/s0.41 L/s

由表4可知，管道1的阻力系數(shù)被較準(zhǔn)確識別，其它管道阻力系數(shù)值基本不變。這是由于管道1為高位水池出水管，其管道流量遠(yuǎn)大于其它管道，處于高負(fù)荷水力狀態(tài)，導(dǎo)致監(jiān)測值對管道1阻力系數(shù)敏感度遠(yuǎn)大于其它管道。上述結(jié)論可通過表1、表2的雅克比矩陣進(jìn)行說明。雅克比矩陣又稱靈敏度矩陣，反映了監(jiān)測值對參數(shù)的靈敏程度，其中元素值越大表明對應(yīng)參數(shù)對模型計(jì)算精度影響越大且越容易被識別，反之亦然。例如：表1中的第一列代表了各節(jié)點(diǎn)水壓對管道1的阻力系數(shù)靈敏度，其中各值比其它各列的值均大了一個(gè)數(shù)量級以上，表明監(jiān)測值對管道1阻力系數(shù)靈敏度遠(yuǎn)大于其它管道，即管網(wǎng)水力模型精度主要取決于管道1的阻力系數(shù)，且管道1阻力系數(shù)更容易識別。

由表5可知，參數(shù)識別前，節(jié)點(diǎn)水壓平均絕對誤差為0.48 m，管道流量平均絕對誤差為0.91 L/s；參數(shù)識別后，節(jié)點(diǎn)水壓平均絕對誤差為0.1 m，管道流量平均絕對誤差為0.41 L/s，模型計(jì)算誤差整體上明顯減小，這表明利用所提出算法識別管網(wǎng)阻力系數(shù)能提高模型計(jì)算精度。此外，節(jié)點(diǎn)1的水壓與真實(shí)值的差異僅為0.01 m，遠(yuǎn)小于監(jiān)測誤差值0.11 m，這表明應(yīng)用所提出的加權(quán)方法能有效防止參數(shù)過擬合。再者，根據(jù)表4可知，整個(gè)參數(shù)識別過程僅需要三次迭代，表明所提出算法計(jì)算效率高。

4.2案例2

為進(jìn)一步驗(yàn)證算法可行性，本案例利用某實(shí)際大型供水管網(wǎng)測試算法。管網(wǎng)基本情況如圖3所示，其中包括43個(gè)節(jié)點(diǎn)、62根管道。假設(shè)在節(jié)點(diǎn)4、6、10、17、30、36及38上設(shè)置水壓監(jiān)測點(diǎn)，監(jiān)測誤差服從N（0， 0.32）正態(tài)分布，且有e4=0.29 m、e6=-0.41 m、e10=0.27 m、e17=0.06 m、e30=-0.18 m、e36=-0.43 m、e38=0.15 m。根據(jù)管材及管齡估計(jì)的管道阻力系數(shù)經(jīng)驗(yàn)值為100（海澄威廉系數(shù)），實(shí)際管道阻力系數(shù)服從N（100，102）的正態(tài)分布。

圖3舉例管網(wǎng)2

Fig.3Example pipe network2應(yīng)用所提出算法識別管網(wǎng)阻力系數(shù)，限于篇幅原因，不對結(jié)果進(jìn)行詳細(xì)列舉。總體而言，與案例1類似，流量較大的主供水管道4、6、11、42、47、51的阻力系數(shù)被較準(zhǔn)確識別，與靈敏度分析結(jié)果一致。經(jīng)管道阻力系數(shù)識別，模型節(jié)點(diǎn)水壓平均絕對誤差由0.76 m降低到0.11 m，最大節(jié)點(diǎn)水壓計(jì)算誤差由1.5 m降低到0.48 m，模型計(jì)算精度有較大改善，這表明所提出算法可用于實(shí)際大型管網(wǎng)參數(shù)識別。此外，整個(gè)參數(shù)識別過程僅需要三次迭代，且6、36節(jié)點(diǎn)水壓計(jì)算誤差小于監(jiān)測值隨機(jī)誤差，表明所提出算法計(jì)算效率高，能有效避免參數(shù)過擬合問題。

5結(jié)論

供水管網(wǎng)阻力系數(shù)識別是欠定的非線性優(yōu)化問題，目前大多數(shù)研究通過管道分組這一參數(shù)化方法將欠定問題轉(zhuǎn)換為超定，并采用遺傳算法或其它類似隨機(jī)搜索算法求解。本文提出了基于先驗(yàn)信息的供水管網(wǎng)阻力系數(shù)識別算法，其根據(jù)管材、管齡等先驗(yàn)信息估計(jì)管道阻力系數(shù)，將估計(jì)值引入目標(biāo)函數(shù)，采用高斯牛頓法進(jìn)行求解。與現(xiàn)有方法相比，所提出算法無需對管道分組，利用先驗(yàn)信息將欠定問題轉(zhuǎn)換為超定，避免了管道分組不唯一的問題。再者，推導(dǎo)了供水管網(wǎng)阻力系數(shù)雅克比矩陣解析式用于搜索向量構(gòu)造，提高了參數(shù)識別計(jì)算效率。最后，采用小型管網(wǎng)闡明了雅克比矩陣計(jì)算及搜索向量構(gòu)造，利用實(shí)際大型管網(wǎng)對算法進(jìn)行測試。

結(jié)果表明，所提出算法通過3次迭代就能獲得最終識別結(jié)果，計(jì)算效率高且能避免參數(shù)過擬合問題。通過分析還發(fā)現(xiàn)，管網(wǎng)水力模型計(jì)算精度主要取決于管網(wǎng)中供水主管管道阻力系數(shù)，通過識別這些管道阻力系數(shù)，保持其它管道阻力系數(shù)不變，不失為一種可行的識別方法。值得說明的是，參數(shù)識別結(jié)果及識別后模型計(jì)算精度與監(jiān)測點(diǎn)數(shù)量、布置位置及數(shù)據(jù)采集時(shí)管網(wǎng)運(yùn)行狀態(tài)密切相關(guān)，通常應(yīng)用高負(fù)荷運(yùn)行狀態(tài)下的水壓監(jiān)測數(shù)據(jù)能得到更準(zhǔn)確的管道阻力系數(shù)校核值。通過優(yōu)化監(jiān)測點(diǎn)布置，采集消火栓放水試驗(yàn)時(shí)監(jiān)測值能改善識別結(jié)果、提高模型計(jì)算精度。鑒于監(jiān)測點(diǎn)布置本身是一個(gè)復(fù)雜的優(yōu)化問題，其超出了本文的研究范圍，在此不進(jìn)行深入探討。在工程實(shí)踐中，Walski [16]建議可將監(jiān)測點(diǎn)布置在用水量較大的節(jié)點(diǎn)及管網(wǎng)外圍（遠(yuǎn)離水源）的節(jié)點(diǎn)。

參考文獻(xiàn)：

[1] 袁一星， 張志軍. 供水管網(wǎng)校核模型參數(shù)估計(jì)與求解方法的研究[J]. 給水排水， 2005， 31（9）：105111.

YUAN Y X， ZHANG Z J. Study on estimation and approach of parameters of calibration model for water distribution network model[J]. Water & Waste Water， 2005， 31（9）：105111.（in Chinese）

[2] 王卓然. 給水管網(wǎng)管道阻力系數(shù)校正試驗(yàn)研究[D]. 哈爾濱工業(yè)大學(xué)， 2012.

WANG Z R. experimental study on calibration of pipeline roughness in water distribution networks[D]. Harbin Institute of Technology， 2012.（in Chinese）

[3] 詹書俊， 陶濤. 基于NSGA—II的供水管網(wǎng)模型校核[J]. 給水排水， 2013， 39（3）：158160.

ZHAN S J， TAO R. Model checking of water supply network based on NSGA II[J]. Water & Waste Water， 2013， 39（3）：158160.（in Chinese）

[4] DINI M， TABESH M. A new method for simultaneous calibration of demand pattern and hazenwilliams coefficients in water distribution systems[J]. Water Resources Management， 2014， 28（7）：20212034.

[5] 劉永鑫， 鄒平華， 馬月璇. 基于遺傳算法的供水管網(wǎng)阻力系數(shù)辨識[J]. 中國給水排水， 2014（23）：113116.

LIU Y X， ZOU P H， MA Y X. Identification of resistance coefficient of water supply network based on genetic algorithm[J]. China Water & Wastewater， 2014（23）：113116.（in Chinese）

[6] 信昆侖， 詹書俊， 陶濤，等. 基于靈敏度分析的供水管網(wǎng)模型多目標(biāo)校核[J]. 同濟(jì)大學(xué)學(xué)報(bào)：自然科學(xué)版， 2014， 42（5）：736739.

XIN K L， ZHAN S J， TAO T. Multi—objeetive Calibration of hydraulic model in water distribution network based on sensitivity analysis[J]. Journal of Tongji University（Natural Science） ， 2014， 42（5）：736739.（in Chinese）

[7] KANG D， LANSEY K， KANG D， et al. Demand and roughness estimation in water distribution systems[J]. Journal of Water Resources Planning & Management， 2010， 137（1）：2030.

[8] WU Z Y， CLARK C. Evolving effective hydraulic model for municipal water systems.[J]. Water Resources Management， 2009， 23（1）：117136.

[9] MALLICK K N， AHMED I， TICKLE K S， et al. Determining pipe groupings for water distribution networks[J]. Journal of Water Resources Planning & Management， 2002， 128（128）：130139.

[10] LANSEY K E， ElSHORBAGY W， AHMED I，et al. Calibration assessment and data collection for water distributionnetworks[J].Journal of Hydraulic Engineering， 2001， 127（4）： 270279.

[11] KANG D， LANSEY K. Realtime demand estimation and confidence limit analysis for water distribution systems[J]. Journal of Hydraulic Engineering， 2009， 135（10）： 825837.

[12] PEREZ R， PUIG V， PASCUL J， et al. Pressure sensor distribution for leak detection in Barcelona water distribution network[J]. Water science and technology： water supply， 2009， 9（6）： 715.

[13] MNDEZ M， ARAYA J A， SNCHEZ L D. Automated parameter optimization of a water distribution system[J]. Journal of Hydroinformatics， 2013， 15（1）：7185.

[14] OSTFELD A， SALOMONS E， ORMSBEE L， et al. Battle of the water calibration networks[J]. Journal of Water Resources Planning & Management， 2012， 138（5）：523532.

[15] ORMSBEE L E， LINGIREDDY S. 1997. Calibration of hydraulic network models [J]. Journal AWWA， 89（2）： 4250.

[16] WALSKI T M. Technique for calibrating network models. 1983 [J]. Journal of Water Resources Planning & Management， 109（4）：360372.