基于自適應(yīng)MOEA/D 的供水管網(wǎng)多目標(biāo)優(yōu)化

董 娜

（內(nèi)蒙古電子信息職業(yè)技術(shù)學(xué)院電子與自動化學(xué)院，內(nèi)蒙古 呼和浩特 010070）

0 引 言

城市供水系統(tǒng)是整座城市的生命線，完善的供水管網(wǎng)能夠保障人民群眾日常生活的正常運(yùn)轉(zhuǎn)。在城市供水系統(tǒng)中，管網(wǎng)優(yōu)化設(shè)計(jì)一直是研究熱點(diǎn)問題之一［1］。在管網(wǎng)設(shè)計(jì)時(shí)，要求在滿足供水需求量和水壓等約束條件下，期望能盡量減少節(jié)點(diǎn)富余水頭總和，以降低管網(wǎng)設(shè)計(jì)的故障率，增強(qiáng)其可靠性。同時(shí)，又期望盡可能降低管網(wǎng)建造總成本，以提高管網(wǎng)建造的經(jīng)濟(jì)性。顯然，總造價(jià)和故障率是一對具有沖突特性的目標(biāo)函數(shù)。隨著總造價(jià)的提高，故障率呈下降趨勢；反之，隨著總造價(jià)的降低，故障率呈上升趨勢。因此，供水管網(wǎng)設(shè)計(jì)顯然是一個(gè)多目標(biāo)優(yōu)化問題［2］。傳統(tǒng)的單目標(biāo)優(yōu)化方法需要利用權(quán)重系數(shù)將多目標(biāo)優(yōu)化問題轉(zhuǎn)化為單目標(biāo)優(yōu)化問題。但是權(quán)重系數(shù)極難確定，而且對算法結(jié)果影響巨大。由于多目標(biāo)進(jìn)化算法一次運(yùn)行能夠獲取多個(gè)Pareto最優(yōu)解，適合于解決管網(wǎng)多目標(biāo)優(yōu)化問題。研究者將其應(yīng)用在城市管網(wǎng)優(yōu)化設(shè)計(jì)當(dāng)中，為市政工程的施工建設(shè)提供有益參考［3，4］。

基于分解的多目標(biāo)進(jìn)化算法（MOEA/D）采用分解方法將一個(gè)多目標(biāo)優(yōu)化問題轉(zhuǎn)化為一組單目標(biāo)優(yōu)化子問題，并以協(xié)作的方式處理它們。MOEA/D公認(rèn)為是多目標(biāo)進(jìn)化領(lǐng)域的主流算法之一［5］。MOEA/D 在實(shí)際工程領(lǐng)域具有廣泛的應(yīng)用，如無人機(jī)航跡規(guī)劃［6］、能源計(jì)劃優(yōu)化［7］、網(wǎng)絡(luò)接入控制［8］、樓宇負(fù)荷優(yōu)化［9］、天線優(yōu)化［10］、飛行器氣動外形優(yōu)化［11］等。MOEA/D在解決一些復(fù)雜多目標(biāo)優(yōu)化問題時(shí)，其難以發(fā)現(xiàn)目標(biāo)空間中Pareto 前沿上的邊界解，而與這些邊界解對應(yīng)的決策空間中的輸入矢量有可能是期望的Pareto最優(yōu)解［12］。

為了取得收斂性和多樣性的最佳平衡，研究者提出多種改進(jìn)方法。例如：LI 等人將差分進(jìn)化算子（DE）引入到MOEA/D中，形成MOEA/D-DE 算法［13］。ZHANG 等人設(shè)計(jì)一種基于精英群體的自適應(yīng)權(quán)重向量調(diào)整策略，并集成到MOEA/D 中形成MOEA/D-AWA 算法。MOEA/D-AWA 通過引入精英群體，以在目標(biāo)空間稀疏區(qū)域添加新的子問題，從而能夠有效解決復(fù)雜MOP問題［14］。LI等人設(shè)計(jì)一種自適應(yīng)交配選擇策略，提出一種簡單有效的穩(wěn)定匹配模型來協(xié)調(diào)選擇過程［15］。此外，如自適應(yīng)交配選擇機(jī)制、鄰域大小集合、可變鄰域策略等均已成功集成到MOEA/D中來解決實(shí)際工程問題。

然而，實(shí)際工程中大多數(shù)多目標(biāo)優(yōu)化問題具有復(fù)雜的Pareto 最優(yōu)前沿，比如長尾巴、尖峰或不連續(xù)。以供水管網(wǎng)優(yōu)化問題來說，總造價(jià)和故障率不一定呈現(xiàn)出標(biāo)準(zhǔn)的Pareto 前沿。我們總是期望，當(dāng)總造價(jià)緩慢增加時(shí)，故障率將會緩慢降低；當(dāng)總造價(jià)緩慢降低時(shí)，故障率將會緩慢增加。然而，現(xiàn)實(shí)情況卻是，在某個(gè)臨界點(diǎn)后，總造價(jià)的些許降低可能導(dǎo)致故障率的極大增加，總造價(jià)的極大增加也可能無法使故障率有些許降低。這就出現(xiàn)了所謂的長尾巴、尖峰等復(fù)雜Pareto 前沿的情況。針對這些復(fù)雜優(yōu)化問題，其中一個(gè)目標(biāo)的細(xì)微變換會導(dǎo)致另一個(gè)目標(biāo)較大的變換，從而導(dǎo)致Pareto 最優(yōu)解集中在Pareto 前沿的中部區(qū)域，邊界區(qū)域很難獲取到Pareto 最優(yōu)解。顯然，傳統(tǒng)MOEA/D難以有效完成這類問題的尋優(yōu)。

主要針對供水管網(wǎng)多目標(biāo)優(yōu)化呈現(xiàn)的Pareto解集前沿面存在尖峰、長尾和不連續(xù)等問題，提出一種帶有兩階段策略和小生境引導(dǎo)策略的MOEA/D 算法（MOEA/D-TPN），旨在解決Pareto 解的多樣性和分布性的不平衡問題。然后利用MOEA/DTPN 算法求解城市供水管網(wǎng)多目標(biāo)優(yōu)化問題，獲取雙環(huán)管網(wǎng)的最優(yōu)解集，為城市供水管網(wǎng)的工程設(shè)計(jì)提供可行參考方案。

1 城市供水管網(wǎng)多目標(biāo)優(yōu)化問題

管網(wǎng)形式主要分環(huán)狀和樹狀，樹狀主要用于邊緣地區(qū)，其他地區(qū)一般布置環(huán)狀管網(wǎng)以保證用水安全。Alperovits 與Shamir 提出的雙環(huán)管網(wǎng)是公認(rèn)的管網(wǎng)優(yōu)化模型，如圖1 所示。它包含了1個(gè)水源、7個(gè)節(jié)點(diǎn)和8個(gè)段長1 000 m 的管段，海曾威廉系數(shù)設(shè)定為130［16］。表1 給出了各節(jié)點(diǎn)的節(jié)點(diǎn)高程、最低自由水頭和需水量等參數(shù)。

圖1 雙環(huán)管網(wǎng)結(jié)構(gòu)圖Fig.1 Structure diagram of double ring pipe network

根據(jù)雙環(huán)管網(wǎng)和單元管段的結(jié)構(gòu)，可得到如下的供水管網(wǎng)經(jīng)濟(jì)性目標(biāo)函數(shù)，即管網(wǎng)的總造價(jià)［16］：

式中：P為管網(wǎng)每年大修與折舊的費(fèi)用在整個(gè)管網(wǎng)造價(jià)中的比例；E=［j（1+j）n］/［（1+j）n-1］為基建投資效果系數(shù)；j為銀行利率；n為資金回收期；np為管網(wǎng)的管段數(shù)；Di、li分別為第i段管道的直徑和長度；a+bDiλ為管網(wǎng)建造費(fèi)用，a、b和λ分別為管道造價(jià)公式的系數(shù)及指數(shù)；H0為水泵靜揚(yáng)程；hi為管段i的水頭損失；Qp為進(jìn)入管網(wǎng)的總流量，K=（24×356×10 000γσ）/（102η）=8 600γ σ/η；σ為泵站電價(jià)；γ為管網(wǎng)使用年限內(nèi)管網(wǎng)供水不均勻所占比；η為泵站效率。

由于節(jié)點(diǎn)富余水頭越大，管網(wǎng)中水壓越高，事故發(fā)生的概率就越大，因此可靠性目標(biāo)函數(shù)定義為節(jié)點(diǎn)富余水頭總和，即管網(wǎng)的故障率［16］：

式中：Hi為節(jié)點(diǎn)i處的自由水頭：Himin為節(jié)點(diǎn)i處要求的最小水壓；I為管網(wǎng)節(jié)點(diǎn)總數(shù)；Hi-Himin為節(jié)點(diǎn)i的富余水頭。

2 MOEA/D原理

MOEA/D利用分解方法將多目標(biāo)優(yōu)化問題轉(zhuǎn)換為一系列單目標(biāo)優(yōu)化子問題，并構(gòu)造問題導(dǎo)向的鄰域關(guān)系，然后以協(xié)同方式同時(shí)處理所有子問題，最終得到完整的Pareto 前沿［5］。常用的分解方法有權(quán)重聚合法、切比雪夫法和基于懲罰的邊界交叉法［5］。由于權(quán)重聚合法不能很好地處理真實(shí)Pareto 前沿為凹狀的問題，基于懲罰的邊界交叉法的性能受懲罰因子的影響較大，因此以切比雪夫法為例進(jìn)行研究。在該方法中，第i個(gè)子問題定義為如下形式：

其中，z*=min｛fi（x）|x∈Ω｝，i=1，2，…，m，λ=（λ1，λ2，…，λm）T為權(quán)重向量，滿足λi＞=0，且有

3 MOEA/D-TPN算法

考慮到公式（3）是z*作為參考點(diǎn)，當(dāng)其處理復(fù)雜凸優(yōu)化問題時(shí)，算法尋找到的Pareto 最優(yōu)解會向Pareto 前沿中部區(qū)域聚集。對于該問題，如果采用znad作為參考點(diǎn)，算法將能夠?qū)ふ业竭吔鐓^(qū)域的Pareto解。上述思想如圖2所示。

圖2 Pareto解集分布情況Fig.2 The distribution of Pareto solution set.

基于此，提出一種自適應(yīng)MOEA/D 算法，其將進(jìn)化過程分成兩個(gè)階段。在第一階段利用z*作為參考點(diǎn)，在第二階段利用znad作為參考點(diǎn)。同時(shí)為了增加算法的多樣性，設(shè)計(jì)了一套小生境引導(dǎo)策略，以選擇不同的交配個(gè)體。

3.1 兩階段策略

第一階段采用公式（3）將多目標(biāo)優(yōu)化問題分解成子問題進(jìn)行協(xié)同優(yōu)化。第二階段采用公式（4）將凸多目標(biāo)優(yōu)化問題分解成子問題進(jìn)行協(xié)調(diào)優(yōu)化。

其中，znad=max｛fi（x）|x∈Ω｝，i=1，2，…，m。

然而，在處理多目標(biāo)優(yōu)化問題時(shí)，一般無法事先知道問題的凸凹性。因此何時(shí)選取公式（3）和公式（4）對問題進(jìn)行求解存在困難。本文所提兩階段策略，在最大迭代步數(shù)的Mr%，使用公式（3）處理多目標(biāo)優(yōu)化問題。在第一階段結(jié)束時(shí)，使用基于擁擠的方法來評估MOEA/D 是否得到了一組均勻解［12］。該方法能夠估算出Pareto 前沿中間區(qū)域和極端區(qū)域的解的密度。如果中間區(qū)域的解比極端區(qū)域的解更密集，則說明多目標(biāo)優(yōu)化問題可能是凸的?；诖?，在第二階段，采用公式（4）處理多目標(biāo)優(yōu)化問題。

下面闡述基于擁擠的評估方法的原理［12］。由于在MOEA/D 權(quán)重向量λ滿足λ1+λ2+ … +λm= 1，根據(jù)幾何平均不等式定理，始終成立，因此λ位于超平面λ1+λ2+ … +λm= 1 的中心。換句話說，越大，權(quán)重向量就越接近中心權(quán)重向量。因此，將MOEA/D 的n個(gè)權(quán)重向量集分為兩個(gè)子集，即中間子集Wm和極端子集We：

其中，Ωλ是MOEA/D 運(yùn)用單純形格子法生成的n個(gè)權(quán)重向量集。公式（5）中設(shè)置系數(shù)為0.5，能夠取得Wm和We的最佳分類。

然后，根據(jù)權(quán)重向量子集Wm和We，將種群分為兩個(gè)子種群，分別稱為中間區(qū)域子種群Pm和極端區(qū)域子種群Pe。兩個(gè)子種群的擁擠程度采用公式（7）和（8）進(jìn)行估計(jì)：

式中：Dmid和Dext分別是中間子種群和極端子種群的擁擠度值；γ(i)表示子種群中第i個(gè)解的擁擠度，定義如下：

式中：d（i，j）為解i與j之間的距離；B（i）和T分別為解i鄰域內(nèi)的個(gè)體和鄰域大小。

3.2 小生境引導(dǎo)策略

傳統(tǒng)MOEA/D 在鄰域交配和更新時(shí)會產(chǎn)生許多類似的解。盡管在MOEA/D-DE 中使用DE 算子取代了模擬二進(jìn)制交叉算子，以克服該問題。但是當(dāng)Pareto前沿具有不連續(xù)區(qū)域時(shí)，這種缺陷會導(dǎo)致種群多樣性急劇降低。

參考文獻(xiàn)［12］提出一種具有自適應(yīng)選擇能力的小生境引導(dǎo)策略。在該策略中，首先計(jì)算出每個(gè)個(gè)體在其相鄰個(gè)體上的生態(tài)位計(jì)數(shù)。每個(gè)個(gè)體i的生態(tài)位計(jì)數(shù)是通過對相鄰個(gè)體求和獲得的一個(gè)共享函數(shù)，用nc（i）來表示：

式中：dij是個(gè)體i和j之間的距離；sh（dij）是測量個(gè)體i和j之間相似程度的共享函數(shù)，定義為：

式中：σshare是一個(gè)預(yù)先定義的生態(tài)位半徑；α是一個(gè)常數(shù)，稱為共享水平。

如果個(gè)體的生態(tài)位計(jì)數(shù)超過了給定的閾值，這意味著個(gè)體與它的鄰接個(gè)體相似，所以最好選擇鄰域之外的個(gè)體作為交配的父代種群。因此，自適應(yīng)更新定義如下：

其中，β是與生態(tài)位半徑σshare密切相關(guān)的閾值，即β主要由σshare決定。由于nc（i）的最大值為T，將β設(shè)為T/2。P和P′產(chǎn)生方法如下：

其中，當(dāng)rand1和rand2在范圍［0，1］內(nèi)獨(dú)立生成隨機(jī)數(shù)時(shí)，P與MOEA/D-DE 中的P保持相同。P′的定義意味著如果個(gè)體i的生態(tài)位計(jì)數(shù)達(dá)到閾值β，它有50%的機(jī)會選擇其周圍的個(gè)體作為交配父代種群。

3.3 MOEA/D-TPN算法框架

其中，是的第k個(gè)分量，xrk1、xrk2和xrk3是從種群中隨機(jī)選擇的3 個(gè)不同個(gè)體。CR和F是兩個(gè)控制參數(shù)。采用多項(xiàng)式突變策略從產(chǎn)生一個(gè)解y=（y1，…，yn）：

其中，rand是［0，1］區(qū)間隨機(jī)數(shù)，η和Pm分別是分布指數(shù)和突變率，bk和ak分別是決策變量的上下界。

4 實(shí)驗(yàn)分析

4.1 參數(shù)設(shè)置

MOEA/D-TPN 的參數(shù)設(shè)置如下：種群規(guī)模N=100，鄰域規(guī)模T=20，nr=2，最大迭代次數(shù)maxIteration=100；從鄰域內(nèi)選取解的概率δ=0.9；差分進(jìn)化的雜交率CR=0.5，尺度因子F=0.5；多項(xiàng)式變異的變異率pm=1/n，分布性指數(shù)η=20。

4.2 基準(zhǔn)測試

為了評價(jià)MOEA/D-TPN 的有效性，采用6 個(gè)具有長尾巴、尖峰等復(fù)雜POF 的MOP （F1-F6）作為測試案例［12］，決策變量個(gè)數(shù)均設(shè)置為20，取值范圍在［0，1］n，具體信息見表2。

表2 具有復(fù)雜POF的MOPTab.2 MOP with complex POF

4.4 實(shí)驗(yàn)結(jié)果

圖3 給出了MOEA/D-TPN 的逼近結(jié)果，并與MOEA/D-DE和MOEA/D-STM 進(jìn)行了對比分析?？梢钥闯觯瑢τ贔1 測試問題，MOEA/D-DE 和MOEA/D-STM 獲取的解會收斂到POF 的中間區(qū)域，難以發(fā)現(xiàn)邊界解，MOEA/D-TPN 容易找到邊界解。對于F2 測試問題，MOEA/D-DE、MOEA/D-STM 和MOEA/D-TPN三個(gè)算法獲取到解的精度基本相同，但是MOEA/D-TPN 的收斂速度最快。F3、F4 和F5 三個(gè)測試函數(shù)極其復(fù)雜，具有長尾和尖峰區(qū)域，MOEA/D-DE 僅在POF 中部區(qū)域上找到一部分解，MOEA/D-TPN 在邊界區(qū)域也能尋找到一部分解。對于更復(fù)雜的三目標(biāo)優(yōu)化問題F6，MOEA/D-TPN 也表現(xiàn)出了更優(yōu)越的尋優(yōu)性能。

圖3 3種算法的逼近結(jié)果Fig.3 Approximation results of three algorithms

圖4給出了3種對比算法30次獨(dú)立運(yùn)行取得最小IGD 值時(shí)的演化曲線?？梢钥闯觯瑢τ贔1-F6 測試函數(shù)，與MOEA/D-DE和MOEA/D-STM 相比，MOEA/D-TPN 的進(jìn)化曲線均表現(xiàn)出較快的收斂速度和較好的逼近精度。

圖4 3種算法的IGD進(jìn)化曲線Fig.4 IGD evolution curve of three algorithms

圖5 給出了3 種對比算法30 次獨(dú)立運(yùn)行繪制的IGD 盒形圖?？梢钥闯?，對于F1、F3、F4 和F5 四個(gè)測試函數(shù)，MOEA/DTPN的IGD總體處于較低范圍，分布性好，尋優(yōu)性能突出。對于F2 和F6，MOEA/D-TPN 與MOEA/D-DE、MOEA/D-STM 的性能相當(dāng)，略占優(yōu)勢。

圖5 3種算法的IGD盒形圖Fig.5 IGD box diagram of three algorithms

4.4 供水管網(wǎng)優(yōu)化

下面利用所設(shè)計(jì)的MOEA/D-TPN 算法對雙環(huán)管網(wǎng)的造價(jià)（經(jīng)濟(jì)性）和節(jié)點(diǎn)富余水頭總和（可靠性）兩個(gè)目標(biāo)進(jìn)行優(yōu)化。MOEA/D-TPN 的參數(shù)設(shè)置同上。圖6 給出了MOEA/D-TPN 算法獲取到的Pareto 前沿。每一個(gè)Pareto 解對應(yīng)的決策空間變量是［x1，x2，x3，x4，x5，x6，x7，x8］，分別表示雙環(huán)管網(wǎng)的8個(gè)管徑。

圖6 管網(wǎng)的Pareto前沿Fig.6 Pareto front of pipe network

本文設(shè)計(jì)一種基于模糊隸屬函數(shù)的智能決策方法，利用該方法能夠在Pareto 前沿上決策出一個(gè)折衷解，從而獲得8 個(gè)管徑的一個(gè)最佳組合。首先利用模糊隸屬函數(shù)定義第i個(gè)目標(biāo)函數(shù)fi的Pareto解xk的滿意度，如下：

其中，fmaxi和fmini分別代表第i個(gè)目標(biāo)函數(shù)fi的最大值和最小值。然后， 可以計(jì)算出xk的規(guī)范化隸屬度函數(shù)μk，如下：

其中，m為目標(biāo)函數(shù)個(gè)數(shù)，|S|是算法獲得的Pareto 解集的元素個(gè)數(shù)。本文選取μk為最小值時(shí)的解，即可圖6中的折衷解。

從圖6中可以看出，MOEA/D-TPN能夠獲取更多的邊界解。此時(shí)，從邊界中進(jìn)行取值能夠取得經(jīng)濟(jì)性和可靠性的最佳平衡。

從圖6 的Pareto 前沿中選取25 個(gè)解，列出決策變量和目標(biāo)函數(shù)的取值，如表3所示。可以看出，管網(wǎng)造價(jià)和節(jié)點(diǎn)富余水頭總和之間存在沖突特性。隨著管網(wǎng)造價(jià)的逐步升高，節(jié)點(diǎn)富余水頭總和越來越低，表明總造價(jià)增大的同時(shí)，故障率逐步降低，可靠性逐步增大；反之總造價(jià)減少的話，故障率逐步抬高，可靠性逐步降低。這在管網(wǎng)施工時(shí)就會產(chǎn)生沖突，施工方一方面期望能盡量降低總造價(jià)，另一方面也期望故障率盡可能地低，可靠性盡可能地高。利用MOEA/D-TPN 決策出的最優(yōu)解將能夠在保證足夠可靠性的前提下盡量降低成本，為城市供水管網(wǎng)優(yōu)化設(shè)計(jì)提供有益參考。

表3 算法獲取的Pareto解Tab.3 Pareto solutions obtained by algorithm

5 結(jié) 語

本文提出一種自適應(yīng)MOEA/D 算法，命名為MOEA/DTPN，利用其解決城市供水管網(wǎng)多目標(biāo)優(yōu)化設(shè)計(jì)問題。為了平衡MOEA/D 解的多樣性和分布性，設(shè)計(jì)了兩階段策略和小生境引導(dǎo)策略。利用6 個(gè)具有復(fù)雜Pareto 前沿的測試函數(shù)進(jìn)行了算法性能驗(yàn)證，并與MOEA/D-DE、MOEA/D-STM 等多個(gè)先進(jìn)算法進(jìn)行了比較，在4 個(gè)測試函數(shù)中MOEA/D-TPN 的尋優(yōu)結(jié)果最好，在其余兩個(gè)測試函數(shù)中性能表現(xiàn)相當(dāng)。根據(jù)我國供水系統(tǒng)的運(yùn)行情況，利用城市供水雙環(huán)官網(wǎng)模型進(jìn)行了仿真測試，確定了總造價(jià)和故障率兩個(gè)沖突的目標(biāo)函數(shù)，運(yùn)用MOEA/D-TPN算法對管網(wǎng)進(jìn)行多目標(biāo)優(yōu)化設(shè)計(jì)，求解出雙環(huán)管網(wǎng)的最優(yōu)解集。下一步工作是利用本文設(shè)計(jì)的MOEA/D-TPN 解決大型復(fù)雜管網(wǎng)的多目標(biāo)優(yōu)化問題。