河道下游存在梯級水庫的潰壩洪水演進模擬研究

賀成山,陳情情,劉 波,吳 昊,郭魯楠

（安徽理工大學,安徽 淮南 232001）

隨著社會的不斷發展,水利設施不斷完善,水庫大壩成為當今水利工程中不可或缺的一部分。我國每年因水庫大壩潰決導致的財產損失與人員傷亡時有發生,興修水庫大壩以及研究潰壩洪水演進規律迫在眉睫。曹洪建等人采用naoe-FOAM-SJTU 求解器模擬潰壩洪水對下游直立方柱沖擊過程中潰壩波的劇烈變形,為預測波浪對結構物的沖擊提供科學依據[1]。牟迪等人對比分析了DNS 和RAS 兩種模型在潰壩波對下游結構物作用的不同結果,二者差別較小[2]。丁偉業等人采用VOF 法與NLk-ε 紊流模型耦合的方式對下游存在不同幾何形狀障礙物和90°彎角的潰壩進行了模擬,得出河床地形變化、潰口發展等多種潰壩水流的影響因素[3]。Pankaj Kumara 等人將SPH 和FVM方法在OpenFOAM中進行了耦合,實現了一種基于相位變化和相分數α 變化的界面檢測算法以及任意時間步長的SPH 區域和FVM區域的識別[4]。羅富強等人采用OpenFOAM 中的多相流求解器interFoam對二維潰壩問題進行了自適應網格策略研究,并通過三維潰壩問題的驗證得出三維自適應網格策略的有效性[5]。Y. Wang 等人在OpenFOAM中建立了一個預測溢洪道射流狀態的模型,采用流體體積法捕捉動態自由曲面,得出湍流對溢洪道下游的流動特性起著重要的作用[6]。Andrea Colagrossi 等人提出了一種處理二維界面流的光滑粒子流體力學(SPH)方法。用此方法討論了密度比變化對潰壩問題的影響,討論了空氣夾帶對負荷的影響。結果表明,該方法穩定、易于處理具有界面斷裂和空氣夾帶的各種空氣- 水流動[7]。

鑒于物理模型實驗的費時費力,本研究采用數值模擬的方法對下游存在梯級水庫時上游水庫潰壩洪水演進過程進行研究。利用CFD 開源軟件OpenFOAM建立數值模型進行模擬,通過ParaView軟件進行數據后處理分析。

1 數學模型

由于潰壩波是在重力作用下發生運動,氣液兩相流的分層尤為明顯,利用流體體積法（VOF）能夠較好的模擬潰壩水流的自由液面,完成對自由液面的捕捉。VOF 涉及到的流體體積函數為：

式中：α 為體積分數,并滿足對流輸運方程：

本研究使用氣液兩相流模型進行數值模擬,假設兩種流體均為不可壓縮流體,控制方程分別如下。

連續性方程

雷諾時均方程

式中：μeff分別表示流體密度及流體有效粘度；u 表示流體速度。

2 模型驗證

為了驗證模型的有效性,本研究采用荷蘭海洋研究所(MARIN)構建的物理模型進行驗證,該模型是尺寸為3.22 m×1 m×1 m 的長方體水箱,該模型頂部為開放式,四周與底部均為墻體。模型右側部分設有0.55 m 高的長方體水柱,水柱代表上游水庫,通過水柱坍塌模擬大壩潰決,模型下游設有箱體障礙物,箱體障礙物兩側留有大小等同的空隙。圖1 為數值模擬（左）和物理模型試驗（右）的結果對比。

圖1 不同時刻的潰壩波自由液面

通過對比分析可知,模擬結果與試驗結果基本一致。在物理模型試驗和數值模擬中都可以看到水柱坍塌后自由液面的形狀,水流在0.5 s 撞擊箱型障礙物的瞬間,潰壩波形狀是相同的,撞擊形成的水舌開始上升,其余水流從箱體障礙物兩側空隙流入后方；在0.84 s 時一部分水流從箱型障礙物兩側流出后撞擊在模型左側墻體上形成反射波,此時水流撞擊箱型障礙物的水舌向前彎曲,撞擊墻體產生的反射波沿著左側墻體爬升,反射波與障礙物前的潰壩波在1.02 s 于障礙物左側空中相遇；反射波與上游的潰壩波碰撞后摻雜著部分空氣在1.7 s 共同向上游行進,模擬結果與試驗結果保持一致。

3 數值模擬與結果分析

本研究在三維潰壩的基礎上研究下游存在梯級水庫時上游潰壩水流的洪水演進規律,為此利用OpenFOAM建立了數值模型,模型尺寸圖見圖2。該模型上游設有0.9 m 高的水庫,模型底部設有一定坡度,假設在距離模型上游墻體1 m 處設置閘門,通過閘門迅速開啟實現瞬時潰壩。下游設有梯級水庫,梯級水庫初始時刻設有一定水深,兩座大壩尺寸一致。

圖2 模型三維尺寸圖

該模型是利用OpenFOAM 中自帶的blockMesh進行網格點的布置以及網格劃分,模型共計704 000個網格,計算時間步長為0.001 s,每隔50 個時間步長輸出一次結果,計算總時間為3 s。根據setFields 定義非均勻初始場,選用多相流中的interFoam 求解器對潰壩水流在下游庫區的流動過程進行模擬,利用VOF 法捕捉控制方程,interFoam 能夠很好的模擬潰壩水流中氣液兩相流的運動情況。圖3 為不同時刻潰壩水流在下游梯級水庫的洪水演進過程。

圖3 不同時刻的自由液面

結合不同時刻的自由液面分析上游水庫潰壩洪水在下游梯級水庫中的演進過程。t=0 s 時下游兩個水庫均設有一定初始水位,模型保持靜止狀態；上游水庫潰壩后水流進入下游第一個水庫并與水庫中原有水體發生碰撞從而形成壅水,直到水位高于壩頂位置后水流開始溢出,t=0.5 s 時下游第一個壩發生漫頂現象,漫頂過程中壩前水位達到最高；t=0.75 s 時由于第一個壩漫頂水流與第二個水庫中原有水體發生碰撞形成壅水現象,第二個壩的壩前水深達到最大值；t=0.9 s 時第二個壩也出現漫頂現象；t=1.1 s 時水流漫過第二個壩頂向下游行進；t=1.7 s 時由于水體與模型右側墻壁碰撞形成反射波與上游來水混合,能量不斷消耗,最終模型整體趨于平衡。

根據圖4 中梯級水庫大壩上監測點得到的數據可知,在0.5 s 時,上游水庫大壩瞬時全潰后水流到達第一個壩前并形成漫頂,壩前水深達到最大,壩頂處壓強和流速急劇上升達到最大值；在0.75 s 時,潰壩洪水漫過第一個壩流入第二個水庫,水庫水位持續上升,形成壅水,此時第二個壩有漫頂趨勢,而壩頂處的壓強和流速逐漸增大,由于第一個壩體的攔截作用,上游洪水的水量減小,第二個壩體壩頂處的水位、流速和壓強的變化趨勢減緩,此時第一個壩前水流由于壩體的攔截作用而形成反射波向上游行進,上游的潰壩波與反射波碰撞,經過一段時間后共同向下游行進,反復幾次后第一個水庫中的水體歸于平靜,此過程中第一個水庫壩頂處的壓強和流速不斷發生變化；在0.9 s 時,第二個壩出現漫頂現象,上游不斷來流的情況下,水體發生擾動,壩頂處的水位和壓強達到最大值；在1.1 s 時,部分水體漫過第二個壩流向下游河道,由于上游來水與壩體漫頂的原因水庫中的水不斷混合碰撞,水體與大量空氣摻雜在一起后形成反射波向上游行進,到達一定高度后動能轉化為勢能,在勢能的作用下再次向下游行進；在1.7 s 時,模型下游河道中的水體撞擊右側墻體后反射回上游,隨著動能和勢能的不斷轉換與消耗,最終整個模型中的水體流速減緩達到平衡。

圖4 梯級水庫壩頂的壓強和流速變化

4 結論

本研究利用CFD 開源軟件OpenFOAM中的流體體積法(VOF)以及標準k-ε 湍流模型模擬下游存在梯級水庫時的三維潰壩問題。采用有限體積法對控制方程離散,利用壓力隱式算子分割法（PISO）對方程組進行求解。根據荷蘭海洋研究所建立的物理模型試驗來驗證模型的有效性和可靠性。研究結果表明,該模型能夠較好的模擬下游存在梯級水庫時,潰壩洪水的流場、水面線等演進過程,為探究潰壩洪水在下游庫區演進問題提供了科學依據。同時,該模型也存在一定問題,實際水庫大壩潰決時并不是瞬時全潰,河床形態、河道走向、梯級水庫之間的距離以及不同的河道坡度都會影響潰壩洪水的演進過程。