人類活動干擾對灞河流域設計洪水的影響研究

慕登睿， 袁衛(wèi)寧， 呂繼強， 羅平平， 范 磊

(1.長安大學 環(huán)境科學與工程學院， 陜西 西安 710061；2.長安大學 旱區(qū)地下水文與生態(tài)效應教育部重點實驗室， 陜西 西安710061)

1 研究背景

我國北方地區(qū)水資源天然稟賦不足，目前水資源緊缺與水環(huán)境惡化問題已限制地區(qū)經(jīng)濟與社會發(fā)展。自20世紀90年代起，區(qū)域城市化進程加快，進一步加劇了水資源供需矛盾。較多學者研究認為[1-3]影響水資源變化最直接因素是氣候條件和土地利用情況[4-6]。當前，半干旱地區(qū)人類活動對流域水文過程的影響日益顯著，增加了徑流變化的復雜性和不確定性[7]。國內(nèi)外研究典型人類活動對流域徑流的影響已逐漸由定性分析轉(zhuǎn)變?yōu)槎坑嬎鉡8-9]。目前研究成果多以洪水要素變化的統(tǒng)計規(guī)律為主[10]，亟需加強人類活動與水文過程相互作用及人類活動對洪水過程影響機理研究。強烈的人類活動不斷改變城市河流的流域下墊面，對不同區(qū)域不同量級洪水過程影響客觀存在，但影響程度不同[11-12]。人類活動干擾類型增多、強度增大，加之水文循環(huán)過程的復雜性，使得人類活動對洪水過程的影響評價結(jié)果帶有很大的不確定性。這些都是下一步人類活動對徑流影響研究中需要重點解決的科學問題[13]。

灞河屬渭河的一級支流，為秦嶺北麓關(guān)中城市群重要水源地之一，是西安市“八水繞長安”中最大的城市型河流。流域內(nèi)農(nóng)業(yè)發(fā)達，灌溉歷史悠久，建有輞灞渠、李家河水庫等諸多引水工程。近年來，流域中游、下游城市新區(qū)建設不斷加快，灞河流域內(nèi)城鎮(zhèn)化面積增加、水利工程建設增多，使得河流水文過程發(fā)生了巨大的變化。本文以灞河為例，利用MCDRM模型模擬1980-2015年間的6期土地利用類型變化和水庫建設對不同重現(xiàn)期的流域設計暴雨洪水的影響，研究結(jié)果對區(qū)域水資源的可持續(xù)利用、防災減災及水利工程的安全高效運行具有重要實踐參考價值。

2 研究區(qū)概況

灞河發(fā)源于秦嶺北麓的藍田縣，由南向北匯入渭河(東經(jīng)109°00′～109°47′、北緯33°50′～34°27′)，總流域面積1 821 km2[14]。流域多年平均降水集水面積主要分布在降水量比較大的秦嶺山區(qū)，約占流域面積的57%，暴雨洪水匯流時間短，具有暴漲猛落的山區(qū)性河流特點。灞河洪枯流量相差懸殊，據(jù)馬渡王水文站(上游控制性水文站)實測資料統(tǒng)計，灞河多年平均流量為15.32 m3/s，年徑流量為4.96×108m3，7-10月份徑流量占全年的56.5%，12-翌年2月占全年徑流量的6.7%，灞河最大洪峰流量為2 160.0 m3/s(1953年8月12日)。流域土地利用類型以耕地、草地和林地為主，集中在流域上中游，2000年后河流下游水域和城市面積增加，截止2015年城市面占流域總面積11.77%。

李家河水庫樞紐位于灞河一級支流輞川河上，地處藍田縣玉川鄉(xiāng)李家河村，水庫壩址以上流域面積362 km2，設計總庫容5 258×104m3，調(diào)節(jié)庫容4 400×104m3，死庫容608×104m3；校核洪水位881.29 m，設計洪水位880.00 m，汛期限制水位880 m，死水位839 m，供水設計流量3.2 m3/s，屬中型水庫。李家河水庫建設的主要任務是以城鎮(zhèn)供水為主，兼顧發(fā)電。2014年，李家河水庫建成后，除汛期有部分洪水量下泄外，其余入庫水量將被水庫調(diào)蓄利用，僅為下游河道預留0.3 m3/s生態(tài)基流量。

采用1975-2010年灞河流域上中游10個氣象站點的氣象數(shù)據(jù)與1個控制性水文站的徑流數(shù)據(jù)，研究區(qū)地理位置及流域水文、氣象站點如圖1所示。

圖1 灞河流域水文、氣象站點位置圖

3 改進的網(wǎng)格分布式模型MCDRM的建立

日本京都大學防災研究所2003年提出的網(wǎng)格化分布式降雨徑流模型CDRM(grid-Cell Distributed Rainfall Runoff Model)用于研究流域暴雨洪水過程[15-17]。CDRM模型可實現(xiàn)水庫實時調(diào)蓄情況下的洪水過程模擬計算，并利用蒙特·卡羅方法計算流域下墊面條件變化的流域產(chǎn)流及匯流參數(shù)合理區(qū)間優(yōu)化，可避免“異參同效”現(xiàn)象，解決空間蒸發(fā)與土壤含水量模擬總和不變，但存在徑流模擬結(jié)果失真的問題。

3.1 MCDRM模型建立

本文根據(jù)灞河流域人類活動實際情況，構(gòu)建不同的人類活動干擾過程模塊，用以概化不同類型人類活動對水文過程的影響機理，提高模型的適用性。本文提出改進后的分布式MCDRM(Modified Cell Distributed Rainfall Runoff Model)模型，并通過蒙特·卡羅方法與人工調(diào)節(jié)參數(shù)方法結(jié)合優(yōu)化模型參數(shù)，對場次洪水過程進行模擬以降低模擬結(jié)果的不確定性。

模型采用運動波模型來計算每一個網(wǎng)格產(chǎn)生的徑流量。每個網(wǎng)格單元必要的模型參數(shù)包括粗糙值和地形屬性。根據(jù)土地使用的情況，分配等效的粗糙值度數(shù)。通過DEM(Digital Elevation Model)定義地形屬性，例如坡度和坡長。該模型模擬了3種側(cè)向流動機制：毛管流、重力流、壤中流。考慮人類活動干擾的MCDRM模型計算流域產(chǎn)匯流過程原理見圖2[18-19]。

圖2 MCDRM模型計算原理示意圖

在MCDRM的每個網(wǎng)格中，當水深低于不飽和流的水深時，(0≤h≤dm),通過有著不飽和傳導系數(shù)km的達西定律模擬流動過程。在毛管孔隙中計算其沿著坡度下降方向的平均速度vm，建立公式為：

(1)

式中：H為水頭，m；x為網(wǎng)格中的水平距離,m；k為非飽和導水率,m/s。

根據(jù)Gardner[20]提出的非飽和導水率計算公式，其表達式為：

(2)

式中：km為毛管水中的飽和滲透系數(shù),m/s;Se為飽和自由度;β為不飽和指數(shù)常數(shù)。

飽和自由度計算公式為：

(3)

式中:θ為體積含水量,%;θm為最大毛細管含水量,%。

因此,公式(1)可以改寫為下式：

(4)

式中：i為山坡梯度;h為水位，m;dm為土層中的毛細管的深度(dm=Dθm),其中D為土層深度,m。通過dm與vm相乘，得到不飽和土層中的網(wǎng)格側(cè)向流的單位寬度qa(0 ≤h≤dm):

(5)

當水深超過不飽和流動的水深，超過的水流作為飽和地下水流，通過飽和水力傳導系數(shù)為ka的達西定律進行模擬。其重力水的坡度下降方向的平均速度va為：

va=kai

(6)

式中：ka為飽和水力傳導系數(shù),m/h。

在飽和情況下，網(wǎng)格的每單位寬度qb的總出流量為(dm≤h≤da)：

qb=dmkmi+(h-dm)kai

(7)

(8)

則可得：

(9)

由h=dm，公式(9)可寫為：

(10)

對于公式(7)來說：

qb=dmkmi+hkai-dmkai

(11)

由公式(11)可得：

(12)

βkmi=kai

(13)

則得到關(guān)系式：

km=ka/β

(14)

當水深超過土層D的有效孔隙率θa(da=Dθa)水流將為表層流，通過曼寧公式模擬平均速度vs：

(15)

式中：n為河道糙率系數(shù)，其值與土地利用類型有關(guān)，m為5/3。

總出流量為地表與地下水量總和：

(16)

(17)

連續(xù)性方程的每個網(wǎng)格單元的坡面流動示意圖如圖4所示，連續(xù)性方程表達式為：

(18)

式中：t為時間，h;φ為傾斜角;r為凈降雨強度,mm/h，模型參數(shù)之間的關(guān)系包括：dm,da,ka,n，β，Lax-Wendroff有限差異組合被用于求解網(wǎng)格單元中每個節(jié)點上的一維運動波動方程。

圖3 單位寬度的流量(q)與深度(h)之間的關(guān)系

圖4 坡面流動示意圖

3.2 MCDRM模型地形數(shù)據(jù)提取

可以由DEM提取水系和子流域特征代表的流域參數(shù)。在DEM中，可以采用柵格的個數(shù)表示流域平面坐標。據(jù)此，就可以計算流域地表上每個柵格的地表坡度、流向、流量分配系數(shù)、集水面積等。

MCDRM模型地形參數(shù)通過DEM處理得到，需要率定產(chǎn)匯流參數(shù)主要有：植被根系區(qū)缺水量、土壤下滲率呈指數(shù)衰減的速率、土壤剛達到飽和時的傳導度、重力排水的時間滯時參數(shù)、根帶最大蓄水能力、根帶土壤的初始缺水量、地表坡面匯流的有效速度、主河道匯流的有效速度。模型各產(chǎn)匯流參數(shù)取值見表1。

表1 MCDRM模型率定參數(shù)表

本文采用蒙特·卡羅方法與手動調(diào)節(jié)參數(shù)相結(jié)合的方式實現(xiàn)參數(shù)變化區(qū)間確定與優(yōu)化。調(diào)節(jié)參數(shù)時，以模擬結(jié)果的相對誤差總和最小和相關(guān)系數(shù)最大構(gòu)建多目標優(yōu)化方程，在給定參數(shù)變化區(qū)間的初始范圍內(nèi)，優(yōu)化確定合理的參數(shù)初始范圍。結(jié)合優(yōu)化結(jié)果，手動調(diào)整參數(shù)優(yōu)化范圍，繼續(xù)優(yōu)化模擬，直到模擬結(jié)果達到滿意為止。多目標優(yōu)化模擬流程見下圖5。

圖5 多目標參數(shù)優(yōu)化計算流程圖

4 灞河流域土地利用變化分析

采用從中國科學院資源環(huán)境科學數(shù)據(jù)中心下載的精度為1 km、年份為 1980、1995、2000、2005、2010和2015年的6期土地利用影像圖。不同設計標準降雨數(shù)據(jù)采用參考文獻[22-23]，降雨數(shù)據(jù)的時間步長為1 h，開展土地利用影響的降雨洪水過程模擬計算與定量研究。灞河流域不同時期土地利用見圖6，不同時期土地利用占比例統(tǒng)計見表2。由圖6可看出，隨著西安市城市快速擴張，城鎮(zhèn)化建設導致灞河流域城市用地比例在逐漸增加。不同土地利用糙率系數(shù)表，見表3。

表2 灞河流域土地利用變化矩陣 %

表3 不同土地利用類型的糙率系數(shù)

西安“八水潤長安”工程建設，在流域內(nèi)修建人工湖、濕地公園以及橡膠壩蓄水形成的人造水景工程，使得灞河流域水體面積略有增加。區(qū)域內(nèi)以城市化為代表的典型人類活動致使城市用地增加及水體面積增加，改變了流域降雨產(chǎn)匯流下墊面情況，必然會影響流域洪水過程。

圖6 不同時期灞河流域土地利用情況

5 土地利用變化、水庫建設的設計洪水過程響應模擬與量化

灞河流域中下游強烈人類活動改變了流域下墊面條件，直接影響產(chǎn)流過程，同時水利工程調(diào)蓄運行，直接改變河道匯流過程，洪水徑流過程必然發(fā)生變化。本文在研究灞河流域人類活動及下墊面要素的變化程度的基礎(chǔ)上，采用基于網(wǎng)格化的分布式降雨洪水模型，模擬土地利用和水庫建設運行對不同設計標準的設計洪水過程、洪量、洪峰等參數(shù)的影響，初步分析人類活動及下墊面要素的變化對洪水要素的影響程度和典型流域的洪水徑流演變特征，開展人類活動影響下洪水徑流過程變化與定量研究。

5.1 水系和子流域特征代表的流域參數(shù)提取

依據(jù)NASA Shuttle Radar Topographic Mission 90m Digital Elevation Model (SRTM DEM)數(shù)字高程影像提取水系和子流域特征代表的流域參數(shù)。據(jù)此，計算流域地表上每個柵格的地表坡度、流向、流量分配系數(shù)、集水面積等。

5.2 推理公式法與MCDRM模型法驗證

采用MCDRM模型模擬計算設計暴雨重現(xiàn)期為T=10年、T=20年、T=50年和T=100年的設計洪水過程，并依據(jù)《西安市實用水文手冊》推薦的推理公式法計算灞河上游馬渡王水文站設計洪水，驗證模型模擬結(jié)果。驗證和評價合格后，開展灞河全流域設計洪水計算。計算結(jié)果見表4。

5.3 基于MCDRM模型的全流域設計洪水響應模擬

以1980年土地利用為基準，對比2015年土地利用變化與水庫運行后不同設計標準暴雨洪水的洪峰、洪水總量變化，定量分析土地利用變化及水庫運行對流域設計洪水過程響應。不同設計標準暴雨洪水模擬計算結(jié)果見圖7。

表4 馬渡王水文站洪峰流量計算成果匯總表

5.4 MCDRM模型結(jié)果分析

以1980土地利用及對應的設計洪水為基準，比較2015年后土地利用與李家河水庫運行情況下的流域設計洪水變化，定量分析土地利用與李家河水庫運行對流域設計洪水影響程度。土地利用變化與水庫運行對灞河全流域設計洪水定量影響見表5。

圖7 不同設計標準暴雨洪水模擬計算結(jié)果對比圖

設計暴雨重現(xiàn)期/a設計洪峰/(m3·s-1)1980年洪峰/(m3·s-1)2015年土地利用與水庫運行后洪峰/(m3·s-1)洪峰變化比例/%設計洪量/108 m31980年洪量/108 m32015年土地利用與水庫運行后洪峰/108 m3洪量變化比例/%101953.91509.722.732.1971.72821.35202527.41965.022.252.6042.04421.51503625.12789.223.063.1472.47521.351004598.73500.423.883.5432.78621.37

對比水量變化及水文響應比例結(jié)果，表明土地利用與水庫運行綜合影響使得2015年后流域設計洪水洪峰與洪水總量減少20%左右。

6 結(jié)論和建議

本文建立的基于MCDRM的灞河流域水文模擬模型能夠較好地模擬灞河流域在不同重現(xiàn)期下的設計洪水過程，得到以下結(jié)論:

(1)根據(jù)灞河流域人類活動實際情況，構(gòu)建不同的人類活動干擾過程模塊，提高模型的適用性。并通過蒙特·卡羅方法與人工調(diào)節(jié)參數(shù)方法結(jié)合優(yōu)化模型參數(shù)，對場次洪水過程模擬以降低模擬結(jié)果的不確定性。依據(jù)《西安市實用水文手冊》推薦的推理公式法計算結(jié)果驗證MCDRM模型模擬效果。設計暴雨重現(xiàn)期為10、20、50、100 a時，MCDRM模型模擬得到的洪峰流量值與經(jīng)驗公式法計算得到的洪峰流量值分別相差5.6%，3.3%，4.3%，6.2%，平均差值為4.85%，模型可用于人類活動影響的洪水過程模擬計算。

(2)在土地利用變化與水庫建設等人類活動干擾的影響下，1980-2015年期間灞河全流域洪峰與洪量的變化總量削減約20%。較1985年，2015年流域設計標準為10、20、50、100年一遇的洪峰流量減少比例(平均減少22.98%)大于設計洪水徑流總量減少比例(平均減少約21.4%)。其中，100年一遇標準設計洪峰徑流減少最大(23.88%)，20年一遇標準設計洪水徑流總量減少最大(21.51%)。

(3)在半干旱區(qū)流域水資源緊缺、城市洪澇現(xiàn)象日益加重的背景下，本文選擇典型城市河流研究人類活動干擾后設計洪水過程變化，對地區(qū)水資源規(guī)劃、防洪減災尤為重要。同時，也對于今后研究半干旱區(qū)域流域洪水過程響應提供了新的選擇。