基于SA-FOA的水庫防洪優化調度研究

王成民，王 昊，刁艷芳

（1.日照水庫管理運行中心，山東 日照 276800； 2.山東農業大學水利土木工程學院，山東 泰安 271018）

0 引 言

洪水災害是目前發生頻率最高的自然災害，據統計，僅2018年全球洪水災害事件共發生109起，死亡人口和受災人口分別為1 995 人和1 262 萬人，直接經濟損失45億美元［1］，其中，中國、印度、印度尼西亞、美國等國家是洪水災害發生次數和遭受損失最多的國家。雖然在近30年間，全球洪水災害死亡人口及受災規模都呈現持續下降趨勢，但經濟損失卻呈上升趨勢。水庫是在河流上筑壩進而起到防洪效果的一種重要水利工程，且隨著城市供水規模的不斷擴大和農業用水灌溉方式的多樣化，水庫在防洪、供水、發電及生態環境保護等方面發揮了巨大的社會和經濟效益。科學的水庫防洪調度設計是安全經濟建設水庫與實時調度發揮設計效益的關鍵。目前，常見的兩類水庫調度方法包括水庫常規調度和水庫優化調度，后者由于可以通過設定優化目標來制定優化作業方案，因此在水庫防洪調度中廣泛應用。起初，水庫優化調度所采用的方法主要有線性規劃［2］、非線性規劃、動態規劃［3］、逐次優化算法［4］及大系統分解協調原理［5］等，然而，面對由水庫群、蓄滯洪區、湖泊及其他防洪工程組成的復雜防洪系統，上述優化算法表現出收斂效率低、維度低等明顯局限性。

現代計算技術的進步使得應用生物學、物理學和人工智能原理的啟發式智能算法得到了廣泛發展，這些算法可以有效解決上述問題，主要包括遺傳算法［6］、蟻群算法［7］、人工魚群算法［8］及粒子群算法［9］等。由于智能啟發式算法相比傳統的優化算法能夠解決非線性、間斷、不可微和多維的復雜問題，因此在水庫優化調度中得到了廣泛的應用。然而，啟發式智能算法也存在不足，例如遺傳算法易早熟且計算量大；蟻群算法、人工魚群算法過于復雜且計算量大；粒子群算法易陷入局部最優值。上述智能算法因其各自尋優時存在的弊端，使得應用受到限制。果蠅優化算法（FOA）［10］是近年來提出的一種新型智能優化算法，與其他算法相比，它具有原理簡單、參數少、計算量小及易于實現等優點，故受到了國內外眾多學科領域學者的廣泛關注。通過對該算法的大量研究發現，該算法在遇到高維優化問題時，存在過早收斂和易于陷入局部最優等問題，影響了算法的穩定性［11］。為彌補FOA的不足，采用設置變動的搜索距離L和模擬退火算法（Simulated Annealing algorithm，SA）的Metropolis準則［12］對FOA的迭代過程進行擾動，有效避免上述缺陷。

基于優先考慮水庫及下游防護區防洪安全的前提，建立水庫雙目標防洪優化調度模型；然后，提出耦合模擬退火算法的果蠅優化算法（SA-FOA）對模型進行求解；最后，以日照水庫2020年的兩場洪水調度過程為例，驗證了該算法的可行性和合理性。

1 水庫雙目標防洪優化調度模型的建立和求解

1.1 水庫雙目標防洪優化調度模型的建立

水庫自身、庫區和下游防護區是水庫防洪調度的3個主要防護對象，其中，為保證水庫及庫區安全，希望水庫汛期盡可能貼近汛限水位運行，故常以水庫調洪最高水位最低或最大防洪庫容最小作為目標函數；為保障下游防護區安全，希望水庫能夠最大限度地削減洪峰流量，故常以最大削峰準則作為目標函數［13］。本文綜合考慮3個主要防護對象，建立水庫雙目標防洪優化調度模型，以實現調洪最高水位最低和最大限度削減洪峰流量。

1.1.1 目標函數

（1）水庫及庫區防洪安全目標。調洪最高水位越小，水庫及庫區越安全，因此以壩前最高水位作為目標函數，即：

式中：T為調洪歷時；Z(t)為t時刻水庫壩上水位，m。

（2）下游防護區安全目標。經證明，水庫最大削峰準則的目標是水庫出庫流量過程平方和最小，即下泄流量最均勻［14］。利用該準則所構建的目標函數的表達式為［14-16］：

式中：q(t)為t時刻水庫的泄流量，m3/s；R（t）為t時刻水庫到下游防護區之間的區間來水流量，如無區間來水，其取值為0，m3/s。

1.1.2 約束條件

（1）水量平衡方程。

（2）水庫下泄流量約束。

（3）水庫泄流能力約束。

（4）庫容約束。

（5）非負約束：以上變量均為非負。

式中：Q(t)為第t時段水庫的入庫流量，m3/s；V（t）為第t時段水庫的庫容，m3；Δt為時段長，h；Φ[V(t)]為水庫蓄水量為V(t)時的最大下泄能力，m3/s；q(t)min、q(t)max分別為第t時段水庫的泄流量下限、上限，m3/s；V(t)min、V(t)max分別為第t時段水庫的庫容下限、上限，m3。

1.2 模型的求解

1.2.1 FOA計算步驟

FOA由潘文超于2011年提出［17］，是一種通過模擬果蠅種群覓食行為尋求最優解的方法，其尋優過程如下：

（1）初始參數的設置，包括果蠅群里初始位置X0和Y0，最大覓食代數Gmax，搜索距離L，種群規模size。

（2）果蠅個體通過嗅覺隨機搜尋食物，其位置為：

式中：i=1，2，…，size。

（3）計算果蠅個體與原點之間的距離Disti，再計算味道濃度判斷值Si，該值可取為距離的倒數，即：

（4）將Si代入味道濃度判斷函數中，求解果蠅個體的味道濃度Smelli：

（5）找出果蠅種群中味道濃度最佳值Smellb，如果目標函數是求解最小值，則味道濃度最小值的個體為最佳，反之亦然。

（6）保留最佳味道濃度Smellb及其位置坐標，此時果蠅種群利用視覺向該位置飛行。

（7）迭代尋優，重復步驟（2）～（4），判斷當前味道濃度是否優于前一次迭代的最佳味道濃度，并判斷當前迭代次數是否小于Gmax，若是執行步驟（6）；否則，結束迭代輸出最優值。

1.2.2 SA算法

SA算法是Kirkpatrick等人［18］受到固體或熱力學系統的物理退火法的啟發而提出的一種概率優化方法。在該方法中，熱力系統的當前狀態、能量方程和基態類似于優化問題的當前調度解、目標函數和全局最優解。SA算法使用基于概率的Me‐tropolis準則搜索全局最優并跳出局部最優，該準則設置接受劣解的概率［18，19］，其定義如下：

式中：G為迭代次數；pG是第G次迭代中的接受概率；EG是第G次迭代中的目標函數值；Eg是歷史最優目標函數值；KbO是Boltzmann系數；TG是第G次迭代中的退火溫度，該溫度最初設置為較大的值，然后通過以下溫度控制函數減小到較小的值：

式中：α為退火系數，其取值區間為［0.8，1.0］。因此，隨著TG的逐漸減小，劣解的接受概率將接近0。

1.2.3 SA-FOA

FOA中果蠅移動的搜索距離L會影響其尋優的結果，L大有利于跳出局部極值進而尋求全局最優，并且提高收斂效率，但其局部搜索能力將下降；L過小會導致陷入局部極值，且降低收斂效率。在算法迭代的過程中，對L大小的需求不同。在算法迭代初期，需要較大的L以提高收斂速度；在迭代后期，需要較小的L以避免錯過最優值。因此，本文采用隨著迭代次數逐漸減小的搜索距離，公式如下：

再者引入SA算法，將目標函數增量?s=SmellbG+1-SmellbG作為衡量標準，其中，SmellbG為第G次迭代中果蠅種群中味道濃度最佳值，根據Metropolis準則依概率對產生的新個體進行選擇，選取進入下一代的群體。具體方法見1.2.4節的第（6）步。該方法可以很好地確定搜索方向，又能依概率舍棄，跳出局部最優，使算法能夠較快的搜索到最優位置。

1.2.4 SA-FOA確定水庫防洪優化調度方案的步驟

SA-FOA確定水庫防洪優化調度方案的步驟如下：

（1）初始參數的設置，包括最大覓食代數Gmax，初始搜索距離L0，當前覓食代數G=1，種群規模size，空間維度dim=T，即洪水時段數，初始溫度TG，退火系數α。

（2）隨機生成初始種群，即水庫初始泄流過程線，即：

（3）將種群代入味道濃度判斷函數中，求出果蠅個體的味道濃度值，即：

式中：為第i個種群在第G次覓食時的最高水位，m；Zmax和Zmin分別為調度期水庫允許的最高水位和最低水位，m；為第i個種群在第G次覓食時泄流量的平方和分別為第G次覓食時，所有種群泄流量平方和的最大值和最小值。

（4）尋求果蠅種群中味道濃度最佳的個體。本文是求最小值問題，則果蠅種群中味道濃度最小值的個體為最佳，即SmellbG。

（5）采用式（13）計算LG，果蠅群體的位置變為：

再由第（3）和第（4）步驟求出SmellbG+1。

（6）計算味道濃度差值?s=SmellbG+1-SmellbG。依據Metropolis準則，當Δs≤0時，SmellbG+1被接受，保存該果蠅個體；當Δs>0時，依據式（11）的概率pG=exp(-Δs/TG)收斂。若此時的pG>rand，則新的最小值被接受，將(t)替代作為下次尋優的初始位置；否則，保持不變，開始進行下一次尋優。

（7）令G=G+1，如果G≤Gmax，轉入第（3）步繼續進行迭代，否則停止迭代，輸出最優解。

2 實例研究

2.1 日照水庫簡介

日照水庫位于山東省日照市東港區，控制流域面積548 km2，總庫容3.18億m3，興利庫容1.82 億m3，是一座集防洪、供水、發電、灌溉等多種功能于一體的大（II）型水庫，同時也是日照市區重要飲用水水源地。日照水庫保護著下游后村鎮、高興鎮、南湖鎮等7個鄉鎮、街道辦事處90個自然村，9.06萬人民群眾的生命財產安全，及204國道、兗石鐵路、日蘭和沈海高速公路和青連高鐵等交通干線以及臨港大中型企業等防洪重點保護對象。日照水庫防洪標準為100年一遇設計、5 000年一遇校核，汛限水位42.50 m、興利水位43.00 m、防洪高水位43.50 m、設計洪水位43.97 m、校核洪水位46.52 m。汛期洪水調度原則為，當雨前或預計雨后水位低于42.50 m時，不泄洪；當雨前或預計雨后水位42.50～43.50 m時，控泄流量1 000 m3/s；當雨前或預計雨后水位高于43.50 m時，溢洪閘自由敞泄；當水位高于43.97 m時，北放水洞參與泄洪。

2.2 結果對比

收集了日照水庫20200813和20200826兩場洪水資料，包括水庫的水位、庫容、入庫流量和泄流量，時段長為2 h，同時收集了水庫的水位～庫容～泄量關系等數據，滿足建立水庫防洪優化調度模型的要求。參考已有文獻［11，17-19］及對兩場洪水的計算，確定SA-FOA的初始參數為Gmax=2 000，L0=2，size=100，T1=106，α=0.9。采用SA-FOA對日照水庫雙目標防洪優化調度模型進行求解，得到非劣解集中各調度方案的目標函數見表1和表2，同時采用FOA進行求解，兩個算法求得的調度方案非劣解前端見圖1和圖2，兩場洪水部分方案的泄流過程見圖3和圖4。

圖1 20200813洪水調度非劣方案集Fig.1 Non-inferior scheme set of 20200813 flood

圖2 20200826洪水調度非劣方案集Fig.2 Non-inferior scheme set of 20200826 flood

圖3 20200813洪水部分調度方案下泄流量過程Fig.3 Discharges of partial operation schemes of 20200813 flood

圖4 20200826洪水部分調度方案下泄流量過程Fig.4 Discharges of partial operation schemes of 20200826 flood

表1 20200813洪水調度非劣方案集Tab.1 Non-inferior scheme set of 20200813 flood

表2 20200826洪水調度非劣方案集Tab.2 Non-inferior scheme set of 20200826 flood

（1）由表1、2和圖1、2可見，由SA-FOA和FOA求得的調度方案集的最高水位和泄流量平方和近似成反比關系，這是因為兩者代表的目標存在沖突導致的。由圖1和圖2可見，在同一最高水位時，SA-FOA求得的泄流量平方和小于FOA，在同一泄流量平方和時，SA-FOA求得的最高水位小于FOA，由此可見，SA-FOA能夠比FOA求得更優的目標值。

（2）20200821洪水的實測最高水位和泄流量平方和分別為42.74 m和3 477 174，由表1可見，均大于SA-FOA求得方案集的目標值，比方案集中最高水位的最大值和最小值分別大0.17%和0.25%，比方案集中的泄流量平方和的最大值和最小值分別大14.87%和22.35%。20200826洪水的實測最高水位和泄流量平方和分別為42.90 m和1 133 554，由表2可見，同樣均大于SA-FOA求得方案集的目標值，比方案集中最高水位的最大值和最小值分別大0.27%和0.35%，比方案集中的泄流量平方和的最大值和最小值分別大5.63%和71.54%。綜上，SAFOA求得的泄流方案優于實測值。

（3）由圖3可見，20200821洪水的4個優化方案的最大泄流量均小于實測值，且比實測泄流過程提前加大了泄流，這是導致優化方案最高水位低于實測值的原因。由圖4可見，20200826洪水的4個優化方案的最大泄流量均大于實測值，然而其不超過日照水庫的控泄流量1 000 m3/s，故對下游防護區也是安全的。

3 結 論

（1）鑒于FOA存在過早收斂、易于陷入局部最優的不足，本文提出了SA和FOA耦合算法SA-FOA，該算法采用SA的Metropolis準則引入FOA，使其在迭代過程中依概率接受最優解，同時采用逐漸遞減的搜索距離L，并將該算法求解水庫防洪優化調度模型。

（2）由日照水庫兩場洪水調度的實例可以看出，SA-FOA求解的調度方案集的最高水位和泄流量平方和均小于實測值和FOA求得的調度方案集，由此表明SA-FOA求得的優化調度方案集既能夠保證水庫和庫區安全，又能提高下游防護區安全。綜上所述，本文建立的水庫防洪優化調度模型及SA-FOA能夠制定出滿意的洪水調度決策方案，為決策者提供技術支持。