黃國新 謝小華 鄧 蕊
(1.江西省水文局,江西 南昌 330000;2.江西省吉安市水文局,江西 吉安 343000)
江西省位于長江中游南岸,河網密布水系發達,獨特的地形氣候條件造成該地區洪澇災害頻繁,洪水預報對水安全至關重要。但由于部分地區水文站網缺失,水文監測資料缺乏,給洪水預報、水庫調度及安全度汛工作帶來了很大的挑戰。吉安市行政區劃范圍內有大中型水庫47座,其中,無長期水文觀測站點的水庫有36座,占全部大中型水庫數量的76.6%。為便于上述無資料地區進行洪水預報,開展地貌單位線適用性研究很有必要。
傳統的流域單位線從系統水文的角度進行研究,把流域看作一個系統,降雨和出流分別視為輸入和輸出。因此,流域單位線成了黑箱模型的產物。流域地貌瞬時單位線是以流域地形地貌情況以及概率學隨機理論為基礎,具有物理意義的流域匯流隨機模型。流域地貌瞬時單位線以研究區獨特的地形地貌特點,揭示流域基本特性與單位線之間的關系;以研究區水文氣象特性來體現與單位線的非線性時變之間的關系,克服了神經網絡等學習方法推求單位線的一些缺陷,與研究區現實情況更為貼切,更能反映區域特點。自1979年委內瑞拉的水文學者提出該模型以來,近幾十年來得到了迅速發展[1]。尤其是隨著計算機技術和GIS技術的高速發展,同時伴隨數字高程信息化的誕生,為研究地貌瞬時單位線與研究區地形地貌特征之間的相互關系,提供了非常有效的技術支持和保障,更加快速地促進了地貌瞬時單位線在區域洪水預報工作中的實踐與研究。地貌單位線采取數字高程模型(Digital Elevation Model)、地理信息系統(Geographic Information System)、遙感(Remote Sensing)等相關先進科學手段,使我們可以實時獲得研究區的地形地貌,同時充分考慮了人類因素的影響導致的流域時空變化,讓研究更接近流域的實際情況。受流域面積及形狀的限制,地貌單位線一般適用于50~1000km2的流域范圍。本文選取吉安市南車水庫流域為研究對象,基于DEM 數字高程信息,提取流域內河網基本情況,獲取地貌參數,然后建立流域匯流模型并驗證分析[2]。
南車水庫位于贛江水系禾水一級支流牛吼江上,地處吉安市東南部,坐落于泰和縣橋頭鎮大壟坑村,興建于1992年9月,2003年6月建成蓄水,水庫壩上斷面以上集水面積459km2,流域范圍內匯集有中朝水、少江水、六八河、高幣水等河流,覆蓋新江鄉、雙橋鄉和五斗江鄉等鄉鎮,于大壟坑村匯入南車水庫。研究區控制河長57.9km,流域范圍見圖1。

圖1 研究區范圍
南車水庫設計總庫容15380萬m3,興利庫容9520萬m3,防洪庫容415萬m3,死庫容2810萬m3,是一座以防洪、灌溉、發電等功能為主的混凝土面板碾壓堆石壩大(2)型水庫。
1.2.1 DEM數據
數字高程圖是數據分析的基礎,本文數字高程數據來源為地理空間數據云GDEM 30m分辨率數字高程數據。數據獲取后由ArcGIS提取處理[3]。
1.2.2 降雨徑流資料
本研究選取南車水庫流域2017—2018年6場降雨及徑流資料進行模型分析及驗證,經過分析和計算提取相關指標,見表1。

表1 降雨及徑流資料
降雨資料來源于流域范圍內測站上報的實時水雨情數據庫,流量、峰現時間值根據庫容變化量按水量平衡原理倒推入庫流量得到。
徑流系數計算公式如下:
(1)
式中:W為場次洪水實測總量,m3;F為流域面積,km2;P為流域面降雨量,mm。
通過查閱《吉安市水文手冊》,南車水庫流域范圍內降雨歷時一般在35h以內,徑流系數在0.5~0.8。由表1可知,6場次流域降雨歷時均值為22h,徑流系數均值為0.6,所選場次洪水能夠反映流域的實際情況。
1.2.3 研究方法
地貌瞬時單位線是既有物理意義又簡便易行的一種匯流模型,將其引進水文模型中,直接參與水文計算。其基本思路是:設有多個(n個)體積、質量相等且互相沒有聯系、作用的水滴,在很短的時間匯到研究區,水滴在運動軌跡內要經過相應的運移時間才移動到流域出口斷面位置,每個水滴都有自己相應的運移時間TB,但它們卻擁有一致的分布函數。基于區域內的特定時間注入的水量與流出的水量之差等于水量變化量(水量動態平衡原理)及大數法則,初步判斷流域瞬時單位線(IUH)相當于水滴在原來位置到流域出口斷面運移時間的概率密度函數。在天然下墊面條件下,無數下落到流域內的水滴順著自己的軌跡,通過填坑、補洼后匯集到達流域出口斷面。每個水滴途經的線路都有無數的方案,這些結果的概率函數可以采用Strahler以及Horton方法求得。計算全部結果的運行軌跡的隨機運行時間的概率密度函數同此運動軌跡的函數概率相乘,然后求和,即為區域內水滴滯留時間的概率函數,就是本次應用的洪水預報方法——地貌瞬時單位線(GIUH)[4]。
假設某一時間降落到研究區內分布均勻的降水量是由無數微小的水滴構成,且水滴之間僅為很弱的相互影響和作用,則每個水質點的運動都可被看作是一個由轉移概率控制的從某級別河道向同級或高級河道運移的Markov過程,那么該流域的瞬時單位線——u(0,t) 同水滴匯集到流域出口時間的分布密度函數fB(t)相等,得到下式:
u(0,t)=fB(t)
(2)
據此,可以判定該流域匯流的S曲線S(t)與水滴匯集到流域出口時間的概率分布函數FB(t)相等,得到下式:
S(t)=FB(t)
(3)
式(2)、式(3)表明,只要求得水滴匯集到流域出口時間的概率分布密度或概率分布函數二者中的其一,就能夠確定流域的瞬時單位線或S曲線,與傳統流域單位線的分析方式不同的是,地貌單位線主要基于區域地形地貌要素等角度與流域水文之間的相互影響,通過以上方式求得[1]。
基于水滴匯流時間T等于其流路長度L與其速度V之商的基本關系式,即T=L/V,根據概率論的理論推導得出下式:

(4)

(5)
式中:g(l)與φ(v)分別為L與V的分布密度;G(l)與ψ(v)分別為L與V的分布函數;vmax為流域中水滴的最大匯集速度。
式(4)、式(5)將確定流域地貌瞬時單位線或S曲線的問題,轉化為推求水滴向流域出口斷面匯集的流路長度和速度的分布函數或分布密度問題。水滴匯流的流路長度可通過DEM直接提取得到,而由于受地形和地貌的影響,匯流速度在空間分布上具有不均勻性,很難直接確定。根據《美國水土保持局工程手冊》中提出的方法,參照Manning公式的形式,根據DEM單元的坡度計算流速:
V=asb
(6)
式中:s為單元的平均坡度;a和b為經驗參數。
因此通過提取坡度分布律的方式來確定水滴匯流速度的分布律。以此為基礎,將式(5)改為離散形式,可得
(7)
式中:m為離散的坡度值的總數;Θ(s)為坡度的分布函數。
求得FB(t)后,就能夠得到某一Δt時段的地貌單位線:
U(Δt,t)=S(t)-S(t-Δt)=FB(t)-FB(t-Δt)
(8)
式中:U(Δt,t)為計算時段為Δt的地貌單位線。
提取得到的流域分水線(山脊線)是水源的起點,所以其匯流累積量為0。因此,在Arcmap中使用水文模塊對區域的匯流累積量的0值進行提取,并將提取的數據與生成的正地形求交,得到消除了存在于負地形區域中的錯誤的山脊線,獲得流域的分水線。
在Arcgis中創建一個線圖層,以生成的分水線作為參考圖層,結合南車水庫流域及其嵌套流域地貌圖,勾勒出流域邊界,并將邊界圖層轉化為面域圖對DEM圖進行裁剪,得到南車水庫流域的DEM圖,見圖2。

圖2 南車水庫流域邊界及流域DEM圖
考慮到內插函數存在一定的局限性和研究區內某些非常規地形的約束,通常會造成人為內插出的DEM高程面有一些凹陷地區,在對這些區域進行水文分析計算和水滴運移軌跡推算時,容易生成不合理或者錯誤的水力條件,進而獲取錯誤的水流方向信息和匯流時間。所以,在對流域內各水滴進行水力分析和水文計算之前,要優先對提取的原始DEM高程數據進行處理,使得各洼地得以填充,得到處理后的研究區的DEM高程圖[5]。
利用地理信息系統Arcgis10.2,利用基于流向的河網提取方法,得到柵格河網數據,見圖3。

圖3 南車水庫流域河網
地貌單位線預報方案中河流級別的劃分采用斯特拉勒法,它是基于修正的霍頓河流級別劃分方式而得出的。斯特拉勒分級法定義從河源出發的河流為1級河流,同級的兩條河流交匯所形成的河流級數增加1級,不同級別的兩條河流交匯所形成的河流級別為二者中較高者。
選取地理信息系統操作軟件Arcgis10.2,利用已經提取的河網信息和參數,運用斯特拉勒分級法對河網進行分級運算,獲取具備所研究流域相應級別的河網屬性數據,以此為基礎可獲取流域各級別河段的數量、流程、集水面積和其他與地貌單位線推求相關的各項地貌特征值,同時運用霍頓河系定律進行計算就能夠獲取流域內的河數率、河長率和面積率,其表達式分別為
河數率:
(9)
河長率:
(10)
面積率:
(11)

各河道具體數據見表2,劃分的集水區域矢量圖見圖4。

表2 各河道屬性數據

圖4 集水區域矢量圖
根據地貌單位線通用公式,進行瞬時單位線計算。為簡化工作量,本文采用編程方式,將上述計算過程封裝為計算工具。把河道三參數輸入工具,計算得到本流域的地貌瞬時單位線,見圖5,進而得到流域1h地貌瞬時單位線[6]。

圖5 地貌瞬時單位線計算
為便于洪水預報分析,將地貌單位線導入到江西省吉安市水庫防汛調度預報系統中,該系統能夠自動讀取歷史、實時數據,輔助完成對預報方案的率定,預報期可以達到72h。該系統現已經在實際洪水預報工作中應用,目前運行穩定可靠。
對南車水庫流域確定的6場洪水過程進行模擬統計,采用地貌單位線進行流域匯流計算。模擬結果洪量、洪峰、峰現時間的平均相對誤差分別為17%、19%和3h,場次洪水綜合預報合格率50%。圖6~圖11 為實測與模擬洪水過程線比較結果。

圖6 20170628洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。

圖7 20170711洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。

圖8 20180606洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。

圖9 20180612洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。

圖10 20180621洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。

圖11 20180622洪水預報過程注 圖中綠色線條代表實測流量過程,紅色為預報流量過程。
場次洪水評定參照水文情報預報規范,洪峰、洪量預報以預見期內實測變幅的20%作為許可誤差,峰現時間按3h計,預報結果見表3[7]。

表3 洪水預報結果統計

續表
本文根據童冰星等[8]提出的地貌單位線確定方法,以吉安市南車水庫流域為例,利用DEM信息提取了流域的地貌單位線,并將之用于該流域的洪水模擬,通過對6場次洪水預報與實測值比較,峰現時間合格率為83.3%,洪峰流量合格率為66.7%,洪水總量合格率為66.7%,場次洪水綜合預報合格率為50%,取得了較好的模擬結果,進一步驗證了通過流路長度分布,確定地貌單位線的有效性,并將此地貌單位線法應用到吉安市31個無資料水庫洪水預報工作中。
通過研究分析預報存在誤差的原因:一方面水庫庫容采用4段制上報,入庫流量采用水量平衡原理倒推入庫,人為地對實測洪峰值進行了均化處理,預報數值計算步長為1h,影響了洪水預報精度;另一方面,預報流域的基流大小對預報結果影響較大,而對無資料地區的基流獲取尚且存在困難。隨著實際應用和實踐,對于無資料地區的水文預報問題,有望展開更加深入的研究并解決上述問題。