基于流固耦合的閘下淹沒出流流場瞬態分析

趙 勇，鄭圣義，蔡衛江，榮 紅，許 棟

(1.南瑞集團(國網電力科學研究院)有限公司，江蘇南京211106；2.河海大學能源與電氣學院，江蘇南京211100)

水工鋼閘門的安全關系到人民生命、財產安全，而水流誘發的閘門劇烈振動很容易造成閘門破壞[1-2]。閘門開啟泄流時，由于周圍流場運動因素不穩定，閘門上出現脈動荷載，從而引發閘門振動。據不完全統計[3]，我國失事的閘門中有20多起是因水流誘發振動，流激振動應予以高度重視。不同水流特性和結構特性情況下，閘門振動的性質大不一樣。特定條件下，工作閘門若是局部開啟，閘后將發生淹沒水躍[4-5]。淹沒水躍循環往復地沖擊門葉結構，旋滾渦流容易造成結構空蝕以及振動，惡劣情況下閘門會因動力失穩而破壞。通過對國內外文獻的查閱總結[6- 8]，應用有限元數值理論，構建流固耦合有限元模型，分析比較各工況下弧形閘門結構動特性，總結結構動特性在各工況下的變化規律，獲得對結構振動影響較大的開度區間。對于水工鋼閘門設計有一定的參考價值。

現階段國內對弧形閘門淹沒出流下瞬態流速變化、動水壓力變化和湍流動能變化等方面的研究存在著大量空白。本文首次基于ANSYS FLOTRAN多場耦合分析平臺，對水體—弧形閘門耦合模型進行了瞬態求解，研究了不同開度下流場瞬態流速、近壁水體動壓力和湍流動能隨時間變化變化規律。

1 水體—閘門耦合理論基礎與基本方程

水—閘門耦合問題屬于流固耦合問題流體動力作用下有限位移的長期問題，在泄流過程中，水體與閘門面板間有明確的耦合界面，水體位移比較有限，且泄流過程時間較長。水—閘門耦合(流固耦合)基本方程包括流體控制方程、耦合界面節點動壓力方程、等效節點載荷方程和結構控制方程[9]等。

(1)流體控制方程

(1)

(2)耦合界面節點動壓力方程

(2)

式中，Np=f+g+h，其中，f，g，h為三部分節點動壓力矩陣的行數；{Qp}為Sp面處的節點動壓力；{Qn}為耦合界面的節點動壓力；{Qu}為其他存在節點上的動壓力。推導得耦合界面Sn上的節點動壓力的基本方程

(3)

(3)等效節點載荷方程。耦合界面上的節點動壓力{Qn}可分解出耦合界面Snk上第k個單元S個節點處的動壓力，即

(4)

式中，{q}(k)為第k個單元s個節點處的動壓力；[A](k)為位號矩陣；SE為耦合界面上的流體單元個數。附加動壓力(結構坐標系下)的等效節點載荷為

(5)

(4)結構控制方程

(6)

2 模型建立

2.1 工程實例

本文以黃河某水電站胸墻式泄洪閘水體—弧形閘門為研究對象。針對該工作門開啟泄流情況，應用數值模擬軟件ANSYS，建立水體—弧門耦合有限元模型，研究淹沒出流下流場的瞬態特性。某胸墻式閘門為深孔閘門，正常水位下門后有一定深度的水體。水庫正常蓄水位134.00 m，下游平均水位125.30 m，高于溢流堰頂，水頭差8.70 m。閘門底檻高程120.192 m，閘門運行的水力特點是低水頭閘下淹沒出流，要求操作方式為動水啟閉。實際運行中，閘門常局部開啟泄流，閘門泄流時門底后出現淹沒水躍，振動問題較為嚴重。

2.2 有限元模型及邊界條件

計算選取上游段水體為孔口正前方水體，賦予其自由液面以下的相應水壓力，不影響計算結果。整個流體區域長20.6 m，寬8.904 m；門前水體長3.5 m，高4.155 m；門后水體長17.1 m，高11.420 m。工作閘門尺寸8.904 m×5.977 m，面板前沿距離支鉸中心9.0 m，支鉸中心高程126.0 m，門底高程120.192 m，閘門設計水頭13.808 m。坐標原點取在兩支鉸中心線的對稱中心，x軸為水流方向，y軸沿鉛直方向，z軸沿主橫梁軸向[10]。工程實例的流體為水，不考慮溫度，只需給出水的密度和粘度。鋼閘門為Q235鋼，彈性模量E=2.06×105MPa，泊松比μ=0.3，容重γ=78.5 kN/m3。弧形閘門是典型的空間薄壁復雜結構體系，本文首次采用實體單元模擬面板。水體以及弧形閘門各構件模擬單元的選擇為：流體(水體)為FLUID142單元，面板、主梁、縱梁、邊梁、小橫梁、支臂、筋板為等SHELL63單元；支臂斜撐桿為BEAM188單元；面板、圓柱鉸、支座與支鉸間的連接塊為SOLID187單元。網格采取不同區域不同網格密度，在流固耦合交界面附近區域的水體網格較密，遠離交界面區域的水體網格相對稀疏。水體—弧形閘門耦合有限元模型如圖1所示。

圖1 水體-弧形閘門耦合有限元模型(開度0.20)

流體邊界條件：假設底面和側面為無滑移固壁，即x、y、z向流動速度均為零；同時考慮計算簡便，上下游上體上表面設置為固壁邊界并施加滑移條件，確保有水流方向的速度。上游進口為壓力入口，出口為自由出流。

結構邊界條件：閘門在支鉸孔和吊耳孔處受x、y、z方向位移約束，x軸沿水平方向，y軸沿豎直方向，z軸沿支鉸中心連線方向。本文選擇RNGk-ε模型模流態[11]，相比于標準k-ε模型，RNGk-ε模型在ε方程中增添了一項用以反映主流時均應變率，從而對耗散率方程的模擬實現了優化，湍流的各向異性得以表達，非常適合用于模擬流線的彎曲流動。

3 水體—弧形閘門耦合流場瞬態分析

3.1 瞬態流速和動水壓力分析

圖2所示為0.20開度下1.00 s時刻對應的流場三維流速矢量。從圖2可以看出該開度下最大流速出現在閘門的底緣附近，過流斷面中心處流速相對左右側壁稍大一些。為便于觀察，選取各流場對稱面流速矢量，并將閘門底緣附近結果進行適當放大。部分開度下瞬態流速矢量如圖3、4、5、6所示。

圖2 0.2開度下1.00 s流場流速矢量(單位：m/s)

通過對各開度下流場瞬態流速和動水壓力的計算分析可知：各開度下水流在流經弧形閘門底緣時，水流的流態發生急劇變化，由于過流斷面的縮減，水流流速增加，且流經門底時流向發生偏轉，門底后方水域均出現了逆向回流。在計算時間內，各開度下門底附近水流的最大流速隨時間有不同幅度的增加，在0～0.16 s時間內，流速急劇增加，開度為0.05～0.20最為明顯，0.05時瞬態流速相對最高，呈現近乎射流狀態；0.05～0.20時流速變化劇烈，流態相對惡劣。

從圖3、4、5、6可以看出，近壁動水壓力大小分布不再滿足“梯形遞增”的形式，閘門底緣處水體壓力不是最大的，門體附近最大水壓力集中在距離閘門底緣約1/3處。閘門底緣前后斷面流體運動規律符合水力學伯努利能量方程，即流速高處壓力低，流速低處壓力高。在計算時間內，各開度下門前水體(除底緣附近水體)壓力隨時間非線性地遞增，開度越小門前水體壓力增加越多，且離上游耦合面下1/3處越近增加越明顯，最大增加約0.3 MPa，增幅為35.1%(開度0.05)。在計算時間內，各開度下門后水體壓力隨時間非線性地遞減，隨著開度的遞增，門后水體壓力減小的幅度逐漸變小。其中，開度為0.05時水體壓力最大降幅為112.4%，而開度為0.75時水壓最大降幅僅為24.8%。反映了開度越小，門后水體壓力下降越明顯。

圖3 0.05開度下瞬態流速矢量(單位：m/s)

圖4 0.05開度下動水壓力云圖(單位：Pa)

圖5 0.375開度下瞬態流速矢量(單位：m/s)

圖6 0.375開度下動水壓力云圖(單位：Pa)

3.2 閘門底緣近壁水體動壓力分析

為了解閘門底緣處水體壓力瞬態變化情況，在閘門底緣處選取4個近壁水體節點，生成節點動壓力時程變化曲線。門底前緣和門底后緣各取2個節點，底緣處近壁節點選取情況如圖7所示。圖8為各開度下4個節點動壓力時程變化曲線。

圖7 閘門底緣處近壁水體節點選取

從圖7的節點動壓力時程變化可以看出：①各開度下不同時刻，水體節點壓力不同。固定開度下，壓力大小關系為節點1>節點2>節點3>節點4，說明在門底前緣上方壓力高于下方壓力，在門底后緣上方壓力低于下方壓力，門底前緣壓力要高于后緣壓力。②對比各開度下各節點壓力峰值，發現其峰值均隨開度的增加而增加，從0.05開度到0.75開度，各節點的壓力峰值的增幅依次為70.8%、88.1%、117.1%和131.3%，說明閘下泄流時過流面積越小流速越高水壓越低，過流面積越大流速越小水壓越高，符合伯努利方程；同時反映出隨著開度的增加，閘門底緣近壁水體中，門底后方水體壓力變化最顯著。③開度0.05～0.20時，門底后緣節點3和節點4壓力先后出現負值，其大小在-6 kPa左右，說明該開度區間下門底產生負壓。當開度大于0.375時負壓不再出現。反映了小開度下淹沒出流造成的閘門底緣渦流空化現象，使得閘門底緣存在著下吸力。④在計算時間內，各節點壓力值變化方式不同。在0～0.2 s時間內，開度為0.05～0.20時，各節點壓力值出現“陡增陡減”變化；而0.375～0.75開度下，各節點壓力值出現“陡增緩減”的變化，反映出開度在0.05～0.20區間，作用在閘門底緣的脈動壓力變化較劇烈。

圖8 各開度下門底近壁節點壓力時程變化

4 結 論

通過對不同閘門開度下的流場流速和近壁水體動壓力的瞬態計算分析，得到以下結論：①開度在0.05～0.20區間，門底流速急劇增加流態相對劇烈；②閘門底緣處水壓力并非最大，門體附近最大水壓力集中在上游面板距離底緣約1/3處；③門后淹沒水躍對閘門產生較強的脈動壓力，開度在0.05～0.20區間，門底后緣水體節點壓力出現負值，表示閘門底緣存在著下吸力，其大小在-6 kPa左右；④開度在0.05～0.20區間下，各節點壓力值出現“陡增陡減”的變化，反映作用在閘門底緣的脈動壓力變化比較劇烈。