基于分布式水文模型的阿倫河流域降雨-徑流計算

黃金柏， 王 斌， 溫佳偉， 朱士江， 甄自強

(1.揚州大學 水利與能源動力工程學院， 江蘇 揚州 225009； 2.東北農業大學 水利與建筑學院， 黑龍江

哈爾濱 150030； 3.揚州大學 信息工程學院， 江蘇 揚州 225009； 4.三峽大學 水利與環境學院， 湖北 宜昌 443002)

基于分布式水文模型的阿倫河流域降雨-徑流計算

黃金柏1， 王 斌2， 溫佳偉3， 朱士江4， 甄自強1

(1.揚州大學 水利與能源動力工程學院， 江蘇 揚州 225009； 2.東北農業大學 水利與建筑學院， 黑龍江

哈爾濱 150030； 3.揚州大學 信息工程學院， 江蘇 揚州 225009； 4.三峽大學 水利與環境學院， 湖北 宜昌 443002)

摘要：[目的] 為黑龍江省西部半干旱區季節性地表徑流計算提供科學的方法，并推求該區雨季各月的徑流系數。[方法] 基于對阿倫河流域下墊面的實際調查，對土壤垂直剖面模型化，利用GIS構建研究區DEM及河網，以運動波理論的基礎方程式構建降雨—徑流計算方法，以對觀測流量的數值模擬檢驗模型實用性。通過對2012和2013年5—10月降雨—徑流計算結果的分析，分別推求了各月徑流量和徑流系數。[結果] 模型計算精度可以達到誤差基準允許的范圍之內(<0.03)；7，8月的月徑流系數在0.5以上，計算時段內的徑流量分別占各年總降水量的34.2%和34.7%。[結論] 構建的降雨—徑流數值模型，量化了計算時期內各月的徑流系數，適用于對研究區降雨—徑流過程的計算。

關鍵詞：阿倫河流域； 降雨—徑流； 分布式水文模型； 數值計算

當前，在中國的很多地區，水資源已成為制約社會經濟發展和生態環境保護協調發展的一個主要障礙[1-2]。如何實現對地區性水資源的準確評估，使可利用水資源有效保證經濟社會的不斷發展已成為一個亟待解決的課題[3]。近年來，探求可準確計算流域降雨—徑流過程的數值計算方法，從而實現對流域水資源的準確推求是國內外水科學研究領域許多研究者致力的課題。基于地形、地質等流域實際物理條件，結合可描述水文過程的數學方程式構建數值模型，從而實現對流域水資源的準確推求，已成為當前國內外水文水資源研究領域的共識[4-6]。

黑龍江省西部半干旱區是中國重要的糧食生產基地[7]，是中國典型的旱作農業區，旱災一直是該地區典型的自然災害之一[8]。該地區屬于寒溫帶半干旱季風氣候，降雨時空分布不均，冬季嚴寒少雪，春季干旱，秋季雨量集中，不利于農業生產的開展[9-10]。

本研究依托具有黑龍江省西部半干旱區典型氣候、水文及地形特征的阿倫河流域，構建分布式水文模型，在檢驗模型實用性及計算精度的基礎上，對研究區降雨—徑流過程進行計算并對結果進行解析，以期為黑龍江省西部半干旱區季節性水資源的準確評估及開發利用提供基礎方法和數據，從而為該地區水資源有效地支撐糧食生產服務，并為該地區水文模型數值計算平臺的構建提供借鑒方法。

1材料與方法

1.1　研究區概況

阿倫河流域(地理位置為東經 122°04′—124°04′，北緯 47°37′—48°48′)發源于內蒙古自治區境內大興安嶺東南麓，是嫩江右岸一級支流，中下游流經黑龍江省齊齊哈爾市甘南縣梅里斯區匯入嫩江。流域面積6 297 km2，河道長度318 km[11-12]，流域海拔153～198 m，地勢由西北向東南階梯式遞降，其中下游位于松嫩平原西部。地質構造屬于大興安嶺新華夏構造帶，表層為第四系覆蓋層[13]。阿倫河流域處于中高緯度地區，屬寒溫帶大陸性氣候。多年平均氣溫2.0～2.5 ℃，多年平均降水量約為450 mm，年降水分布不均，7，8月降水量占全年降水總量的60%左右。年蒸發量超過1 000 mm，屬于半干旱地區[11-13]。流域冬季寒冷漫長，一般凍土期從當年的11月至次年4月。據鉆孔調查，表層土主要成分是壤土(黑壚土)，其下是亞黏土及黏土層，底層為砂巖層。

1.2　水文觀測

降雨觀測采用翻倒式雨量計(7852M-L10，φ165 mm×H240 mm，地表徑流流量采用水位計 (HM-910-02-309)。徑流觀測斷面(48°1′22″N,S123°36′28″E)的上游集水面積為4 993 km2，占流域總面積的79.3%，自河源到觀測點的主河道長229.25 km。降雨觀測地點為(48°1′31″N, 123°34′7″E)。水文數據的觀測自2012年初開始，由于外界不確定因素的影響，自觀測數據在時間序列上存在不同程度的缺失。基于選取的流域尺度，對降雨—徑流過程的計算需要多點數據支撐。除了自觀測水文數據，其他數據來源為阿倫河流域那吉水文站(48°5′42″N,123°28′8″E)的降雨資料和黑龍江省甘南縣氣象局提供的降雨數據(數據類型為多點平均)。對降雨—徑流過程計算時，基于對多源數據的綜合分析，進行合理的選取及應用。如同期發生全流域降雨時，取多點時間序列雨量的平均值；流域內發生區域性降雨時，結合流域數字高程模型(DEM)和河網上各分布式小流域的具體位置，使降雨與流域空間位置對應，進行計算。觀測的地表徑流水位數據轉換為流量數據的方法為：通過橫斷測量的方法，確定觀測斷面的形狀，利用曼寧平均流速公式，可以將觀測時段內的水位數據，轉換成流量數據[14]，作為基礎數據，與模型計算(計算徑流)結果進行比較。

1.3　流域土壤垂直剖面模型化

基于對研究區地形及基本水文地質條件的調查結果，對土壤垂直剖面的分層情況模型化，即構建自地面開始至地下某一含水層的土壤垂直剖面模型，作為承載雨水降落到地面后在垂直方向的運動載體，使水的運動方式受下墊面實際物理條件的約束。以手工結合簡單機械鉆孔(直徑5 cm)的方式，對流域內多點的土壤垂直剖面分層情況進行調查，調查深度至潛水含水層底面的弱透水層。基于鉆孔調查并結合對流域內現有多口潛水井的井深和水深季節性變動的調查結果，篩選土壤垂直剖面分層結構的特征參數，構建的土壤垂直剖面模型如圖1所示。模型由坡面和河道構成，坡面區間自地面開始至潛水含水層底部(砂巖層表面)被分成兩層，第1層由壤土和亞黏土構成，其厚度為5—10 m，潛水含水層位于第1層的下部。其下是黏土層(第1個弱透水層)，厚度超過20 m，黏土層的下部是第一承壓含水層和砂巖層。河道區間自河床表面至砂巖層表面被開發成1層。坡面區間的潛水含水層自由水面在坡面下端與河道徑流相通，模型上各水流入和流出成分如圖1所示。

1.4　降雨-徑流計算模型構建

1.4.1計算公式運動波模型被廣泛地應用于對降雨—徑流過程的計算[15-16]，該模型以物理性參數描述流域以及水流運動過程，不但可以對地表徑流進行準確計算，而且可以通過達西公式描述和計算滲透以及地下水徑流過程[17-18]，所以，采用運動波理論的基礎方程式構建降雨—徑流計算方法，以坡面為例，給出計算公式如下：

地表徑流連續方程式：

(1)

式中：dt——計算的時間步長為1 s(s)； dx——水流方向上計算的空間步長為1 km(km)；h——水深(m)；q——單寬流量(m2/s)；r——降雨(m/s)；f1——第1層土壤平均滲透速度(m/s)；α——為引入系數，用于表達在不同水深條件下第1層的實際入滲速度。下同。

圖1　土壤垂直剖面模型化示意圖

曼寧平均流速公式：

(2)

式中：n——糙率(s·m-1/3)；R——水力半徑(m)；I——坡度。

在上述設定的時間步長dt和水流方向的空間步長dx條件下，經驗證，數值計算滿足穩定計算的條件，即Umaxdt≤dx(Umax為最大流速，m/s)。

滲透流(地下水)連續方程式(以第1層為例)：

(3)

達西定律：

(4)

第2層用于地下水計算的公式與第1層采用的公式相同，只是由于層號的不同導致部分參數腳標發生變化；河道與坡面采用同一組方程式進行計算，只是由于個別水流入或流出成分的變化導致公式在形式上略有差別。

1.4.2初始條件和邊界條件作為計算的初始條件，計算開始時刻，地表徑流水深在各最末級分布式小流域源點處被設為0，地下水的初始(邊界)水深由計算開始時刻的實際調查結果確定。垂直方向的邊界條件如第1層和第2層的滲透速度分別根據地表徑流和地下水的水深由公式(1)和公式(3)經過計算確定，如公式(1)中的α等于水深(hi)與計算時間步長(dt)的比值(α=hi/dt)，當hi/dt

1.4.3計算參數用于計算的主要參數通過實際調查、實驗以及利用流域數字高程模型(DEM)確定。例如：河道和坡面坡度利用DEM結合水準測量確定；土層厚度以鉆孔調查并參照現有潛水井實際條件確定；表層土壤垂向滲透系數以及土壤孔隙率通過滲透實驗確定。因為同類參數隨河網上各分布式小流域的不同而有所差異，土壤水力學特性參數呈現出時空變動的特性。另外，流域植被條件的季節性變化明顯，如在作物的生長期，地面幾乎被作物覆蓋，而作物的蒸散發量也因生長期內氣象條件和作物本身的生物條件而有所不同，而在非作物生長期，超過50%的面積為裸地。由于地面植被季節性變化明顯，導致蒸散發的準確計算難于實現，在降雨—徑流過程計算時，引入一個損失系數近似地評價蒸散發(表1)。對于第2層土壤的有效孔隙率、垂向滲透系數等難以利用上述方法確定的參數，采用首先給這些未定參數賦予合理的初值，通過模型海量計算考察數值模擬誤差的方法來確定[19]。表1列出了模型主要參數的特征值，即各參數處于相應的量級。

表1　主要計算參數特征值索引

2結果及討論

2.1　模型檢驗

2.1.1數值模擬利用計算機軟件Fortran開發計算程序(數值模型)，通過對觀測地表徑流的數值模擬來檢驗模型的實用性和計算精度。計算區域為地表徑流觀測斷面的上游集水區(面積4 993 km2)，以該區域為依托對降雨—徑流過程的計算結果，在尺度上可以較充分地反應阿倫河流域的地表徑流特征。計算區域的河網基于該區的數字高程模型(DEM)利用GIS中的ArcMap生成(圖2a)，對河網分割后再集中化的各分布式小流域的空間連接關系如圖2b所示。

圖2　計算區域的河網及各小流域空間連接關系

對觀測斷面徑流過程進行模擬，結果(觀測流量和計算結果的比較)如圖3所示。由圖3可知，計算流量與觀測值之間的擬合效果較好，在計算時段內二者之間沒有明顯的差別，計算結果很好地再現了觀測徑流發生的過程。

2.1.2誤差分析采用式(5)作為誤差判斷的基準[20-21]，該基準要求觀測值與計算值之間的誤差小于0.03。

(5)

式中：Er——誤差；n——計算次數；Qo(i)——i時刻的觀測流量(m3/s)；Qc(i)——i時刻的計算流量(m3/s)；Qop——計算時段內的最大觀測流量(m3/s)。

對數值模擬結果進行誤差分析，由誤差計算結果可知，2013年6月11—13日降雨—徑流模擬結果的誤差為0.007(圖3a)；2013年7月15—18日模擬結果的誤差為0.01(圖3b)，其中，7月15日12時至7月16日6時長歷時降雨事件發生期間降雨—徑流過程數值模擬的誤差為0.018；對其他時段的降雨-徑流模擬結果進行隨機誤差計算，結果均在誤差判斷基準允許的范圍之內(<0.03)。由誤差分析的結果可知，本研究開發的分布式降雨—徑流數值模型適用于對阿倫河流域降雨—徑流過程進行計算。

a 20130611—20130613　　　　　　　　　　　　　　　　　　　b 20130715—20130718

圖3降雨-徑流數值模擬結果

2.2　結果及分析

研究流域所在地區自11月至次年4月為凍土期(包括凍融期)，在此其間降水的形式主要為降雪，沒有集中降雨發生。利用2012(降水量522 mm)，2013(降水量651 mm)年觀測的降雨數據，通過模型計算得到的降雨—徑流過程線如圖4所示。計算時段為每年有集中降雨事件發生的5—10月，涵蓋研究區整個雨季(6—9月)。因為計算時間步長為1 s，數據點過密，流量計算結果以“日均值”、降雨數據以1 d序列給出。

根據雨量觀測結果，2012和2013年計算時段內的降雨量分別為476.2和622.9 mm，分別占當年降雨總量的91.2%和95.6%；對計算時段內各月徑流結果進行統計，并推求逐月徑流系數，結果如表2所示。

圖4　降雨-徑流計算結果

項目 2012年5月6月7月8月9月10月2013年5月6月7月8月9月10月降雨量/mm36.6091.20135.2053.80120.8038.6022.10189.70218.3077.6090.7024.50徑流量/107m33.8211.5037.7020.6011.304.162.541.8156.2029.2016.606.58徑流系數0.210.250.560.770.190.220.230.190.520.750.370.54

由表2結果可知，在計算時段內，2012，2013年5，6月徑流系數較小，因為次降雨事件所產生的徑流量受降雨強度、降雨量及降雨歷時等因子以及下墊面條件(如植被、表層土壤含水率)的影響[22]，該時期內降雨事件的次降雨量不大且雨強均值較低，另一方面，下墊面具有較強的入滲能力同時伴隨較強烈的蒸發也是導致該時期徑流系數較小的原因。7，8月的月徑流系數較大，其值超過0.5，其中8月為最大，其主要原因為2012和2013年的7月降雨量都為當年最多，其前期降雨使土壤含水量增高，隨著降雨事件的持續發生，產生了較為集中的徑流量；而8月發生的降雨多為降雨量相對集中的降雨事件，同時，8月陰雨天較多，空氣濕度相對較大導致蒸散發減弱，從而導致徑流量相對于其他月份更為集中。2012年9和10月徑流系數較小，但2013年同期月徑流系數相對較大，其主要原因是2012年此間發生的降雨多為次降雨量不大且雨強較小的降雨事件，而2013年同期發生了降雨量相對集中，雨強較大的次降雨事件，從而導致兩年在該時段內徑流系數有較大的差別。

基于降雨—徑流數值計算結果，可以推求出2012和2013年5—10月的徑流量，分別占當年總降雨量的34.2%和34.7%。2012和2013年5—10月的降雨量均占當年總降水量的90%以上，年內其他時間的降雨(水)量較小，總量低于年總降水量的10%，其間的徑流形式主要為融雪徑流且其占年總徑流量的份額很小，由以上分析可知，2012，2013年的年徑流系數應略高于0.34。

3結 論

(1) 本研究構建的分布式水文數值模型適用于對阿倫河流域降雨—徑流過程進行計算，其計算精度可以達到誤差基準允許的0.03范圍之內。

(2) 基于模型計算結果，分別推求了2012和2013年5—10月的徑流量和徑流系數，兩年7,8月的徑流系數在0.5以上；5—10月的徑流量，分別占當年總降雨量的34.2%和34.7%。

研究結果期待為同一地區季節性地表水資源的準確評估以及水文過程數值計算平臺的構建提供基礎數據。

[參考文獻]

[1]謝家澤,陳志愷.中國水資源[J].地理學報,1990,45(2):210-219.

[2]李中鋒,劉昌明,楊志峰,等.對我國水資源問題的哲學思考[J].科技導報,2002(9):39-43.

[3]段春青,邱林,黃強,等.灌區農業水資源承載力研究[J].西北農林科技大學學報:自然科學版,2005,33(4):135-138.

[4]Jain M K, Singh V P. DEM-based modelling of surface runoff using diffusion wave equation[J]. Journal of Hydrology, 2005,302(1/4):107-126.

[5]Du Jinkang, Xie Shunping, Xu Youpeng, et al. Development and testing of a simple physically-based distributed rainfall-runoff model for storm runoff simulation in humid forested basins [J]. Journal of Hydrology, 2007,336(3/4):334-346.

[6]Chua Lloyd H C, Wong Tommy S W, Wang X H. Information recovery from measured data by linear artificial neural networks: An example from rainfall-runoff modeling[J]. Applied Soft Computing, 2011,11(1):373-381.

[7]趙雨森,魏永霞.黑龍江省西部半干旱區土壤水分入滲規律及其模擬研究[J].灌溉排水學報,2008，27(4):110-112.

[8]姜秋香,付強,王子龍.黑龍江省西部半干旱區土壤水分空間變異性研究[J].水土保持學報,2007,21(5):118-122.

[9]魏永霞,張忠學,王立敏.半干旱區坡耕地抗旱保水技術集成對大豆水分利用效率的影響[J].灌溉排水學報,2007,26(6):73-75,82.

[10]王宇,陳麗華,楊啟紅,等.東北半干旱區主要農田防護林樹種蒸騰速率研究[J].水土保持通報,2008,28(4):48-51.

[11]劉新宇,趙嶺,王立剛,等.阿倫河流域水土保持林土壤抗蝕性研究[J].防護林科技,2000(3):21-23.

[12]王繼常,李利.梅里斯阿倫河流域生物治理工程的探討[J].防護林科技,2014(3):100-101.

[13]張玉峰.甘南縣阿倫河流域水環境質量現狀、變化趨勢及防治對策[J].黑龍江環境通報,2012,36(4):35-26,39.

[14]周方錄,黃金柏,王斌.基于柵格的不規則斷面水深—流量關系曲線確定方法[J].水資源研究,2013(2):109-113.

[15]Yomoto A, Islam M N. Kinematic analysis of flood runoff for a small-scale upland field[J]. Journal of Hydrology, 1992,137(1/4):311-326.

[16]Chua Lloyd H C, Wong Tommy S W, Sriramula L K. Comparison between kinematic wave and artificial neural network models in event based runoff simulation for an overland plane[J]. Journal of Hydrology, 2008,357(3/4):337-348.

[17]Cabral M C, Garrote L, Bras R L, et al. A kinematic model of infiltration and runoff generation in layered and sloped soils[J]. Advances in Water Resources, 1992,15(5):311-324.

[18]Sarkar R, Dutta S. Field investigation and modeling of rapid subsurface storm flow through preferential pathways in a vegetated hillslope of Northeast India[J]. Journal of Hydrologic Engineering, 2012,17(2):333-341.

[19]Huang Jinbai, Hinokidani O, Yasuda H, et al. Study on characteristics of the surface flow of the upstream region of Loess Plateau[C]. Annual Journal of Hydraulic Engineering, JSCE, 2008,52:1-6.

[20]黃金柏,王斌,檜谷治,等.耦合融雪的分布式流域“降雨—徑流”數值模型[J].水科學進展,2012，23(2):194-199.

[21]Huang Jinbai, Hinokidani O, Yasuda H, et al. Effects of the check dam system on water redistribution in the Chinese Loess Plateau[J]. Journal of Hydrologic Engineering, 2013,18(8):929-940.

[22]Huang Jinbai, Wen Jiawei, Hinokidani O, et al. Runoff and water budget of the Liudaogou Catchment at the wind-water erosion crisscross region on the Loess Plateau of China[J]. Environmental Earth Sciences, 2014,72(9):3623-3633.

Rainfall-Runoff Calculation of Alun River Basin Based on Distributed Hydrological Model

HUANG Jinbai1, Wang Bin2, WEN Jiawei3, ZHU Shijiang4, ZHEN Ziqiang1

(1.CollegeofHydraulic,EnergyandPowerEngineering,YangzhouUniversity,Yangzhou,Jiangsu

225009,China; 2.CollegeofWaterConservancyandArchitectureofNortheastAgriculturalUniversity,

Harbin,Heilongjiang150030,China; 3.CollegeofInformationofYangzhouUniversity,Yangzhou,Jiangsu225009,

China; 4.CollegeofHydraulic&EnvironmentalEngineeringofThreeGorgesUniversity,Yichang,Hubei443002China)

Abstract：[Objective] To provide a scientific method for surface runoff calculation and to estimate monthly runoff coefficient in rainy season for the semiarid region of the Western Heilongjiang Province.[Methods] The Alun River basin which flows through the Western Heilongjiang Province was chosen as the study area. Modeling for the vertical soil profile was achieved based on investigating the physical properties of the underlying surface. Digital elevation model(DEM) and river channel network of the study area was generated by GIS-ArcMap. Algorithm on rainfall-runoff was established by kinematic wave equations. Model validation was carried out by numerical simulation for the observed flow. Monthly runoff and runoff coefficient were estimated through analyzing runoff calculation results in the period from May to October in 2012 and 2013.[Results] The model accuracy was within the allowable range of the error criterion(<0.03); Monthly runoff coefficient in July and August were more than 0.5; Runoff in the period from May to October accounted for 34.2% and 34.7 % of annual total precipitation in 2012 and 2013, respectively. [Conclusion] The developed numerical model for rainfall-runoff calculation can be applied to the study area and monthly runoff coefficient in the study period.

Keywords:Alun River basin; rainfall-runoff; distributed hydrological model; numerical calculation

文獻標識碼：B

文章編號：1000-288X(2015)01-0224-06

中圖分類號：TV121.1， P343.1

收稿日期：2014-09-25修回日期：2014-10-08

資助項目：黑龍江省教育廳海外學人科研項目“黑龍江省西部半干旱區降雨—徑流數值解析方法的研究”(1251H017)

第一作者：黃金柏(1974—)，男(漢族)，黑龍江省樺南縣人，博士，副教授，主要從事水文過程數值模型、數字流域、河流工學方面的研究。E-mail:huangjinbai@aliyun.com。