信昆侖，盛希夫，陶 濤，項(xiàng)寧銀

（同濟(jì)大學(xué) 環(huán)境科學(xué)與工程學(xué)院，上海 200092）

近幾年，城市供水管網(wǎng)的突發(fā)污染事故屢有發(fā)生，除水源地污染引起管網(wǎng)飲用水問(wèn)題之外，如管道錯(cuò)接等供水管網(wǎng)內(nèi)部的原因也會(huì)造成水質(zhì)的局部乃至全局惡化.因此，如何快速準(zhǔn)確地定位污染源是飲用水安全保障急需解決的課題.關(guān)于污染源的定位研究，Laird等［1］提出了用非線性規(guī)劃方法來(lái)求解這一問(wèn)題.首先用源頭追蹤算法替代EPANET水質(zhì)模擬計(jì)算引擎從而產(chǎn)生必要的數(shù)據(jù)，然后通過(guò)帶正則項(xiàng)的非線性規(guī)劃方法求出污染源侵入的時(shí)間和地點(diǎn).針對(duì)這一方法解的非唯一性問(wèn)題，Laird等［2］又提出了先用混合整數(shù)二次規(guī)劃的方法確定一個(gè)候選污染源節(jié)點(diǎn)集，然后在此污染源節(jié)點(diǎn)集的基礎(chǔ)上利用前一方法篩選出最有可能是污染源的節(jié)點(diǎn).但這個(gè)方法對(duì)于大中型管網(wǎng)的適用性和此模型與實(shí)際管網(wǎng)相符程度仍然有待研究.Guan等［3］利用EPANET作為內(nèi)嵌的水質(zhì)模擬器再結(jié)合負(fù)梯度法提出了模擬優(yōu)化算法.這一方法在已知管網(wǎng)模型的運(yùn)行模式和水質(zhì)監(jiān)測(cè)數(shù)據(jù)充足可靠的前提下能夠較好地解出污染源侵入的時(shí)間和地點(diǎn).但只考慮了水質(zhì)反應(yīng)為線性的情況，也沒(méi)有給出確定候選污染源節(jié)點(diǎn)的方法.Pries等［4］針對(duì)此問(wèn)題提出了關(guān)系樹-線性規(guī)劃算法.通過(guò)構(gòu)造上下游節(jié)點(diǎn)的污染物濃度關(guān)系樹并形成線性規(guī)則，再利用關(guān)系樹的線性規(guī)則反向求解線性規(guī)劃問(wèn)題，得到污染源節(jié)點(diǎn)位置及污染注入屬性.Cristo等［5］提出了離散最優(yōu)化的方法來(lái)求解污染源的識(shí)別問(wèn)題.主要通過(guò)構(gòu)造水質(zhì)污染矩陣來(lái)確定可能的污染源注入點(diǎn)，再通過(guò)最小化模擬值和監(jiān)測(cè)值的差來(lái)從可能的污染源點(diǎn)中確定污染源的位置.另外，Liu等［6］同樣利用EPANET作為內(nèi)嵌的水質(zhì)模擬計(jì)算引擎用遺傳算法隨機(jī)搜索最優(yōu)解.由于以上方法均需要管網(wǎng)的水質(zhì)監(jiān)測(cè)點(diǎn)可以報(bào)告實(shí)時(shí)的污染物濃度，對(duì)于水質(zhì)監(jiān)測(cè)點(diǎn)僅能判斷是否出現(xiàn)污染物的情況，Yang等［7］提出用貝葉斯方法計(jì)算已檢出污染物的水質(zhì)監(jiān)測(cè)節(jié)點(diǎn)其上游節(jié)點(diǎn)為污染源的概率，然后根據(jù)概率的大小判斷出可能的污染源.

本質(zhì)上來(lái)說(shuō)，供水管網(wǎng)污染源的識(shí)別問(wèn)題是管網(wǎng)水質(zhì)模擬的反演問(wèn)題，從數(shù)學(xué)的角度來(lái)看也是一個(gè)優(yōu)化問(wèn)題.本文基于Guan等［3］利用負(fù)梯度法求解供水管網(wǎng)污染源識(shí)別問(wèn)題，結(jié)合Cristo等［5］提出的候選污染源節(jié)點(diǎn)集合方法，進(jìn)一步考慮供水管網(wǎng)水質(zhì)反應(yīng)的復(fù)雜性，用差分步長(zhǎng)法代替權(quán)重系數(shù)法求解梯度方向，采用EPANET作為內(nèi)嵌的水質(zhì)模擬計(jì)算引擎結(jié)合負(fù)梯度法來(lái)求解污染源侵入的時(shí)間、地點(diǎn)和侵入過(guò)程，并結(jié)合實(shí)例管網(wǎng)進(jìn)行算法求解性能的分析.

1 模型構(gòu)建

在城市供水管網(wǎng)擁有相當(dāng)精度的水力水質(zhì)模型且管網(wǎng)上安裝了若各個(gè)水質(zhì)監(jiān)測(cè)點(diǎn)的情況下，污染源的識(shí)別問(wèn)題可以轉(zhuǎn)化為一個(gè)數(shù)學(xué)上的最優(yōu)化問(wèn)題.本文以模擬時(shí)段內(nèi)各個(gè)水質(zhì)監(jiān)測(cè)點(diǎn)的實(shí)際污染物濃度值與模擬值差的平方和為目標(biāo)函數(shù)，通過(guò)分析水質(zhì)監(jiān)測(cè)點(diǎn)不同時(shí)刻模擬濃度值與實(shí)測(cè)濃度值差來(lái)修正計(jì)算候選污染源在特定時(shí)刻的節(jié)點(diǎn)濃度.待模擬時(shí)段內(nèi)所有的候選污染源節(jié)點(diǎn)濃度更新后，再次進(jìn)行水質(zhì)模擬.當(dāng)模擬結(jié)果達(dá)到預(yù)設(shè)的條件后，輸出候選污染源節(jié)點(diǎn)濃度矩陣.由于求解所得的候選污染源污染物累計(jì)注入量與其作為真實(shí)污染源的可能性成正比，因此可以此作為確定污染源節(jié)點(diǎn)的依據(jù).

首先定義

式中：xi（t）表示候選污染源節(jié)點(diǎn)i在t時(shí)刻的污染物濃度值，對(duì)于k個(gè)時(shí)間段，則xi是一個(gè)行向量，即xi＝［xi（t1），xi（t2），…，xi（tk）］；yj（t）表示水質(zhì)監(jiān)測(cè)點(diǎn)j在t時(shí)刻的模擬污染物濃度值；n表示候選污染源節(jié)點(diǎn)總數(shù)；m表示布置的水質(zhì)監(jiān)測(cè)點(diǎn)總數(shù)；t表示模擬時(shí)刻.

監(jiān)測(cè)點(diǎn)模擬濃度可以表示為

結(jié)合監(jiān)測(cè)點(diǎn)的實(shí)際監(jiān)測(cè)結(jié)果，污染源識(shí)別問(wèn)題的目標(biāo)函數(shù)可以表示為［3］

式中：X＊為污染源識(shí)別問(wèn)題的最優(yōu)解，y′i（t）為水質(zhì)監(jiān)測(cè)點(diǎn)i在t時(shí)刻的實(shí)際濃度值.為了便于求解，這里需作幾點(diǎn)假設(shè)：①污染物的流量相對(duì)于管網(wǎng)水量而言可忽略不計(jì)，不影響管網(wǎng)原來(lái)的水力工況；②污染物在管網(wǎng)節(jié)點(diǎn)處投加，投加點(diǎn)可以是多個(gè)，投加量可以隨時(shí)間變化，但在每一水質(zhì)時(shí)間步長(zhǎng)之內(nèi)是恒定不變的；③污染物侵入管網(wǎng)中每個(gè)節(jié)點(diǎn)的概率是相同的，也就是說(shuō)污染物是隨機(jī)侵入的；④管網(wǎng)中的水質(zhì)監(jiān)測(cè)點(diǎn)能夠監(jiān)測(cè)到任意濃度的污染物.

2 模型求解

由于式（3）表示的最優(yōu)化問(wèn)題目標(biāo)函數(shù)和約束條件均為非線性函數(shù)，且約束條件是由供水管網(wǎng)恒定流方程、連續(xù)性方程、管道水質(zhì)反應(yīng)動(dòng)力學(xué)方程共同決定的復(fù)雜非線性等式約束，因此該問(wèn)題為高度復(fù)雜的非線性規(guī)劃問(wèn)題.考慮到約束條件可采用成熟的EPANET水力水質(zhì)分析模塊進(jìn)行隱式求解，故基于負(fù)梯度法，通過(guò)差分法求目標(biāo)函數(shù)梯度向量對(duì)此非線性規(guī)劃問(wèn)題進(jìn)行求解.求解過(guò)程如下：

步驟1 給定初始近似點(diǎn)X（0）及精度ε＞0，若‖▽f（X（0））‖2≤ε，則X（0）即為近似極小點(diǎn).需要說(shuō)明的是梯度向量的每一個(gè)元素均為目標(biāo)函數(shù)對(duì)X＝［x1，x2，…，xn］相應(yīng)元素的偏導(dǎo)數(shù).該偏導(dǎo)數(shù)通過(guò)差分法近似求解

式中：ΔX是第i個(gè)元素為δ且其余元素為0的列向量；δ為差分步長(zhǎng)，其取值大小應(yīng)該合理確定；f即目標(biāo)函數(shù)（見式（3））.當(dāng)水質(zhì)反應(yīng)為線性關(guān)系時(shí)，水質(zhì)監(jiān)測(cè)點(diǎn)所測(cè)得的污染物濃度值可以表達(dá)成污染物侵入節(jié)點(diǎn)的污染物濃度加權(quán)值［3］，即

偏導(dǎo)數(shù)的求法可以簡(jiǎn)化為

式中：aj，i為水質(zhì)監(jiān)測(cè)點(diǎn)j對(duì)候選污染源節(jié)點(diǎn)i的權(quán)重系數(shù)，即當(dāng)候選污染源節(jié)點(diǎn)i有一個(gè)單位的增量時(shí)，在水質(zhì)監(jiān)測(cè)點(diǎn)j有相應(yīng)的污染物濃度增量.權(quán)重系數(shù)可以通過(guò)EPANET水質(zhì)分析中的源頭追蹤功能得到.

步驟2 若‖▽f（X（0））‖2＞ε，求步長(zhǎng)λ0，并計(jì)算

求最佳步長(zhǎng)本文采用黃金分割法，最佳步長(zhǎng)即為使得目標(biāo)函數(shù)最小的步長(zhǎng).若‖▽f（X（1））‖2≤ε，則輸出X（1），否則轉(zhuǎn)入步驟3.

步驟3 一般地說(shuō)，若‖▽f（X（k））‖2≤ε，則X（k）即為所求的近似解；若‖▽f（X（k））‖2＞ε，則求步長(zhǎng)λk，并確定下一個(gè)近似點(diǎn)

如此迭代下去，直至達(dá)到所設(shè)定精度要求為止.

3 算例研究

以供水管網(wǎng)實(shí)例Net3為例進(jìn)行污染源的識(shí)別問(wèn)題求解.Net3管網(wǎng)模型含節(jié)點(diǎn)數(shù)92個(gè)，水庫(kù)2座，網(wǎng)中水塔3個(gè)，管段數(shù)117個(gè)，水泵2臺(tái)，如圖1所示.圖中虛線范圍內(nèi)的節(jié)點(diǎn)為候選污染源節(jié)點(diǎn).

圖1 管網(wǎng)拓?fù)浣Y(jié)構(gòu)圖Fig.1 Topology layout of the case study network

假定在節(jié)點(diǎn)101和119同時(shí)注入污染物，注入方式為設(shè)置點(diǎn)注入［8］，注入時(shí)間段為3∶00am～10∶00am，污染物注入質(zhì)量濃度變化模式見圖2.管網(wǎng)水質(zhì)模擬時(shí)間總長(zhǎng)24h，水力時(shí)間步長(zhǎng)1h，水質(zhì)時(shí)間步長(zhǎng)5min.在管網(wǎng)中隨機(jī)選取6個(gè)節(jié)點(diǎn)作為水質(zhì)監(jiān)測(cè)點(diǎn)（分別為連接節(jié)點(diǎn)121，151，111，267，205，213），如圖1所示.通過(guò)水質(zhì)模擬，以這些節(jié)點(diǎn)模擬污染物質(zhì)量濃度數(shù)據(jù)作為監(jiān)測(cè)結(jié)果數(shù)據(jù)（見表1）.

表1 監(jiān)測(cè)點(diǎn)模擬污染物質(zhì)量濃度Tab.1 Simulated concentrations at monitoring sites mg·L-1

3.1 確定候選節(jié)點(diǎn)集

候選節(jié)點(diǎn)集是指根據(jù)水質(zhì)監(jiān)測(cè)點(diǎn)對(duì)污染物的監(jiān)測(cè)情況所確定的可能是管網(wǎng)污染事故污染源發(fā)生位置的節(jié)點(diǎn)集合.目前關(guān)于管網(wǎng)突發(fā)污染事故污染源追蹤的研究，仍然以通過(guò)經(jīng)驗(yàn)法確定候選節(jié)點(diǎn)集合為主［3］，候選節(jié)點(diǎn)集合確定的范圍和準(zhǔn)確程度都會(huì)影響到污染源識(shí)別問(wèn)題最優(yōu)化算法的性能.為此，本研究根據(jù)文獻(xiàn)［5，9］提出的最小覆蓋集方法進(jìn)行候選節(jié)點(diǎn)集合的確定.首先求出節(jié)點(diǎn)的監(jiān)測(cè)時(shí)間矩陣，以24h為它們的服務(wù)水平，可以確定出每個(gè)監(jiān)測(cè)點(diǎn)所覆蓋的范圍.通過(guò)表1可知，本次污染事故假定模式下，實(shí)際上僅水質(zhì)監(jiān)測(cè)點(diǎn)151，111，267，205，213能監(jiān)測(cè)到污染物.因此可以把5個(gè)監(jiān)測(cè)點(diǎn)共同覆蓋的節(jié)點(diǎn)作為候選節(jié)點(diǎn)集（如圖1所示），即

3.2 結(jié)果分析

考慮到按照負(fù)梯度方向搜索得到的是局部最優(yōu)解，且絕大多數(shù)候選污染源節(jié)點(diǎn)均不是污染源注入點(diǎn)，可將初始解向量的每個(gè)元素均設(shè)置為0.求解過(guò)程中，差分步長(zhǎng)δ取0.0100，最高迭代次數(shù)取360.圖2中顯示采用差分法計(jì)算梯度向量求解得到的結(jié)果，由于求解的多數(shù)候選點(diǎn)的質(zhì)量濃度值均很小，圖中僅顯示30個(gè)候選污染源節(jié)點(diǎn)中所求解的注入質(zhì)量濃度較大的4個(gè)節(jié)點(diǎn)的注入質(zhì)量濃度隨時(shí)間變化情況.候選污染源節(jié)點(diǎn)101是各時(shí)刻累計(jì)投加質(zhì)量濃度值最大的節(jié)點(diǎn)，初步選定為污染源節(jié)點(diǎn).同時(shí)，候選污染源節(jié)點(diǎn)10從3∶00am到10∶00am基本上以100mg·L-1的質(zhì)量濃度持續(xù)注入.圖1表明，節(jié)點(diǎn)10和節(jié)點(diǎn)101是同一個(gè)管段的起止節(jié)點(diǎn).該現(xiàn)象表明有一個(gè)污染源在以節(jié)點(diǎn)10和節(jié)點(diǎn)101為起止節(jié)點(diǎn)的管道中.事實(shí)上，候選污染源節(jié)點(diǎn)101的注入質(zhì)量濃度變化曲線基本與實(shí)際污染源節(jié)點(diǎn)101的質(zhì)量濃度變化曲線相吻合.此外，候選污染源節(jié)點(diǎn)119的投加質(zhì)量濃度呈拋物線型變化也與實(shí)際的污染源119投加質(zhì)量濃度變化相吻合，可以確定這里也是一個(gè)污染源注入點(diǎn).另外，候選污染源節(jié)點(diǎn)151作為污染源節(jié)點(diǎn)119的鄰近節(jié)點(diǎn)，也有一定質(zhì)量濃度的污染物注入，這就更表明了另一個(gè)污染源就是節(jié)點(diǎn)119或者在其附近.

圖2 污染物注入質(zhì)量濃度變化曲線（差分法）Fig.2 Injection concentrations of contamination source（finite difference method）

3.3 影響求解結(jié)果的因素分析

（1）差分步長(zhǎng)

對(duì)于圖1所顯示的小型算例管網(wǎng)，將迭代次數(shù)均設(shè)定為360次，差分步長(zhǎng)δ依次取值為0.1000，0.0100，0.0010和0.0001時(shí)，計(jì)算結(jié)果的梯度向量模長(zhǎng)值依次為16.52，0.24，0.15和0.33.為了定量地衡量計(jì)算結(jié)果與實(shí)際結(jié)果的擬合程度，計(jì)算污染源節(jié)點(diǎn)101和119各個(gè)時(shí)段假定注入濃度和解出注入濃度差的平方和總和，依次為：1564.09，1914.74，7703.17，64201.58，如表2所示.

由表2可知，差分步長(zhǎng)取0.0100時(shí)對(duì)應(yīng)的計(jì)算結(jié)果（梯度向量模長(zhǎng)為0.24，目標(biāo)函數(shù)值為1914.74，而梯度向量模長(zhǎng)和目標(biāo)函數(shù)值是反映梯度法求解質(zhì)量的兩個(gè)主要指標(biāo)）最好，此時(shí)差分步長(zhǎng)取值正好與EPANET水質(zhì)選項(xiàng)中設(shè)定的水質(zhì)公差值相同.其原因在于水質(zhì)公差表示管網(wǎng)水質(zhì)的最小變化值，也即一部分水體與另一部分水體水質(zhì)的最小差值［9］，這個(gè)值也表征了管網(wǎng)中水質(zhì)監(jiān)測(cè)點(diǎn)的污染物最低監(jiān)測(cè)限度，因此當(dāng)差分步長(zhǎng)取水質(zhì)公差值時(shí)，計(jì)算結(jié)果最好.若差分步長(zhǎng)繼續(xù)減小，由于梯度計(jì)算將逐步超出管網(wǎng)模型的計(jì)算精度范圍，計(jì)算結(jié)果反而變差.

表2 差分步長(zhǎng)求解結(jié)果對(duì)比Tab.2 Comparison of results with differentδs

（2）梯度計(jì)算方法

針對(duì)本文研究的算例管網(wǎng)，分別分析了通過(guò)差分法（見式4）以及線性比例法（見式6）進(jìn)行梯度計(jì)算的求解效果.線性比例法計(jì)算后得到的結(jié)果如圖3所示，其中迭代次數(shù)為400次，收斂時(shí)間為1min，梯度向量的模長(zhǎng)為3230.83；通過(guò)差分法編碼計(jì)算后得到的結(jié)果如圖2所示，其中迭代次數(shù)為400次，收斂時(shí)間為5min，梯度向量的模長(zhǎng)為0.25.

圖3 污染物注入質(zhì)量濃度變化曲線（線性比例法）Fig.3 Injection concentrations of contamination sources with linearly proportional method

圖3表明在候選污染源節(jié)點(diǎn)115，117，119和120處有污染物注入，所定位出的污染源位置在假定污染源注入點(diǎn)附近區(qū)域，但對(duì)于污染源節(jié)點(diǎn)污染物濃度的變化識(shí)別精度較差；而通過(guò)差分法計(jì)算（見圖2）則能夠直接確定污染源注入的時(shí)間、地點(diǎn)和侵入過(guò)程.此外，迭代次數(shù)相同時(shí)，采用線性比例法計(jì)算結(jié)果的梯度向量模長(zhǎng)遠(yuǎn)大于采用差分法計(jì)算的梯度向量模長(zhǎng).因而采用線性比例法計(jì)算的結(jié)果要比用差分步長(zhǎng)編碼計(jì)算的結(jié)果差，但用線性方式編碼計(jì)算的收斂速度快于用差分法計(jì)算的速度.在采用線性比例法進(jìn)行計(jì)算時(shí)，還需要為迭代變量設(shè)置一個(gè)可取值區(qū)間，以保證解的非負(fù)性和所能取的最大值，否則將會(huì)出現(xiàn)無(wú)法收斂的情況.而用差分法計(jì)算時(shí)，僅需要保證解的非負(fù)性，無(wú)需設(shè)置解的最大限值.分析其原因在于，線性比例法計(jì)算時(shí)，由于水質(zhì)監(jiān)測(cè)點(diǎn)對(duì)于某一候選污染源節(jié)點(diǎn)的權(quán)重系數(shù)，是取相應(yīng)時(shí)間步長(zhǎng)時(shí)段內(nèi)的權(quán)重系數(shù)平均值，因此所求得的梯度向量并不能精確地給出目標(biāo)函數(shù)的收斂方向，導(dǎo)致計(jì)算結(jié)果較差.因此，建議在需要快速確定污染源點(diǎn)位時(shí)可采用線性比例法，當(dāng)需要更精確地分析污染物注入時(shí)間、濃度及變化模式時(shí)，采用基于管網(wǎng)水質(zhì)模擬的差分法進(jìn)行計(jì)算.

4 結(jié)論

負(fù)梯度法可有效解決供水管網(wǎng)單點(diǎn)污染物注入模式和多點(diǎn)污染物注入模式下的污染源定位問(wèn)題.相對(duì)于遺傳算法等隨機(jī)搜索算法，沿負(fù)梯度方向搜索具有更高的搜索效率，通過(guò)引入差分算法可有效降低目標(biāo)函數(shù)對(duì)于決策變量求導(dǎo)的難度.而通過(guò)將利用最小覆蓋集確定候選污染源節(jié)點(diǎn)集合的方法與梯度法相結(jié)合，可進(jìn)一步提高對(duì)污染源識(shí)別問(wèn)題的求解能力.此外，通過(guò)不同差分步長(zhǎng)下的計(jì)算結(jié)果表明，當(dāng)其值取水質(zhì)模擬所采用的水質(zhì)公差值時(shí)，算法求解結(jié)果最好.而相對(duì)于用線性比例法求解梯度向量，差分法可獲得更高的求解質(zhì)量，不但可以定位污染源，并能以較高精度識(shí)別污染物注入的時(shí)間和濃度變化模式.

由于污染源識(shí)別問(wèn)題是利用有限個(gè)水質(zhì)監(jiān)測(cè)點(diǎn)獲得的數(shù)據(jù)反演污染源的侵入過(guò)程，因此，水質(zhì)監(jiān)測(cè)點(diǎn)的數(shù)量及空間布置對(duì)于問(wèn)題的求解有直接的影響.因此，在將來(lái)的研究中，本方法應(yīng)進(jìn)一步與水質(zhì)監(jiān)測(cè)點(diǎn)的優(yōu)化布置結(jié)合，以更進(jìn)一步提高算法的求解效率.此外，對(duì)于無(wú)法有效獲取管網(wǎng)中水質(zhì)的精確濃度數(shù)據(jù)的情況，如何對(duì)突發(fā)污染事件進(jìn)行污染源的辨識(shí)，是未來(lái)需要進(jìn)一步研究的問(wèn)題.

［1］Laird C D，Biegler L T，van Bloemen Waanders B G，et al.Contamination source determination for water networks［J］.Water Resources Planning and Management，2005，131（2）：125.

［2］Laird C D，Biegler L T，van Bloemen Waanders B G.Mixedinteger approach for obtaining unique solutions in source inversion of water networks［J］.Water Resources Planning and Management，2006，132（4）：242.

［3］GUAN Jiabao，Aral M M，Maslia M L，et al.Identification of contaminant sources in water distribution systems using simulation optimization method：case study ［J］. Water Resources Planning and Management，2006，132（4）：252.

［4］Pries A，Ostfeld A.Contamination source identification in water systems：a hybrid model trees-linear programming scheme［J］.Water Resources Planning and Management，2006，132（4）：263.

［5］Cristo D C，Leopardi A.Pollution source identification of accidental contamination in water distribution networks［J］.Water Resources Planning and Management，2008，134（2）：197.

［6］LIU Li，Ranji Ranjithan S，Mahinthakumar G.Contamination source identification in water distribution systems using an adaptive dynamic optimization procedure［J］.Water Resources Planning and Management，2011，137（2）：183.

［7］YANG Xueyao，Boccelli D L，De Sanctis A E.A Bayesian approach for probabilistic contaminantion source identification［C］／／World Environmental and Water Resource Congress.［S.l.］：ASCE，2011.

［8］Rossman L A.EPANET2 users’manual［M］.Cincinnati：U.S.Environmental Protection Agency，2000.

［9］王鳳仙.供水管網(wǎng)水質(zhì)管理軟件應(yīng)用和預(yù)警技術(shù)研究［D］.上海：同濟(jì)大學(xué)環(huán)境科學(xué)與工程學(xué)院，2010.WANG Fengxian.Study on application of water quality management software and early-warning technology for water distribution systems［D］.Shanghai：College of Enviromental Science and Engineering of Tongji University，2010.