鄭金德，粟 忠，范玉潔，許景璇，李勇軍，蘇 釗，吳樹誠，方榮杰

(1.桂林理工大學，廣西環(huán)境污染控制理論與技術重點實驗室，廣西 桂林 541004；2.桂林理工大學，巖溶地區(qū)水污染控制與用水安全保障協(xié)同創(chuàng)新中心，廣西 桂林 541004；3.桂林水文中心，廣西 桂林 541001)

1 引言

桂林市以其獨特的人文和山水自然風光而聞名于世，已被評為中外著名的國家級風景和旅游城市。但是桂林市每年汛期的暴雨洪水幾乎都會給這座秀美的城市帶來不同程度的洪澇災害，而且在枯水期漓江的水量又嚴重不足，給桂林市的城鄉(xiāng)飲用水、漓江旅游用水等造成一定的困擾。隨著旅游產(chǎn)業(yè)以及經(jīng)濟的發(fā)展，桂林市對防汛的要求也越來越高。桂林市急需準確的雨情、水情信息，為桂林市的防汛工作提供相應的參考依據(jù)。

目前，關于漓江流域河段流量演變規(guī)律的研究相對較少。馬斯京根法主要應用于河段流量演算中，由McCarthy于1938年第一次提出。該方法從提出至今，一直被國內(nèi)外廣大從事水文工作的人員廣泛使用。馬斯京根模型的主要參數(shù)為K、x，對演算參數(shù)求解的方法大多采用人工試算法。在運用馬斯京根模型進行流量演算的過程中，發(fā)現(xiàn)演算出來的結果精度不是很高，為此許多科研工作者對馬斯京根法進行了改進。在對K、x求解的方面，JyS.Wu等提出了改進圖解法， O，Donnell提出了一種直接三參數(shù)馬斯京根法，運用矩陣計算方法直接推求K、x的值，且矩陣法較圖解法推求出來的結果更合理[1,2]。最小二乘法[3,4]是基于數(shù)學層面的計算方法，其概念清晰，在水文統(tǒng)計過程中被普遍應用。例如，董曉華等[5]采用最小二乘法，對馬斯京根法河道流量演算精度的影響做了探究，并提出了采用最小二乘法有助于提高河道流量計算的精度；衛(wèi)小蓉[6]、薛夢楠[7]等也基于最小二乘法對不同河道的流量進行模擬，二者均取得了較好的演算結果。

本文針對漓江流域上游三街水文站至桂林水文站段干支流，選取三街水文站、靈川水文站和桂林水文站的逐日實測流量數(shù)據(jù)，在有約束條件下推導出流量演進系數(shù)的最優(yōu)估計值，客觀地建立適用于漓江流域上游的馬斯京根流量演算模型，并探究該河段的流量演進規(guī)律，進而為漓江流域防洪工作地進行及水資源的合理調(diào)度提供相應的理論依據(jù)。

2 材料與方法

2.1 流域概況

漓江流域上游桂林水文站以上的集水區(qū)總面積為2762 km2，平均海拔約150 m，河段主要干流長105 km，流經(jīng)桂林市區(qū)的河段長4.93 km，整個流域的地勢自西北向東南傾斜[8]。且桂林水文站以上流域內(nèi)雨水充足，水資源豐富；流域多呈羽狀，支流短，河道類型為峽谷式河槽，河床的坡降大，漓江流域上游存在多個暴雨區(qū)。桂林市氣候?qū)賮啛釒Ъ撅L氣候，多年平均降雨量為1889.4 mm，汛期大多集中在3～8 月份[9]，汛期降雨量約占年降水量的75%，在汛期內(nèi)漓江流域降水較多，頻繁發(fā)生水澇災害。桂林市區(qū)三街水文站至桂林水文站主要河流基本特征見表1，河流及水文站網(wǎng)分布見圖1。

表1 漓江桂林水文站以上主要河流基本特征

2.2 數(shù)據(jù)來源

本次模型探究所使用的資料由桂林水文中心所提供，主要研究河段為漓江上游的三街水文站至桂林水文站，該河段內(nèi)甘棠江為漓江的一級支流，甘棠江上布設有靈川水文站，其它支流并未布設有相應的水文測站。三街水文站處于漓江上游河段的干流上，該水文站以上的流域內(nèi)建有斧子口水庫、川江水庫和小溶江水庫，各水庫所處位置見圖1，且流域以上支流眾多，地形地勢復雜。進行流量演進的前提是演算河道上下游的流量資料要相互對應，三街水文站為新建的水文站，所具有的流量監(jiān)測數(shù)據(jù)年份較少，并且由前文可知桂林市的汛期主要集中在每年的3～8月份。因此，本次研究所使用的數(shù)據(jù)為三街水文站、靈川水文站和桂林水文站2019、2020年3～8月份的逐日實測流量數(shù)據(jù)。以2019年3～8月份的逐日實測流量數(shù)據(jù)進行模型參數(shù)的率定，以2020年3～8月份的逐日實測流量數(shù)據(jù)對模型進行了驗證。

圖1 桂林市區(qū)河流及水文站網(wǎng)分布

2.3 馬斯京根流量演算模型

馬斯京根模型主要是一種在河道流量(洪水)演算中被廣泛應用的模型，它利用水量平衡方程式代替連續(xù)方程，以槽蓄方程代替復雜的水動力方程，極大地簡化了計算過程[10]。本文依據(jù)漓江流域上游河段干支流洪水匯集的特性，采用先合后演[11]的流量演算方法，將干、支流斷面上三街水文站與靈川水文站的實測流量相互疊加在一起，以此為基礎作為漓江流域上游河段的總入流量，按照無支流河段的演算方法進行連續(xù)馬斯京根演算。

2.3.1 馬斯京根方程

在第1、2時段內(nèi)對水量平衡方程和槽蓄方程[12,13]進行差分求解，可推導出流量演算方程式：

Q2=C0I2+C1I1+C2Q1

(1)

式(1)中：

(2)

(3)

(4)

C0+C1+C2=1.0

(5)

因此，在使用馬斯京根方程時，只要確定了參數(shù)K、x與△t的值，就可根據(jù)式(2)、(3)、(4)分別推導得出C0、C1、C2的值，然后再根據(jù)流量演算方程式(1)演算出下斷面的最終流量過程Q(t)。

顯然，確定模型演算參數(shù)x、K的值，成為使用馬斯京根流量演算模型的關鍵。模型參數(shù)x、K值的確定，一般根據(jù)馬斯京根法W～Q`呈單一直線的基本假設對x、K的值進行確定，傳統(tǒng)的確定x、K值的方法有試錯法、最小二乘法、最小面積法等。利用傳統(tǒng)的方法推求x、K值是使W～Q`呈單一直線的最優(yōu)估計值，但是根據(jù)所推導出的x、K值計算的流量演算系數(shù)未必能使演算出流量與實測出流量的擬合誤差最小[14]。而基于最小二乘法，直接對流量演算系數(shù)C0、C1、C2的最優(yōu)估值進行推導，最后再根據(jù)相應公式進行x、K值的反算，再利用所推求出的系數(shù)進行河道流量演算，則可取得較好的演算結果。

2.3.2 模型參數(shù)優(yōu)選依據(jù)

模型參數(shù)的優(yōu)選以演算出流量與實測出流量的離差平方和最小為依據(jù)，構造目標函數(shù)[15]表達式(6)與約束方程表達式(7)：

(6)

(7)

式(6)中S為離差平方和，n為實測流量點數(shù)。

由式(2)、(3)、(4)可知C0、C1、C2必須滿足式(5)的要求。因此可以根據(jù)拉格朗日乘數(shù)法構造一個非線性規(guī)劃問題，由式(6)、式(7)構成等式約束方程，再構造拉格朗日函數(shù)，將等式所需要約束的極值問題轉化為一個無線性約束的極值問題，其極值函數(shù)為式(8)：

L(λ,C0,C1,C2)=

(8)

式(8)中，λ為拉格朗日乘數(shù)。式(8)達到最小的C0、C1、C2值即為式(7)的解，根據(jù)極值存在的必要條件，則有式(9)、(10)：

(9)

(10)

經(jīng)處理化簡后得式(11)、(12)、(13)：

(11)

(12)

(13)

聯(lián)立式(5)與式(11)、(12)、(13)即可求出流量演算系數(shù)C0、C1、C2的最優(yōu)估值。然后根據(jù)式(14)、(15)即可反求出x、K的值。

(14)

(15)

2.4 模型的驗證方法

本文主要以確定性系數(shù)[16]和相關系數(shù)[17]作為模型預測精度的評價標準，確定性系數(shù)DC按式(16)計算，相關系數(shù)R按式(17)計算：

(16)

(17)

3 結果與分析

3.1 馬斯京根流量演算系數(shù)的推求

選取三街水文站、靈川水文站至桂林水文站2019年3～8月份的逐日實測流量資料，代入式(8)對馬斯京根模型的流量演算系數(shù)C0、C1、C2進行推求，然后根據(jù)式(14)與式(15)推求出K、x的值；在對傳統(tǒng)的馬斯京根模型的流量演算系數(shù)進行推求，兩者所推求出的流量演算系數(shù)值見表2。把推求得到的C0、C1、C2值代入式(1)進行模型率定期的流量演算，演算出來的流量值與實測流量值成果見圖2。

表2 河段流量演算系數(shù)值

圖2 演算流量與實測流量對比(率定期)

由圖2可知，利用上文所介紹的方法與傳統(tǒng)的方法對馬斯京根法的流量演算系數(shù)C0、C1、C2進行最優(yōu)取值，然后在進行河道流量演算，得出實測流量值與演算流量值具有大致相同的變化趨勢，演算出的流量值與水文站實測的流量值擬合較好。

3.2 模型檢驗

把由2019年3～8月份的逐日實測流量數(shù)據(jù)推求出的河段流量演算系數(shù)值，代入到三街水文站至桂林水文2020年3～8月份的逐日流量數(shù)據(jù)進行模擬，模擬結果見圖3。從圖2、圖3可以看出，無論是在模型的率定期或者驗證期，所演算得出的流量值與實際監(jiān)測得到的流量值之間是具有較大偏差的。在流量較小的階段，演算值與實測值之間可以達到很好的擬合效果，而在流量較大的階段，演算值與實測值之間盡管存在著大體一致的變化趨勢，但是仍然存在著一定的偏差。產(chǎn)生誤差的主要原因是，在進行流量演算的過程中并未加入桃花江與南溪河的流量進行運算，而由表1 可知桃花江的多年平均流量為12.2 m3/s，如若加上桃花江與南溪河的流量進行演算，則可取得更好的結果。

由圖3和表3可以看出，把由2019 年的流量數(shù)據(jù)確定出來的馬斯京根演算系數(shù)C0、C1、C2值，代入到三街水文站2020 年的流量監(jiān)測數(shù)據(jù)進行演算，演算出的結果與桂林水文站實測的流量結果擬合程度較好，二者都達到了乙級標準，可為洪水預報工作提供相應的參考依據(jù)。

由表3可以明顯看出，基于改進的馬斯京根法所推求出來的流量值，不論是在率定期還是驗證期，其確定性系數(shù)與相關系數(shù)都比由馬斯京根法推求出來的高，但是對比二者的平均相對誤差值，改進的馬斯京根法卻比馬斯京根法的高。從圖2、圖3可以看出不論是在模型的率定期還是驗證期，改進的馬斯京根法與實測流量值之間具有更好的擬合效果，曲線的變化趨勢更加一致。改進的馬斯京根法是基于最小二乘法，直接推導出馬斯京根模型的河道流量演算系數(shù)，其計算過程簡單，取得的擬合效果更好。而由馬斯京根法推求的流量演算系數(shù)，計算過程繁瑣，演算系數(shù)的推求受計算者的影響較大，所取得的演算結果往往與實測的結果擬合效果較差。因此，改進的馬斯京根法比馬斯京根法更適合應用于漓江上游的流量演算。

圖3 演算流量與實測流量對比(驗證期)

表3 模型精度評定結果

4 結論

(1)基于改進的馬斯京根流量演算法，推求出了適合漓江上游流域的流量演算系數(shù)C0、C1、C2(0.37；0.62；0.01)的值。改進的馬斯京根流量演算法較傳統(tǒng)的馬斯京根法取得更好的擬合效果，在模型驗證期內(nèi)，改進的馬斯京根法所推求出的結果確定性系數(shù)達到了0.812，相關系數(shù)達到了0.940，平均相對誤差為-0.195，達到了乙級精度?？蔀楣鹆质械暮樗A報作業(yè)提供相應的參考依據(jù)。

(2)改進的馬斯京根法較傳統(tǒng)的馬斯京根法，減少了繁瑣的參數(shù)試算過程，提高了參數(shù)的選取效率，同時也增強了演算結果的可靠性。

(3)本次模型探究所使用的數(shù)據(jù)是水文測站每日監(jiān)測的流量數(shù)據(jù)，其精度可能會存在一定的問題，且使用的數(shù)據(jù)資料年份較少，會給參數(shù)的確定帶來一定的影響，再加上有一些支流上并未布設水文測站，所以在今后的工作中應及時掌握流域的來水情況，及時地對模型參數(shù)進行相應的修正。