基于MIKE21 FM模型的洞庭湖區(qū)平原城市洪水演進(jìn)模擬

朱世云， 于永強(qiáng)， 俞芳琴， 劉 俊， 游志康

(1.河海大學(xué) 水文水資源學(xué)院， 江蘇 南京 210098； 2.中國(guó)電建集團(tuán)中南勘測(cè)設(shè)計(jì)研究院有限公司，湖南 長(zhǎng)沙 410014； 3.南京市浦口區(qū)水利局，江蘇 南京 211800)

1 研究背景

城市洪澇災(zāi)害作為城市常見(jiàn)自然災(zāi)害之一，不但阻礙城市的經(jīng)濟(jì)發(fā)展，還會(huì)造成重大人員傷亡和社會(huì)動(dòng)蕩[1]。我國(guó)是洪澇災(zāi)害頻發(fā)的國(guó)家,除沙漠、極端干旱地區(qū)和高寒地區(qū)外,我國(guó)大約2/3的國(guó)土面積都存在著不同程度和不同類型的洪澇災(zāi)害。國(guó)外對(duì)洪澇問(wèn)題做出許多積極探索，總結(jié)出豐厚的理論基礎(chǔ)。20世紀(jì)50年代，對(duì)密西西比河和俄亥俄河流域建立了水動(dòng)力學(xué)模型[2-4]，Herman等[5]利用情景模型對(duì)潰堤洪水進(jìn)行分析，Jonkman等[6-7]、Pistrika等[8]采用基于情景模擬的方法給出了定量的潰堤洪水影響結(jié)果，Segura[9]在進(jìn)行洪災(zāi)影響分析時(shí)利用GIS疊加危險(xiǎn)性、暴露性和脆弱性因子估算出概率背景下的洪災(zāi)損失。國(guó)內(nèi)洪澇分析方法主要采用了水文學(xué)方法和水力學(xué)方法，水文學(xué)主要通過(guò)數(shù)理統(tǒng)計(jì)得出洪澇災(zāi)與相關(guān)因素的相關(guān)規(guī)律對(duì)洪澇進(jìn)行分析，水力學(xué)方法主要通過(guò)對(duì)圣維南方程的推導(dǎo)和迭代進(jìn)行洪澇模擬分析[10-12]。

地處長(zhǎng)江中游平原的洞庭湖區(qū)是承納湘、資、沅、澧四水和吞吐長(zhǎng)江的洪道型調(diào)蓄湖泊，經(jīng)歷了數(shù)次湖盆擴(kuò)大和縮小的過(guò)程[13]。近年來(lái)，隨著長(zhǎng)江四口帶來(lái)的泥沙淤積和人工圍湖，導(dǎo)致湖盆抬升，人水爭(zhēng)地矛盾顯著。這一系列問(wèn)題造成洞庭湖區(qū)城市洪澇災(zāi)害頻發(fā)，洪澇損失日益增大。同時(shí)湖區(qū)城市地形平坦多位于湖區(qū)防洪高水位之下，依靠修建堤防營(yíng)造防洪包圍圈來(lái)保證水安全，堤防失事后果將是巨大經(jīng)濟(jì)損失和人員傷亡，城市將陷入癱瘓，對(duì)洪澇無(wú)法完全消除地區(qū)開(kāi)展研究是復(fù)雜且有必要的。

根據(jù)洞庭湖區(qū)平原城市獨(dú)特自然地理特征，提出洞庭湖區(qū)平原城市洪澇模擬方法，制定不同情景分析城市洪水方案，對(duì)湖區(qū)平原城市洪水問(wèn)題進(jìn)行動(dòng)態(tài)分析，可為城市洪澇災(zāi)害防治提供決策協(xié)助，指導(dǎo)洞庭湖區(qū)平原城市水利建設(shè)和規(guī)劃。對(duì)洞庭湖區(qū)平原城市的洪澇防治策略和工程規(guī)劃布局具有積極意義，對(duì)我國(guó)自然災(zāi)害研究和城市安全體系構(gòu)建做出積極探索。

2 研究方法

城市洪水分析原則上應(yīng)采用水力學(xué)法，尤其在因河道泛濫而導(dǎo)致洪水的平原城市要求采用水力學(xué)方法。水動(dòng)力學(xué)模型是以圣維南方程組為控制方程，具有明確物理意義，在計(jì)算過(guò)程中可以考慮城市中阻水建筑物和水利工程對(duì)洪水演進(jìn)的影響[14]，并且可以很好地模擬設(shè)防城市潰堤漫堤等極端水災(zāi)害情況。其中MIKE21 FM模型是基于數(shù)值解的二維淺水方程，對(duì)恒定體積的N-S方程進(jìn)行積分，對(duì)質(zhì)量、動(dòng)量、溫度、密度等都守恒，可以模擬大部分豎向均勻的二維自由水面流態(tài)[15]。

二維模型中二維水動(dòng)力學(xué)模型的控制方程如下所示：

連續(xù)方程：

(1)

動(dòng)量方程：

(2)

(3)

式中：t為時(shí)間，s；U、V分別為x、y方向的單寬流量，m2/s；H、h分別為水位和水深，m；n為糙率系數(shù)；g為重力加速度，m/s2;q為源匯項(xiàng)。

本文針對(duì)研究區(qū)域采用MIKE21 FM模型，選擇堤防設(shè)計(jì)水位與歷史最高水位作為邊界條件，擬定外洪量級(jí)，根據(jù)堤防險(xiǎn)工段排查確定潰口位置，分析不同模擬情景下的洪水淹沒(méi)特征。

3 實(shí)例分析

3.1 研究區(qū)域概況及資料來(lái)源

岳陽(yáng)市位于長(zhǎng)江中下游兩湖平原的洞庭湖平原，處于湖南省東北部，環(huán)抱洞庭，瀕臨長(zhǎng)江，地勢(shì)東北高、西南低，呈東西走向。城區(qū)位于洞庭湖出口與長(zhǎng)江交匯帶，多為平原地形，外圍有低矮丘陵夾雜，高程分布較均勻。岳陽(yáng)以上流域集水面積約130×104km2，其中長(zhǎng)江占104×104km2，洞庭湖約26×104km2。屬?gòu)闹衼啛釒虮眮啛釒н^(guò)渡濕潤(rùn)大陸性季風(fēng)氣候，溫暖濕潤(rùn)，四季分明，季節(jié)性強(qiáng)。湖區(qū)氣候均一，山地氣候懸殊。多年平均降水量為1 439.1 mm，春夏季多于秋冬季，東部多于西部，年內(nèi)分布和年際分布都不均勻，春夏季降雨量占年均降水量70%，年平均降雨最多時(shí)發(fā)生在1954年(2 191.4 mm)，降雨最少時(shí)發(fā)生在2011年(945.7 mm)。岳陽(yáng)市是典型的洞庭湖區(qū)平原城市，地理位置突出，河湖關(guān)系復(fù)雜，洪澇災(zāi)害頻繁。

本文收集研究區(qū)域內(nèi)1∶500和1∶10000矢量地圖及行政區(qū)劃圖等地形資料，以及鄰近主要控制站(岳陽(yáng)、七里山、蓮花塘等基本站)1954、1998、1996年及其他大水年實(shí)測(cè)洪水資料。

3.2 研究區(qū)劃分

研究區(qū)域中東風(fēng)湖大堤和南湖大堤內(nèi)側(cè)地勢(shì)較低，人口稠密，且有內(nèi)滑坡、風(fēng)浪沖刷、散浸和蟻患等險(xiǎn)工段，選定東風(fēng)湖潰口和南津港潰口進(jìn)行分析，潰口位置見(jiàn)圖1。外洪研究范圍根據(jù)潰口位置及地形條件，劃分為東風(fēng)湖計(jì)算分區(qū)和南湖計(jì)算分區(qū)，面積分別為18.15 km2和96.77 km2。

3.3 模型邊界條件設(shè)置

建立岳陽(yáng)市主城區(qū)二維水力學(xué)模型，潰口采用水位過(guò)程控制。考慮到汛期洞庭湖水位較高，水量較大，潰口后對(duì)洞庭湖水位影響甚微，水位過(guò)程考慮采用定值。

堤防標(biāo)準(zhǔn)是由流域防洪標(biāo)準(zhǔn)制定，當(dāng)外湖水位達(dá)到堤防設(shè)計(jì)水位時(shí)，可能發(fā)生潰堤險(xiǎn)情。在此基礎(chǔ)上，岳陽(yáng)市西臨洞庭湖出口處遭遇歷時(shí)最高洪水位時(shí)，則是工程最不利情景，此時(shí)外洪致災(zāi)性歷時(shí)最大。因此為使情景模擬方案更符合實(shí)際情況，對(duì)岳陽(yáng)市臨洞庭湖洪水水量級(jí)設(shè)置為堤防設(shè)計(jì)水位和歷史最高水位對(duì)應(yīng)不同流量量級(jí)兩種情景。

圖1 潰口位置示意圖

根據(jù)岳陽(yáng)、七里山等水文站和流量站實(shí)測(cè)資料，參照歷史洪水系列，岳陽(yáng)市外湖洪水堤防設(shè)計(jì)和歷史最高標(biāo)準(zhǔn)下的邊界條件分別采用城陵磯七里山站1954年最高水位32.61 m、1998年最高水位34.0 m。

洞庭湖區(qū)堤垸潰口一般在100～500 m之間，由潰口內(nèi)外地形條件、設(shè)計(jì)水位等分析潰口可能發(fā)展寬度，擬定潰口流量量級(jí)范圍，最終根據(jù)擬定流量量級(jí)經(jīng)模型試算確定潰口寬度。外洪潰口洪水量級(jí)設(shè)置見(jiàn)表1。

表1 岳陽(yáng)市外洪洪水量級(jí) m3/s

根據(jù)洞庭湖區(qū)堤垸歷史潰決調(diào)查資料，采用瞬間潰的形式，潰口參數(shù)設(shè)置見(jiàn)表2。

3.4 模型構(gòu)建

模型網(wǎng)格剖分采用不規(guī)則三角網(wǎng)格，局部地形變化較大的區(qū)域進(jìn)行局部加密。共剖得網(wǎng)格18 785個(gè)，網(wǎng)格平均面積0.003 8 km2，最大網(wǎng)格面0.008 9 km2。為保證模型精度和可靠度，采用高階方法，時(shí)間步長(zhǎng)在0.01～30 s之間，臨界CFL數(shù)為0.8。

岳陽(yáng)市內(nèi)的南湖、東風(fēng)湖水草、石塊數(shù)量較少，湖面開(kāi)闊，淤塞雜草叢生的情況較為少見(jiàn)，糙率n在0.025～0.028之間，陸地面積不透水較高，下墊面硬化程度高，多為灰色建筑物，糙率范圍0.014～0.016，城市中的綠地范圍糙率取值在0.018～0.021，洪水演進(jìn)模型糙率選取標(biāo)準(zhǔn)見(jiàn)表3，糙率場(chǎng)取值情況見(jiàn)圖2。

表2 潰口參數(shù)設(shè)置

表3 不同下墊面洪水模型糙率選取標(biāo)準(zhǔn)

采用MIKE21中水工建筑物堰來(lái)展現(xiàn)外洪潰口特性，采用寬頂堰類型，無(wú)閥控制流向。南津港潰口處堰底高程27 m，根據(jù)不同流量級(jí)堰寬寬度在97～312 m之間；東風(fēng)湖潰口處堰底高程27.68 m，根據(jù)不同流量級(jí)堰寬寬度在115～405 m之間，設(shè)置的潰口堰位置見(jiàn)圖3。

本文采用正壓模型，參考溫度恒為10℃，參考鹽度恒為32PSU。渦黏系數(shù)選用Smagorinsky公式表達(dá)，Smagorinsky系數(shù)選用恒定值0.28。根據(jù)東洞庭湖洪水來(lái)源，考慮東風(fēng)湖潰口和南津港潰口，堤防設(shè)計(jì)水位和歷史最高水位標(biāo)準(zhǔn)設(shè)2000、4000、6000、8000 m3/s流量級(jí)，洪水演進(jìn)模型共設(shè)置8種模擬情景，見(jiàn)表4。

3.5 模型成果

由模型運(yùn)行結(jié)果可以得出岳陽(yáng)市發(fā)生高量級(jí)洪水時(shí)，洪水演進(jìn)不同時(shí)段的淹沒(méi)情況、洪水流態(tài)和進(jìn)洪量。岳陽(yáng)市洪水演進(jìn)模擬進(jìn)洪量情況見(jiàn)表5，各潰口各方案進(jìn)洪流量過(guò)程見(jiàn)圖4～6。

表4 岳陽(yáng)市8種洪水演進(jìn)模型模擬情景方案

表5 不同方案條件下不同進(jìn)洪歷時(shí)的進(jìn)洪量模擬結(jié)果 108 m3

圖2洪水模型糙率取值示意圖圖3洪水演進(jìn)模型潰口位置

圖4南津港潰口情景進(jìn)洪流量變化情況圖5東風(fēng)湖潰口情景進(jìn)洪流量變化情況

圖6 不同模擬情景進(jìn)洪量變化情況

不同潰口，進(jìn)洪量不同，東風(fēng)湖潰口堤防設(shè)計(jì)水位標(biāo)準(zhǔn)和歷史最高水位標(biāo)準(zhǔn)最大進(jìn)洪量分別為0.352×108m3和0.442×108m3，南湖則分別為0.876×108m3和1.049×108m3，二者進(jìn)洪量相差約一倍；不同洪水標(biāo)準(zhǔn)，進(jìn)洪量不同，以南津港潰口為例，歷史最高水位標(biāo)準(zhǔn)比堤防設(shè)計(jì)水位標(biāo)準(zhǔn)進(jìn)洪量多2 700×104m3。潰口流量量級(jí)越大，洪水演進(jìn)越快，以南津港潰口為例，堤防設(shè)計(jì)水位標(biāo)準(zhǔn)2 000 m3/s流量級(jí)，14 h后進(jìn)洪流量達(dá)到平衡，而歷史最高水位標(biāo)準(zhǔn)8 000 m3/s流量級(jí)，5 h后進(jìn)洪流量即達(dá)到平衡，二者相差9 h。東風(fēng)湖和南津港潰口最大淹沒(méi)范圍見(jiàn)圖7。

分別以東風(fēng)湖潰口2 000 m3/s流量級(jí)和南津港潰口4 000 m3/s流量級(jí)為例，不同進(jìn)洪歷時(shí)淹沒(méi)范圍見(jiàn)圖8和圖9。

岳陽(yáng)市中心城區(qū)地勢(shì)相對(duì)較高，東風(fēng)湖和南津港潰口對(duì)中心城區(qū)影響不大，主要影響南湖及東風(fēng)湖周邊地勢(shì)較低矮的區(qū)域。

不同方案，洪水淹沒(méi)范圍不同。潰口均位于內(nèi)湖臨洞庭湖出口堤段，潰堤洪水先入內(nèi)湖，臨近區(qū)域地勢(shì)相對(duì)較低，影響范圍相對(duì)較小，南湖潰口和東風(fēng)湖潰口相比，淹沒(méi)范圍相差9 km2，主要原因是南湖調(diào)蓄區(qū)范圍本身較大。

圖7 不同潰口最大淹沒(méi)范圍

圖8 東風(fēng)湖潰口2 000 m3/s洪水方案淹沒(méi)范圍變化

圖9 南津港潰口4 000 m3/s洪水方案淹沒(méi)范圍變化

4 結(jié) 論

本文采用MIKE21 FM模型分別計(jì)算研究區(qū)域8個(gè)洪水模擬方案，得出的主要結(jié)論如下：

(1)研究區(qū)歷史最高水位情景下淹沒(méi)范圍較大，堤防設(shè)計(jì)水位情景下淹沒(méi)歷時(shí)更長(zhǎng)。水位越高，洪水演進(jìn)速度越快，淹沒(méi)深度也越大。外湖潰堤洪水影響范圍有限，建議開(kāi)展城區(qū)內(nèi)澇風(fēng)險(xiǎn)研究。

(2)研究結(jié)果合理可靠，可為湖區(qū)平原城市防汛救災(zāi)工作提供參考。表明MIKE21 FM模型在湖區(qū)平原城市的洪水?dāng)?shù)值模擬中具有一定應(yīng)用價(jià)值。

參考文獻(xiàn)：

[1] 古葒歡.基于社區(qū)尺度的上海市自然災(zāi)害社會(huì)脆弱性評(píng)估[D]. 上海：上海師范大學(xué), 2016.

[2] TODOROVIC P, ZELENHASIC E. A stochastic model for flood analysis[J]. Water Resources Research, 2010,6(6):1641-1648.

[3] ZELENHASIC E. On the extreme streamflow drought analysis[J].Water Resources Management, 2002,16(2):105-132.

[4] TODOROVIC P. On some problems involving random number of random variables[J]. Annals of Mathematical Statistics, 1970,41(3):1059-1063.

[5] HERMAN J， LOUAT T， HUANG Q， et al. Autoimmune and inflammatory disorder therapy: WO， US9301961[P]. 2016.

[6] JONKMAN S N， MAASKANT B， BOYD E， et al. Loss of life caused by the flooding of New Orleans after Hurricane Katrina: analysis of the relationship between flood characteristics and mortality.[J]. Risk Analysis， 2009， 29(5):676-698.

[7] JONKMAN S N. Global perspectives on loss of human life caused by floods[J].Natural Hazards,2005,34(2):151-175.

[8] PISTRIKA A K, JONKMAN S N. Damage to residential buildings due to flooding of New Orleans after hurricane Katrinal[J].Natural Hazards,2010,54(2):413-434.

[9] SEGURA G, BADILLA E, OBANDO L. Susceptibilidad al deslizamiento en el corredor siquirres-turrialba[J]. Revista Geológica De América Central, 2011(45):101-121.

[10] 毛德華.洪災(zāi)綜合風(fēng)險(xiǎn)分析的理論方法與應(yīng)用研究[M].北京：中國(guó)水利水電出版社，2009.

[11] 劉華振，劉 俊，左 君，等.馬斯京根法在黃河吳堡-龍門(mén)區(qū)間洪水演算中的應(yīng)用[J].水電能源科學(xué)，2012，30(6):53-55.

[12] 程曉陶. 我國(guó)推進(jìn)洪水風(fēng)險(xiǎn)圖編制工作基本思路的探討[J]. 中國(guó)水利, 2005(17):11-13.

[13] 王月容. 洞庭湖退田還湖區(qū)錢(qián)糧湖垸景觀格局、土壤質(zhì)量與土地承載力研究[D]. 武漢：華中農(nóng)業(yè)大學(xué), 2010.

[14] 薛文宇. 城市暴雨積水及街道洪水模擬模型研究[D]. 天津：天津大學(xué), 2015.

[15] 房克照, 尹 晶, 孫家文,等. 基于二維淺水方程的滑坡體興波數(shù)值模型[J]. 水科學(xué)進(jìn)展, 2017, 28(1):96-105.