高 建，張運鑫，2

（1.河北工程大學水利水電學院，河北 邯鄲 056038；2.河北省水資源高效利用工程技術(shù)研究中心，河北 邯鄲 056038）

0 引 言

中國是世界上13 個貧水國家之一。2021 年我國用水總量為5 920.2 億m3，其中農(nóng)業(yè)用水量3 644.3 億m3，占總用水量的61.5%[1]。水資源短缺問題日益成為制約我國農(nóng)業(yè)和經(jīng)濟社會發(fā)展的重要因素。因此，優(yōu)化灌區(qū)輸配水渠系運行管理對提升灌區(qū)渠道輸配水效率，改善灌區(qū)灌溉供用水狀況具有重要意義。

明渠是被廣泛應用于灌區(qū)的輸配水方式[2]。灌區(qū)輸配水渠系一般由干、支、斗、農(nóng)、毛五級渠道構(gòu)成，其中支、斗兩級渠道是灌區(qū)渠系配水由續(xù)灌轉(zhuǎn)變?yōu)檩喒嗟年P(guān)鍵銜接部分，對實現(xiàn)灌區(qū)渠系優(yōu)化配水和提高渠系水利用系數(shù)方面具有重要作用。

有學者對支、斗渠合理配水調(diào)度進行了研究，Suryavanshi和Reddy提出將輪灌分組問題轉(zhuǎn)化為數(shù)學問題，引入了流管概念，將0-1規(guī)劃模型中的二元決策變量Xij看作第i個流管的第j個分水口，它的值為0或1，0表示該分水口關(guān)閉，1表示該分水口打開[3]，通過這種數(shù)學規(guī)劃的方式獲得最優(yōu)配水調(diào)度方案。呂宏興等改進了Suryavanshi 的模型，通過進行引水時間的均一化處理，解決了渠道配水時需要多次調(diào)節(jié)進水閘的問題[4]。

目前，大多數(shù)灌區(qū)仍憑借以往管理經(jīng)驗來確定灌溉渠系的水量和流量分配，即采用經(jīng)驗配水法。盡管采用經(jīng)驗配水方法可以進行灌區(qū)配水，但經(jīng)常出現(xiàn)渠系輸配水時間長、小流量配水以及水資源浪費等問題。經(jīng)驗配水方法無法實現(xiàn)精準配水，也難以對供需變化做出及時響應[5]。隨著水資源供需矛盾日益加劇，對提升灌區(qū)供需水管理水平的要求也越來越高[6]。為了實現(xiàn)渠系的高效優(yōu)化配水，學者們開始研究將智能優(yōu)化算法引入到渠系配水中，如遺傳算法[6-12]、粒子群算法[12-15]、模擬退火算法[16]、灰狼優(yōu)化算法[17]等，通過對這些算法應用獲得渠系最優(yōu)調(diào)度方案。馬孝義等建立了下級渠道引水流量不相等情況下的渠道優(yōu)化配水模型，通過遺傳算法獲得優(yōu)化配水方案，輸水損失更低[9]；劉照等以渠道輸水損失最小與輪灌組間引水持續(xù)時間差異值最小為目標建立多目標優(yōu)化配水模型，通過雙層粒子群算法求解優(yōu)化方案，不但可減少渠道棄水量，同時統(tǒng)一了輪灌組閘門關(guān)閉時間，減少了閘門啟閉次數(shù)[13]；劉葉等以渠系輸水損失最小為目標構(gòu)建兩級灌溉渠系優(yōu)化配水模型，通過遺傳算法確定的方案，不但能滿足灌溉要求，同時減少渠系配水時間以及閘門啟閉次數(shù)[10]；徐淑琴等以灌區(qū)渠道輸水損失最小及干渠水位變化均勻為優(yōu)化目標，運用帶精英策略非支配排序遺傳算法確定運行方案，在滿足優(yōu)化目標的前提下，可以縮短各生育期的配水時間，且流量變化幅度較小[11]。

針對優(yōu)化離散空間約束問題，Kennedy 和Eberhart 在傳統(tǒng)粒子群算法(PSO)基礎上改進并提出離散二進制粒子群算法(BPSO)[18,19]。本文在已有研究基礎上，綜合考慮離散二進制粒子群算法(BPSO)和遺傳算法(GA)的優(yōu)缺點[20]，研究提出了混合二進制粒子群算法（GA-BPSO），用于灌區(qū)支、斗兩級渠系優(yōu)化配水研究。用二元決策變量Xij描述每一個支渠向斗渠配水的分水口配水運行狀態(tài)，其中賦值0表示該分水口關(guān)閉；賦值為1 則表示該分水口開啟。應用MATLAB 對BPSO 算法和GABPSO算法進行模擬計算，并通過應用案例進行檢驗分析。

1 支斗渠系優(yōu)化配水模型

當支渠向斗渠配水時，往往將斗渠分為幾個組，按照“組間續(xù)灌，組內(nèi)輪灌”的模式運行，支渠向斗渠配水灌溉如圖1所示。

圖1 支渠向斗渠配水示意圖Fig.1 Irrigation schematic from branch canal to lateral canal

設支渠的凈流量為Q凈，各斗渠的流量為q，斗渠渠道水利用系數(shù)為η，并用向上取整公式ceil 和向下取整公式floor 進行取整，灌溉編組數(shù)M計算公式如下：

本研究提出的模型以支渠向斗渠輸配水運行時間最短作為目標，模型的決策變量主要包括灌溉組數(shù)、分水口數(shù)和每個分水口的連續(xù)引水時間。約束條件包括引水流量和任意分水口開啟次數(shù)等，優(yōu)化模型如下：

式中：Tb表示支渠運行時間，s；Tl表示M 個灌溉斗渠編組中運行時間最長的灌溉斗渠編組運行時間，s；M 為灌溉編組個數(shù)；N 為被配水斗渠條數(shù)；ti為第i 組灌溉斗渠編組運行時間，s，其中i=1,2,…,M；tj表示第j條斗渠配水運行時間，s，其中j=1,2,…,N；Sj表示第j 條斗渠所需的灌水量，m3；qj表示第j 條斗渠配水流量，m3/s；ηj表示第j條斗渠渠道水利效率；Xij表示第i個灌溉編組中第j個分水口的狀態(tài)，Xij=0 為i個灌溉編組中第j個分水口是關(guān)閉狀態(tài)，即第i個灌溉編組內(nèi)的第j個分水口沒有配水計劃；Xij=1 為i 個灌溉編組中第j 個分水口是打開狀態(tài)，即第i個灌溉編組內(nèi)的第j個分水口有配水計劃。

約束條件如下：

設第i 灌溉編組供水時間為ti，則其總供水時間不應大于灌溉周期T，即：

任意分水口在灌溉周期T內(nèi)只能開啟一次，即：

決策變量取值約束：

各灌溉斗渠編組凈流量之和不大于對應支渠來水凈流量，即：

式中：Qi表示第i個灌溉斗渠編組的凈流量，m3/s。

2 混合二進制粒子群算法

2.1 粒子群算法原理（PSO）

傳統(tǒng)粒子群算法（Particle Swarm Optimization，PSO）[18]是在連續(xù)區(qū)域內(nèi)對所求問題進行優(yōu)化求解的一種算法，粒子群中的每一個粒子個體均表示一種可行的解決方案，其算法原理是通過粒子之間的簡單行為，經(jīng)過交互得出最優(yōu)解，具有較快的收斂效率和局部搜索能力但容易陷入局部最優(yōu)解的特點。

粒子群算法的基本求解計算表達式如下：

式中：v(i)(t)表示第i 個粒子在第t 次的移動速度；x(i)(t)表示第i個粒子在第t次粒子所在的位置；x(i)(t+1)表示第i個粒子在第t+1次所在的位置；c1、c2為粒子的學習因子；rand1、rand2為分布在[0,1]的隨機數(shù)；p(i)best為第i個粒子歷史最優(yōu)解；gbest為當前歷史最優(yōu)解。

進一步分析可將式(10)中第i個粒子在第t+1次的移動的速度可分解為三部分：w v(i)(t)為第一部分，稱為“慣性”部分，這一部分反映了粒子原有的運動趨勢；c1rand1[p(i)best-x(i)(t)]為第二部分，稱為“個人認知”部分，這一部分就反映了粒子向自身粒子最優(yōu)解移動的趨勢；c2rand2[gbest- x(i)(t)]為第三部分，稱為“社會認知”部分，這一部分反映了粒子向全部粒子中的最優(yōu)解移動的趨勢。粒子移動過程圖和粒子群算法基本流程圖分別如圖2和圖3所示。

圖2 PSO流程圖Fig.2 Flow chart of PSO

圖3 PSO粒子移動過程圖Fig.3 PSO particle movement process diagram

典型粒子群算法步驟如下：

（1）按照需求對粒子進行初始化，包括粒子速度和粒子位置；

（2）構(gòu)建適應度函數(shù)，對每個已經(jīng)生成的粒子按照適應度函數(shù)計算相應的適應度；

（3）根據(jù)已經(jīng)計算出的適應度評價出個體最佳粒子以及全局最佳粒子；

（4）根據(jù)公式（10）和公式（11）更新每個粒子的速度和位置；

（5）重新計算更新過后的每個粒子的適應度；

（6）重新計算過后的適應度與個體最佳粒子和全局最佳粒子進行比較，重新選出個體最佳粒子和全局最佳粒子；

（7）依照結(jié)束條件進行判斷，滿足結(jié)束條件結(jié)束算法；若不滿足則重復執(zhí)行上述步驟4）、5）、6）直到滿足條件。

2.2 離散二進制粒子群算法原理（BPSO）

Kennedy 和Eberhart 為解決離散型或二進制類型的問題從而改進了常規(guī)的粒子群算法，提出了離散二進制粒子群算法（Discrete Binary Particle Swarm Optimization Algorithm，BPSO）[19]。該算法是采用二進制編碼0或1對粒子中某一維度的狀態(tài)進行表示，速度更新公式同基本粒子群算法的速度更新公式相同，但沒有使用粒子群算法位置更新公式，取而代之使用sigmoid函數(shù)將速度值映射到[0,1]區(qū)間，以此表示取0或1 的概率。通過rand 隨機數(shù)與sigmoid 函數(shù)映射之后的值進行比較，得到位置的具體值[21]。

下式分別為sigmoid函數(shù)以及位置更新公式：

2.3 遺傳算法（GA）

遺傳算法（Genetic Algorithm，GA）借鑒了達爾文的進化論和孟德爾的遺傳學說，是通過模擬生物在自然環(huán)境中的遺傳和進化的過程而形成的一種自適應全局優(yōu)化搜索算法，具有全局搜索能力強，局部搜索能力較弱的特點[20]。圖4 為遺傳算法流程圖。

圖4 GA流程圖Fig.4 Flow chart of GA

其算法流程大致如下：

（1）按照編碼需求（二進制編碼，實數(shù)編碼等）對種群進行初始化；

（2）構(gòu)建適應度函數(shù)，對每個已經(jīng)生成的個體按照適應度函數(shù)計算相應的適應度；

（3）根據(jù)需求對已經(jīng)生成的個體進行選擇，一般選擇適應度較高的個體為父代；

（4）按照設定好的交叉概率及交叉規(guī)則，將選擇出來的父代進行隨機配對，進行交叉運算，生成新的個體；

（5）按照設定好的變異概率及變異規(guī)則，將選擇出來的父代進行變異操作，生成新的個體；

（6）重新計算適應度并更新種群；

（7）依照結(jié)束條件進行判斷，滿足結(jié)束條件結(jié)束算法；若不滿足重復執(zhí)行上述步驟（3）、（4）、（5）、（6）直到滿足條件。

2.4 混合二進制粒子群算法（GA-BPSO）

根據(jù)GA 和BPSO 的特點，在BPSO 的基礎上融入GA 的方法，融合兩種算法的優(yōu)點并彌補相應不足，得到一種整體性能更優(yōu)的混合二進制粒子群算法(GA-BPSO)。圖5為GA-BPSO基本流程圖。

圖5 GA-BPSO流程圖Fig.5 Flow chart of GA-BPSO

其算法流程大致如下：

（1）按照需求對粒子進行初始化，包括粒子速度和粒子位置；

（2）對每個已經(jīng)生成的粒子按照適應度函數(shù)計算相應的適應度；

（3）根據(jù)已經(jīng)計算出的適應度評價出個體最佳粒子以及全局最佳粒子；

（4）根據(jù)公式（10）更粒子新速度，根據(jù)公式（12）、（13）更新粒子的位置；

（5）計算更新后的粒子適應度；

（6）重新選出個體最佳粒子和全局最佳粒子；

（7）按照需求選擇出一定數(shù)量的個體作為父代，并進行交叉的相關(guān)操作；

（8）按照需求選擇出一定數(shù)量的個體作為父代，并進行變異的相關(guān)操作；

（9）對新生成粒子進行適應度計算，根據(jù)結(jié)果更新粒子并重新選出個體最佳粒子和全局最佳粒子；

（10）依照結(jié)束條件進行判斷，滿足結(jié)束條件結(jié)束算法；若不滿不足重復執(zhí)行步驟（4）、（5）、（6）、（7）、（8）、（9）直到滿足條件。

3 實例應用

某灌區(qū)一條支渠下有7 條斗渠，支渠設計流量2.5 m3/s，渠道水利用系數(shù)為0.92，各斗渠設計流量、控制面積以及渠道水利用系數(shù)如表1 所示。支渠灌溉周期為7 d，本次灌水定額確定為900 m3/hm2，優(yōu)化配水目標為渠系輸配水時間最短。

表1 各斗渠參數(shù)Tab.1 Parameters of each lateral canal

3.1 案例分析

根據(jù)該灌區(qū)支渠、斗渠相關(guān)信息以及公式(1)確定灌溉編組數(shù)M=4 最佳。根據(jù)相關(guān)約束可知每條斗渠只能在一個灌溉編組中運行，結(jié)合算法可將每一條斗渠看作一個1×4 的列向量，均由0 或1 組成。其中0 代表該斗渠在某灌溉編組內(nèi)的配水口為關(guān)閉狀態(tài)，即該灌溉編組內(nèi)該支渠沒有配水計劃；1代表該斗渠在某灌溉編組內(nèi)的配水口為打開狀態(tài)，即該灌溉編組內(nèi)該斗渠有配水計劃。如[1,0,0,0]T表示該斗渠在第一個灌溉組中運行，而不在其余灌溉組中運行，示意如圖6(a)所示。

圖6 算法分析圖Fig.6 Algorithm analysis graph

按照上述設想，如果將灌溉編組和渠道成生一個只包含0或1 的二維數(shù)組，初始化各灌溉斗渠編組分水口狀態(tài)如表2 所示，則當前配水狀態(tài)為7條斗渠均在第一個灌溉編組中，其余灌溉編組無斗渠配水計劃，按此設定方式會降低算法的運行效率。所以對BPSO 的粒子初始化做出改變，在生成粒子時，生成一個7 維的向量，每一個維度可看成每一條斗渠變化狀態(tài)，其中：假設這個粒子初始化時為[0,0,0,0,0,0,0]，當粒子某一維度發(fā)生變化時則對應斗渠分水口在灌溉編組中工作狀態(tài)就會發(fā)生變，從第i個灌溉編組變換到第i+1個，具體如圖6(b)所示。當粒子由[0,0,1,0,0,0,0]變?yōu)閇0,0,0,0,0,0,0]則第三條斗渠會從第三個灌溉編中有配水計劃變?yōu)樵诘谒膫€灌溉編組中有配水計劃。當某斗渠在最后一個灌溉編組中有配水計劃時下一次變換就會變換到在第一個灌溉編組中有配水計劃；如果粒子某一維度未發(fā)生改變則對應的斗渠分水口在灌溉編組中的狀態(tài)就保持不變。

表2 初始化各灌溉斗渠編組分水口狀態(tài)Tab.2 Initialize the status of each irrigation lateral grouping manifold

3.2 算法應用

本文在PC 端使用MATLAB 軟件求解算法。將渠系配水時間作為適應度函數(shù)值，最佳配水方案將具有最小的適應度函數(shù)值。算法參數(shù)確定：種群Z=20，迭代次數(shù)N=100，慣性權(quán)重W=0.9，個體學習因子C1=2，社會學習因子C2=2，交叉概率Pc=0.7，變異概率Pm=0.05。按照公式（1）確定灌溉斗渠分組數(shù)M=4，根據(jù)各灌溉斗渠編組分水口狀態(tài)，按照公式（4）計算出每個灌溉斗渠分組運行的時間ti，再根據(jù)公式（3）得出所有ti運行時間最長的灌溉斗渠分組運行的時間Tl，將Tl作為算法的適應度函數(shù)值，使其盡可能小以達到算法求解要求。

3.3 結(jié)果分析

針對上述兩級渠系輸配水問題案例，可得到如下分析：

（1）兩種算法得到的配水編組的方式都可以將渠系配水時間降到最低，具體方案如表3和表4所示。

表3 BPSO和GA-BPSO計算最優(yōu)配水方案結(jié)果（1）Tab.3 Results of BPSO and GA-BPSO calculations of optimal water distribution schemes (1)

表4 BPSO和GA-BPSO計算最優(yōu)配水方案結(jié)果（2）Tab.4 Results of BPSO and GA-BPSO calculations of optimal water distribution schemes (2)

（2）通過表5 以及圖7 可知，在參數(shù)相同的情況下BPSO在大概21 代左右能得到最優(yōu)解，而GA-BPSO 則在大概12 代左右得到最優(yōu)解，且兩者最優(yōu)解相同?？梢?，針對BPSO 的不足之處結(jié)合GA 之后的新算法GA-BPSO 的整體性能優(yōu)于BPSO。

表5 BPSO和GA-BPSO迭代終止時的性能對比Tab.5 Performance comparison between BPSO and GA-BPSO at iterative termination

圖7 算法結(jié)果對比Fig.7 Comparison of algorithm results

（3）由于在BPSO 中包含速度和位置兩個變量，在與GA結(jié)合時交叉、變異操作會導致速度、位置均發(fā)生相應改變，而速度變化會影響下一次迭代時位置發(fā)生改變的概率，所以在模擬過程中會出現(xiàn)個別GA-BPSO效果不如BPSO的情況。

4 結(jié) 論

本文針對灌區(qū)支斗兩級渠系輸配水優(yōu)化問題，根據(jù)BPSO和GA的優(yōu)缺點，在BPSO基礎上融合GA的思想，構(gòu)建一種更為高效的GA-BPSO。通過GA-BPSO 對案例問題進行優(yōu)化，模擬結(jié)果表明GA-BPSO 達到最優(yōu)解的平均收斂代數(shù)要比BPSO達到最優(yōu)解的平均收斂代數(shù)要低，說明結(jié)合后的GA-BPSO 的整體性能相比BPSO 要更好。證實了支斗兩級渠系在“組間續(xù)灌，組內(nèi)輪灌”的輸配水模式下GA-BPSO 具有更好的可行性和高效性，可為灌區(qū)用水調(diào)度提供一種較優(yōu)的解決方案。由于遺傳操作對速度的改變可能會對下一次迭代時位置的改變產(chǎn)生影響，在未來算法優(yōu)化的方向上可以進一步研究GABPSO 運行過程中的速度更新方式及更高效的交叉、變異方式。