村鎮(zhèn)易澇區(qū)動態(tài)洪水風(fēng)險分析模型研究

牛 帥，崔信民

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

隨著我國城鎮(zhèn)化進(jìn)程的加快，村鎮(zhèn)地區(qū)社會經(jīng)濟(jì)水平發(fā)展迅速，洪水災(zāi)害造成的損失日益嚴(yán)重。在工程與非工程措施并重的防洪減災(zāi)新思路下，對村鎮(zhèn)地區(qū)進(jìn)行動態(tài)洪水風(fēng)險分析，可以為村鎮(zhèn)地區(qū)的水利規(guī)劃、防洪減災(zāi)等提供關(guān)鍵信息參考，具有非常重要的現(xiàn)實意義。

近年來，國內(nèi)外學(xué)者對村鎮(zhèn)地區(qū)的防洪標(biāo)準(zhǔn)、洪水風(fēng)險、洪災(zāi)損失等內(nèi)容進(jìn)行了廣泛的研究[1,2]，并取得了較大進(jìn)展。王高丹等[3]利用多源數(shù)據(jù)通過GIS方法、層次分析法進(jìn)行了海南島中部村鎮(zhèn)地區(qū)山洪災(zāi)害風(fēng)險分析。申邵洪等[4]采用一維、二維耦合水動力學(xué)模型開展了保康縣城鎮(zhèn)區(qū)域的洪水影響及淹沒分析研究。練繼建等[5]采用二維水動力模型分析了萬泉河下游農(nóng)村地區(qū)洪水淹沒情況,繪制了脆弱性分布圖。由于洪水來源的不確定性，傳統(tǒng)的靜態(tài)洪水風(fēng)險圖不能滿足實際需要，在實際防汛中需要進(jìn)行動態(tài)洪水風(fēng)險分析[6,7]。吳濱濱等[8]研發(fā)了大陸澤、寧晉泊洪水風(fēng)險實時分析系統(tǒng)。寇嘉瑋等[9]結(jié)合水動力學(xué)模型和WebGIS技術(shù),構(gòu)建了洪澤湖周邊滯洪區(qū)暴雨內(nèi)澇預(yù)報預(yù)警系統(tǒng)。

本文根據(jù)村鎮(zhèn)易澇區(qū)的暴雨內(nèi)澇特點，采用一、二維水動力學(xué)模型構(gòu)建了易澇區(qū)暴雨內(nèi)澇耦合計算模型，模型可以根據(jù)未來預(yù)測降雨過程和預(yù)測潮位過程進(jìn)行易澇區(qū)內(nèi)澇風(fēng)險的動態(tài)分析，可以為村鎮(zhèn)地區(qū)的實際防汛工作提供及時有效的決策參考信息。

1 研究區(qū)域概況

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

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

2 易澇區(qū)洪水分析模型

2.1 河網(wǎng)一維水流模型

采用一維水動力學(xué)模型對易澇區(qū)內(nèi)河網(wǎng)的水流運動進(jìn)行模擬，描述河流一維水流運動的圣維南方程組包括連續(xù)性方程和動量方程：

(1)

(2)

式中：Q為流量；Z為斷面水位；A為過水?dāng)嗝婷娣e；q為旁側(cè)入流量；g為重力加速度；K為流量模數(shù)；x為空間坐標(biāo)；t為時間坐標(biāo)。

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

2.2 地表二維水動力學(xué)模型

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

(3)

(4)

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

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

圖2 河網(wǎng)概化與計算網(wǎng)格Fig.2 Generalization of river network and computation grid

2.3 模型耦合

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

(5)

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

(6)

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

2.4 模型驗證

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

圖3 調(diào)查點地理位置圖Fig.3 Geographical position of investigation point

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

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

3 積澇風(fēng)險動態(tài)分析

3.1 模型封裝與接口調(diào)用

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

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

圖4 模塊調(diào)用流程圖Fig.4 Chart of module transfer flow

3.2 設(shè)計頻率靜態(tài)洪水風(fēng)險圖

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

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

4 結(jié) 語

本文構(gòu)建了萬城鎮(zhèn)易澇區(qū)一維、二維耦合水動力學(xué)模型，并采用歷史洪水資料進(jìn)行了模型驗證，可以有效地模擬暴雨條件下的萬城鎮(zhèn)易澇區(qū)洪水風(fēng)險。面對實際防汛需要，通過降雨、潮位邊界條件的動態(tài)設(shè)定構(gòu)建了萬城鎮(zhèn)易澇區(qū)動態(tài)洪水風(fēng)險分析模型?？梢詾榇彐?zhèn)地區(qū)的洪水預(yù)測預(yù)警、防洪減災(zāi)等提供關(guān)鍵技術(shù)支撐，具有良好的實際應(yīng)用價值。