基于MIKE21水動力模型計算感潮河段閘下水位流量關系

王新強

(1.福建省水利水電勘測設計研究院，福州 350001；2.福建省水動力與水工程重點實驗室，福州 350001)

1 概述

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

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

2 MIKE21 水動力模型

2.1 控制方程

MIKE21 水動力數學模型基于Navier-Stokes方程[11]。

具體控制方程如下：

(1)

(2)

(3)

式中:

t——時間；

x,y——笛卡爾坐標系坐標；

η——水位；

d——靜止水深；

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

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

f——Coriolis系數，f=2ωsinφ，φ表示所在地緯度，ω是指地球自轉的角速度；

g——重力加速度；

ρ——水的密度；

sxx、sxy、syy——輻射應力的分量；

S——源匯項；

Txx、Txy、Tyx、Tyy——水平粘滯應力；

us,vs——源匯項的不同方向的水流速度。

2.2 求解方法

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

(4)

式中：

U——守恒變量的向量；

F——通量向量函數；

S——源項向量。

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

(5)

式中：

I——無粘性通量；

V——粘性通量。

3 模型建立與驗證

3.1 模型建立

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

圖1 木蘭陂—江陰島網格分布示意

3.2 模型相關參數

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

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

3.3 模型驗證

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

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

圖2 大潮潮位計算值與實測值對比

圖3 大潮流速計算值與實測值對比

圖4 小潮潮位計算值與實測值對比

圖5 小潮流速計算值與實測值對比

由驗證結果可知，3個測點的潮位與流速計算值與實測值吻合性較好，符合相關規程的要求，比較準確地反映了木蘭溪下游河道的水流運動規律，因此認為模型參數的選取是合理的，可以將此模型用于計算分析木蘭溪寧海閘閘閘下水位流量關系。

4 寧海閘水位流量關系數值計算

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

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

4.1 計算條件

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

4.2 閘下潮位合理性分析

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

圖7 現狀條件下東甲站實際潮位與計算潮位過程

圖8 現狀條件下閘下多年平均潮位過程

4.3 閘下水位流量關系計算成果

本次利用二維水動力數學模型仿真計算水閘下泄各級流量分別與閘下低潮位、高潮位遭遇情況時寧海閘閘下水位流量關系，具體情況如下。

1) 泄流量遇低潮工況

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

表1 寧海閘泄流量遇低潮工況時，閘下水位流量關系比較

圖9 寧海閘泄流量遭遇低潮時的水位流量關系

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

2) 泄流量遇高潮工況

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

圖10 寧海閘泄流量遭遇高潮時的水位流量關系

表2 寧海閘泄流遇高潮工況時，閘下水位流量關系

5 結語

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