村鎮易澇區動態洪水風險分析模型研究

牛 帥，崔信民

(南京水利科學研究院水文水資源研究所，南京 210029)

隨著我國城鎮化進程的加快，村鎮地區社會經濟水平發展迅速，洪水災害造成的損失日益嚴重。在工程與非工程措施并重的防洪減災新思路下，對村鎮地區進行動態洪水風險分析，可以為村鎮地區的水利規劃、防洪減災等提供關鍵信息參考，具有非常重要的現實意義。

近年來，國內外學者對村鎮地區的防洪標準、洪水風險、洪災損失等內容進行了廣泛的研究[1,2]，并取得了較大進展。王高丹等[3]利用多源數據通過GIS方法、層次分析法進行了海南島中部村鎮地區山洪災害風險分析。申邵洪等[4]采用一維、二維耦合水動力學模型開展了?？悼h城鎮區域的洪水影響及淹沒分析研究。練繼建等[5]采用二維水動力模型分析了萬泉河下游農村地區洪水淹沒情況,繪制了脆弱性分布圖。由于洪水來源的不確定性，傳統的靜態洪水風險圖不能滿足實際需要，在實際防汛中需要進行動態洪水風險分析[6,7]。吳濱濱等[8]研發了大陸澤、寧晉泊洪水風險實時分析系統?？芗维|等[9]結合水動力學模型和WebGIS技術,構建了洪澤湖周邊滯洪區暴雨內澇預報預警系統。

本文根據村鎮易澇區的暴雨內澇特點，采用一、二維水動力學模型構建了易澇區暴雨內澇耦合計算模型，模型可以根據未來預測降雨過程和預測潮位過程進行易澇區內澇風險的動態分析，可以為村鎮地區的實際防汛工作提供及時有效的決策參考信息。

1 研究區域概況

萬城鎮易澇區位于海南省萬寧市，瀕臨南海，地理位置如圖1所示。萬城鎮轄居10個社區，32個村委會，2個鎮辦農場，擁有東星工業開發區和烏場港經濟開發區。易澇區處于海南省暴雨高發區，暴雨多、強度大。易澇區內河網水系較多，河道下游連接小海，小海洪潮受潮汐通道束縛影響，水位易于高漲，形成頂托，排水不暢。根據近年來統計數據，年平均受淹造成不同程度損失達7～8次之多(下游農村區受淹次數居多)。

圖1 易澇區位置示意圖Fig.1 Geographical position of waterlogged area

2 易澇區洪水分析模型

2.1 河網一維水流模型

采用一維水動力學模型對易澇區內河網的水流運動進行模擬，描述河流一維水流運動的圣維南方程組包括連續性方程和動量方程：

(1)

(2)

式中：Q為流量；Z為斷面水位；A為過水斷面面積；q為旁側入流量；g為重力加速度；K為流量模數；x為空間坐標；t為時間坐標。

采用廣泛使用的Pressman格式對方程組進行數值離散求解[10]。該格式采用隱式離散格式，具有較好的穩定性，能夠滿足模擬計算的要求。根據河道地形資料，對易澇區內的銅鼓溪、面前溪、太陽河舊河等河道進行一維建模。共構建概化斷面750個，概化河段733個，河道汊點9個，河段平均間距約70 m。河網模型上游邊界取至太陽河左岸堤防處，考慮到易澇區內河道下游受到小海潮位頂托，下游邊界采用小海潮位邊界條件。

2.2 地表二維水動力學模型

易澇區內地表的洪水演進模擬采用二維水動力學模型，可以計算得到洪水淹沒水深、到達時間、流速分布等信息。若不考慮科氏力和風力的影響，引入靜水壓力假設，則二維淺水方程的守恒形式為：

(3)

(4)

式中：U為守恒變量；E、F分別為x和y方向通量；S為源項；S0x、S0y分別為x和y方向的底坡；Sfx、Sfy分別為x和y方向上的摩阻比降；h為水深；g為重力加速度；u、v分別為x和y方向的垂線平均流速。

基于非結構網格采用有限體積法對二維淺水方程進行數值離散求解[11,12]，采用Roe格式的近似Riemann解來求解界面通量，對底坡源項特征分解以保證和諧性。根據易澇區地形特征、水系分布等資料，確定萬城鎮易澇區建模范圍面積約59.5 km2。采用非結構網格對易澇區區域進行網格剖分，見圖2。網格單元共17 325 個，網格平均邊長55 m。為了避免高程點數據失準導致局部高程凹陷造成的淹沒水深的奇異解，對萬城鎮易澇區的數字高程模型進行數據校核，采用距離加權法對網格單元進行高程插值。根據萬城鎮下墊面的土地利用分類情況(水田、草地、林地等)資料，結合糙率參數取值表對網格單元的糙率參數進行賦值。土壤下滲和蒸發分別按照下滲能力和蒸發能力在降雨中進行扣除后得到凈雨過程，作為暴雨邊界條件，采用直接降雨法進行地表二維模型積澇計算。

圖2 河網概化與計算網格Fig.2 Generalization of river network and computation grid

2.3 模型耦合

為了實現河道水流通過河道堤防漫溢進入易澇區及易澇區內洪水回流河道，對一、二維水動力學模型進行耦合計算[13]，采用下面的寬頂堰流公式計算堤防漫溢水量實現一、二維模型水量交換，同時滿足水量守恒條件。

(5)

h1=max(Zu,Zd)-Ab；h2=min(Zu,Zd)-Zb

(6)

式中：q為通過堤防的單寬流量；Zu、Zd分別為堤防處河道內外的水位；Zb為河道堤防高程。

2.4 模型驗證

采用2010年10月實測暴雨洪水資料對易澇區洪水風險分析模型進行驗證分析。2010年9月30日至10月7日受到臺風影響，萬寧市出現持續強降雨，20個觀測站累計降雨量均超過800 mm，其中3個超過1 000 mm，過程降雨為歷史罕見。由于在實際中很難得到暴雨造成的整個淹沒區域內所有積澇點的淹沒水深數據，考慮到暴雨內澇時的積水一般積蓄在地勢低洼地帶，所以通常選擇一些代表性低洼處的最大淹沒水深作為模型計算的驗證資料。通過實地調研分析收集到5個地勢低洼處的代表性積澇點的最大淹沒水深數據，調查點位置見圖3。

圖3 調查點地理位置圖Fig.3 Geographical position of investigation point

將本次實測降雨過程作為降雨輸入邊界條件，采用模型進行模擬計算，將實際調查點的計算最大淹沒水深與實測最大淹沒水深進行比較驗證，見表1。從比較結果可知，計算淹沒水深與實測淹沒水深基本一致，滿足模型驗證要求，說明模型是合理可靠的。

表1 實測水深與計算水深比較 m

3 積澇風險動態分析

3.1 模型封裝與接口調用

為了實現對預測降雨、潮位等突發的洪水風險因子下的易澇區積澇風險動態分析，將降雨邊界條件和下游潮位邊界條件允許外部訪問對兩個邊界條件過程進行動態修改，然后調用模型進行洪水淹沒風險計算。既可以針對特定水文頻率制作靜態洪水風險圖，也可以針對未來預測降雨過程、潮位過程進行易澇區積澇風險動態分析。

采用計算機程序語言對計算模型進行封裝處理。設計模型與外部系統交互的接口共3個：①模型輸入接口，用于供外部系統輸入任意的預測降雨過程、潮位過程等模型計算邊界條件，作為模型計算的動態風險因子；②模型計算接口，用于供外部系統顯式調用模型進行洪水淹沒方案的模擬計算；③模型輸出接口，用于模型計算結果(淹沒水深、流速等)輸出供外部系統展示洪水演進、查詢淹沒水深等。各個模塊調用流程圖見圖4。

圖4 模塊調用流程圖Fig.4 Chart of module transfer flow

3.2 設計頻率靜態洪水風險圖

根據《廣東省暴雨參數等值線圖》及《廣東省暴雨徑流查算圖表》中與海南省相關的水文設計分析成果，推求計算萬城鎮易澇區的5年一遇、10年一遇、20年一遇、50年一遇24 h設計暴雨過程。模型計算邊界條件中，以設計暴雨過程經過扣損后的凈雨過程作為二維模型的輸入邊界條件，以小海的多年平均高潮位作為河網一維模型的下游潮位邊界條件。對5年一遇、10年一遇、20年一遇、50年一遇設計頻率下的萬城鎮易澇區積澇風險進行計算，統計得到各個設計頻率暴雨下的最大淹沒水深分布圖，如圖5所示。

圖5 不同設計降雨下的最大淹沒水深分布圖Fig.5 Maximum submerged water depth distribution map under different design rainfall

4 結 語

本文構建了萬城鎮易澇區一維、二維耦合水動力學模型，并采用歷史洪水資料進行了模型驗證，可以有效地模擬暴雨條件下的萬城鎮易澇區洪水風險。面對實際防汛需要，通過降雨、潮位邊界條件的動態設定構建了萬城鎮易澇區動態洪水風險分析模型。可以為村鎮地區的洪水預測預警、防洪減災等提供關鍵技術支撐，具有良好的實際應用價值。