基于區間兩階段魯棒優化模型的灌區水資源優化配置

陳紅光 王瓊雅 李曉寧 王中君 李晨陽

(1.東北農業大學水利與建筑學院， 哈爾濱 150030； 2.東北農業大學工程學院， 哈爾濱 150030)

0 引言

隨著社會經濟和人口的快速增長，水資源短缺問題日益突出。隨著糧食需求不斷增加，農業作為用水大戶，水資源供需問題不斷加重[1-2]。合理配置農業用水、提高農業用水效率對緩解水資源供需矛盾和用水壓力具有十分重要的意義[3-5]。灌區作為農業用水主體，其水資源的合理優化配置對提高水資源利用率、促進農業水資源的可持續開發利用、權衡灌區水資源經濟和風險之間的關系具有重要作用。

近年來，水資源優化配置得到了廣泛發展，出現了多種優化配置方法。婁帥等[6]基于免疫遺傳算法構建多階段群決策優化模型，解決漳河流域水資源優化配置問題。孫冬營等[7]運用模糊聯盟合作博弈方法將水資源進行兩次分配，得到流域水資源的合理配置方案。陳述等[8]將粒子群人工蜂群算法運用到灌區渠-塘-田的水資源優化配置中。上述研究基于確定的水資源系統進行配置，但在實際中，灌區灌溉系統受到氣象條件、作物需水量以及政策變化等因素影響，是一個復雜的動態配水系統[9]。傳統方法在處理系統中的不確定因素時，存在一定的局限性。因此，為了解決這些包含不確定因素的問題，提出了不確定優化方法，其中兩階段隨機規劃模型得到了廣泛應用[10]。HUANG等[11]將不確定性優化和兩階段隨機規劃方法應用到水資源系統管理中；LI等[12]將區間兩階段機會約束模型運用到水資源規劃中；ZHANG等[13-14]運用改進后的區間兩階段隨機規劃算法，對三江平原水資源進行優化配置。區間兩階段優化算法在大量學者[15-17]的努力下，逐漸得到完善，用以處理水資源配置過程中的不確定性問題。區間兩階段隨機效用模型處理不確定因素十分有效，但卻忽略了系統風險問題，模型結果不具有絕對可行性。而魯棒優化方法則能夠在規劃過程中對風險進行有效規避，并權衡系統中可變隨機值和追索成本的關系[18-20]。

本文針對灌區水資源配置模型中存在風險的特點，構建多水源、多作物配水模型，以兩階段線性規劃為基礎，建立區間兩階段魯棒優化模型(Interval-parameter two-stage robust stochastic programming model，ITRM)，采用概率密度函數、離散區間表示系統的不確定性，用魯棒系數表示系統的風險。在不同來水情況下，對各個灌區的不同作物進行水資源的優化配置，以期為灌區管理者提供風險可控、成本最優的配水方案。

1 模型建立與求解

1.1 模型建立

灌區多水源優化配置是一個復雜的系統過程，供水情況受來水量等因素影響，存在一定的隨機特點[21]。且多水源聯合優化調度目的是使水資源得到合理的分配，滿足灌區用水需求，同時盡量減少用水成本。因此，本文以灌區農作物需水量為決策變量，引入用水成本、缺水懲罰系數，并分兩個階段確定灌區水資源的最優配置。在保證灌區水資源承載力的前提下，以正常水平年作物需水量為依據，確定作物預先配水目標，并將其作為第1階段的決策變量；由于受灌區初始儲水量、天然來水量、水源蒸發量等因素的影響，供水量可能小于供水目標，這就需要決策者進行調整，減少灌區供水量或者外調水源進行補水，而減少供水量會影響作物產量，外調水源則會增加用水成本，都會產生經濟懲罰，為降低用水成本，需要對第1階段的配水量進行調整，將缺水量作為第2階段的決策因素。故兩階段隨機規劃模型為

(1)

式中Cij——水源i向作物j的輸水成本，元/m3

Wij——水源i向作物j的預先目標配水量，m3

E——隨機變量的期望值

Dij——水源i未滿足作物j預先目標配水量時的缺水懲罰系數，元/m3

Sij——水源i未滿足作物j預先目標配水量時的缺水量，m3

i——水源，i為1、2分別表示地表水和地下水

j——作物，j為1、2、3分別表示水稻、玉米、大豆

f——系統總成本，元

(2)

式中Sijh——來水量水平為h時，i水源向j作物供水，未滿足預先目標配水量的缺水量，m3

已有的研究較少考慮系統風險問題，無法保證模型最優解的絕對可行性。魯棒優化方法將風險體現在函數當中，對風險進行有效測評，并在規劃過程中規避風險，平衡系統的成本和風險的關系，可以有效增加模型求解的可行性和系統的穩定性[22]。因此，本文在兩階段隨機規劃模型的基礎上引入魯棒優化方法，建立兩階段隨機魯棒優化模型。模型表述為

(3)

式中ρ——魯棒系數

(4)

式中θih——松弛變量

θih主要作用是保證模型的穩定性和非負性。作物預先目標配水量(需水量)Wij具有不確定性，且作物的價格和產量變動導致懲罰系數Dij也不確定。為表示不確定性，引入區間參數來表示不確定性參數。“+”表示參數的上限，“-”表示參數的下限，建立區間兩階段魯棒優化模型(ITRM)

(5)

(1)需水量約束

(6)

式中Wijmax——作物j正常生長的最大需水量

Wijmin——作物j正常生長的最小需水量

(2)水源可用水量約束

地表水約束為

(7)

地下水約束為

(8)

式中Qij——灌區初期水源i的儲水量，m3

Qim——灌區末期水源i的蓄水量，m3

不錯，《罕哈冉惠傳》中出現了不少反映佛教思想的內容。如學者們指出的，當哈冉惠和他的兩個兄弟遭蟒古斯暗害，誤食了蟒古斯投下的毒藥后，哈冉惠變成了長有九十五顆頭的大黑蟒古斯，他的胞弟烏蘭岱·莫日根變成了一尊石人，烏蘭岱的坐騎變成了一尊石馬，他的義兄弟吉爾吉斯·賽因·貝托爾則變成了一頭黃色的野豬。是兩位公主派出了維蘭·索龍嘎的兒子去向菩薩求救，菩薩將三件寶物——萬能的金套繩、能起死回生的仙丹、智慧的金盤賜給了小童，小童用此三寶使三位英雄得救，恢復了原貌。菩薩，為眾人所知，當然屬于重要的佛的形象，然而當菩薩向小童授予三件寶物之前，卻說出下面一段頗令人費解的話：

qih——水源i的來水量，m3

Qsi——灌區水源i的蒸發量，m3

Qimin——水源i最低蓄水量，m3

(3)追索變量約束

追索變量約束為

(9)

(4)水源最大供水量約束

水源最大供水量約束為

(10)

式中Wimax——水源i最大可供水量，m3

(5)非負約束

非負約束為

(11)

1.2 模型求解

將ITRM分為兩個子模型進行求解，對應的目標函數下限子模型為

(12)

約束條件為

(13)

(14)

約束條件為

(15)

(16)

(17)

zij=zijopt(?i,j)

(18)

最優配水量為

(19)

式中Aij——水源i向作物j的最優配水量

其算法流程如圖1所示。

圖1 區間兩階段魯棒優化方法算法流程圖Fig.1 Algorithm flow chart of interval-parameter two-stage robust stochastic programming method

2 案例研究

2.1 區域概況

三江平原地區牡丹江灌區位于黑龍江省東南部(131°13′～133°37′E，44°57′～46°12′N)，耕地面積為3.07×109m2，所轄12個現代化農場。灌區內土壤肥沃，水資源豐富，主要農作物為水稻、大豆、玉米，以井灌為主。然而，近些年隨著作物種植面積的不斷增加，作物需水量增大，導致灌區用水量壓力增大，水資源危機加重。實際配置水資源量如果小于農作物的需水量,供水不足將造成農作物減產，增大懲罰成本，若滿足農作物需水量，則又增加用水成本和系統供水風險。因此，本文運用ITRM有效的平衡系統成本和風險，解決灌區水資源分配問題。

2.2 參數確定

選取水稻、大豆、玉米3種作物為研究對象，根據《牡丹江統計年鑒》[24]、《黑龍江墾區統計年鑒》[25]等相關資料可以得到3種作物的灌區種植面積、單位面積灌溉用水區間值，作物灌溉面積比例參照文獻[27]，水稻、大豆、玉米分別為100%、10%、10%，根據文獻[24-25]和當地多年統計數據，得到作物單位面積最小灌溉量和最大灌溉量(表1)。以灌區2011年3種作物種植面積數據為已知條件，將作物單位面積最小(最大)用水量乘以作物種植面積、作物灌溉面積比例，得到作物最小(最大)需水量，查閱文獻[24-25]得到灌區水量分配的相關經濟系數(表2)。根據《黑龍江省農墾水利志》[26](1987—2015年)對灌區多年天然降水量的統計結果，可知1992、2002、2008、2011年灌區的降水量為400～450 mm，1987、1990、2009、2013年降水量為600 mm以上，其他年份降水量在450～600 mm之間。把來水量水平分為高、中、低3個水平，由統計結果可知，中流量年份比高降流量和低流量年份多，且高流量和低流量出現概率相近，故假設預測年份的降水量水平高、中、低出現的概率為0.2、0.6、0.2。查閱文獻[27]獲得地下水和地表水的供水比例為4∶1。其中地下水的主要補給來源為側向徑流補給[28]，2011年全年牡丹江地區地表水過境水量為2.482×1010m3，但其提水量2.3×107m3，僅占過境水量的0.093%。根據《牡丹江統計年鑒》[24]、《黑龍江墾區統計年鑒》[25]以及多年數據統計分析，得到灌區在不同來水量水平下的可供水量區間值(表3)，灌區地表水最大可供水量為6.2×108m3，最小蓄水量為1.1×108m3，地下水最大可供水量為9.7×108m3，最小儲蓄量為2.55×108m3。

表1 灌區作物種植面積及單位面積灌溉用水量Tab.1 Irrigation crops area and water per unit in irrigation district

表2 灌區作物最大最小需水量、預先配水量和經濟系數Tab.2 Irrigation corps maximum and minimum water distribution, target value of crops water distribution in advance and economic parameters

表3 不同來水量水平下農作物灌溉水源可用水量Tab.3 Water available for crops irrigation under different water inflow levels

2.3 模型計算結果與分析

圖2為灌區用于水稻、玉米、大豆的最優目標配水量，其中地表水最優目標配水量為6.39×107m3，地下水為2.696 8×108m3，在中水平流量下灌區可用的地表水、地下水資源量依次為[2.25×108,2.42×108] m3、[3.25×108，3.47×108] m3，此時配水量充足，不存在缺水。而在低流量情況下，灌區的地下水可用水量為[2.58×108,2.76×108] m3，此時可分配水量幾乎達到上限，作物供水存在缺失。可能的原因是:①水資源不合理開發，過度利用地下水。②缺少合理的規劃和利用。

最優目標配水量由模型第1階段得到，最優分配水量則由式(19)確定。魯棒優化方法旨在第2階段中對成本期望進行追索，并對系統的風險進行評估。表4中給出了ρ取不同值時，對應的缺水量和最優配水量。當ρ=0時，模型為普通區間兩階段隨機規劃模型，代表了決策者對風險持中立態度，不考慮成本可變性。由表4可知：①對于水稻，地表水的決策變量為0，當ρ=0時，在低、中、高3種來水水平下，其缺水量分別為[3.68×106，5.68×106] m3、[1.53×106，3.53×106] m3、[2.10×105，1.21×106] m3，因此，最優配水量為[4.419×107，4.619×107] m3、[4.634×107，4.834×107] m3、[4.866×107，4.966×107] m3。由表3可知，水稻的地表水最小需水量為[4.321×107，4.444×107] m3,可滿足水稻的供水需求。對于玉米，地表水的決策變量為1，當ρ=0時不產生缺水量，最優配水量為6.03×106m3。在供水量充足的情況下，系統更偏好避免懲罰風險，所以將目標配水量上限作為最優目標配水量。②隨著ρ的增大，缺水量不斷增加。例如，在低水平流量下，當ρ=0.4，1，2時，水稻的地下水缺水量為[2.179×107，2.788×107] m3、[3.305×107，3.877×107] m3、[3.628×107，4.356×107] m3，最優配水量為[1.716 1×108，1.767 0×108] m3、[1.607 2×108，1.664 4×108] m3、[1.559 3×108，1.632 1×108] m3，玉米的地表水缺水量為0 m3、[8.3×105，9.7×105] m3、[1.02×106，1.34×106] m3，最優配水量為6.03×106m3、[5.06×106，5.20×106] m3、[4.69×106，5.01×106] m3。低水平時系統來水量較少，通過減少配水量來調節系統的最優配水量，即，魯棒系數越大，系統的穩定性越強。模型增加了作物缺水量，降低了作物的配水量，從而使系統的穩定性增加。③當ρ逐漸增大時，作物的缺水量隨之增大，但當ρ≥3之后，缺水量幾乎不再變化。例如，當ρ=2,3,5時，在中水平流量下，大豆的地表水缺水量為[8.80×105，1.19×106] m3、[1.32×106，2.58×106] m3、[1.36×106，2.60×106] m3，最優配水量為[6.81×106，7.12×106] m3、[5.42×106，6.68×106] m3、[5.46×106，6.64×106] m3，玉米的地下水缺水量為[6.23×106，6.88×106] m3、[7.22×106，8.67×106] m3、[7.22×106，8.67×106] m3，最優配水量為[1.723×107，1.788×107] m3、[1.544×107，1.689×107] m3、[1.544×107，1.689×107] m3，缺水量的增加會使系統穩定性增加，但為了保證作物正常生長的最小需水量，缺水量不再增大，此時系統趨于穩定。

表4 不同ρ取值的作物缺水量和最優配置水量Tab.4 Crops water shortage and water optimal allocation of different values of ρ

續表4

表5 決策變量值Tab.5 Value of decision variable

圖2 作物最優目標配水量Fig.2 Optimized allocation targets for different crops

最優系統成本區間如圖3、4所示，系統成本隨著魯棒系數的變化而變化，呈遞增趨勢。如圖3所示，當ρ=0時，最優系統總成本為[1.104 32×109,2.049 95×109]元，當ρ=1時，模型最優系統總成本為[1.331 55×109，2.235 76×109]元。當ρ=5時，最優系統總成本為[1.943 77×109，2.657 69×109]

元。在低流量水平下(圖4a)，最優系統成本在[6.835 5×108，9.286 7×108]元(ρ=5)和[3.693 9×108，6.962 5×108]元(ρ=0)之間變化。在中流量水平下(圖4b)，最優系統成本在[6.496 7×108，9.130 5×108]元(ρ=5)和[3.683 3×108，6.864 5×108]元(ρ=0)之間變化。在高流量水平下(圖4c)，系統成本在[6.105 5×108，8.160 2×108]元(ρ=5)和[3.666 1×108，6.672 5×108]元(ρ=0)之間變化。對比圖3、4，隨著水資源最優分配量的變化，系統成本呈現一定的變化規律：①魯棒系數增加，引起系統成本增大，當ρ≥3之后，成本幾乎不變，說明系統已經趨于穩定。②隨著魯棒系數增加，成本的上下限差值變小，系統的穩定性增加，經濟性和穩定性得到了較好的平衡。③較高的成本對應較高的缺水水平。當可用水量較高時，決策者可利用的水資源量也會較多，如果實際分配量較少，則會產生較高的系統風險和較多的懲罰成本；相反，如果可用水量較少，則決策者需要減少實際供水量，采取保守決策，降低系統風險增加穩定性。決策者可以根據系統分析結果，針對灌區實際情況，制定風險和經濟相協調的水資源配置策略。

圖3 最優系統總成本Fig.3 Optimized net system cost

圖4 不同流量水平下最優系統成本Fig.4 Optimized net system cost under different water inflow levels

3 結論

(1)針對水資源分配過程中存在風險的問題，將魯棒優化方法與兩階段規劃方法耦合，建立了區間兩階段魯棒優化模型，以三江平原牡丹江灌區農業水資源配置為例進行了研究。模型結果表明，系統總成本隨著魯棒系數的變化有一定的規律。當模型不考慮系統風險時，即ρ=0時，系統成本在[1.104 32×109,2.049 95×109]元之間變化，隨著魯棒系數的增大，模型的缺水量增加，使得系統的穩定性增強、成本增加，當ρ=0.4、1、2時，系統成本分別在[1.263 25×109，2.185 67×109]元、[1.331 55×109,2.235 76×109]元、[1.608 79×109,2.415 52×109]元之間變化，但當ρ=3、5時，系統缺水量不再增大，此時系統達到穩定狀態，成本在[1.903 27×109,2.634 75×109]元、[1.943 77×109,2.657 69×109]元之間變化。由此可知，在不同來水量水平下，通過增加魯棒系數增加模型結果的可行性，對灌區的用水成本、系統的穩定性和系統的風險三者之間進行了充分的權衡，使配置方案更具有實際操作性和靈活性。

(2)與傳統區間兩階段隨機規劃模型相比，區間兩階段魯棒優化方法不但可以有效地解決不確定條件下的隨機問題和區間問題，魯棒優化方法的加入可以捕獲規劃過程中產生的風險問題，避免優化結果出現高風險狀態，彌補了模型存在風險的缺陷，增強了系統的穩定性。將區間兩階段魯棒優化方法應用到灌區水資源優化配置中，驗證了模型的應用性和有效性。模型結果表明，通過魯棒系數的變化，生成一系列對應的不同風險水平、不同情境的可行性方案，顯示了系統經濟性和系統穩定性之間的權衡。