黃 馗，穆振宇，艾學(xué)山，3，4，陳森林，3，郭佳俊，丁 杰

（1.廣西電網(wǎng)電力調(diào)度中心，廣西 530013；2.武漢大學(xué)水利水電學(xué)院，武漢 430072；3.武漢大學(xué)水資源與水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，武漢 430072；4.海綿城市建設(shè)水系統(tǒng)科學(xué)湖北省重點(diǎn)實(shí)驗(yàn)室（武漢大學(xué)），武漢 430072）

0 引 言

在水庫(kù)調(diào)度過(guò)程中，入庫(kù)流量計(jì)算是最重要的基礎(chǔ)工作之一，水庫(kù)洪水預(yù)報(bào)調(diào)度方案的編制、水庫(kù)調(diào)度圖的編制、水庫(kù)調(diào)度經(jīng)濟(jì)評(píng)價(jià)以及水庫(kù)洪水調(diào)節(jié)演算等都需要入庫(kù)流量資料［1］。

水庫(kù)入庫(kù)流量一般采用水量平衡方法進(jìn)行反推計(jì)算，目前大部分水庫(kù)均采用靜庫(kù)容方法計(jì)算，利用水位庫(kù)容關(guān)系曲線計(jì)算開(kāi)始和終止時(shí)刻的水位對(duì)應(yīng)的水庫(kù)庫(kù)容，由此得到該時(shí)段蓄水量的變化，同時(shí)，計(jì)算該時(shí)段水庫(kù)總的出庫(kù)水量，利用水量平衡原理，得到水庫(kù)在該時(shí)段內(nèi)的總?cè)霂?kù)水量，再與時(shí)段長(zhǎng)相除，得到該計(jì)算時(shí)段內(nèi)的平均入庫(kù)流量［2］。該方法體現(xiàn)了水量平衡原理、流量連續(xù)性原理等主要規(guī)律，但由于計(jì)算過(guò)程中出庫(kù)流量與時(shí)段庫(kù)容變化量是依據(jù)庫(kù)容曲線、泄流特性曲線及機(jī)組耗流量特性曲線等計(jì)算得到，故會(huì)因水位波動(dòng)、各曲線擬合誤差、時(shí)段長(zhǎng)度等因素影響，使得最終得到的入庫(kù)流量結(jié)果具有不可忽略的偏差，而在梯級(jí)水電站水庫(kù)調(diào)度過(guò)程中，該偏差會(huì)隨著計(jì)算的延續(xù)而不斷擴(kuò)大，進(jìn)而演變成梯級(jí)水電站之間的水量不平衡問(wèn)題，比如可能出現(xiàn)兩種極端情況：區(qū)間徑流系數(shù)大于1.0或?yàn)樨?fù)值。

另一方面，由于調(diào)峰電廠需要隨時(shí)調(diào)整負(fù)荷，其變動(dòng)會(huì)引起發(fā)電流量的變動(dòng)，同時(shí)由于測(cè)量誤差和曲線誤差等原因，通過(guò)插值等方式得到的時(shí)段庫(kù)容變化量也不準(zhǔn)確，二者共同作用使得短歷時(shí)內(nèi)計(jì)算所得到的入庫(kù)流量數(shù)據(jù)波動(dòng)很大，在其曲線上表現(xiàn)為明顯的“鋸齒”現(xiàn)象［3］，極不利于后續(xù)調(diào)度過(guò)程的實(shí)施，然而，在計(jì)算時(shí)段較長(zhǎng)時(shí)，鋸齒現(xiàn)象表現(xiàn)的不明顯，這一是因?yàn)殡S著時(shí)間段的增加，各種影響入庫(kù)反推的誤差會(huì)相互抵消［4］，二是在長(zhǎng)時(shí)段計(jì)算中發(fā)電流量等因素的變動(dòng)減小，誤差會(huì)被均化。

為減少入庫(kù)流量數(shù)據(jù)中的鋸齒現(xiàn)象，唐海華等［1］提出對(duì)庫(kù)區(qū)河道進(jìn)行分段優(yōu)選，并根據(jù)各分段靜庫(kù)容曲線和代表水位計(jì)算入庫(kù)流量，使得其波動(dòng)幅度顯著減??；武煒等［4］從入庫(kù)流量反推方式出發(fā)，提出使用五點(diǎn)三次平滑算法對(duì)入庫(kù)流量進(jìn)行反推從而有效降低入庫(kù)流量計(jì)算過(guò)程中的干擾。另一方面，梯級(jí)水電站之間的水量不平衡問(wèn)題是由入庫(kù)流量反推中的誤差演變而來(lái)，采取物理或數(shù)學(xué)模型進(jìn)行入庫(kù)流量預(yù)報(bào)可以有效避免此誤差的發(fā)生，基于此思想，劉曉陽(yáng)等［5］基于降雨徑流模型對(duì)梅山水庫(kù)入庫(kù)流量進(jìn)行模擬，體現(xiàn)出降雨在徑流模擬和預(yù)報(bào)中的重要性；康玲等［6］基于Morlet 小波對(duì)丹江口水庫(kù)入庫(kù)流量周期性進(jìn)行分析并揭示了入庫(kù)流量的多尺度震蕩特性和周期性變化規(guī)律；劉章君等［7］構(gòu)建Copula-BFS 模型并以三峽水庫(kù)汛期入庫(kù)流量概率預(yù)報(bào)為例對(duì)模型進(jìn)行檢驗(yàn)，為入庫(kù)流量預(yù)報(bào)提供了一種新型途徑；李茂貴等［8］建立了針對(duì)酒埠江水庫(kù)入庫(kù)流量的逐步回歸預(yù)報(bào)模型；Herbert 等［9］利用深度學(xué)習(xí)模式進(jìn)行長(zhǎng)時(shí)期的入庫(kù)流量預(yù)報(bào)，有效提高了入庫(kù)流量預(yù)報(bào)的精度；Zhang等［10］基于多種模型開(kāi)發(fā)出一種機(jī)器學(xué)習(xí)算法并用于水庫(kù)入庫(kù)流量的預(yù)報(bào)，取得了一定成果。然而，上述方法雖然在一定程度上解決了出庫(kù)流量鋸齒現(xiàn)象以及誤差演變?yōu)樗坎黄胶鈫?wèn)題，但是，在短時(shí)段預(yù)報(bào)中上述方法的預(yù)報(bào)精度還不能有效滿足調(diào)度要求，同時(shí)，運(yùn)用數(shù)學(xué)模型或算法進(jìn)行預(yù)報(bào)難以反映流域站點(diǎn)間隱含的物理機(jī)制，說(shuō)服力不足。

由此可見(jiàn)，流域梯級(jí)水電站之間水流演進(jìn)計(jì)算存在兩個(gè)方面的問(wèn)題，一是因庫(kù)容曲線、泄流特性曲線及機(jī)組耗流量特性曲線等不夠精確，計(jì)算的出庫(kù)流量和時(shí)段庫(kù)容變化量不可避免的存在誤差，從而導(dǎo)致梯級(jí)水電站之間的水量不平衡問(wèn)題（如區(qū)間徑流系數(shù)大于1.0 或?yàn)樨?fù)值），二是短時(shí)段計(jì)算所得到的入庫(kù)流量過(guò)程鋸齒現(xiàn)象嚴(yán)重。針對(duì)這兩個(gè)問(wèn)題，本文提出一種受上游水庫(kù)調(diào)控影響的日時(shí)段入庫(kù)流量預(yù)報(bào)方法，將入庫(kù)流量分為上游水庫(kù)出庫(kù)流量與區(qū)間匯流兩部分，以日為時(shí)段長(zhǎng)，應(yīng)用API 模型計(jì)算區(qū)間產(chǎn)匯流過(guò)程，并在上游水庫(kù)出流的河道演算中引入修正系數(shù)處理區(qū)間來(lái)水或引水問(wèn)題，通過(guò)郁江流域西津電站的實(shí)際應(yīng)用進(jìn)行檢驗(yàn)。

1 日時(shí)段入庫(kù)流量預(yù)報(bào)方法

日時(shí)段入庫(kù)流量預(yù)報(bào)方法的基本思路是：以日為時(shí)段長(zhǎng)，以需要預(yù)報(bào)的水電站上游水庫(kù)的出庫(kù)流量、區(qū)間降雨量和水電站控制面積為基礎(chǔ)，以水庫(kù)反推的入流數(shù)據(jù)作為擬合標(biāo)準(zhǔn)進(jìn)行水庫(kù)日平均入庫(kù)流量的模擬計(jì)算，具體計(jì)算步驟如圖1所示。

圖1日時(shí)段入庫(kù)流量預(yù)報(bào)方法計(jì)算流程圖Fig.1 Calculation flow chart of daily inflow forecast method

首先，根據(jù)區(qū)間降雨量等數(shù)據(jù)進(jìn)行區(qū)間產(chǎn)流計(jì)算，同時(shí)對(duì)上游出庫(kù)流量進(jìn)行修正；其次，對(duì)區(qū)間產(chǎn)流進(jìn)行匯流計(jì)算，對(duì)修正后的上游出庫(kù)流量過(guò)程進(jìn)行河道匯流計(jì)算；最后，將這兩個(gè)過(guò)程疊加作為下游水庫(kù)的日時(shí)段入庫(kù)流量過(guò)程。

1.1 區(qū)間產(chǎn)流計(jì)算方法

區(qū)間流域產(chǎn)流計(jì)算應(yīng)用API 模型［11］，流域前期影響雨量計(jì)算公式為：

式中：Pa，t+1為t + 1時(shí)段前期影響雨量；K為消退系數(shù)；Pa，t為t時(shí)段前期影響雨量；Pa，M為流域最大前期影響雨量；Pt為t 時(shí)段降雨量。在長(zhǎng)時(shí)段計(jì)算中初始前期影響雨量可設(shè)置為0。

假設(shè)時(shí)段徑流系數(shù)βt與Pa，t、Pt的關(guān)系為：

式中：α 為指數(shù)參數(shù)，表示由降雨轉(zhuǎn)化為徑流的轉(zhuǎn)換率，作為參數(shù)需要率定。

則t時(shí)段降雨所形成的徑流深Rt為：

則t時(shí)段由降雨形成的總流量Qt為：

式中：A為該計(jì)算水電站所控制流域面積；As為計(jì)算水電站上游電站控制流域面積。

1.2 上游水庫(kù)出庫(kù)流量的修正方法

上游出庫(kù)流量數(shù)據(jù)是由出力、水位、閘門開(kāi)度等數(shù)據(jù)計(jì)算，在計(jì)算過(guò)程中經(jīng)常會(huì)用到曲線插值計(jì)算，并且往往將一個(gè)時(shí)間點(diǎn)測(cè)得的瞬時(shí)值作為該時(shí)段的平均值用于計(jì)算，由于曲線的精度問(wèn)題以及計(jì)算方法所產(chǎn)生的誤差，導(dǎo)致上游出庫(kù)流量不太準(zhǔn)確，同時(shí)水電站之間還存在區(qū)間來(lái)水或引水問(wèn)題，故以下游入庫(kù)流量數(shù)據(jù)為標(biāo)準(zhǔn)，設(shè)置修正系數(shù)aQc 對(duì)上游水庫(kù)的出庫(kù)流量進(jìn)行修正，即假設(shè)t時(shí)段能夠到達(dá)下游的實(shí)際流量：

式中：Qc，t為t時(shí)段上游出庫(kù)流量。

1.3 基于匯流系數(shù)法的日時(shí)段匯流計(jì)算方法

上游水庫(kù)出庫(kù)流量的河道匯流演算和區(qū)間匯流演算，可以應(yīng)用趙人俊教授于1962年推導(dǎo)出的馬斯京根分段連續(xù)演算匯流曲線［12］，即n個(gè)單元河段單位入流的出流過(guò)程為：

其中：

如果將C0、C1、n 作為參數(shù)（其中n 為整數(shù)，C0、C1為0 到1 之間的小數(shù)），應(yīng)用式（6）和式（7）即可得到匯流系數(shù)過(guò)程φ0，n～φm，n，如圖2所示，為了保證水量平衡，m的取值應(yīng)滿足：

圖2 匯流系數(shù)過(guò)程示意圖Fig.2 Schematic diagram of confluence coefficient process

由于區(qū)間匯流與河道匯流機(jī)制不同，需要分別率定區(qū)間匯流系數(shù)（φ0，n～φm，n）以及河道匯流系數(shù)

1.4 入庫(kù)流量合成計(jì)算方法

t時(shí)段由區(qū)間降雨形成的入庫(kù)流量ΔQR，t：

t時(shí)段由上游出庫(kù)流量形成的入庫(kù)流量ΔQc，t：

則t時(shí)段計(jì)算水庫(kù)的日平均入庫(kù)流量QAr，t為：

1.5 參數(shù)率定與檢驗(yàn)

該入庫(kù)流量預(yù)報(bào)方法需要率定的參數(shù)見(jiàn)表1和表2。

表1 產(chǎn)流及流量修正參數(shù)Tab.1 Parameter table of runoff generation and runoff conversion

表2 匯流參數(shù)Tab.2 Parameter table of runoff concentration

入庫(kù)流量數(shù)據(jù)的不準(zhǔn)確體現(xiàn)在與實(shí)際流量之間的誤差以及鋸齒現(xiàn)象，但隨著計(jì)算時(shí)間推移，誤差能在一定程度上相互抵消，故在長(zhǎng)計(jì)算期內(nèi)，歷史入庫(kù)流量數(shù)據(jù)仍能代表站點(diǎn)的實(shí)際入流過(guò)程以及實(shí)際入庫(kù)總水量，故本文以歷史入庫(kù)流量數(shù)據(jù)作為擬合標(biāo)準(zhǔn)。根據(jù)資料情況劃分率定期與檢驗(yàn)期，考慮到各個(gè)參數(shù)之間共同對(duì)結(jié)果造成影響，所以本文采用分步率定，步驟如下：

（1）將計(jì)算期內(nèi)修正后上游出庫(kù)總水量與區(qū)間產(chǎn)流總水量之和與實(shí)際入庫(kù)總水量的差值絕對(duì)值最小作為目標(biāo)函數(shù)，用于率定產(chǎn)流及出庫(kù)流量修正參數(shù)（見(jiàn)表1），其形式為：

式中：Wsy(t)指t 時(shí)段內(nèi)經(jīng)過(guò)修正后的上游出庫(kù)水量；Wcl(t)指t時(shí)段內(nèi)降雨產(chǎn)流所形成的水量；Wrk(t)指t 時(shí)段內(nèi)的實(shí)際入庫(kù)水量。

（2）將計(jì)算期內(nèi)最后得到的各時(shí)段模擬入庫(kù)流量與各時(shí)段實(shí)測(cè)入庫(kù)流量的納什效率系數(shù)作為擬合標(biāo)準(zhǔn)，用于率定匯流參數(shù)（見(jiàn)表2），其形式為：

式中：Qm(t)表示t時(shí)段內(nèi)的模擬入庫(kù)流量；Qs(t)表示t時(shí)段內(nèi)的實(shí)測(cè)入庫(kù)流量表示計(jì)算期內(nèi)實(shí)測(cè)平均入庫(kù)流量。

2 實(shí)例應(yīng)用

2.1 實(shí)例概況

郁江位于云南省東部和廣西壯族自治區(qū)南部，是西江流域一級(jí)支流，發(fā)源于云南省文山壯族苗族自治州廣南縣，于廣西桂平市注入西江潯江段，郁江流域水庫(kù)群有12 座主要電站，其中西津電站位于廣西橫縣，是一座以發(fā)電為主，兼顧航運(yùn)的水電站，在南寧市下游約100 km，其上游有屬右江支流的金雞灘水電站和屬左江支流的山秀水電站，具體關(guān)系見(jiàn)圖3。

圖3 西津水電站上下游關(guān)系示意圖Fig.3 Relationship between upstream and downstream of Xijin hydropower

2.2 參數(shù)率定及統(tǒng)計(jì)結(jié)果

選取西津水庫(kù)2013-01-01-2020-01-01 的日降雨量、入庫(kù)流量以及上游兩座水庫(kù)的出庫(kù)流量，以2013-01-01-2017-12-31 作為率定期，2018-01-01-2020-01-01 作為檢驗(yàn)期，以納什效率系數(shù)作為模型的評(píng)判標(biāo)準(zhǔn)，由于西津上游存在兩個(gè)電站，故率定兩套流量修正參數(shù)以及河道匯流參數(shù)。

本文采用遺傳算法［13］根據(jù)節(jié)1.5目標(biāo)函數(shù)及率定步驟進(jìn)行參數(shù)率定，產(chǎn)匯流及流量修正系數(shù)參數(shù)結(jié)果見(jiàn)表3、4 以及圖4、5及圖6。

圖4 區(qū)間流域匯流系數(shù)Fig.4 Concentration coefficient of interval watershed

圖5 金雞灘水電站出庫(kù)河道匯流系數(shù)Fig.5 Concentration coefficient of Jinjitan hydropower station outflow channel

圖6 山秀水電站出庫(kù)河道匯流系數(shù)Fig.6 Concentration coefficient of Shanxiu hydropower station outflow channel

由率定參數(shù)計(jì)算西津水電站各時(shí)期及逐年納什效率系數(shù)，結(jié)果見(jiàn)表5 和圖7 所示。2013年和2014年的徑流模擬結(jié)果見(jiàn)圖8和圖9所示。

圖7 各年納什效率系數(shù)Fig.7 NSE of each year

圖8 西津水庫(kù)2013年徑流模擬結(jié)果Fig.8 Simulation results of runoff in Xijin Reservoir in 2013

圖9 西津水庫(kù)2014年徑流模擬結(jié)果Fig.9 Simulation results of runoff in Xijin Reservoir in 2014

2.3 成果分析

本方法通過(guò)對(duì)上游出庫(kù)流量進(jìn)行修正，并基于降雨進(jìn)行區(qū)間產(chǎn)流計(jì)算，最后將這兩部分采取匯流計(jì)算得到模擬入庫(kù)流量數(shù)據(jù)，計(jì)算步驟簡(jiǎn)潔且符合水文學(xué)基本原理，率定得出的降雨產(chǎn)流參數(shù)以及匯流系數(shù)很好地代表了該流域區(qū)間及河道的匯流特征，由于2016年短時(shí)強(qiáng)降雨較多，使得該年的模擬精度相較其他時(shí)段略有下降，但在率定期本方法模擬計(jì)算納什效率系數(shù)達(dá)到了0.943，檢驗(yàn)期納什效率系數(shù)達(dá)到了0.914，取得了很好的計(jì)算效果，同時(shí)本方法有效解決了梯級(jí)水電站之間的水量不平衡問(wèn)題，形成的流量過(guò)程也避免了傳統(tǒng)方式所產(chǎn)生的鋸齒現(xiàn)象。

表4 匯流參數(shù)率定表Tab.4 Calibration table of runoff concentration parameter

表5日尺度納什效率系數(shù)Tab.5 NSE in each period of daily scale

2.4 參數(shù)敏感性分析

根據(jù)各參數(shù)特性，以率定參數(shù)為基礎(chǔ)，采取固定其他參數(shù)，在一定范圍內(nèi)變動(dòng)某一參數(shù)方式進(jìn)行敏感性分析，參數(shù)變動(dòng)范圍見(jiàn)表6、7，敏感性分析結(jié)果見(jiàn)圖10、11和圖12。

圖10 納什效率系數(shù)與消退系數(shù)K關(guān)系Fig.10 Relationship of regression coefficient K and NSE

圖11 納什效率系數(shù)與匯流參數(shù)n關(guān)系Fig.11 Relationship of runoff concentration coefficient n and NSE

圖12 納什效率系數(shù)與產(chǎn)流參數(shù)、出庫(kù)修正系數(shù)及部分匯流參數(shù)關(guān)系Fig.12 Relationship of runoff generation，outflow revised，part of runoff concentration coefficient and NSE

表6 消退系數(shù)及部分匯流參數(shù)取值Tab.6 Values of regression coefficient and partial concentration parameters

通過(guò)參數(shù)敏感性分析研究可得到以下結(jié)論：

（1）由于實(shí)際資料中西津水庫(kù)上游出庫(kù)水量之和小于該站入庫(kù)水量之和，缺失的水量由降雨過(guò)程來(lái)彌補(bǔ)。消退系數(shù)K最為敏感，其對(duì)流域降雨產(chǎn)流過(guò)程具有重要影響，納什效率系數(shù)隨前期最大影響雨量Pa，M、產(chǎn)流指數(shù)α變動(dòng)所產(chǎn)生的變幅也較大。

（2）山秀站年平均流量要高于金雞灘站年平均流量，故山秀水電站出庫(kù)流量修正系數(shù)aQc*、山秀水電站出庫(kù)河道匯流參數(shù)n*變動(dòng)所產(chǎn)生的變幅要比金雞灘站對(duì)應(yīng)參數(shù)變動(dòng)所產(chǎn)生的變幅稍大。

表7 最大前期雨量、產(chǎn)流指數(shù)、出庫(kù)修正系數(shù)以及部分匯流參數(shù)取值Tab.7 Values of PaM，α，and partial concentration parameters

（3）以納什效率系數(shù)為標(biāo)準(zhǔn)，相較于上游出庫(kù)部分，降雨產(chǎn)流及匯流部分對(duì)整個(gè)模擬過(guò)程有更強(qiáng)的影響。

3 結(jié) 論

通過(guò)理論分析及在廣西郁江流域西津水庫(kù)的實(shí)例計(jì)算，得到如下結(jié)論：

（1）本文提出的日時(shí)段入庫(kù)流量預(yù)報(bào)方法以水量平衡理論、API模型以及馬斯京根匯流系數(shù)法為理論基礎(chǔ)，以上游出庫(kù)流量和區(qū)間降雨產(chǎn)流作為入庫(kù)流量?jī)蓚€(gè)組成部分，并通過(guò)出庫(kù)流量修正解決上游出庫(kù)流量數(shù)據(jù)不準(zhǔn)確以及區(qū)間來(lái)水或取水情況不明確等實(shí)際問(wèn)題，符合水文學(xué)基本原理，且對(duì)現(xiàn)有資料運(yùn)用充分，有效解決了梯級(jí)水電站之間水量不平衡問(wèn)題以及傳統(tǒng)方式所產(chǎn)生的流量過(guò)程鋸齒現(xiàn)象。

（2）參數(shù)敏感性分析結(jié)果表明，消退系數(shù)K對(duì)于納什效率系數(shù)值有著重要影響，為參數(shù)選取與率定提供了必要指導(dǎo)。

（3）本方法所需率定參數(shù)少，率定速度快，在實(shí)例應(yīng)用中，率定期模擬計(jì)算納什效率系數(shù)達(dá)到了0.943，檢驗(yàn)期納什效率系數(shù)達(dá)到了0.914，實(shí)現(xiàn)了很好的模擬效果，對(duì)水電站水文預(yù)報(bào)以及水庫(kù)調(diào)度方案的制定與實(shí)施具有積極意義和重要參考價(jià)值。