二維水動(dòng)力學(xué)模型在蓄滯洪區(qū)洪水演進(jìn)模擬分析中的應(yīng)用

李允軍，徐 青，周元斌，吉慶偉，向小華

（1.南瑞集團(tuán)有限公司，江蘇 南京 210003；2.江蘇省駱運(yùn)水利工程管理處，江蘇 宿遷 223800；2.河海大學(xué)，江蘇 南京 210003）

1 概述

洪水時(shí)空遷移、波動(dòng)特性復(fù)雜，加之近年來(lái)極端性氣候頻發(fā)，水文情勢(shì)多變，影響了洪水調(diào)度決策的科學(xué)性。開(kāi)展洪水演進(jìn)規(guī)律研究，對(duì)于蓄滯洪區(qū)的科學(xué)運(yùn)用和防洪減災(zāi)具有重要的指導(dǎo)意義。

國(guó)內(nèi)外洪水淹沒(méi)分析方法經(jīng)歷了從根據(jù)歷史洪災(zāi)調(diào)查資料勾畫(huà)洪水范圍，到采用水文學(xué)和水力學(xué)方法建立數(shù)學(xué)模型模擬洪水演進(jìn)過(guò)程獲得豐富的洪水水力特征值的過(guò)程。二維水動(dòng)力學(xué)模型用于蓄滯洪區(qū)、堤防保護(hù)區(qū)等洪水風(fēng)險(xiǎn)分析始于20世紀(jì)末[1]，結(jié)合當(dāng)時(shí)的最新技術(shù)，中國(guó)水利水電科學(xué)研究院研發(fā)了基于二維不規(guī)則網(wǎng)格和有限體積法的二維水動(dòng)力學(xué)模型軟件。采用二維數(shù)學(xué)模型能夠提供更加詳細(xì)的水情信息，目前國(guó)內(nèi)外越來(lái)越注重采用水力學(xué)方法進(jìn)行數(shù)值模擬。

2 二維水動(dòng)力學(xué)模型

2.1 基本方程

河流平面二維水動(dòng)力學(xué)模型基本方程[2]為：

式中：h為水深，g為重力加速度，zb為河底高程，n為糙率系數(shù)，x、y為縱向坐標(biāo)和橫向坐標(biāo)，u、v為垂向平均流速的縱向分量和橫向分量。

令qu=hu、qv=hv，將上述方程寫(xiě)成矢量形式[3]如下：

且Ut+F（U）x+G（U）y=S。

定義雅可比矩陣AF、AG：

由式（3）將式（2）寫(xiě)成非守恒形式即：

分析其特征結(jié)構(gòu)，計(jì)算特征值和相應(yīng)的右特征向量分別為：

特征向量

通常情況下（即h≠0時(shí)）雅可比矩陣AF、AG的特征值 SF1、SF2、SF3、SG1、SG2、SG3兩兩不相等，說(shuō)明式（2）、式（4）為雙曲型方程組。

2.2 模型的數(shù)值解

模型的數(shù)值解法包括有限差法（FDM）、特征法（MOC）、有限元法（FEM）、有限體積法（FVM）等，F(xiàn)VM將計(jì)算域劃分成若干規(guī)則或不規(guī)則形狀的單元或控制體。在計(jì)算出通過(guò)每個(gè)控制體邊界沿法向輸入（出）的流量和動(dòng)量通量后，對(duì)每個(gè)控制體分部進(jìn)行水量和動(dòng)量的平衡計(jì)算，便得到計(jì)算時(shí)段末各控制體平均水深和流速。因此，F(xiàn)VM正是對(duì)于推導(dǎo)原始微分方程所用控制體的回歸，與FDM和FEM的數(shù)值逼近相比，其物理意義更直接明晰[4-5]，其控制體示意圖如圖1。

圖1 控制體示意圖

在圖1中的控制體上建立積分方程，并進(jìn)行離散，可得控制體的離散形式：

式中：ΔVi為控制體i的面積，N為控制體的界面數(shù)為源項(xiàng)S在控制體i上的平均值，δlij為控制體i、j之間的界面長(zhǎng)度。

2.3 定解條件

2.3.1 初始條件

二維模型給定各計(jì)算網(wǎng)格點(diǎn)上水位、流速初值：

2.3.2 邊界條件

固壁邊界為：

開(kāi)邊界給定水位過(guò)程線(xiàn)：

2.3.3 模型計(jì)算的條件設(shè)定

（1）初始條件設(shè)定

網(wǎng)格產(chǎn)生以后，必須給各個(gè)計(jì)算單元賦予初始狀態(tài)。初始條件可以是網(wǎng)格結(jié)點(diǎn)的初始水面高程或水深，或者x、y方向上的初始流速。初始條件的規(guī)定，一是根據(jù)問(wèn)題的物理需求，如靜水或均勻流；二是根據(jù)部分地點(diǎn)的觀(guān)測(cè)數(shù)據(jù)，所缺空間分布由內(nèi)插估計(jì)。常設(shè)初始流動(dòng)為已達(dá)平衡態(tài)的恒定流，初始條件的誤差隨著時(shí)間會(huì)很快衰減[6]。

（2）邊界條件設(shè)定

邊界可分為2類(lèi)：一是陸邊界（閉邊界），是實(shí)際存在的，是水域與陸地或器壁的交界面；二是水邊界（開(kāi)邊界），是人為規(guī)定的，是截取的一部分水體所形成的有界計(jì)算域。與初始條件相比，邊界條件對(duì)數(shù)值計(jì)算的結(jié)果影響很大。邊界處理的2個(gè)基本要求是：使計(jì)算問(wèn)題在數(shù)學(xué)上適定，在物理上合理，盡量不影響內(nèi)點(diǎn)數(shù)值解的精度和穩(wěn)定性，內(nèi)點(diǎn)與邊界點(diǎn)的格式不一致和開(kāi)邊界對(duì)輸入波的虛假反射都是誤差擾動(dòng)源。對(duì)于陸邊界，一般使用無(wú)滑移條件來(lái)設(shè)定，即認(rèn)為水深在邊界的法向方向沒(méi)有變化，而水流速度在邊界的法向方向?qū)?shù)為零。

（3）動(dòng)邊界處理

動(dòng)邊界是水平計(jì)算域中有水與無(wú)水區(qū)域的界線(xiàn)。水陸邊界的外移是由于內(nèi)側(cè)水位高于外側(cè)地面，而內(nèi)縮則由于內(nèi)側(cè)水位低于同側(cè)地面。在邊界附近水深通常較小，同時(shí)邊界處存在法向流速，不同于一般的陸地邊界。

（4）參數(shù)的選擇

糙率是反映下墊面地表粗糙程度的參數(shù)，是水力模型中的重要參數(shù)。模型控制參數(shù)主要有演算時(shí)段控制參數(shù)和結(jié)果輸出控制參數(shù)。演進(jìn)時(shí)段控制參數(shù)包括起始計(jì)算時(shí)間、終止計(jì)算時(shí)間和計(jì)算時(shí)間步長(zhǎng)。結(jié)果輸出控制參數(shù)是結(jié)果文件輸出的時(shí)間步長(zhǎng)間隔數(shù)[7]。

3 黃墩湖滯洪區(qū)概況

黃墩湖滯洪區(qū)列入了國(guó)家蓄滯洪區(qū)名錄，是沂沭泗流域防洪體系重要組成部分。黃墩湖滯洪區(qū)位于駱馬湖西側(cè)，2009年國(guó)務(wù)院批復(fù)水利部《全國(guó)蓄滯洪區(qū)建設(shè)與管理規(guī)劃》明確，以邳洪河大堤（中運(yùn)河、駱馬湖二線(xiàn)大堤）房亭河南堤、徐洪河、廢黃河高灘地為界，地跨徐州、宿遷2市，分屬邳州、睢寧、宿豫3縣（市、區(qū)），總面積230 km2。蓄滯洪區(qū)內(nèi)有黃墩湖滯洪閘、爆破口門(mén)（勝利段、曹甸段）、防洪大堤等防洪建筑物。

按國(guó)家防汛防旱總指揮部《關(guān)于沂沭泗河洪水調(diào)度方案的批復(fù)》規(guī)定：“如預(yù)報(bào)駱馬湖水位超過(guò)26.0 m，當(dāng)駱馬湖水位達(dá)到25.5 m時(shí)，啟用黃墩湖滯洪區(qū)滯洪，確保宿遷大控制安全”。黃墩湖的滯洪方式為，開(kāi)啟進(jìn)洪閘和爆破口門(mén)相結(jié)合。當(dāng)駱馬湖水位超過(guò)24.5 m并預(yù)報(bào)繼續(xù)上漲時(shí)，退守宿遷大控制。同時(shí)，做好黃墩湖滯洪準(zhǔn)備，如預(yù)報(bào)駱馬湖水位超過(guò)26.0 m，當(dāng)水位達(dá)到25.5 m時(shí)，啟用黃墩湖滯洪。

4 應(yīng)用分析

根據(jù)黃墩湖滯洪區(qū)不同分洪組合得到邊界洪水過(guò)程，作為黃墩湖滯洪區(qū)二維水流數(shù)學(xué)模型的輸入邊界條件，根據(jù)黃墩湖滯洪閘、勝利口門(mén)、曹甸口門(mén)的不同開(kāi)啟組合，可制定出不同的洪水計(jì)算方案，采用模型進(jìn)行模擬計(jì)算得到方案結(jié)果。對(duì)傳播，由于蓄滯洪區(qū)南部地勢(shì)比較高，所以洪水主要集中在中部低洼地帶。上述分析表明，流場(chǎng)計(jì)算結(jié)果符合水流運(yùn)動(dòng)的基本規(guī)律。

（3）洪水特征要素分析

圖2 洪水特征要素圖

對(duì)模型計(jì)算結(jié)果進(jìn)行洪水演進(jìn)分析計(jì)算，統(tǒng)計(jì)得到最大淹沒(méi)范圍、最大淹沒(méi)水深、洪水到達(dá)時(shí)間、洪水淹沒(méi)歷時(shí)等洪水特征要素，如圖2所示。模型計(jì)算結(jié)果進(jìn)行洪水演進(jìn)分析計(jì)算，得到可能最大淹沒(méi)范圍、淹沒(méi)水深等洪水風(fēng)險(xiǎn)計(jì)算結(jié)果。選取黃墩湖滯洪閘啟用、勝利口門(mén)不啟用、曹甸口門(mén)不啟用的計(jì)算方案進(jìn)行洪水計(jì)算成果進(jìn)行合理性分析，包括水量守恒、局部流場(chǎng)、淹沒(méi)過(guò)程分析等。

（1）水量平衡分析

由于黃墩湖滯洪閘啟用后，駱馬湖洪水進(jìn)入蓄滯洪區(qū)，隨著閘門(mén)泄洪過(guò)程，蓄滯洪區(qū)內(nèi)的蓄水總量加上駱馬湖的當(dāng)前水量應(yīng)該等于滯洪閘啟用時(shí)駱馬湖初始水位下的蓄水庫(kù)容。

（2）局部流場(chǎng)與淹沒(méi)過(guò)程分析

選擇黃墩湖滯洪閘下游附近的局部流場(chǎng)進(jìn)行分析，結(jié)合蓄滯洪區(qū)的地形和典型時(shí)刻流場(chǎng)分布圖，由下面的洪水淹沒(méi)流場(chǎng)與水深分布圖可知，滯洪閘下游處流速較大，距離閘下游越近流速越大。駱馬湖洪水由滯洪閘進(jìn)入蓄滯洪區(qū)后，在水流慣性和地形地勢(shì)的雙重影響下向周?chē)屯莸貛U(kuò)散

5 結(jié)語(yǔ)

在黃墩湖滯洪區(qū)建立二維淺水流洪水演進(jìn)模擬模型，采用基于有限體積法的數(shù)值求解，準(zhǔn)確反映了蓄滯洪區(qū)洪水波運(yùn)動(dòng)傳播特性，有效捕捉了洪水間斷波，準(zhǔn)確描述洪水傳播特性和傳播時(shí)間。

針對(duì)黃墩湖蓄滯洪區(qū)分洪問(wèn)題采用滯洪閘、勝利口門(mén)、曹甸口門(mén)來(lái)模擬不同組合分流方式下洪水演進(jìn)過(guò)程，并對(duì)模型計(jì)算成果進(jìn)行了合理性分析，可知：

（1）基于FVM格式的二維水動(dòng)力學(xué)模型計(jì)算的河道水位、流速和斷面實(shí)測(cè)流態(tài)具有較好的吻合度，計(jì)算精度較高，適用于河道二維流態(tài)的全程模擬。

（2）整個(gè)洪水模擬過(guò)程，水量是守恒的，統(tǒng)計(jì)洪水淹沒(méi)歷時(shí)、淹沒(méi)范圍、淹沒(méi)水深等致災(zāi)參數(shù)均比較合理。

（3）此法為蓄滯洪區(qū)的洪水演進(jìn)模擬分析提供了一條行之有效的途徑。