基于MIKE21水動力模型計(jì)算感潮河段閘下水位流量關(guān)系

王新強(qiáng)

(1.福建省水利水電勘測設(shè)計(jì)研究院，福州 350001；2.福建省水動力與水工程重點(diǎn)實(shí)驗(yàn)室，福州 350001)

1 概述

閘下水位流量關(guān)系計(jì)算作為水利水電工程設(shè)計(jì)中的關(guān)鍵性基礎(chǔ)工作，具有非常重要的作用。在水利水電工程規(guī)劃設(shè)計(jì)中，遇到設(shè)計(jì)斷面附近無實(shí)測水位、流量或者實(shí)測資料不足的情況時(shí)，常常采用多種方法在水文調(diào)查和臨時(shí)測流的基礎(chǔ)上擬定水位流量關(guān)系[1-3]。目前，針對水位流量關(guān)系的研究較多，宋運(yùn)凱等采用正交多項(xiàng)式和冪函數(shù)法計(jì)算得到春河伊春站水位流量關(guān)系[4]；梅立庚利用水力學(xué)法，采用一維恒定非均勻漸變流計(jì)算得到斷面水位流量關(guān)系[5]；楊曉華等在研究自適應(yīng)加速遺傳算法求解的基礎(chǔ)上得到水位流量關(guān)系[6]。李玉榮等通過對比三峽水庫建庫前后主要水文站實(shí)測資料，研究了三峽水庫建庫前后荊江低水水位流量關(guān)系[7]；張英駿等用密集的數(shù)據(jù)替代數(shù)學(xué)公式的模式,由計(jì)算機(jī)對收集的盡量多和新的實(shí)測水位流量數(shù)據(jù)兩兩直線內(nèi)插后加權(quán)平均,直接得到推流流量結(jié)果,即無線推流方法[8]；程雪源等利用Excel擬合求解河流水位-流量關(guān)系曲線，提高了水流-流量關(guān)系曲線擬合的質(zhì)量以及精度、效率等[9]；王浩以感潮河段水位預(yù)報(bào)作為研究對象，從河道洪水演算水文學(xué)的方法出發(fā)，綜合考慮感潮河道洪水運(yùn)動的基本特征，河道水位流量轉(zhuǎn)換的基本方程等，研究出一種適用于感潮河段的水位預(yù)報(bào)方法[10]。

一般而言，水位流量關(guān)系指的是河道中某一橫斷面的流量與同時(shí)點(diǎn)水位之間的對應(yīng)關(guān)系，一般以水位為縱坐標(biāo)、流量為橫坐標(biāo)的水位流量關(guān)系曲線來表示，也可以選取適當(dāng)?shù)臄?shù)學(xué)方程式，并根據(jù)曲線和數(shù)學(xué)方程式列出便于查讀的水位流量關(guān)系表[12]，傳統(tǒng)水位流量關(guān)系擬合方法計(jì)算比較方便，但是對于水流條件較為復(fù)雜的河道計(jì)算精度不高，尤其是河口感潮段河流，同時(shí)受徑流和潮流共同作用，水流往復(fù)運(yùn)動，閘下水位流量關(guān)系不易準(zhǔn)確確定，本文根據(jù)木蘭溪實(shí)測地形以及水文條件，利用MIKE21 軟件建立木蘭溪下游二維水動力數(shù)學(xué)模型，能夠根據(jù)設(shè)計(jì)需要，模擬計(jì)算感潮段水閘下泄不同流量分別與高、低潮位遭遇的情況，得到閘下水位流量關(guān)系。

2 MIKE21 水動力模型

2.1 控制方程

MIKE21 水動力數(shù)學(xué)模型基于Navier-Stokes方程[11]。

具體控制方程如下：

(1)

(2)

(3)

式中:

t——時(shí)間；

x,y——笛卡爾坐標(biāo)系坐標(biāo)；

η——水位；

d——靜止水深；

h——總水深，h=η+d；

u,v——x,y方向上的速度分量；

f——Coriolis系數(shù)，f=2ωsinφ，φ表示所在地緯度，ω是指地球自轉(zhuǎn)的角速度；

g——重力加速度；

ρ——水的密度；

sxx、sxy、syy——輻射應(yīng)力的分量；

S——源匯項(xiàng)；

Txx、Txy、Tyx、Tyy——水平粘滯應(yīng)力；

us,vs——源匯項(xiàng)的不同方向的水流速度。

2.2 求解方法

本次水動力數(shù)值模擬計(jì)算所用模型在空間上進(jìn)行有限體積法離散，該方法將空間劃分為不重復(fù)的控制單元，能夠保證水動力模型中水量和動量守恒，模型在時(shí)間上的離散方法為顯性歐拉法[12]。 二維淺水方程積分形式一般可以寫成如下：

(4)

式中：

U——守恒變量的向量；

F——通量向量函數(shù)；

S——源項(xiàng)向量。

在Cartesian 坐標(biāo)中，可以將上面的二維淺水方程表示成如下形式：

(5)

式中：

I——無粘性通量；

V——粘性通量。

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

3.1 模型建立

為了能夠使模型計(jì)算結(jié)果更加合理、準(zhǔn)確，本次水動力數(shù)學(xué)模型采用非結(jié)構(gòu)網(wǎng)格對模型進(jìn)行劃分，模型的上邊界定于木蘭陂，并且考慮興化灣江陰島位置的潮位過程不受木蘭溪洪水流量的影響，將數(shù)學(xué)模型建模下邊界建至外海江陰島，基于2007年實(shí)測地形進(jìn)行驗(yàn)證計(jì)算，基于2010年1月實(shí)測地形進(jìn)行工程計(jì)算，河道邊界為防洪堤線(如圖1所示)。

圖1 木蘭陂—江陰島網(wǎng)格分布示意

3.2 模型相關(guān)參數(shù)

本次水動力數(shù)值模擬計(jì)算采用的參數(shù)主要有河道干濕邊界、河道糙率、渦粘系數(shù)以及時(shí)間步長等。

模擬中設(shè)定完全濕單元與完全干單元的水深閾值分別為0.1 m與0.005 m，即水深大于0.1的單元為完全濕單元，水深小于0.005 m的單元不參加計(jì)算，水深介于兩者之間的單元為半濕潤單元。河道糙率在數(shù)值模擬中具有十分重要的作用，本模型根據(jù)實(shí)測水流資料率定、調(diào)試模型糙率，最終確定模型計(jì)算河段的糙率取值區(qū)間在0.028～0.035之間，模型中渦粘系數(shù)取值0.28 m2/s，時(shí)間步長取為0.01～30 s。

3.3 模型驗(yàn)證

為了使模型驗(yàn)證更加準(zhǔn)確，分別采用2007年8月14—15日大潮及8月7—8日小潮兩次水文測驗(yàn)資料進(jìn)行驗(yàn)證計(jì)算，計(jì)算時(shí)下游出口邊界由三江口實(shí)測潮位過程控制，河道內(nèi)共有洋埕閘、寧海橋下和三江口3個水文測量斷面。

大潮以及小潮時(shí)，沿程三江口測點(diǎn)、寧海橋下測點(diǎn)以及洋埕閘測點(diǎn)潮位、流速過程計(jì)算值與實(shí)測值的對比見圖2～圖5所示。

圖2 大潮潮位計(jì)算值與實(shí)測值對比

圖3 大潮流速計(jì)算值與實(shí)測值對比

圖4 小潮潮位計(jì)算值與實(shí)測值對比

圖5 小潮流速計(jì)算值與實(shí)測值對比

由驗(yàn)證結(jié)果可知，3個測點(diǎn)的潮位與流速計(jì)算值與實(shí)測值吻合性較好，符合相關(guān)規(guī)程的要求，比較準(zhǔn)確地反映了木蘭溪下游河道的水流運(yùn)動規(guī)律，因此認(rèn)為模型參數(shù)的選取是合理的，可以將此模型用于計(jì)算分析木蘭溪寧海閘閘閘下水位流量關(guān)系。

4 寧海閘水位流量關(guān)系數(shù)值計(jì)算

木蘭溪位于福建省東部沿海莆田市境內(nèi)，下游為濱海南北洋平原(又稱興化平原)，木蘭溪流域面積為 1 732 km2，河長為105 km，系閩中最大河流，木蘭陂以下34.1 km河道為木蘭溪下游感潮河段，潮流作用明顯。擬建寧海閘水利工程處于木蘭溪下游河口段，距離寧海橋下約800 m(木蘭溪流域水系和擬建閘位置如圖6所示)。水閘樞紐工程由閘室、門庫、弧門支座、消能工、上下游接線護(hù)岸、控制房、管理房等建筑物組成，寧海閘閘室由河中央常規(guī)閘加兩側(cè)平面立軸式弧形鋼閘門組成，其中立軸式弧門共布置2孔，單孔凈寬為90.0 m，閘室底板頂高程為-0.50 m。其中常規(guī)閘共布置3孔，單孔凈寬為25.0 m，閘室底板頂高程為-2.50 m。寧海閘設(shè)計(jì)洪水標(biāo)準(zhǔn)為100年一遇(P=1%)，校核洪水標(biāo)準(zhǔn)為300年一遇(P=0.33%)，設(shè)計(jì)防潮標(biāo)準(zhǔn)為50年一遇高潮位(P=2%)。閘上正常蓄水位為2.80～4.30 m，主汛期閘上的限制蓄水位為2.80 m。

圖6 木蘭溪流域水系及擬建寧海閘位置示意

4.1 計(jì)算條件

本次寧海閘閘下水位流量關(guān)系數(shù)值模擬計(jì)算中，考慮水閘泄流分別遭遇低潮和高潮，通過利用MIKE21水動力模型仿真計(jì)算閘下水位流量關(guān)系，主要用于寧海閘閘下消能防沖的設(shè)計(jì)以及復(fù)核寧海閘的泄流能力。模型上邊界為水閘各級下泄流量，下邊界條件基于多年平均潮位過程進(jìn)行數(shù)值模擬計(jì)算。

4.2 閘下潮位合理性分析

擬建寧海閘位于驗(yàn)證點(diǎn)寧海橋下游側(cè)，通過驗(yàn)證結(jié)果可知實(shí)測值與計(jì)算值吻合良好(見圖2)，另外計(jì)算現(xiàn)狀條件下外海邊界為多年平均潮位過程時(shí)，東甲潮位站潮位過程線，并與東甲站實(shí)際潮位過程線對比(如圖7所示)，結(jié)果顯示計(jì)算結(jié)果良好，所以，利用數(shù)值模擬計(jì)算得到閘下的多年平均潮位過程線是可行的(見圖8)。由過程線可知，河道整體潮位分布形態(tài)正常，計(jì)算漲潮時(shí)東甲最高潮位為5.04 m，閘下最高潮位為5.03 m，仿真計(jì)算成果合理。

圖7 現(xiàn)狀條件下東甲站實(shí)際潮位與計(jì)算潮位過程

圖8 現(xiàn)狀條件下閘下多年平均潮位過程

4.3 閘下水位流量關(guān)系計(jì)算成果

本次利用二維水動力數(shù)學(xué)模型仿真計(jì)算水閘下泄各級流量分別與閘下低潮位、高潮位遭遇情況時(shí)寧海閘閘下水位流量關(guān)系，具體情況如下。

1) 泄流量遇低潮工況

利用MIKE21水動力模型仿真計(jì)算寧海閘下泄各級流量分別與低潮位遭遇工況下的閘下水位流量關(guān)系，并統(tǒng)計(jì)計(jì)算結(jié)果(見圖9、表1)。

表1 寧海閘泄流量遇低潮工況時(shí)，閘下水位流量關(guān)系比較

圖9 寧海閘泄流量遭遇低潮時(shí)的水位流量關(guān)系

通過對比分析本次仿真計(jì)算成果與水文復(fù)核計(jì)算成果，可以看出，在相同流量下，水動力仿真計(jì)算得到的水位值稍低，兩種方法計(jì)算得到水位流量關(guān)系整體趨勢一致，可據(jù)此進(jìn)行水閘消能防沖設(shè)計(jì)。

2) 泄流量遇高潮工況

利用MIKE21水動力數(shù)學(xué)模型計(jì)算水閘下泄各級流量分別匹配閘下多年平均潮位過程，計(jì)算得到寧海閘閘下水位流量關(guān)系。由于水文計(jì)算成果是考慮下游水位為平均高潮位2.19 m，故本次數(shù)模計(jì)算提取相同條件下的閘下水位與之對比。結(jié)果顯示，在流量較小時(shí)，相同流量情況下，水文公式計(jì)算得到的水位值略小，隨著水閘下泄流量的增加，相同流量下，水文公式計(jì)算的水位略高，兩者結(jié)果相差較小，趨勢較為一致，說明水動力仿真計(jì)算能夠較準(zhǔn)確的模擬下游潮位與徑流共同作用下閘下水位流量關(guān)系，因此，對于洪水工況的水閘泄流能力問題，可依據(jù)水動力數(shù)模仿真計(jì)算成果進(jìn)行復(fù)核(如圖10、表2)。

圖10 寧海閘泄流量遭遇高潮時(shí)的水位流量關(guān)系

表2 寧海閘泄流遇高潮工況時(shí)，閘下水位流量關(guān)系

5 結(jié)語

通過實(shí)測地形、水文資料建立MIKE21水動力數(shù)學(xué)模型，在驗(yàn)證合理的基礎(chǔ)上，利用該模型仿真計(jì)算寧海閘下泄不同量級流量分別遭遇低潮、高潮的工況，得到閘下水位流量關(guān)系。結(jié)果表明利用水動力數(shù)學(xué)模型仿真計(jì)算方法能夠較準(zhǔn)確的得到徑流和潮流共同作用下感潮河段閘下水位流量關(guān)系，該方法相對傳統(tǒng)水位流量關(guān)系計(jì)算方法更加仿真、準(zhǔn)確，并且簡明、直觀、適應(yīng)多種較復(fù)雜條件，計(jì)算結(jié)果可用于水閘消能防沖設(shè)計(jì)和復(fù)核水閘泄流能力。對于受潮流作用明顯的閘下水位流量關(guān)系計(jì)算具有一定的參考價(jià)值，并且該模型建立后也可以進(jìn)行河道水面線以及防洪影響分析計(jì)算等。