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

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

（安徽理工大學(xué),安徽 淮南 232001）

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

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

1 數(shù)學(xué)模型

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

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

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

連續(xù)性方程

雷諾時均方程

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

2 模型驗證

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

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

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

3 數(shù)值模擬與結(jié)果分析

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

圖2 模型三維尺寸圖

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

圖3 不同時刻的自由液面

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

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

圖4 梯級水庫壩頂?shù)膲簭姾土魉僮兓?/p>

4 結(jié)論

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