王 貝，朱 迪，何錫君，陳 東，徐長江，周研來，陳 華

（1. 浙江省水文管理中心，浙江 杭州 310009； 2. 武漢大學(xué) 水資源與水電工程科學(xué)國家重點(diǎn)實(shí)驗(yàn)室， 湖北 武漢 430072；3. 開化縣水利局，浙江 開化 324300； 4. 長江水利委員會水文局，湖北 武漢 430010）

0 引 言

水庫群的建設(shè)與投運(yùn)，一方面通過對天然徑流的調(diào)節(jié)實(shí)現(xiàn)興利除害，另一方面也改變了河流原有的水文情勢，可能對河流生態(tài)系統(tǒng)的動態(tài)平衡造成影響［1-3］。因此，為維護(hù)河流生態(tài)健康，堅持人與自然和諧共生的生態(tài)理念，開展水庫群生態(tài)調(diào)度的研究是非常有必要的。

近年來，國內(nèi)外學(xué)者圍繞水庫群生態(tài)調(diào)度開展了一列研究，在生態(tài)調(diào)度模型構(gòu)建和優(yōu)化求解方面取得了有意義的成果。在建模方面，Babel等［4］考慮工業(yè)、農(nóng)業(yè)、生活及生態(tài)等不同部門需水要求，構(gòu)建了以需水滿足度最大、用水效益最大的優(yōu)化配置模型；Yang等［5］采用IHA∕RVA 法確定漢江下游襄陽斷面的生態(tài)流量過程，構(gòu)建了考慮河道水文情勢的丹江口水庫生態(tài)調(diào)度模型；龍凡等［6］采用年內(nèi)展布法和改進(jìn)的流量歷時曲線法計算生態(tài)流量過程，設(shè)置多種生態(tài)流量約束方案，構(gòu)建溪洛渡—向家壩生態(tài)調(diào)度模型，分析不同生態(tài)流量約束下的水電站發(fā)電變化特征；王昊等［7］基于Vague 集理論構(gòu)建了瀾滄江景洪水電站生態(tài)調(diào)度評價模型，優(yōu)化原有水電站調(diào)度方案，增加電站下游生態(tài)效益。在求解算法方面，發(fā)展了一批以傳統(tǒng)優(yōu)化算法和智能優(yōu)化算法為代表的求解算法。例如，高仕春等［8］采用改進(jìn)大規(guī)模線性規(guī)劃方法優(yōu)化沙潁河水庫群、閘壩等水工程生態(tài)調(diào)度過程，提高沙潁河流域關(guān)鍵斷面生態(tài)用水的保證率；陳悅云等［9］以贛江流域?yàn)槔瑯?gòu)建了考慮發(fā)電、供水以及生態(tài)的多目標(biāo)水庫群，采用改進(jìn)的多目標(biāo)粒子群算法進(jìn)行優(yōu)化求解，剖析了不用優(yōu)化目標(biāo)間的相互關(guān)系；白濤等［10］建立水電站發(fā)電和減水河段生態(tài)保證程度的多目標(biāo)生態(tài)調(diào)度模型，采用非支配遺傳算法（Non-dominated Sorting Genetic Algorithm-Ⅱ， NSGA-Ⅱ）進(jìn)行求解，量化了減水河段興利和生態(tài)目標(biāo)的脅迫關(guān)系。

目前研究對象多數(shù)集中于大流域平原地區(qū)，對于山區(qū)性中小流域的生態(tài)調(diào)度研究較為少見。馬金溪流域位于衢江上游，是錢塘江南源之一，是典型的山區(qū)性中小流域，具有“九山半水半分田”的特征［11］。流域所在的開化縣為全國生態(tài)示范區(qū)、浙江生態(tài)縣，是全國17個具有全球意義的山地保護(hù)區(qū)之一［12］。目前，馬金溪流域目前建有中型水庫2 座，在建大型水庫1 座，隨著水庫群的建設(shè)運(yùn)行，將會對流域生態(tài)環(huán)境和區(qū)域經(jīng)濟(jì)社會發(fā)展產(chǎn)生影響。因此，如何充分發(fā)揮水庫群水資源調(diào)配作用，貫徹開化“生態(tài)立縣”的發(fā)展戰(zhàn)略，協(xié)同生態(tài)和興利目標(biāo)的均衡優(yōu)化，是值得進(jìn)一步研究的問題。

以馬金溪流域已建的2 座中型水庫和在建的1 座大型水庫為研究實(shí)例，考慮水庫綜合利用功能及下游控制站點(diǎn)的生態(tài)需求，構(gòu)建面向生態(tài)、發(fā)電和供水需求的馬金溪流域水庫群多目標(biāo)生態(tài)調(diào)度模型，采用NSGA-Ⅱ算法進(jìn)行優(yōu)化求解，研究不同頻率來水條件下生態(tài)、發(fā)電和供水目標(biāo)間的關(guān)系，為馬金溪流域水庫群的生態(tài)調(diào)度決策提供參考。

1 研究區(qū)域概化及數(shù)據(jù)

馬金溪河道全長91.2 km，流域面積1 011.3 km2，主要支流有何田溪、村頭溪、中村溪、池淮溪、龍山溪、馬尪溪等6條支流，流域內(nèi)大中型水庫有齊溪、茅崗和開化水庫（在建），流域水系概化圖如圖1所示。齊溪水庫位于馬金溪上游，是一座以防洪、發(fā)電為主，結(jié)合灌溉、養(yǎng)漁等綜合利用的中型水庫，總庫容4 517 萬m3。茅崗水庫位于馬金溪支流中村溪上游，是一座以灌溉為主，結(jié)合防洪、發(fā)電等綜合利用的中型水庫，總庫容1 125 萬m3。開化水庫位于齊溪水庫下游，是一座防洪、供水和改善流域生態(tài)環(huán)境為主的大型水庫，總庫容1.84 億m3。根據(jù)齊溪、茅崗和開化水庫的初步設(shè)計報告書，以及《開化縣小水電下游河道生態(tài)核算報告書》等資料，統(tǒng)計水庫特征參數(shù)情況如表1所示。

表1 水庫特征參數(shù)Tab.1 The characteristics of reservoirs

圖1 馬金溪流域水系概化圖Fig.1 Generalization diagram of the Majinxi basin

《開化縣小水電下游河道生態(tài)核算報告書》基于偏工程不利原則，選擇文圖水文站1976年、1981年以及1968年分別作為平水年（P=50%）、枯水年（P=75%）以及特枯水年（P=95%）的典型年。本文以馬金溪下游文圖水文站為生態(tài)保障斷面，將文圖站年徑流系列排頻，依據(jù)《開化縣小水電下游河道生態(tài)核算報告書》的計算成果，分別選取上述來水頻率典型年，計算時段為月，水庫入庫和支流來水基于同期數(shù)據(jù)按同倍比放大得到。此外，根據(jù)《開化縣開化水庫工程水資源論證報告》等資料，確定用水區(qū)生活、農(nóng)業(yè)等部門的需水量、流域綜合耗水率等數(shù)據(jù)。

2 模型構(gòu)建

基于圖1所示的馬金溪流域水系概化圖，本文面向生態(tài)、發(fā)電和供水等要求，構(gòu)建馬金溪流域水庫群多目標(biāo)生態(tài)調(diào)度模型，相應(yīng)的數(shù)學(xué)表達(dá)式如下。

2.1 目標(biāo)函數(shù)

（1）生態(tài)目標(biāo)。一般來說，天然水文情勢下的河流生物多樣性和生態(tài)系統(tǒng)完整性最好［13］，因此，本文以文圖水文站調(diào)度后的流量與天然流量偏差最小作為生態(tài)目標(biāo)，選擇修正全年流量偏差函數(shù)（Amended Annual Proportional Flow Deviation，AAPFD）作為目標(biāo)函數(shù)表達(dá)式［14］，具體如下：

式中：Qd（t）表示調(diào)度后文圖水文站第t時段流量；QN（t）表示文圖水文站第t時段天然流量；QˉN表示文圖水文站調(diào)度期內(nèi)天然流量平均值；T表示調(diào)度時長。

（2）發(fā)電目標(biāo)。調(diào)度期內(nèi)水庫水電站總發(fā)電量最大，表達(dá)式如下：

式中：km表示第m個水庫水電站出力系數(shù)；Qfd，m（t）表示第m個水庫水電站第t時段發(fā)電流量；Hm（t）表示第m個水庫水電站第t時段水頭；Δt表示計算時段長。

（3）供水目標(biāo)。調(diào)度期內(nèi)流域生活、農(nóng)業(yè)缺水量最小，表達(dá)式如下：

式中：Wxu，j（t）和Wqu，j（t）分別表示第j個取水區(qū)第t時段需水和取水量；J表示取水區(qū)數(shù)量。

從生態(tài)、發(fā)電和供水目標(biāo)函數(shù)上看，生態(tài)目標(biāo)與發(fā)電、供水目標(biāo)存在競爭關(guān)系，發(fā)電和供水目標(biāo)的改善，會改變河流的天然水文情勢，進(jìn)而增大AAPFD的數(shù)值，影響生態(tài)目標(biāo)效益。

2.2 約束條件

（1）水庫水量平衡。

式中：Vm（t）表示第m個水庫第t時段庫容；Qmin（t）和Qmout（t）分別表示第m個水庫第t時段入庫和出庫流量；Wqu，m（t）表示第m個水庫第t時段取水量。

（2）水位限制約束。

式中：Zm（t）表示第m個水庫第t時段水位；Zmmax（t）和Zmmin（t）分別表示調(diào)度期內(nèi)第m個水庫第t時段水位的上下限。

（3）水庫起調(diào)水位約束。

式中：Zmstart表示第m個水庫起調(diào)水位值，一般取為水庫的正常高水位。

（4）水庫出庫流量約束

式中：Qmmax（t）和Qmmin（t）分別表示第m個水庫出庫流量的上下限。

（5）水電站出力限制約束。

式中：Nm（t）表示第m個水庫水電站第t時段出力；Nmma（xt）和（t）分別表示調(diào)度期內(nèi)第m個水庫水電站第t時段出力的上下限。

（6）取水節(jié)點(diǎn)約束。

式中：Rj（t）表示第j個取水區(qū)第t時段取水后的河道流量；Qsy，j（t）和Qqj，j（t）分別表示第j個取水區(qū)第t時段上游和區(qū)間來水。

（7）綜合用水節(jié)點(diǎn)約束：

式中：WL，j（t）和WA，j（t）分別表示第j個取水區(qū)第t時段的生活用水、農(nóng)業(yè)用水量；WL，jxu（t）和WA，jxu（t）分別表示第j個取水區(qū)第t時段的生活、農(nóng)業(yè)需水量。

（8）退水量計算約束：

式中：Wtui，j（t）表示第j個取水區(qū)第t時段的退水量；α表示流域綜合耗水系數(shù)。

（9）非負(fù)約束。各種變量均為非負(fù)值。

3 求解方法

研究構(gòu)建的模型綜合考慮流域水庫群生態(tài)、發(fā)電和供水等需求，是一個多目標(biāo)優(yōu)化調(diào)度問題，水庫優(yōu)化調(diào)度常用的線性規(guī)劃（Linear Programming， LP）和動態(tài)規(guī)劃（Dynamic Programming， DP）等算法不再適用［9］。近年來，隨著多目標(biāo)智能算法的興起，諸如多目標(biāo)粒子群算法（Multiple Objective Particle Swarm Optimization， MOPSO）［15］、多目標(biāo)差分進(jìn)化算法［16］、NSGA-Ⅱ算法［17］等被廣泛用于求解此類水庫多目標(biāo)優(yōu)化調(diào)度問題［18］。其中，Deb 等［19］提出的NSGA-Ⅱ算法，具有搜索速度快，收斂性好等優(yōu)點(diǎn)。因此，本文采用NSGA-Ⅱ算法對馬金溪水庫群多目標(biāo)生態(tài)調(diào)度模型進(jìn)行優(yōu)化求解，其求解流程如下。

（1）設(shè)置NSGA-Ⅱ算法的相關(guān)參數(shù)，包括：最大迭代次數(shù)為K、交叉和變異的概率分別為p1和p2，種群規(guī)模NP等。

（2）令迭代次數(shù)k=1，隨機(jī)初始化一個父代種群Pk，計算種群內(nèi)不同個體的生態(tài)、發(fā)電和供水目標(biāo)的適應(yīng)度值，并將所有個體按非支配關(guān)系排序，計算其擁擠度。

（3）采用選擇、交叉、變異算子產(chǎn)生子代種群Qk。

（4）將父代種群與子代種群合并組成Rk，種群規(guī)模為2NP；對Rk進(jìn)行非支配排序，產(chǎn)生非支配集并計算擁擠度，并生成新的父代種群Pk+1。

（5）通過選擇、交叉、變異算子產(chǎn)生子代種群Qk+1，迭代次數(shù)k=k+1。

（6）重復(fù)步驟（4）、（5），直至迭代次數(shù)達(dá)到K，獲取面向生態(tài)、發(fā)電和供水等多目標(biāo)優(yōu)化調(diào)度的非劣解集，及其相應(yīng)的水庫群調(diào)度方案。

為有效處理水量平衡、出庫流量、取用水計算等復(fù)雜約束，基于馬金溪流域水資源系統(tǒng)概化情況（如圖1（b）所示），本文以各水庫各時段出庫流量和用水區(qū)取水量為決策變量，即：u=｛Q1ou（t1），…，Q1ou（t12），Q2ou（t1），…，Q2ou（t12），Q3ou（t1），…，Q3out（12），Wqu，（11），…，Wqu，（112），Wqu，（21），…，Wqu，（212），Wqu，（31），…，Wqu，（312），Wqu，（41），…，Wqu，（412），｝，共計84個。參考相關(guān)文獻(xiàn)［10］、［17］的計算參數(shù)設(shè)置，經(jīng)計算驗(yàn)證，本文中，K=500，p1=0.8，p2=0.2，NP=200。

4 計算結(jié)果及分析

4.1 生態(tài)、發(fā)電、供水多目標(biāo)關(guān)系分析

基于來水、需水、水庫群等相關(guān)資料，采用NSGA-Ⅱ算法對上述多目標(biāo)生態(tài)調(diào)度模型進(jìn)行求解，得到P=50%、P=75%以及P=95%等不同頻率來水下的馬金溪流域水庫群生態(tài)、發(fā)電和供水等目標(biāo)的非劣解集，如圖2 所示。不同頻率來水下，生態(tài)、發(fā)電和供水目標(biāo)之間的相關(guān)關(guān)系圖如圖3至圖5所示。

圖2 不同頻率來水非劣解集分布結(jié)果Fig.2 The distribution of non-inferior solutions under different inflow frequencies

圖3 生態(tài)目標(biāo)和發(fā)電目標(biāo)之間關(guān)系Fig.3 The relationship between the objectives of ecology and power generation

對比圖2（a）～（c）可知，在不同頻率來水條件下，隨著流域來水增加，文圖水文站的AAPFD 從0.28～0.39（P=95%）減小至0.14～0.27（P=50%），水庫水電站發(fā)電量從5 313～6 388 萬kWh（P=95%）增加至8 521～9 793 萬kWh（P=50%），而流域供水的缺水量從0～805 萬m3（P=95%）隨之減小至0～105 萬m3（P=50%），3 個目標(biāo)均有所改善。圖3 反映的是在不同頻率來水條件下，供水缺水量取不同數(shù)值時，生態(tài)目標(biāo)和發(fā)電目標(biāo)之間的關(guān)系。從圖3 可知，當(dāng)供水缺水量為定值時，AAPFD 會隨發(fā)電量的增加而增加，即生態(tài)效益會隨發(fā)電目標(biāo)的改善而減??；當(dāng)供水缺水量增加時，生態(tài)目標(biāo)和發(fā)電目標(biāo)之間的點(diǎn)據(jù)呈現(xiàn)上移趨勢，表現(xiàn)為發(fā)電效益隨著生態(tài)和供水效益的減少而增加。圖4反映的是在不同頻率來水條件下，發(fā)電量取不同數(shù)值時，生態(tài)目標(biāo)和供水目標(biāo)之間的關(guān)系。從圖4可知，當(dāng)發(fā)電量固定時，隨著供水缺水量的減少，AAPFD 的數(shù)值沒有太大幅度的變化，說明生態(tài)目標(biāo)和供水目標(biāo)之間的競爭關(guān)系并不明顯；當(dāng)發(fā)電目標(biāo)改善時，生態(tài)目標(biāo)和供水目標(biāo)之間的點(diǎn)據(jù)呈現(xiàn)上移趨勢，表現(xiàn)為供水效益隨發(fā)電效益的增加而減小。圖5反映的是在不同頻率來水條件下，AAPED 取不同數(shù)值時，發(fā)電目標(biāo)和供水目標(biāo)之間的關(guān)系。從圖5 可知，供水目標(biāo)和發(fā)電目標(biāo)間的競爭關(guān)系較為明顯，強(qiáng)于圖4中的生態(tài)目標(biāo)和供水目標(biāo)之間的關(guān)系。

圖4 生態(tài)目標(biāo)和供水目標(biāo)之間關(guān)系Fig.4 The relationship between the objectives of ecology and water supply

圖5 發(fā)電目標(biāo)和供水目標(biāo)之間關(guān)系Fig.5 The relationship between the objectives of power generation and water supply

總體而言，馬金溪流域生態(tài)、發(fā)電和供水目標(biāo)之間呈現(xiàn)競爭關(guān)系，改善一個目標(biāo)需要以犧牲其余兩者中至少一個目標(biāo)為代價。其中，生態(tài)目標(biāo)與發(fā)電目標(biāo)、供水目標(biāo)與發(fā)電目標(biāo)間的競爭關(guān)系較為明顯，而生態(tài)目標(biāo)和供水目標(biāo)間的競爭關(guān)系相比之下較弱。

4.2 典型非劣解集下游斷面流量過程分析

分別選擇APPFD 值最?。ū硎旧鷳B(tài)目標(biāo)最優(yōu)，記為方案1）、總發(fā)電量最大（表示發(fā)電目標(biāo)最優(yōu)，記為方案2）以及總?cè)彼孔钌伲ū硎竟┧繕?biāo)最優(yōu)，記為方案3）所對應(yīng)的調(diào)度方案，計算不同來水頻率下各方案的目標(biāo)值，如表2 所示。各方案下文圖水文站的流量過程如圖6所示。

表2 不同來水頻率下選取方案的目標(biāo)值Tab.2 The objective values of the selected schemes under different inflow frequencies

圖6 不同來水頻率文圖水文站調(diào)度前后流量過程Fig.6 The streamflow of the Wentu station before and after operation under different inflow frequencies

從表2 可知，對于P=50%和P=95%頻率來水，方案1 計算的發(fā)電量要明顯小于方案2，發(fā)電量分表減少1 272 萬和1 075萬kWh，但缺水量可以削減至0，而對于P=50%頻率來水，方案1的發(fā)電量與方案2區(qū)別較小，但AAPFD 值和缺水量明顯更小，生態(tài)和供水目標(biāo)更優(yōu)。對于不同頻率來水情況，方案3 可以將缺水量減少到0，且發(fā)電和生態(tài)目標(biāo)分別相較方案1和方案2有所改善，是一種偏向折中的方案。

從圖6可知，在不同來水頻率下，對于徑流占比較大的月份（豐水期），按方案1調(diào)度后的流量低于天然徑流，但高于其他方案；對于徑流占比較小的月份（枯水期），按方案1調(diào)度后的流量與天然徑流接近，低于其他方案。而方案2 調(diào)度后的流量過程與方案1 相反，即，豐水期方案2 調(diào)度后的流量最小，低于同時期的天然徑流和其他方案調(diào)度結(jié)果；枯水期方案2 調(diào)度后的流量過程要高于同時期的天然徑流和其他方案調(diào)度結(jié)果。方案3調(diào)度后的流量過程則是部分月份與方案1 相似，部分月份與方案2 相似。具體而言，以圖6（c）為例，在P=95%來水條件下，5月和7月的徑流占比較大，1-2月以及8-11月的徑流占比較?。环桨? 在5 月和7 月調(diào)度后的流量相較天然徑流和其他方案偏低，1-2 月以及8-11 月調(diào)度后的流量過程則高于天然徑流和其他方案，而方案1 的調(diào)度結(jié)果則與之相反。另外，方案3 在1-3月、6 月的調(diào)度流量過程表現(xiàn)出與方案1 的相似性，在5 月、9 月和10 月的調(diào)度流量過程又與方案2 較為接近。由此可以看出，方案1 和方案2 調(diào)度后的流量過程呈現(xiàn)相異性，反映出生態(tài)目標(biāo)與發(fā)電目標(biāo)競爭性較強(qiáng)，而供水目標(biāo)與前兩者的競爭性相對較弱。

5 結(jié) 論

以馬金溪流域2 座建成的中型水庫及1 座在建的大型水庫為研究對象，以下游文圖水文站AAPFD 值最小、流域水庫水電站群總發(fā)電量最大以及供水缺水量最少作為優(yōu)化目標(biāo)，構(gòu)建了面向生態(tài)、發(fā)電和供水的馬金溪流域水庫群多目標(biāo)生態(tài)調(diào)度模型，并采用NSGA-Ⅱ算法進(jìn)行求解，分析了P=50%、P=75%以及P=90%等不同頻率來水條件下各目標(biāo)間的關(guān)系，得到結(jié)論如下。

（1）在不同頻率來水條件下，隨著流域來水的增加，文圖水文站AAPFD值減小，水電站總發(fā)電量增加，供水缺水量減少，生態(tài)、發(fā)電和供水目標(biāo)都能有所改善。

（2）在同一頻率來水條件下，發(fā)電目標(biāo)與生態(tài)、供水目標(biāo)呈現(xiàn)出競爭關(guān)系，其中，發(fā)電與生態(tài)目標(biāo)之間的競爭性較強(qiáng)，與供水目標(biāo)之間的競爭性較弱。

（3）分別選取AAPFD 值最小、水電站總發(fā)電量最大以及供水缺水量最少對應(yīng)的非劣解集作為方案1～3，其中，相較其它方案，方案1豐水期調(diào)度后的文圖水文站流量更大，枯水期流量更小，方案2則與之相反；而從目標(biāo)值的計算和文圖站流量過程上看，方案3是一種更偏向折中的方案。

研究解析了山區(qū)性中小流域水庫群生態(tài)效益和興利效益之間的相關(guān)關(guān)系，研究結(jié)果可為保障馬金溪流域生態(tài)流量和增加水庫群綜合效益，提供科學(xué)指導(dǎo)和參考依據(jù)。但本文主要依據(jù)歷史典型年資料，且生態(tài)目標(biāo)采用AAPFD 計算。因此，未來研究將考慮結(jié)合水文預(yù)報信息，采用多種生態(tài)流量計算方法核算生態(tài)流量，開展保障山區(qū)性中小流域生態(tài)流量的預(yù)報調(diào)度，以進(jìn)一步驗(yàn)證模型方法的可行性。