山丘區(qū)小流域洪水風(fēng)險(xiǎn)模型研究及應(yīng)用

李曉英，談亞琦，蘇志偉，田佳樂(lè)

(1.河海大學(xué) 水利水電學(xué)院，南京 210098； 2.杭州誠(chéng)禹水利科技有限公司，杭州 310008)

為預(yù)防和減輕洪水帶來(lái)的危害，抵御洪水侵襲，研究洪水演進(jìn)的規(guī)律顯得十分重要。隨著計(jì)算機(jī)科學(xué)的不斷進(jìn)步，越來(lái)越多的洪水演進(jìn)模型得到開(kāi)發(fā)研究，成為防洪減災(zāi)理論體系的組成部分和分析工具。

洪水演進(jìn)數(shù)值模擬是研究洪水演進(jìn)規(guī)律的有效方法之一。Caleffi等人[1]采用二維淺水方程對(duì)意大利托切河進(jìn)行洪水演進(jìn)模擬。Bladé E等人[2]基于有限體積法分別建立一維、二維水動(dòng)力學(xué)模型，并提出基于數(shù)值通量的耦合方法，加強(qiáng)一維～二維模型之間的動(dòng)量傳遞，在位于西班牙的埃布羅河進(jìn)行驗(yàn)證；李傳奇等人[3]利用構(gòu)建的一維二維水動(dòng)力耦合模型，針對(duì)不同典型降雨過(guò)程，對(duì)模型進(jìn)行驗(yàn)證，并最終用于模擬濟(jì)南市不同重現(xiàn)期下洪水淹沒(méi)情況；苑希民等人[4]建立了一、二維水動(dòng)力耦合模型，其中一維模型采用Preissmann格式離散，二維利用Roe格式離散，應(yīng)用模型于黃河青銅峽河西灌區(qū)潰堤洪水的模擬，較為真實(shí)地反映了洪水在計(jì)算區(qū)域的演進(jìn)過(guò)程與淹沒(méi)范圍；吳天蛟等人[5]以MIKE11為基礎(chǔ)，區(qū)間入流采用分布式模型，對(duì)三峽庫(kù)區(qū)進(jìn)行洪水演進(jìn)模擬；謝作濤等人[6]從求解Holly-Preissmann格式出發(fā)，結(jié)合有限差分法建立一維洪水演進(jìn)模型；槐文信等[7]建立一、二維水動(dòng)力學(xué)模型對(duì)渭河下游河道及洪泛區(qū)洪水進(jìn)行數(shù)值仿真模擬；付成威等人[8]建立了一、二維實(shí)時(shí)動(dòng)態(tài)耦合模型，并利用濕水深和干水深理論，改進(jìn)傳統(tǒng)洪水演進(jìn)模型，將模型在谷堆圩蓄滯洪區(qū)進(jìn)行了驗(yàn)證，結(jié)果良好。

基于上述研究，本文應(yīng)用圣維南方程組作為一維模型的基本方程；應(yīng)用水流運(yùn)動(dòng)的能量守恒原理和質(zhì)量守恒推導(dǎo)得到平面二維非恒定流數(shù)值模型的基本方程組，進(jìn)而得到二維水動(dòng)力學(xué)模型。通過(guò)側(cè)向連接的方式建立一維、二維耦合水動(dòng)力學(xué)模型，利用該模型于太平溪山丘區(qū)小流域，對(duì)該流域作P=10%設(shè)計(jì)暴雨條件下的洪水演進(jìn)模擬分析，為該區(qū)域防洪風(fēng)險(xiǎn)管理提供數(shù)據(jù)參考。

1 原理介紹[9-12]

1.1 一維洪水演進(jìn)模型原理

以圣維南方程組為一維水動(dòng)力學(xué)模型采用的基本方程。方程表示如下：

式中：Q為流量，m3/s；A為斷面面積，m2；B為河道水面寬度，m；α為動(dòng)量校正系數(shù)；g為重力加速度,m/s2。

1.2 二維洪水演進(jìn)模型原理

根據(jù)水流運(yùn)動(dòng)的質(zhì)量守恒以及能量守恒原理，可以推導(dǎo)出平面二維非恒定流數(shù)值模型的基本方程組，該方程組由連續(xù)性方程以及X、Y方向的動(dòng)量方程組成：

式中：h為實(shí)際水深,m；η為水面高程,m；u、v分別為X方向和Y方向的流速,m/s；S為源匯項(xiàng),m2/s；ρ0為水的密度,kg/m3；f為科式力系數(shù)；τsx、τsy分別為X方向和Y方向風(fēng)切應(yīng)力；τbx、τby分別為X方向和Y方向底部切應(yīng)力；Tih為黏滯項(xiàng)、紊動(dòng)項(xiàng)和擴(kuò)散項(xiàng)的總和。

1.3 一維、二維模型耦合方式及模型求解過(guò)程

本文以側(cè)向連解的方式完成一維洪水演進(jìn)模型與二維洪水演進(jìn)模型的耦合，耦合模型求解過(guò)程見(jiàn)圖1。

圖1 耦合模型求解過(guò)程Fig.1 Coupling model solving process

2 模擬實(shí)例[13-14]

2.1 區(qū)域概況

太平溪流域?qū)儆谏酱ㄠl(xiāng)在安吉縣的南部，該流域是苕溪流域的源頭支流，安吉縣境內(nèi)區(qū)域受饅頭山水庫(kù)控制，控制范圍內(nèi)的集雨面積約為27.6 km2，河流長(zhǎng)約為11.9 km。該流域所屬區(qū)域包括安吉縣山川鄉(xiāng)大里村、九畝村、船村3個(gè)行政村，該區(qū)域在汛期來(lái)臨時(shí)，天然來(lái)水量十分巨大且流速迅急，給該流域及其下游防護(hù)區(qū)安全造成極大的威脅。

2.2 太平溪流域流水網(wǎng)生成

在特定的流域，徑流匯聚的多少與水流累積值成正比。因此，可用水流累積個(gè)數(shù)來(lái)定義河流的起始點(diǎn)。在河網(wǎng)生成的第一步中，必須先設(shè)置水流累積矩陣?yán)鄯e值的閾值，認(rèn)為在水流累積個(gè)數(shù)小于該值的網(wǎng)格上不能產(chǎn)生足夠的徑流而形成水道，而大于該值的網(wǎng)格上能形成水道。由此，水流累積個(gè)數(shù)大于該值的網(wǎng)格所組成的流水網(wǎng)，即為一個(gè)柵格形式的河網(wǎng)，可稱之為模擬河網(wǎng)。

當(dāng)采用閾值來(lái)確定河網(wǎng)起始點(diǎn)時(shí)，由于流域地形、下墊面產(chǎn)流機(jī)制以及流域所處區(qū)域氣候條件等因素的不同，其閾值也不同。本文基于DEM數(shù)據(jù)，生成了太平溪流域的模擬河網(wǎng)；采用自然流域分塊的方法生成流域子單元，將太平溪流域分為6個(gè)子單元(圖2)，其中子單元3、4、5位于太平溪的干流；子單元6位于上游；子單元1和2分別為區(qū)間入流。

圖2 太平溪流域模擬河網(wǎng)和劃分子單元示意圖Fig.2 A schematic diagram of simulating river network and dividing sub - unit in Taipingxi basin

2.3 太平溪流域網(wǎng)格劃分

本文基于無(wú)結(jié)構(gòu)網(wǎng)格來(lái)生成太平溪流域坡面二維洪水演進(jìn)的計(jì)算網(wǎng)格。

二維洪水演進(jìn)模型的邊界可以分為開(kāi)邊界與陸地邊界，開(kāi)邊界就是指流量過(guò)程、水位過(guò)程等人為設(shè)定一部分水體所形成的有界計(jì)算域；而陸地邊界是指水流運(yùn)動(dòng)的邊界，一般使用無(wú)滑移邊界來(lái)設(shè)定，認(rèn)為水流不能越過(guò)陸地邊界。陸地邊界是實(shí)際存在的邊界。由于流域邊界為山脊分水嶺，一般水流不能越過(guò)分水嶺，故本文以太平溪流域邊界作為生成網(wǎng)格的陸地邊界。此外，考慮到二維水動(dòng)力模型需要與一維水動(dòng)力模型進(jìn)行耦合，故本文在生成網(wǎng)格時(shí)將太平溪干流區(qū)域也當(dāng)作單獨(dú)的邊界進(jìn)行處理。生成網(wǎng)格時(shí)，每個(gè)網(wǎng)格的最大面積限制在2 000 km2以內(nèi)。生成網(wǎng)格后將太平溪內(nèi)的高程信息插值至各計(jì)算網(wǎng)格，結(jié)果見(jiàn)圖3。

圖3 太平溪流域二維洪水演進(jìn)網(wǎng)格示意圖Fig.3 Schematic diagram of two-dimensional flood evolution grid in Taipingxi basin

3 耦合模型構(gòu)建

3.1 邊界條件

二維洪水演進(jìn)模型的邊界條件包括兩類：①水陸邊界條件：此類邊界條件為圖4中二維洪水演進(jìn)模擬范圍的邊界，認(rèn)為水流無(wú)法越過(guò)該邊界，在水陸邊界法向處的流量及流速均為零。②流量/水位邊界條件：此類邊界條件為太平溪干流向流域坡面漫灘的過(guò)程。本文采用側(cè)向連接的方式(圖4)將一維太平溪干流與二維太平溪流域坡面進(jìn)行耦合，利用堰流計(jì)算通過(guò)側(cè)向連接的水流:

式中：q為交換水量，m3/s；W為寬度，m，一般取單元格和河道相連的邊長(zhǎng)；C為堰流系數(shù)；k為堰指數(shù)；Hus為堰上游水位,m；Hds為堰下游水位,m；Hw為堰頂高程,m。

圖4 側(cè)向連接示意圖Fig.4 Side connection diagram

3.2 動(dòng)邊界處理

由于洪水在流域坡面上運(yùn)動(dòng)的特殊性，計(jì)算區(qū)域中會(huì)存在隨水位漲落而變化的動(dòng)邊界, 為保證模型計(jì)算的連續(xù)性, 采用干/濕處理技術(shù)來(lái)對(duì)動(dòng)邊界進(jìn)行處理。其中，干水深設(shè)置為0.005 m，濕水深設(shè)置為0.07 m，淹沒(méi)水深設(shè)置為0.05 m。即當(dāng)計(jì)算區(qū)域水深小于0.005 m時(shí)，該計(jì)算區(qū)域記為“干”，不參加計(jì)算；當(dāng)水深大于0.07 m時(shí),該計(jì)算區(qū)域記為“濕”,重新參加計(jì)算。

3.3 初始條件及模型參數(shù)

在太平溪干流未發(fā)生漫灘時(shí)流域坡面上的積水較淺，本文假定此時(shí)區(qū)域內(nèi)是無(wú)水的，認(rèn)為初始水位為流域坡面的地面高程值。

由于缺少實(shí)測(cè)資料，各類地形、區(qū)域取值因地制宜，但大致保持在一定范圍內(nèi)。本文依據(jù)相關(guān)經(jīng)驗(yàn)及相關(guān)文獻(xiàn)各糙率值，確定為如下：一、二維的河道取值為0.03，村莊取值為0.07，水田取值為0.05，樹(shù)叢取值為0.05，道路以及空地取值為0.035。

4 太平溪流域洪水淹沒(méi)風(fēng)險(xiǎn)分析

選取P=10%頻率下的50 h的設(shè)計(jì)暴雨，基于本文構(gòu)建的耦合洪水演進(jìn)模型對(duì)太平溪流域的洪水淹沒(méi)風(fēng)險(xiǎn)進(jìn)行分析。圖5為P=10%設(shè)計(jì)暴雨條件下不同計(jì)算單元的出流過(guò)程。

圖5 太平溪流域P=10%設(shè)計(jì)洪水模擬過(guò)程 (子單元1、2、6)Fig.5 P=10% design flood simulation process in Taipingxi basin(subunit 1、2、6)

根據(jù)以上邊界條件模擬了洪水太平溪干流及坡面上的演進(jìn)過(guò)程。圖6為太平溪干流在P=10%時(shí)最高水位的沿程分布情況，最高水位出現(xiàn)在第20 h，也即邊界入流達(dá)到峰值的時(shí)刻。從圖6可以看出，太平溪干流大多數(shù)斷面所面臨的風(fēng)險(xiǎn)情況較小，水位均在干流主槽內(nèi)，未發(fā)生漫灘的情況。但仍有3個(gè)斷面發(fā)生了漫灘的風(fēng)險(xiǎn)事件，這3個(gè)斷面位置見(jiàn)圖7，其計(jì)算起點(diǎn)距分別為1 913.4、3 044.9和3 676.7 m。

圖6 太平溪干流P=10%最高水位分布Fig.6 P=10% maximum water level distribution in Taipingxi mainstream

圖7 太平溪干流P=10%發(fā)生風(fēng)險(xiǎn)事件斷面示意圖Fig.7 P=10% risk profile schematic in the Taipingxi basin mainstream

圖8～圖11分別為設(shè)計(jì)暴雨時(shí)間內(nèi)35 h的太平溪流域坡面淹沒(méi)的范圍及深度示意圖。從這些圖可以看出，由于漫灘事件的發(fā)生，從設(shè)計(jì)暴雨的20 h起洪水由太平溪干流進(jìn)入流域坡面進(jìn)行傳播，發(fā)生淹沒(méi)的區(qū)域位于風(fēng)險(xiǎn)斷面2，并造成太平溪流域船村的風(fēng)險(xiǎn)事件。但在P=10%設(shè)計(jì)暴雨條件下，洪水淹沒(méi)的范圍相對(duì)有限，主要集中在干流附近，并且主要向下游方向進(jìn)行演進(jìn)，未對(duì)船村造成過(guò)大的影響。從淹沒(méi)發(fā)生的時(shí)間過(guò)程來(lái)看，發(fā)生漫灘事件后2 h(圖9)的淹沒(méi)深度達(dá)到最大值，隨后隨著洪水的擴(kuò)散與歸槽，淹沒(méi)深度逐漸減小，發(fā)生漫灘事件后10 h后的淹沒(méi)高程已顯著降低。

圖8 太平溪流域淹沒(méi)范圍及深度示意圖(19 h)Fig.8 Flooding scope and depth schematic in Taipingxi basin(19 h)

圖9 太平溪流域淹沒(méi)范圍及深度示意圖(22 h)Fig.9 Flooding scope and depth schematic in Taipingxi basin(22 h)

圖10 太平溪流域淹沒(méi)范圍及深度示意圖(25 h)Fig.10 Flooding scope and depth schematic in Taipingxi basin(25 h)

圖11 太平溪流域淹沒(méi)范圍及深度示意圖(35 h)Fig.11 Flooding scope and depth schematic in Taipingxi basin(35 h)

通過(guò)對(duì)比設(shè)計(jì)暴雨過(guò)程、徑流過(guò)程以及淹沒(méi)過(guò)程發(fā)現(xiàn)，由于太平溪流域的匯流時(shí)間較短，故暴雨中心到洪峰的時(shí)間以及漫灘風(fēng)險(xiǎn)發(fā)生的時(shí)間均較短，在P=10%設(shè)計(jì)暴雨條件下的時(shí)間僅為1 h。為此，需要建立洪水預(yù)報(bào)警報(bào)系統(tǒng)。通過(guò)將實(shí)時(shí)降雨、水文等數(shù)據(jù)輸入模型得到預(yù)報(bào)信息，在必要的情況下，啟用警報(bào)系統(tǒng)，為抗洪救災(zāi)、居民撤離等提供寶貴的準(zhǔn)備時(shí)間，減少經(jīng)濟(jì)損失。

從工程措施來(lái)看，通過(guò)提升3個(gè)風(fēng)險(xiǎn)斷面的堤防可以防止漫灘事件的發(fā)生。具體來(lái)說(shuō)，風(fēng)險(xiǎn)斷面1加高5 cm、風(fēng)險(xiǎn)斷面2加高15 cm、風(fēng)險(xiǎn)斷面3加高20 cm。通過(guò)數(shù)值實(shí)驗(yàn)分析發(fā)現(xiàn)，在加高以上3個(gè)斷面的高程后，太平溪流域在P=10%設(shè)計(jì)暴雨條件下無(wú)風(fēng)險(xiǎn)事件發(fā)生。

5 結(jié) 語(yǔ)

1) 基于DEM數(shù)據(jù)生成太平溪流域的模擬河網(wǎng)，基于無(wú)結(jié)構(gòu)網(wǎng)格生成流域坡面二維洪水演進(jìn)的計(jì)算網(wǎng)格，以側(cè)向連解方式完成一維洪水演進(jìn)模型與二維洪水演進(jìn)模型的耦合。

2) 利用建立的耦合模型對(duì)太平溪流域P=10%設(shè)計(jì)暴雨條件下的洪水淹沒(méi)進(jìn)行分析，并根據(jù)分析結(jié)果，為山丘區(qū)太平溪流域提出風(fēng)險(xiǎn)應(yīng)對(duì)對(duì)策及工程和非工程措施。

3) 采用基于數(shù)學(xué)模型技術(shù)的洪水過(guò)程研究，可以準(zhǔn)確地分析流域的洪水過(guò)程特征和洪水風(fēng)險(xiǎn)因子，利用耦合模型進(jìn)行洪水風(fēng)險(xiǎn)分析研究將成為洪水風(fēng)險(xiǎn)管理的一項(xiàng)重要的技術(shù)支撐。