基于多目標粒子群算法的停泵水錘防護優化

劉亞萌，蔣 勁，李 婷，劉 承，陳英強

(武漢大學動力與機械學院 水力機械過渡過程教育部重點實驗室，武漢 430072)

水錘現象一直以來是泵站及輸水工程的常見問題。事故停泵造成的壓力波動會造成管道產生巨大的水錘壓力和汽化現象，嚴重時將導致管道破裂、設備損壞等事故，威脅泵站的運行安全[1,2]。通常需要增加必要的水錘防護措施，以保證泵站系統的安全運行，但由此會帶來工程投資的增加，因此為了解決安全和經濟之間的矛盾，對其進行優化計算是十分必要的[3]。

目前實際工程和科學研究的對象大多是多目標優化問題，與單目標優化不同的是，多目標優化涉及的參數較多，目標量不唯一，導致在算法選擇和計算工作量上有較大的難度[3-5]。本文計算了國內某輸水工程的停泵過渡過程，以水泵出口閥門的關閥規律和單向調壓塔的直徑、初始液位作為設計變量，以系統最大、最小壓力以及單向調壓塔的有效體積作為目標函數，建立數學模型并采用多目標粒子群算法對該工程進行優化計算，得出了有效可靠且經濟的水錘防護方案。

1 泵站水錘數值計算

國內某輸水工程，全長2.4 km，其中引水隧洞長為347 m，直徑為2.2 m，泵站布置4臺臥式單級、雙吸離心泵，總裝機容量4×2 240 kW，水泵額定轉速1 480 r/min，轉動慣量為223.24 kg/m2，設計揚程為164.8 m，設計流量2.90 m3/s，整個輸水管線距離短、揚程高、流量大且管線布置起伏較大，駝峰點較多，具體的管線布置圖如圖1所示。

圖1 輸水工程管線縱剖面圖Fig.1 Longitudinal section of pipeline

1.1 水錘計算特征線法

泵站水錘的計算是針對整個水泵抽水裝置，包括管道內點及與管道連接的泵裝置中的各部分。在水錘計算中，對于管道系統內點的計算是求解水錘基本方程，即由運動方程和連續性方程組成的雙曲型偏微分方程組[6,7]。采用特征線方法將該偏微分方程組離散化，沿特征線方向將它轉換為水錘全微分方程：

(2)

由公式(1)和(2)進行有限差分近似，可得編入計算機程序的相容性方程為：

(4)

式中：B=a/gA；R=fΔx/2gDA2。

1.2 單向調壓塔模型及工作原理

單向調壓塔與管道中間裝有止回閥，只允許塔中水單向流入主干管中，其構造原理圖如圖2所示。水泵正常運行時，由于單向調壓塔內水位低于管中正常工作水頭，止回閥處于關閉狀態，向主管道補給水體的流量為零，即Qp3=0；發生事故停泵后，當管道中的壓力低于單向調壓塔預先設定的壓力值時，即Hp

圖2 單向調壓塔構造原理圖Fig.2 Schematic diagram of one way Surge Tower

由連續性原理可知：

(5)

若任意時刻單向調壓塔內水位HP3，可由此計算時段Δt初的水位H3和流出水體的體積對應的水深，求得：

HP3=H3-0.5Δt(QP3+Q3)/Ast

(6)

主管道的相容性特征線方程為：

QP1=(Cp-Hp)/B

QP2=(Cp-Hm)/B

(7)

由公式(5)～(7)，可解出QP3的表達式為：

QP3=g(CdAp)2-{0.5 (B+Δt/Ast)+{0.25 (B+Δt/Ast)2-

[Cp+Cm+Δt/AstQ3]-2 (H3+ZT)/[g(CdAp)]2}1/2}

(8)

聯立公式(7)、(8)求得：

(9)

式中：Qp1、Qp2分別為Δt流進、流出調壓塔的流量；Qp3、Q3分別為Δt時段末、初由調壓塔流向主管道的流量；Hp3、H3分別為Δt時段末、初調壓塔的水位；Ast為單向調壓塔的面積；Hp為Δt時段末的壓力水頭。

1.3 無防護措施事故停泵計算結果

泵站4臺機組同時正常運行時系統最大壓力出現在水泵出口，壓力值為197.5 m；最不利的運行工況即4臺機組同時事故停泵，無防護措施兩階段關閉閥門(4 s快關至0.3，40 s全關)，最大壓力值為415.5 m，達到穩態運行最大壓力的2.1倍，引水管及沿線多處出現-8 m的壓力，圖3為正常運行和事故停泵關閥的管道壓力包絡線。實踐證明不加防護措施僅控制泵出口的閥門不能改善水錘壓力過大和負壓情況，因此設置單向調壓塔進行防護，但是單向調壓塔的尺寸以及閥門的關閉規律對系統的影響很大，為使整個系統安全運行并降低工程的投資成本，將對調壓塔尺寸和閥門的關閉規律進行優化。

圖3 正常運行和事故停泵關閥的管道壓力包絡線Fig.3 Pipe pressure envelope for normal operation and shutdown of pump

2 優化模型建立

本文采用的多目標粒子群優化算法(Multi-Objective Particle Swarm Optimization,MOPSO)對停泵水錘計算進行優化，該算法將原來只能用于單目標上的粒子群算法(PSO)經過改進應用于多目標優化問題中,具體算法如下：

設粒子群算法中粒子本身找到的最優解為個體極值(pBest)，整個群體所經歷過的最好位置為全局極值(gBest)，粒子 的速度和位置將按如下公式進行更新[11,12]：

Vi=wVi+c1rand() (pBest[i]-Xi)+

c2Rand() (pBest[g]-Xi)

(10)

Xi=Xi+Vid

(11)

式中：c1，c2稱為學習因子，均為常數；rand()和Rand()是[0，1]上的隨機數；w為慣性權重。

解空間的粒子根據式(10)和式(11)不斷調整自己的位置和速度。一個具有n個設計變量，m個目標的多目標優化可以描述為：

miny=F(x)=[f1(x),f2(x),…,fm(x)]

(12)

s.t.gi(x)≤0

(13)

式中：x=(x1,x2,…,xm)∈X?Rn為n維決策向量，X為n維決策空間；y=(y1,y2,…,xn)∈Y?Rm為m維的目標變量；Y為m目標空間；gi(x)≤0是系統約束。

通常情況下，各個目標函數之間可能相互制約甚至是矛盾沖突的，難以使每個目標函數同時達到最優，而多目標優化的最終目的是在各個目標之間進行協調權衡和折衷處理，尋找到相對較優解[13,14]。多目標粒子群優化算法可按照如下公式對傳統粒子群算法采用變異策略：

vm=2 (r3-1)βVmax

(14)

xdi(t)=xdi(t)+vm

(15)

式中：vm是變異值；β∈[0,1]為變異系數，可調節變異程度；r3為在[0,1]范圍隨機變化的值；xdi表示第i個粒子的隨機選中的第d維。

經過改進后，多目標粒子群算法可提高解的多樣性，增加粒子的全局搜索能力，且由變異產生的優異粒子能夠對其他解產生吸引，以此增強粒子對局部最優的逃逸能力，避免其陷入局部最優[7,15]。通常優化計算的過程流程圖如圖4所示：首先需要建立模型，明確工程中的設計變量、優化目標和約束條件，進而確定設計變量的范圍和參數初值，在優化算法的驅動下進行迭代計算，若目標值的計算結果不滿足工程要求，需要重新給定初值進行優化，若滿足可輸出結果。

圖4 優化流程圖Fig.4 Optimization flow chart

針對本文計算的系統，在泵站其他參數均不變的前提下，僅對約束變量進行優化，選取系統中對水錘防護影響較大的參數進行約束，其中設計變量為：單向調壓塔的直徑、初始液位、閥門快關時間、快關百分比、慢關時間；目標變量為：系統最大壓力、出水管最小壓力、引水管最小壓力和單向調壓塔的有效體積，約束的條件需滿足《泵站設計規范》以及本工程的設計要求，據此設計變量、目標量和約束的設定具體如表1和表2所示。

表1 設計變量和約束條件的設定Tab.1 Design variables and constraints

3 優化計算

表2 預期優化目標值Tab.2 expected optimization target value

注：D、H分別為單向調壓塔的直徑和初始液位。

參照表1和表2列出的設計變量、目標變量和約束條件，在模型優化計算前初步試算了20種方案，并提取出其中5種有代表性的結果，如表3所示。

從方案1、2、3可看出，在不改變調壓塔尺寸的情況下，增加快關速率和慢關時間，系統最大水錘壓力明顯減小且負壓減輕；對比方案3、4、5可得，關閥方案相同時，改變單向調壓塔的直徑和初始液位，調壓塔的有效體積越大，負壓情況改善明顯，但同時工程造價也就越高。這5種方案的計算結果雖能改善最大壓力和負壓，但均沒有達到工程要求，因此以方案3的參數作為初始值，并結合表1的約束和工程要求，設定各參數的范圍進行優化計算，其中：快關時間在4～6 s內、慢關時間在40～60 s內，慢關百分比在0.1～0.2內(快關百分比=1-慢關百分比，故不再單獨設定)、單向調壓塔的直徑和初始液位分別在3～5 m和3～7 m內。

表3 試算方案的參數設定及結果Tab.3 Parameter setting and results

采用多目標粒子群算法進行優化后，選取快關時間和調壓塔初始液位2個設計變量的優化歷史，如圖5所示，圖中兩條直線相交的點即為最優的取值。在優化過程中，可清晰地看出初始設計變量的取值遍布了整個變量可行域，說明優化算法中的粒子具有全局搜索的能力，隨后逐漸縮小取值范圍，沿較優的結果方向進行計算，經過600次的迭代后，得到了優化目標的全局最優解。結果顯示在第392次迭代時獲得最優規律：(水泵在第10 s發生事故停泵，故快關和慢關時間需減去10 s)快關時間為4.54 s，慢關百分比為0.10(即4.54 s關閥至0.1)，慢關時間為59.96 s，單向調壓塔的初始液位6.19 m，直徑3.18 m，在該工況下系統最大壓力為247 m，出水管最小壓力-1.53 m，引水管2點處的最小壓力-0.18 m，且單向調壓塔有效體積為49.16 m3相對較小，滿足運行可靠和經濟的要求。

圖6和圖7是最大壓力和最小壓力與調壓塔尺寸及關閥參數的三維圖形。

圖6中最大壓力值與左圖的關閥參數在XZ、YZ軸上規律性較強，相比而言難以提取最大壓力與調壓塔尺寸之間的關系，說明在該泵站系統中慢關時間、關閥速率對最大壓力的貢獻率較大。分析原因可知，停泵后管道中的水流先由正向流動逐漸減小到零流速，在到達零流速之前快關閥門至某一角度，會造成很大的局部阻力，而閥門此時尚未全部關閉，因此水流在重力作用下倒流會使升壓減小；再逐漸慢關閥門，使流速變化的增量減小，由于壓力的升高與流速的變化成正比，從而可使管道的升壓限制在允許的范圍之內。因此合理控制關閥速率和慢關時間，可達到減小水錘壓力升高的防護效果。

同理由圖7可知，單向調壓塔的直徑和初始液位對負壓的改善較為明顯，而最小壓力與關閥規律的曲面圖起伏較大，凸起點較多規律不明確。單向調壓塔須具有足夠的容積，以滿足主管道內產生負壓及汽穴空腔時所需的補給水量，保證塔內水體不被泄空；其次由于塔內水體是靠重力作用向管道補水，故調壓塔的初始液位需滿足一定的要求，如果水壓不夠，向管道補水不及時則難以達到破壞水柱分離的效果，因此單向調壓塔的直徑和初始液位對負壓有著較大的影響。

圖5 快關時間和調壓塔初始液位的優化過程Fig.5 The optimization process of the quick closing time and the initial level of the surge tank

圖6 最大壓力與關閥參數和調壓塔尺寸的三維圖Fig.6 3D diagram of maximum pressure between closing valve parameters and size of surge tank

圖7 最小壓力與單向調壓塔尺寸和關閥參數的三維圖Fig.7 3D diagram of the minimum pressure between closing valve parameters and size of surge tank

4 結果對比

將初始方案和優化計算后的結果進行對比，具體見表4，可知優化后調壓塔的有效體積由87.96 m3降至49.16 m3，減小了44.1%；最大壓力值由276.71 m降至247 m，最小壓力值雖有初始的-0.61 m降為-1.53 m，但依然控制在了-2 m以內，符合《泵站設計規范》以及本工程的設計要求，這也說明了多目標優化的解的特點，即該解對一個或幾個目標函數不能進一步改進，且其他目標值不至于劣化，即為尋求到的最優解。

由于優化方案的設計變量取值較為精確，實際操作中難以實現，故提出了建議值。建議方案與優化計算的方案相比，在關閥規律和調壓塔參數上做了調整，更容易泵站的實際操作，且調整后建議方案的各個目標量的結果均符合設計要求。

表4 計算結果對比 Tab.4 Comparison of calculation results

圖8和圖9分別為初始方案和優化方案的壓力包絡線，以及二者調壓塔液位的變化曲線對比圖，由圖8可明顯看出水錘壓力有較大程度的減小，最大壓力由穩態運行壓力的1.4倍降為1.25倍。圖9中優化后單向調壓塔的初始液位由初始方案的7 m降為6.19 m，初始液位的降低將改變調壓塔的有效體積，進而設計調壓塔的實際尺寸時也可適當減小，降低了工程投資，達到了經濟、可靠的目的。

圖8 初始方案和優化方案的管道壓力包絡線對比圖Fig.8 Comparison of pipeline pressure envelope with initial scheme and optimization scheme

圖9 初始方案和優化方案的調壓塔液位變化圖Fig.9 Comparison of liquid level of pressure regulating tower in initial scheme and optimization scheme

5 結 語

本文以國內某輸水工程為例，采用多目標粒子群算法進行優化計算，將閥門關閉規律和單向調壓塔初始液位和直徑作為設計變量，控制調壓塔體積和最大、最小壓力，尋找到最優的方案，且得到以下結論。

(1)采用多目標粒子群算法對本工程進行優化計算后，尋找到了折衷最優解，單向調壓塔體積由87.96 m3降至49.16 m3，最大水錘壓力由穩態壓力的1.4倍降至1.25倍。

(2)本文結合約束條件和試算方案中的較優值，給定了設計變量的范圍，從而通過精細化計算得到更優的結果，為優化計算減少了工作量，提高了運算效率。

(3)多目標粒子群算法的優化結果與初始方案比較，各目標值均有較大的改善，能夠在工程要求和經濟投資中尋找到較優解，說明該優化算法能夠有效地計算水錘問題，且為其他工程實際中的多目標優化問題提供了解決思路。

□

[1] 王 娜. 泵站壓力管道的水錘研究[J]. 水利規劃與設計, 2016,(1):85-88.

[2] 王 圃, 許蘭森, 王 穎. 高揚程多起伏輸水管道水錘模擬及其防護[J]. 給水排水, 2014,(6):112-114.

[3] 肖曉偉, 肖 迪, 林錦國,等. 多目標優化問題的研究概述[J]. 計算機應用研究, 2011,28(3):805-808.

[4] 謝 濤, 陳火旺. 多目標優化與決策問題的演化算法[J]. 中國工程科學, 2002,4(2):59-68.

[5] DEB K,PRATAP A,AAGRWALI S,et al. A fastand elitist multiobjective genetic algorithm: NSGA-II[J]. IEEE Trans on Evolutionary Computation,2002,6(2):182-197.

[6] 蘭 剛, 蔣 勁, 李東東,等. 空氣罐在長距離輸水管線水錘防護中的應用[J]. 水電站機電技術, 2013,(4):42-44.

[7] 楊祖強, 羅 希, 王文文. 基于改進型多目標遺傳算法的含氣水錘防護方案優化研究[J]. 水資源與水工程學報, 2014,25(1):200-204.

[8] 劉梅清, 劉志勇, 蔣 勁. 基于遺傳算法的單向調壓塔尺寸優化研究[J]. 中國給水排水, 2008,24(23):56-60.

[9] 徐艷艷. 長距離高揚程多起伏輸水管道采用箱式雙向調壓塔等措施的水錘防護研究[D]. 西安：長安大學, 2008.

[10] 劉光臨, 劉志勇, 王聽權,等. 單向調壓塔水錘防護特性的研究[J]. 給水排水, 2002,28(2):82-85.

[11] SIERRA M R, COELLO C A C. Improving PSO-based multi-objective optimization using crowding, multation and Edominance[C]∥ Lecture Notes in Computer Scence. Proc of the 3rd International Conference on Evolutionary Multi-Criterion Optimization,2005:505-519.

[12] 劉衍民, 牛 奔, 趙慶禎. 多目標優化問題的粒子群算法仿真研究[J]. 計算機應用研究, 2011,28(2):458-460.

[13] TSAI S J,SUN T Y,LIU C C. An improved multi-objective particle swarm optimizer for problems[J].Expert Systems with Applications,2010,37:872-886.

[14] 畢榮山, 楊 霞, 譚心舜,等. 基于動態Pareto解集的微粒群優化算法及其在多目標規劃中的應用[J]. 計算機工程與應用, 2004,40(32):85-88.

[15] 劉衍民, 牛 奔, 趙慶禎. 多目標優化問題的粒子群算法仿真研究[J]. 計算機應用研究, 2011,28(2):458-460.