黃國新 謝小華 鄧 蕊
(1.江西省水文局,江西 南昌 330000;2.江西省吉安市水文局,江西 吉安 343000)
江西省位于長江中游南岸,河網(wǎng)密布水系發(fā)達,獨特的地形氣候條件造成該地區(qū)洪澇災(zāi)害頻繁,洪水預(yù)報對水安全至關(guān)重要。但由于部分地區(qū)水文站網(wǎng)缺失,水文監(jiān)測資料缺乏,給洪水預(yù)報、水庫調(diào)度及安全度汛工作帶來了很大的挑戰(zhàn)。吉安市行政區(qū)劃范圍內(nèi)有大中型水庫47座,其中,無長期水文觀測站點的水庫有36座,占全部大中型水庫數(shù)量的76.6%。為便于上述無資料地區(qū)進行洪水預(yù)報,開展地貌單位線適用性研究很有必要。
傳統(tǒng)的流域單位線從系統(tǒng)水文的角度進行研究,把流域看作一個系統(tǒng),降雨和出流分別視為輸入和輸出。因此,流域單位線成了黑箱模型的產(chǎn)物。流域地貌瞬時單位線是以流域地形地貌情況以及概率學(xué)隨機理論為基礎(chǔ),具有物理意義的流域匯流隨機模型。流域地貌瞬時單位線以研究區(qū)獨特的地形地貌特點,揭示流域基本特性與單位線之間的關(guān)系;以研究區(qū)水文氣象特性來體現(xiàn)與單位線的非線性時變之間的關(guān)系,克服了神經(jīng)網(wǎng)絡(luò)等學(xué)習(xí)方法推求單位線的一些缺陷,與研究區(qū)現(xiàn)實情況更為貼切,更能反映區(qū)域特點。自1979年委內(nèi)瑞拉的水文學(xué)者提出該模型以來,近幾十年來得到了迅速發(fā)展[1]。尤其是隨著計算機技術(shù)和GIS技術(shù)的高速發(fā)展,同時伴隨數(shù)字高程信息化的誕生,為研究地貌瞬時單位線與研究區(qū)地形地貌特征之間的相互關(guān)系,提供了非常有效的技術(shù)支持和保障,更加快速地促進了地貌瞬時單位線在區(qū)域洪水預(yù)報工作中的實踐與研究。地貌單位線采取數(shù)字高程模型(Digital Elevation Model)、地理信息系統(tǒng)(Geographic Information System)、遙感(Remote Sensing)等相關(guān)先進科學(xué)手段,使我們可以實時獲得研究區(qū)的地形地貌,同時充分考慮了人類因素的影響導(dǎo)致的流域時空變化,讓研究更接近流域的實際情況。受流域面積及形狀的限制,地貌單位線一般適用于50~1000km2的流域范圍。本文選取吉安市南車水庫流域為研究對象,基于DEM 數(shù)字高程信息,提取流域內(nèi)河網(wǎng)基本情況,獲取地貌參數(shù),然后建立流域匯流模型并驗證分析[2]。
南車水庫位于贛江水系禾水一級支流牛吼江上,地處吉安市東南部,坐落于泰和縣橋頭鎮(zhèn)大壟坑村,興建于1992年9月,2003年6月建成蓄水,水庫壩上斷面以上集水面積459km2,流域范圍內(nèi)匯集有中朝水、少江水、六八河、高幣水等河流,覆蓋新江鄉(xiāng)、雙橋鄉(xiāng)和五斗江鄉(xiāng)等鄉(xiāng)鎮(zhèn),于大壟坑村匯入南車水庫。研究區(qū)控制河長57.9km,流域范圍見圖1。
圖1 研究區(qū)范圍
南車水庫設(shè)計總庫容15380萬m3,興利庫容9520萬m3,防洪庫容415萬m3,死庫容2810萬m3,是一座以防洪、灌溉、發(fā)電等功能為主的混凝土面板碾壓堆石壩大(2)型水庫。
1.2.1 DEM數(shù)據(jù)
數(shù)字高程圖是數(shù)據(jù)分析的基礎(chǔ),本文數(shù)字高程數(shù)據(jù)來源為地理空間數(shù)據(jù)云GDEM 30m分辨率數(shù)字高程數(shù)據(jù)。數(shù)據(jù)獲取后由ArcGIS提取處理[3]。
1.2.2 降雨徑流資料
本研究選取南車水庫流域2017—2018年6場降雨及徑流資料進行模型分析及驗證,經(jīng)過分析和計算提取相關(guān)指標,見表1。
表1 降雨及徑流資料
降雨資料來源于流域范圍內(nèi)測站上報的實時水雨情數(shù)據(jù)庫,流量、峰現(xiàn)時間值根據(jù)庫容變化量按水量平衡原理倒推入庫流量得到。
徑流系數(shù)計算公式如下:
(1)
式中:W為場次洪水實測總量,m3;F為流域面積,km2;P為流域面降雨量,mm。
通過查閱《吉安市水文手冊》,南車水庫流域范圍內(nèi)降雨歷時一般在35h以內(nèi),徑流系數(shù)在0.5~0.8。由表1可知,6場次流域降雨歷時均值為22h,徑流系數(shù)均值為0.6,所選場次洪水能夠反映流域的實際情況。
1.2.3 研究方法
地貌瞬時單位線是既有物理意義又簡便易行的一種匯流模型,將其引進水文模型中,直接參與水文計算。其基本思路是:設(shè)有多個(n個)體積、質(zhì)量相等且互相沒有聯(lián)系、作用的水滴,在很短的時間匯到研究區(qū),水滴在運動軌跡內(nèi)要經(jīng)過相應(yīng)的運移時間才移動到流域出口斷面位置,每個水滴都有自己相應(yīng)的運移時間TB,但它們卻擁有一致的分布函數(shù)?;趨^(qū)域內(nèi)的特定時間注入的水量與流出的水量之差等于水量變化量(水量動態(tài)平衡原理)及大數(shù)法則,初步判斷流域瞬時單位線(IUH)相當于水滴在原來位置到流域出口斷面運移時間的概率密度函數(shù)。在天然下墊面條件下,無數(shù)下落到流域內(nèi)的水滴順著自己的軌跡,通過填坑、補洼后匯集到達流域出口斷面。每個水滴途經(jīng)的線路都有無數(shù)的方案,這些結(jié)果的概率函數(shù)可以采用Strahler以及Horton方法求得。計算全部結(jié)果的運行軌跡的隨機運行時間的概率密度函數(shù)同此運動軌跡的函數(shù)概率相乘,然后求和,即為區(qū)域內(nèi)水滴滯留時間的概率函數(shù),就是本次應(yīng)用的洪水預(yù)報方法——地貌瞬時單位線(GIUH)[4]。
假設(shè)某一時間降落到研究區(qū)內(nèi)分布均勻的降水量是由無數(shù)微小的水滴構(gòu)成,且水滴之間僅為很弱的相互影響和作用,則每個水質(zhì)點的運動都可被看作是一個由轉(zhuǎn)移概率控制的從某級別河道向同級或高級河道運移的Markov過程,那么該流域的瞬時單位線——u(0,t) 同水滴匯集到流域出口時間的分布密度函數(shù)fB(t)相等,得到下式:
u(0,t)=fB(t)
(2)
據(jù)此,可以判定該流域匯流的S曲線S(t)與水滴匯集到流域出口時間的概率分布函數(shù)FB(t)相等,得到下式:
S(t)=FB(t)
(3)
式(2)、式(3)表明,只要求得水滴匯集到流域出口時間的概率分布密度或概率分布函數(shù)二者中的其一,就能夠確定流域的瞬時單位線或S曲線,與傳統(tǒng)流域單位線的分析方式不同的是,地貌單位線主要基于區(qū)域地形地貌要素等角度與流域水文之間的相互影響,通過以上方式求得[1]。
基于水滴匯流時間T等于其流路長度L與其速度V之商的基本關(guān)系式,即T=L/V,根據(jù)概率論的理論推導(dǎo)得出下式:
(4)
(5)
式中:g(l)與φ(v)分別為L與V的分布密度;G(l)與ψ(v)分別為L與V的分布函數(shù);vmax為流域中水滴的最大匯集速度。
式(4)、式(5)將確定流域地貌瞬時單位線或S曲線的問題,轉(zhuǎn)化為推求水滴向流域出口斷面匯集的流路長度和速度的分布函數(shù)或分布密度問題。水滴匯流的流路長度可通過DEM直接提取得到,而由于受地形和地貌的影響,匯流速度在空間分布上具有不均勻性,很難直接確定。根據(jù)《美國水土保持局工程手冊》中提出的方法,參照Manning公式的形式,根據(jù)DEM單元的坡度計算流速:
V=asb
(6)
式中:s為單元的平均坡度;a和b為經(jīng)驗參數(shù)。
因此通過提取坡度分布律的方式來確定水滴匯流速度的分布律。以此為基礎(chǔ),將式(5)改為離散形式,可得
(7)
式中:m為離散的坡度值的總數(shù);Θ(s)為坡度的分布函數(shù)。
求得FB(t)后,就能夠得到某一Δt時段的地貌單位線:
U(Δt,t)=S(t)-S(t-Δt)=FB(t)-FB(t-Δt)
(8)
式中:U(Δt,t)為計算時段為Δt的地貌單位線。
提取得到的流域分水線(山脊線)是水源的起點,所以其匯流累積量為0。因此,在Arcmap中使用水文模塊對區(qū)域的匯流累積量的0值進行提取,并將提取的數(shù)據(jù)與生成的正地形求交,得到消除了存在于負地形區(qū)域中的錯誤的山脊線,獲得流域的分水線。
在Arcgis中創(chuàng)建一個線圖層,以生成的分水線作為參考圖層,結(jié)合南車水庫流域及其嵌套流域地貌圖,勾勒出流域邊界,并將邊界圖層轉(zhuǎn)化為面域圖對DEM圖進行裁剪,得到南車水庫流域的DEM圖,見圖2。
圖2 南車水庫流域邊界及流域DEM圖
考慮到內(nèi)插函數(shù)存在一定的局限性和研究區(qū)內(nèi)某些非常規(guī)地形的約束,通常會造成人為內(nèi)插出的DEM高程面有一些凹陷地區(qū),在對這些區(qū)域進行水文分析計算和水滴運移軌跡推算時,容易生成不合理或者錯誤的水力條件,進而獲取錯誤的水流方向信息和匯流時間。所以,在對流域內(nèi)各水滴進行水力分析和水文計算之前,要優(yōu)先對提取的原始DEM高程數(shù)據(jù)進行處理,使得各洼地得以填充,得到處理后的研究區(qū)的DEM高程圖[5]。
利用地理信息系統(tǒng)Arcgis10.2,利用基于流向的河網(wǎng)提取方法,得到柵格河網(wǎng)數(shù)據(jù),見圖3。
圖3 南車水庫流域河網(wǎng)
地貌單位線預(yù)報方案中河流級別的劃分采用斯特拉勒法,它是基于修正的霍頓河流級別劃分方式而得出的。斯特拉勒分級法定義從河源出發(fā)的河流為1級河流,同級的兩條河流交匯所形成的河流級數(shù)增加1級,不同級別的兩條河流交匯所形成的河流級別為二者中較高者。
選取地理信息系統(tǒng)操作軟件Arcgis10.2,利用已經(jīng)提取的河網(wǎng)信息和參數(shù),運用斯特拉勒分級法對河網(wǎng)進行分級運算,獲取具備所研究流域相應(yīng)級別的河網(wǎng)屬性數(shù)據(jù),以此為基礎(chǔ)可獲取流域各級別河段的數(shù)量、流程、集水面積和其他與地貌單位線推求相關(guān)的各項地貌特征值,同時運用霍頓河系定律進行計算就能夠獲取流域內(nèi)的河數(shù)率、河長率和面積率,其表達式分別為
河數(shù)率:
(9)
河長率:
(10)
面積率:
(11)
各河道具體數(shù)據(jù)見表2,劃分的集水區(qū)域矢量圖見圖4。
表2 各河道屬性數(shù)據(jù)
圖4 集水區(qū)域矢量圖
根據(jù)地貌單位線通用公式,進行瞬時單位線計算。為簡化工作量,本文采用編程方式,將上述計算過程封裝為計算工具。把河道三參數(shù)輸入工具,計算得到本流域的地貌瞬時單位線,見圖5,進而得到流域1h地貌瞬時單位線[6]。
圖5 地貌瞬時單位線計算
為便于洪水預(yù)報分析,將地貌單位線導(dǎo)入到江西省吉安市水庫防汛調(diào)度預(yù)報系統(tǒng)中,該系統(tǒng)能夠自動讀取歷史、實時數(shù)據(jù),輔助完成對預(yù)報方案的率定,預(yù)報期可以達到72h。該系統(tǒng)現(xiàn)已經(jīng)在實際洪水預(yù)報工作中應(yīng)用,目前運行穩(wěn)定可靠。
對南車水庫流域確定的6場洪水過程進行模擬統(tǒng)計,采用地貌單位線進行流域匯流計算。模擬結(jié)果洪量、洪峰、峰現(xiàn)時間的平均相對誤差分別為17%、19%和3h,場次洪水綜合預(yù)報合格率50%。圖6~圖11 為實測與模擬洪水過程線比較結(jié)果。
圖6 20170628洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
圖7 20170711洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
圖8 20180606洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
圖9 20180612洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
圖10 20180621洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
圖11 20180622洪水預(yù)報過程注 圖中綠色線條代表實測流量過程,紅色為預(yù)報流量過程。
場次洪水評定參照水文情報預(yù)報規(guī)范,洪峰、洪量預(yù)報以預(yù)見期內(nèi)實測變幅的20%作為許可誤差,峰現(xiàn)時間按3h計,預(yù)報結(jié)果見表3[7]。
表3 洪水預(yù)報結(jié)果統(tǒng)計
續(xù)表
本文根據(jù)童冰星等[8]提出的地貌單位線確定方法,以吉安市南車水庫流域為例,利用DEM信息提取了流域的地貌單位線,并將之用于該流域的洪水模擬,通過對6場次洪水預(yù)報與實測值比較,峰現(xiàn)時間合格率為83.3%,洪峰流量合格率為66.7%,洪水總量合格率為66.7%,場次洪水綜合預(yù)報合格率為50%,取得了較好的模擬結(jié)果,進一步驗證了通過流路長度分布,確定地貌單位線的有效性,并將此地貌單位線法應(yīng)用到吉安市31個無資料水庫洪水預(yù)報工作中。
通過研究分析預(yù)報存在誤差的原因:一方面水庫庫容采用4段制上報,入庫流量采用水量平衡原理倒推入庫,人為地對實測洪峰值進行了均化處理,預(yù)報數(shù)值計算步長為1h,影響了洪水預(yù)報精度;另一方面,預(yù)報流域的基流大小對預(yù)報結(jié)果影響較大,而對無資料地區(qū)的基流獲取尚且存在困難。隨著實際應(yīng)用和實踐,對于無資料地區(qū)的水文預(yù)報問題,有望展開更加深入的研究并解決上述問題。