MIKE21模型在精河防洪保護(hù)區(qū)避洪轉(zhuǎn)移中的應(yīng)用研究

韓 嶺，蓋永崗，李榮容，雷 鳴

(黃河勘測(cè)規(guī)劃設(shè)計(jì)有限公司，鄭州 450003)

0 引 言

新疆北部地區(qū)年平均降雨量相對(duì)較小，但由于其獨(dú)特的地形地貌特征，暴雨融雪型洪水也較為突發(fā)，對(duì)于山前平原河流，這種突發(fā)的洪水對(duì)其影響更為顯著，嚴(yán)重影響河流兩岸防洪保護(hù)區(qū)內(nèi)的居民生命財(cái)產(chǎn)安全。當(dāng)洪水發(fā)生時(shí)，如何保護(hù)災(zāi)區(qū)居民生命安全，緊急疏散轉(zhuǎn)移受災(zāi)群眾，保障轉(zhuǎn)移群眾基本生活，實(shí)現(xiàn)快速、高效轉(zhuǎn)移是防洪應(yīng)急預(yù)案體系建設(shè)中急需解決的問題[1-3]。隨著防洪理念由控制洪水向洪水管理的轉(zhuǎn)變，洪水?dāng)?shù)值模擬已成為洪水風(fēng)險(xiǎn)分析的一個(gè)重要手段，是洪水管理的科學(xué)依據(jù)[4-6]。

目前比較流行的水動(dòng)力數(shù)值模擬軟件有荷蘭的Delft 3D模型、英國的InfoWorks RS模型、丹麥的MIKE模型及其他自主開發(fā)的軟件模型[7]，雖然這些軟件在功能實(shí)現(xiàn)、計(jì)算算法和建模方法上各有不同，但都得到了廣泛的使用和工程驗(yàn)證，具有很高的可信性。其中MIKE模型因其強(qiáng)大的前、后處理功能，是目前世界范圍內(nèi)應(yīng)用最為廣泛的水力學(xué)模型之一[8-12]。本文擬采用MIKE 21 模型，對(duì)精河防洪保護(hù)區(qū)內(nèi)的避洪轉(zhuǎn)移進(jìn)行分析研究。

1 研究區(qū)概況

精河流域位于新疆維吾爾自治區(qū)西北部，總體地勢(shì)向北逐漸降低，精河發(fā)源于婆羅科努山北坡，經(jīng)出山口后進(jìn)入平原區(qū)，河槽逐漸變?yōu)閷挏\，河流坡降逐漸減小。本次研究區(qū)域位于出山口至芒丁鄉(xiāng)之間，屬河道中下游，是精河流域主要的工農(nóng)業(yè)生產(chǎn)區(qū)，河道兩岸居民點(diǎn)密布。該河段河道長(zhǎng)約30 km，河道基本順直，河床寬度200～400 m，縱坡1.48%～0.622%，自上而下河道逐漸展寬。研究區(qū)域內(nèi)的防洪保護(hù)區(qū)范圍主要覆蓋精河縣的托里鄉(xiāng)、芒丁鄉(xiāng)以及精河鎮(zhèn)，總?cè)丝诩s6萬人，耕地面積約2.33 萬hm2，農(nóng)林牧漁業(yè)總產(chǎn)值13.5億元。研究區(qū)行政區(qū)劃圖見圖1。

圖1 研究區(qū)行政區(qū)劃圖Fig.1 Study regional administrative divisions

2 模型原理

精河河道及防洪保護(hù)區(qū)采用MIKE21 FM (非結(jié)構(gòu)化網(wǎng)格)水動(dòng)力模塊建模，其原理是基于數(shù)值解的二維淺水方程，可以模擬各種作用力下產(chǎn)生的水位和水流變化及任何忽略分層的二維自由表面流。在平面上采用非結(jié)構(gòu)化網(wǎng)格，采用的數(shù)值方法是單元中心的有限體積法；控制方程離散時(shí)，結(jié)果變量u、v位于單元中心，跨邊界通量垂直于單元邊；有限體積法中法向通量通過在沿外法向建立單元水力模型并求解一維黎曼問題而得到；采用顯式時(shí)間積分。

2.1 控制方程

Mike21 FM二維非恒定流計(jì)算模塊的原理基于二維不可壓縮流體雷諾平均應(yīng)力方程，服從布辛涅斯克假設(shè)和靜水壓力假設(shè)。

h=η+d

(1)

二維非恒定淺水方程組為：

(2)

(3)

(4)

(5)

2.2 數(shù)值解法

(1)空間離散。計(jì)算區(qū)域的空間離散是用有限體積法，將該連續(xù)統(tǒng)一體細(xì)分為不重疊的單元，單元可以是三角形或者四邊形。

在笛卡爾坐標(biāo)系中，二維淺水方程組可以寫為：

(6)

式中：U為守恒型物理矢量；F為通量矢量；S為源項(xiàng)，上標(biāo)I、V分別為無黏性的和黏性通量。對(duì)方程(6)第i個(gè)單元積分，并運(yùn)用Gauss原理重寫可得出：

(7)

式中：Ai為單元Ωi的面積；Fi為單元的邊界；ds為沿著邊界的積分變量。使用單點(diǎn)求積分來計(jì)算面積的積分，該求積分點(diǎn)位于單元的質(zhì)點(diǎn)，同時(shí)使用中點(diǎn)求積法來計(jì)算邊界積分，方程(7)可以寫為：

(8)

式中：Ui和Si分別為第i個(gè)單元的U和S的平均值，并位于單元中心；NS是單元的邊界數(shù)；ΔΓj為第j個(gè)單元的長(zhǎng)度。

一階解法和二階解法都可以用于空間離散求解，為了避免數(shù)值振蕩，模型使用了二階TVD格式。

(2)時(shí)間積分。對(duì)于二維模擬，淺水方程的求解有兩種方法：一種是低階方法，另一種是高階方法。低階方法即低階顯式的Euler方法：

Un+1=Un+ΔtG(Un)

(9)

式中：Δt為時(shí)間步長(zhǎng)。高階的方法為使用了二階的Runge Kutta方法：

(10)

3 模型建立與驗(yàn)證

3.1 模型構(gòu)建

精河防洪護(hù)區(qū)建模范圍南至下天吉水庫壩下，北至精河縣北郊，西至大莊子村，東至東莊村，河道長(zhǎng)約30 km，建模范圍189 km2。二維水動(dòng)力學(xué)模型建模范圍見圖1。

整個(gè)防洪保護(hù)區(qū)網(wǎng)格采用不規(guī)則三角形網(wǎng)格進(jìn)行剖分，為保證模型計(jì)算精度，①在高山邊界、精河主河道、精伊霍鐵路和G30連霍高速路沿線兩側(cè)網(wǎng)格適當(dāng)加密；②在河道內(nèi)任一橫截面上的網(wǎng)格數(shù)量不少于4個(gè)來確保河道內(nèi)有足夠模擬河道水流演進(jìn)的實(shí)際計(jì)算點(diǎn)。模型范圍內(nèi)共剖分網(wǎng)格數(shù)為83 956個(gè)，平均網(wǎng)格面積為0.002 km2。網(wǎng)格剖分結(jié)果見圖2。網(wǎng)格剖分后，采用1∶10 000比例尺的DEM高程數(shù)據(jù)進(jìn)行網(wǎng)格高程插值，得到二維模型的地形云圖3。

圖2 模型網(wǎng)格剖分圖 Fig.2 Mesh generation of model

圖3 模型地形云圖Fig.3 Model terrain cloud of model

3.2 邊界條件

精河防洪保護(hù)區(qū)二維水動(dòng)力學(xué)模型的上邊界為流量邊界，設(shè)置于下天吉水庫壩下精河干流河道處。下天吉水庫壩址以上天然設(shè)計(jì)洪水根據(jù)精河山口水文站實(shí)測(cè)長(zhǎng)系列資料推求，選擇精河山口水文站1999年7月19日0時(shí)-7月26日0時(shí)的7日洪水過程，采用同頻率放大法計(jì)算的天然設(shè)計(jì)洪水過程線，按照下天吉水庫一期工程擬定的防洪運(yùn)用方式對(duì)不同量級(jí)設(shè)計(jì)洪水過程進(jìn)行調(diào)蓄計(jì)算，得到百年一遇洪水的洪水過程線，見圖4。

圖4 精河百年一遇洪水過程線(1999年)Fig.4 Hundred-year flood hydrograph of Jing River

模型下邊界為水位邊界，設(shè)置于精河縣城北郊，通過模型試算找出模型下邊界線中可能的出流位置，下邊界的斷面水位流量關(guān)系線根據(jù)實(shí)測(cè)斷面用曼寧公式進(jìn)行推求，計(jì)算中所采用的比降系根據(jù)實(shí)測(cè)該河段河底高程和河段長(zhǎng)計(jì)算得出，作為水面比降，斷面糙率參考精河山口水文站及相關(guān)設(shè)計(jì)成果，并結(jié)合現(xiàn)場(chǎng)查勘情況，取0.028～0.035。水位流量關(guān)系線見圖5。

圖5 模型下邊界水位流量關(guān)系Fig.5 Relationship between water level and water level of model Lower Boundary

3.3 模型參數(shù)選取與驗(yàn)證

(1)模型參數(shù)選取。在精河防洪保護(hù)區(qū)二維水動(dòng)力學(xué)模型中，模型計(jì)算主要是用來模擬洪水在淹沒區(qū)的演進(jìn)過程，提取淹沒水深、流速和洪水到達(dá)時(shí)間等風(fēng)險(xiǎn)要素，故對(duì)于計(jì)算結(jié)果無明顯影響的計(jì)算參數(shù)均采用默認(rèn)值，需要設(shè)置的主要參數(shù)有：計(jì)算時(shí)間和步長(zhǎng)、糙率、干濕邊界等。

經(jīng)反復(fù)調(diào)試，本次計(jì)算區(qū)二維模型的最大時(shí)間步長(zhǎng)設(shè)定為30 s，最小時(shí)間步長(zhǎng)設(shè)定為0.01 s，模型根據(jù)網(wǎng)格質(zhì)量、進(jìn)洪情況及區(qū)域地形復(fù)雜程度會(huì)自動(dòng)調(diào)整計(jì)算時(shí)間步長(zhǎng)。

在二維水動(dòng)力模型中，糙率是影響模型計(jì)算的另一個(gè)比較重要的參數(shù)。精河風(fēng)險(xiǎn)圖保護(hù)區(qū)網(wǎng)格是由面積小于0.005 km2的不規(guī)則三角形網(wǎng)格構(gòu)成，在MIKE21中可以設(shè)定每個(gè)網(wǎng)格不同的糙率值，也可全區(qū)統(tǒng)一設(shè)定為一個(gè)固定的糙率值。

在本次計(jì)算中，利用本區(qū)域土地利用分類圖、遙感衛(wèi)星影像圖，并結(jié)合現(xiàn)場(chǎng)調(diào)查，考慮分區(qū)內(nèi)的地形、地貌、植被狀況，根據(jù)《洪水風(fēng)險(xiǎn)圖編制技術(shù)細(xì)則》[13]、《水力計(jì)算手冊(cè)》[14]等資料，對(duì)不同下墊面賦予不同的糙率值，精河防洪保護(hù)區(qū)內(nèi)地面插值后的糙率值云圖見圖6。

圖6 模型地面糙率值云圖Fig.6 Surface roughness map of model

干濕邊界是MIKE21水動(dòng)力學(xué)模型中為避免模型計(jì)算出現(xiàn)不穩(wěn)定性和不收斂而設(shè)定的參數(shù)，當(dāng)某一網(wǎng)格單元的水深小于濕水深時(shí)，在此單元上的水流計(jì)算會(huì)被相應(yīng)調(diào)整，而當(dāng)水深小于干水深時(shí)，會(huì)被凍結(jié)而不參與計(jì)算。干水深、浸沒水深和濕水深分別取0.005、0.05和0.1 m。

(2)模型參數(shù)驗(yàn)證。根據(jù)精河山口水文站、精河渡槽斷面水位流量關(guān)系線，通過在二維模型中設(shè)置水位、流量監(jiān)控?cái)嗝妫玫娇刂茢嗝嫠涣髁筷P(guān)系線，將模型計(jì)算結(jié)果與設(shè)計(jì)報(bào)告成果進(jìn)行比較，見圖7和圖8。

圖7 山口水文站斷面水位流量關(guān)系Fig.7 Relationship between water level and water level in the section of Shankou Hydrological Station

圖8 精河渡槽斷面水位流量關(guān)系Fig.8 Relationship between water level and water level in the section of Jinghe River aqueduct

可見，精河山口水文站、精河渡槽斷面二維模型計(jì)算結(jié)果與設(shè)計(jì)報(bào)告成果非常接近。根據(jù)相關(guān)設(shè)計(jì)成果報(bào)告，精河渡槽斷面堤防設(shè)計(jì)防流量為118 m3/s，相應(yīng)設(shè)計(jì)水位為403.76 m。二維模型計(jì)算得到118 m3/s相應(yīng)水位為403.769 m，二者的絕對(duì)誤差為9 cm。從水位流量關(guān)系和堤防設(shè)計(jì)水位成果比較結(jié)果來看，本次率定的參數(shù)可以較好的模擬該河段的洪水演進(jìn)過程。設(shè)計(jì)報(bào)告水位和本次計(jì)算水位相差比較小，可以認(rèn)為模型搭建合理。

4 模型計(jì)算及避洪轉(zhuǎn)移分析

4.1 模型計(jì)算分析

當(dāng)精河山口水文站發(fā)生100年一遇洪水時(shí)，經(jīng)下天吉水庫調(diào)蓄后洪峰流量為332 m3/s。洪水進(jìn)入保護(hù)區(qū)后，根據(jù)地形演進(jìn)，演進(jìn)至精伊霍鐵路處，由于河道過流能力較大，洪水均在河道內(nèi)演進(jìn)，精伊霍鐵路不起阻水作用；演進(jìn)至精河渡槽處，洪水開始發(fā)散并產(chǎn)生漫溢，主河道左岸部分地區(qū)受淹；演進(jìn)至G30連霍高速處，部分洪水被高速阻滯，沿高速向東西兩方向演進(jìn)，并經(jīng)由連霍高速大橋及高速上的涵洞繼續(xù)向北演進(jìn)；約30 h后洪水淹沒面積基本達(dá)到最大，此后洪水流量及各淹沒區(qū)域的淹沒水深、淹沒面積逐漸趨于穩(wěn)定，洪水演進(jìn)不同時(shí)間分布見圖9。

圖9 洪水演進(jìn)不同時(shí)間分布Fig.9 Flood routing distribution at different time

4.2 避洪轉(zhuǎn)移分析

精河防洪保護(hù)區(qū)劃分的危險(xiǎn)區(qū)面積小于1 000 km2，根據(jù)《避洪轉(zhuǎn)移圖編制技術(shù)要求》[15]，避洪轉(zhuǎn)移單元按行政村考慮；同時(shí)滿足水深<1.0 m、流速<0.5 m/s，且具有可容納該區(qū)域人口的安全場(chǎng)所和設(shè)施的，釆取就地安置方式，不滿足上述條件的區(qū)域采取轉(zhuǎn)移安置方式。

根據(jù)模型計(jì)算結(jié)果，當(dāng)發(fā)生100年一遇洪水時(shí)，研究區(qū)域內(nèi)水深、流速及到達(dá)時(shí)間情況如圖10、圖11、圖12。

圖10 水深圖Fig.10 Distribution map of water depth

圖11 流速圖Fig.11 Distribution map of flow velocity

圖12 到達(dá)時(shí)間Fig.12 Time of arrival map

經(jīng)對(duì)淹沒區(qū)域內(nèi)水深、流速及到達(dá)時(shí)間的分析，根據(jù)《避洪轉(zhuǎn)移圖編制技術(shù)要求》，可知：有避洪轉(zhuǎn)移任務(wù)的行政村(社區(qū))為：芒丁鄉(xiāng)的巴西莊子村、河西村、五棵數(shù)村、夾巴溝村和精河鎮(zhèn)的鎮(zhèn)綠園社區(qū)、濱河社區(qū)及平社區(qū)，涉及人口4 200人。

但，這些區(qū)域淹沒水深均小于1 m，其中絕大部分區(qū)域淹沒水深小于0.5 m，淹沒水深在0.5～1.0 m的區(qū)域僅涉及五棵數(shù)村沿河的少部分區(qū)域，且需進(jìn)行避洪轉(zhuǎn)移的區(qū)域洪水流速均小于0.5 m/s。

因此，精河防洪保護(hù)區(qū)受洪水影響的人口均符合就地安置條件，不需進(jìn)行轉(zhuǎn)移安置，精河防洪保護(hù)區(qū)就地安置人口見表1。

表1 精河防洪保護(hù)區(qū)就地安置人口表Tab.1 Population locally settlement in Jing River Flood control protected area

通過典型調(diào)查得知，精河防洪保護(hù)區(qū)內(nèi)受洪水影響的區(qū)域集中在精河縣城附近，經(jīng)濟(jì)社會(huì)發(fā)展程度相對(duì)較高，居民房屋的類型以混凝土框架結(jié)構(gòu)、磚混結(jié)構(gòu)和磚木結(jié)構(gòu)為主，房屋質(zhì)量相對(duì)較好；根據(jù)安置區(qū)劃定的就近原則、安全性、通達(dá)性等原則，可對(duì)受影響人員進(jìn)行就地安置，農(nóng)村居民可直接安置于自己房屋的屋頂，城鎮(zhèn)居民可直接安置與樓房的二層或以上樓層，安置場(chǎng)所是比較可靠的。

5 結(jié) 語

防洪保護(hù)區(qū)避洪轉(zhuǎn)移分析主要是在前期洪水模擬的基礎(chǔ)上，結(jié)合現(xiàn)場(chǎng)調(diào)研情況，進(jìn)行危險(xiǎn)區(qū)及轉(zhuǎn)移單元的確定、轉(zhuǎn)移方式的選擇及人口分析，確定避洪轉(zhuǎn)移安置點(diǎn)[16]。本文通過MIKE 21建立了精河河道及防洪保護(hù)區(qū)內(nèi)的整體二維水動(dòng)力學(xué)模型，合理模擬了該河段遭遇百年一遇洪水時(shí)的風(fēng)險(xiǎn)情況，表明MIKE 21模型對(duì)類似山前平原河段的洪水模擬具有良好的適應(yīng)性；同時(shí)在模型計(jì)算結(jié)果的基礎(chǔ)上，分析了防洪保護(hù)區(qū)內(nèi)的淹沒水深、流速和洪水到達(dá)時(shí)間等洪水風(fēng)險(xiǎn)要素，根據(jù)分析結(jié)果確定了防洪保護(hù)區(qū)內(nèi)需轉(zhuǎn)移單元、轉(zhuǎn)移人口和可安置區(qū)域等信息，為防汛應(yīng)急指揮部門在實(shí)際避洪轉(zhuǎn)移安置工作中提供了輔助決策和技術(shù)支撐作用。

□