線性求解器在中長期水量調度中的應用研究

趙亞威，方洪斌

（1.黃河勘測規劃設計研究院有限公司，河南 鄭州 450003；2.水利部黃河流域水治理與水安全重點實驗室（籌），河南 鄭州 450003）

1 引 言

中長期水量調度是在中長期時間尺度下（年計劃月調度或月計劃日調度），協調當地地表水、地下水、非常規水、外調水等多個水源供水量和用戶需水量的一項水資源優化配置工作，是實現水資源系統水量供需平衡和水資源系統綜合效益最大化的重要手段。隨著用水結構、供水來源的多樣化發展，區域水資源系統的復雜性不斷提升，使得中長期水量調度模型變量和約束增多、規模變大，求解難度隨之增加，尋求快速、精確的模型求解方法成為目前中長期水量調度的重要研究方向之一。

國內外許多學者對中長期水量調度模型求解算法進行了研究，主要分為模擬算法和優化算法。模擬算法［1－4］是一種經驗方法，主要依靠調度經驗將模型的一些復雜約束條件概化為對應的調度規則，進而降低模型的求解難度，實際應用較為廣泛，然而這種方法缺乏最優化理論支撐，獲取的解往往不是最優解；優化算法種類較多，主要分為以線性規劃、動態規劃等為代表的傳統優化算法［5－8］和以遺傳算法、粒子群算法等為代表的智能優化算法［9－10］。傳統優化算法比智能優化算法求解效率低，例如線性規劃在求解大規模復雜模型時計算量較大、計算時間較長，但線性規劃全局尋優的優良特性驅使國內外學者對提升其求解效率的研究不斷深入。隨著計算機硬件和軟件技術的不斷突破，計算機整體性能逐步提升，并行計算成為加速線性規劃計算效率的重要方式，例如目前比較流行的CPU并行、GPU并行可以實現多核、多線程并行計算，能夠大幅提升計算效率。針對線性規劃求解效率低的問題，許多學者和軟件廠商基于并行技術開發了線性規劃軟件。線性求解器是專門用于求解線性規劃問題的軟件或程序集（大部分求解器還可求解非線性規劃問題，如二次規劃），目前國際上流行的線性求解器主要有Cplex、Gurobi、Xpress、Lingo／Lindo等，免費開源求解器主要有Lpsolver、Scip、Glpk、Matlab等；國內線性求解器發展起步較晚，目前應用較為廣泛的主要有Cmip、Leaves等。經過大量實例驗證，線性求解器在求解大規模復雜線性模型時顯現出優秀的計算性能。

本文將線性求解器運用于中長期水量調度模型求解中，并將其與遺傳算法進行對比，以期探索模型的快速求解技術，為中長期水量調度方案編制提供參考。

2 中長期水量調度模型

2.1 水資源系統網絡拓撲

中長期水量調度模型一般基于網絡拓撲構建，網絡拓撲由節點、連線構成，節點表示輸入水量、輸出水量、蓄水量等，連線表示水流方向。對于水資源系統，節點可對應水庫、配水樞紐、控制斷面、水廠、用水戶等，其中水庫有蓄水庫容、配水樞紐有調節池，均具有一定調蓄能力，故蓄水量可變化；連線可對應河流、渠道、隧洞、管線等，在進行中長期水量調度時，一般將連線蓄水量設為恒定值或0。節點和連線遵循水量平衡原則，構成供需水網絡。節點、連線水量平衡示意圖分別見圖1、圖2。圖1中，R上游、R區間、R退水、R地表、R地下、R其他分別為節點i的上游相鄰節點出流水量、區間水量、退水量、當地地表水量、地下水量、其他水量（再生水、收集雨水、淡化水等），均為節點i的輸入水量；C下游、C退水、C損失、C生活、C工業、C農業、C生態分別為節點i的下游相鄰節點入流水量、退水量、損失水量、生活供水量、工業供水量、農業供水量、生態供水量，均為節點i的輸出水量；S時段初、S時段末分別為節點i時段初、末蓄水量。圖2中，R節點i為連線ij上游節點i的出流水量（輸入水量）；C損失、C節點j分別為連線ij的損失水量、下游節點j的入流水量（輸出水量）；連線ij的箭頭表示水流方向。

圖1 節點水量平衡示意

圖2 連線水量平衡示意

2.2 基于網絡拓撲的水量調度模型

網絡拓撲是描述水資源系統各單元水量供需關系的重要工具，依據網絡拓撲節點和連線的上下游關系、水量平衡約束、其他約束條件以及調度目標可構建中長期水量調度模型。

2.2.1 目標函數

以水資源系統總缺水量最小為優化目標，目標函數表達式為

式中：I為節點總數；J、j為節點上用水戶總數、序號，j取值為1、2、3、4，分別表示生活、工業、農業、生態用水戶；T為調度時段總數；Qijt、Cijt分別為節點i用戶j時段t的需水量和實際供水量。

2.2.2 約束條件

（1）節點水量平衡。表達式為

其中

式中：Ni、n分別為節點i上游相鄰節點總數、序號；Mi、m分別為節點i下游相鄰節點總數、序號；Rsyint、Rqint、Rtsint及ain、bin、fin分別為節點i時段t的上游第n個相鄰節點出流水量、區間入流水量、退水量及對應的沿程損耗系數；Cxyimt、Ctsimt、Cgsit、Cdcit分別為節點i時段t的下游第m個相鄰節點出流水量、退水量、總供水量、調出到外系統水量；Vi（t－1）、Vit分別為節點i時段t初、末蓄水量；cij為節點i用戶j的退水系數；dim為節點i與下游相鄰節點m之間的退水量占節點i總退水量的比例；ei為節點i上供水管網損耗系數；K、k分別為節點上水源數量、序號，其中k取值為1、2、3、4，分別表示當地地表水、地下水、非常規水和過境水（上游相鄰節點出流水量、區間水量、退水量之和）；Dijkt為時段t水源k向節點i用戶j的實際供水量。

（2）連線水量平衡。若網絡拓撲中節點i與節點p直接相連，水流由節點i流向節點p，則連線ip滿足如下水量平衡方程：

式中：Cxyivt為時段t從節點i流向節點p的出流水量；Rsyput為時段t由節點i流入節點p的水量；v為節點p在節點i下游相鄰節點中的序號，1≤v≤Mi；u為節點i在節點p上游相鄰節點中的序號，1≤u≤Ni。

（3）節點出流上限、下限約束。約束條件為

式中：、分別為節點i下游相鄰的第m個節點時段t的出流量下限、上限，為節點的過流能力，一般取0，若當前節點存在生態流量約束，則取生態流量。

（4）節點蓄水量上限、下限約束。約束條件為

（5）節點各水源向各用戶的實際供水量上限、下限約束。約束條件為

式中：、分別為時段t節點i上水源k向用戶j的實際供水量下限、上限。

3 基于線性求解器的模型求解

本文所構建的中長期水量調度模型中不含非線性表達式，為標準的線性規劃模型，可采用線性求解器進行求解。線性求解器的種類較多，但主要的數學原理、求解方法、操作流程有相似之處。數學原理方面，大多數求解器是基于單純形法、分支定界法等傳統的線性規劃方法；求解方法方面，大多數求解器是輸入標準化線性規劃模型的目標函數系數向量、約束條件系數矩陣和常數向量、目標優化方式Max／Min、決策變量個數、整數變量個數等信息，然后進行求解，最后得到最優決策向量和目標值；操作流程方面，使用線性求解器求解模型一般需要編寫程序，主要有集成開發環境（IDE）和基于應用程序編程接口（API）兩種方式。IDE的優點是程序編寫方式與數學模型的書寫習慣較為統一，程序結構簡單；缺點是不同求解器的程序語法往往不同，在使用多種求解器時需要熟悉多種IDE程序語法，另外IDE一般不能與其他語言進行程序混編，無法集成到C／C＋＋、Java、Python等其他應用程序開發系統中。API優點是支持C／C＋＋、Java、Python等當前主流編程語言，程序語法易于學習，調用方便，便于進行多語言混編，可用于大型應用系統的優化計算和模塊開發；缺點是很多求解器不提供API，有的僅支持部分主流編程語言，且API功能相較于IDE相對單薄。

本文采用C＃與Cplex、Lpsolver線性求解器API混合編程實現模型求解，步驟如下。

（1）模型前處理。將模型進行標準化處理，包括將目標函數轉化為“Min”形式以及將所有不等式約束條件轉化為“≤”形式。

（2）模型關鍵參數提取。提取中長期水量調度數學模型中的關鍵信息，包括目標函數中決策變量系數向量C、約束條件中決策變量系數矩陣A和常數向量b、決策變量上下限、決策變量為整數或非整數等信息，形成矩陣數據或向量數據。

（3）模型求解。將步驟（2）提取的所有數學模型矩陣數據和向量數據，整理成Cplex或Lpsolver的API標準化輸入形式，調用求解器中的線性規劃方法對模型進行求解。

（4）結果輸出。按照Cplex或Lpsolver輸出的規則，將計算得到的最優結果解析為方便進行方案分析的數據格式。

小兒柴桂退熱顆粒的UPLC指紋圖譜及聚類、主成分分析…………………………………………………… 林 源等（4）：474

4 實例分析

以北方某城市2025年規劃水平年中長期水量調度為例進行分析。該市境內主要河流有5條，分別命名為河流A、B、C、D、E。該市現狀供水水源主要由當地地表水、地下水、非常規水和從A河引水組成。2019年該市從A河引水量8.7億m3，占該市總用水量的38.3%。從A河引水進入該市水資源系統后，主要由位于河流B的水庫1和位于河流D的水庫2進行調節，供給該市生活、工業、農業、生態環境需水。該市水資源系統概化見圖3，網絡拓撲概化見圖4。采用Cplex、Lpsolver、GA（遺傳算法）對模型進行求解。

圖3 某北方城市水資源系統概化

圖4 某北方城市水資源系統網絡拓撲概化

4.1 不同方案優化結果分析

該市每年各月按照一定比例（見表1）從河流A的控制水庫引水，按照河流A來水量豐枯變化，分別設置從A河年引水量4億、5億、6億、7億、8億m35種方案，優化結果見表2。

表1 某北方城市各月從河流A的引水量比例

表2 從A河引水不同方案優化結果

從優化結果來看，各方案Cplex、Lpsolver優化得到的缺水量和缺水率均小于GA的優化結果，原因是Cplex、Lpsolver是基于線性規劃的全局優化算法，而GA則存在易陷入局部最優的缺點，無法保證收斂到全局最優解。從計算時間來看，采用Cplex、Lpsolver的各方案計算時間均小于0.1 s，而采用GA的計算時間最短需978.89 s，最長則需1 188.83 s，Cplex、Lpsolver計算時間遠遠少于GA，原因是Cplex、Lpsolver采用了并行計算，充分發揮了計算機性能，極大提升了計算效率，而GA不僅存在繁瑣的選擇、交叉、變異和迭代等過程，且未采用并行計算，計算效率較低。

4.2 不同規模模型優化結果分析

為了進一步檢驗不同算法的性能，改變模型規模進行對比分析。以圖4模型19個節點（大規模）為基礎，將節點數精簡為8個（小規模），分別采用Cplex、Lpsolver、GA進行優化求解，結果見表3。

表3 小規模模型從A河引水不同方案優化結果

通過對比表2、表3可知，由于Cplex、Lpsolver求解效率較高，因此模型規模大小對其計算時間影響不大；小規模模型和大規模模型不同A河引水方案下，GA的平均計算時長分別為111.70、1 107.84 s，GA限于自身算法的復雜性和未采用并行機制，隨著模型規模的增大，計算時間增加。

4.3 從A河引水時機分析

由表2、表3可知，從A河引水8億m3時，仍不能滿足總水量需求，原因是各月從A河引水比例與該市需水比例變化趨勢不完全一致，而該市境內主要依靠水庫1、水庫2來調節從A河引水分配過程，這兩座水庫興利庫容分別為0.651億m3和0.684億m3，相較于從A河引水量4億～8億m3來說均很小，調蓄能力非常有限，無法進行時空調節。為解決這一問題，可以從增加水庫興利庫容或者改變引水時機兩方面入手，但該市近期沒有新建水庫或增加現有水庫興利庫容的規劃，因此只能從改變引水時機角度考慮改善缺水狀況。

在實際調度中，該市年內各月從A河的引水量為固定值，可表示為

式中：qt、q0t分別為第t月該市從A河的引水量、固定引水量。

對原模型進行改進，將“各月從A河引水量為固定值”等式約束修改為“各月從A河引水量之和為固定值”等式約束，計算公式為

式中：Q為從A河年引水總量。Q與q0t的關系為

以16個節點模型為例，以從A河引水量4億m3為起點、8億m3為終點，等值離散為若干個方案，采用Cplex算法對不同引水量下不考慮引水時機（原模型）和考慮引水時機（改進模型）兩種模型進行求解，得到從A河引水量與缺水量的關系，見圖5。

圖5 從A河引水量與缺水量的關系

（1）不同引水方案下，改進模型的缺水量均小于原模型，且隨著從A河引水量的增加，改進模型較原模型缺水量逐漸減小，說明改變引水時機能夠減小從A河引水量，節約系統供水成本。

（2）當從A河引水量為7.65億m3時，改進模型的缺水量為0，該引水量下水庫1、水庫2的水位變化過程見圖6和圖7，可以看出，原模型水庫1僅在2月、9月發揮了調蓄作用，水庫2僅在5月、8月、9月發揮了調蓄作用，而改進模型的水庫1和水庫2的蓄泄變化過程更加明顯，說明改進模型充分發揮了水庫的調蓄作用。原模型、改進模型從A河引水過程見圖8，可知：7月、9—12月當地水量較富裕，改進模型比原模型減小了引水量，而在1—4月、6月加大了從A河引水量，彌補需水缺口；5月需水量較大，而改進模型從A河引水量小于原模型，原因是改進模型1—4月從A河引水量較大，增加了水庫1、水庫2的蓄水量（5月初水庫1、水庫2蓄水量較大），5月從A河引水，水庫1和水庫2同時向用戶供水，故改進模型5月從A河引水量小于原模型；8月當地水量非常富裕，原模型從A河引水量較小，而改進模型引水量較大，原因是8月當地水主要來自節點N1002、N1004、N1006、N2001、N2002、N3002、N4001、N5002、N5003，除N2001（水 庫1）、N4001（水庫2）外，其他節點當地水量（為3.35億m3，占8月當地水總量的92.76%）均無法存入水庫1和水庫2，不能被調蓄，因此需要加大從A河引水量，來滿足9—12月用戶需求。從圖6、圖7可以看出，8月末水庫1和水庫2蓄水量明顯增加，且9—12月為開閘供水狀態，水庫水位不斷降低，直至水庫放空。綜上所述，與原模型相比，改進模型充分發揮了當地水庫的調蓄作用。

圖6 水庫1水位變化過程

圖7 水庫2水位變化過程

圖8 原模型、改進模型從A河引水過程

5 結 語

本文以北方某城市水資源系統為例，建立中長期水量調度模型，采用線性求解器對模型進行求解，并與遺傳算法進行了對比。主要結論如下：①線性求解器Cplex、Lpsolver編程簡單，能夠獲得全局最優解，求解精度高、速度快；②當模型規模增大時，遺傳算法計算效率明顯下降，而Cplex、Lpsolver計算效率變化不明顯；③通過改變引水時機能夠充分發揮現有水庫的調蓄能力，增加系統總供水量。

線性求解器應用的前提是模型必須為線性模型，因此如何實現含非線性約束的中長期水量調度模型的線性化處理，是下一步的研究重點。