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

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

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

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

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

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

式中：α 為體積分?jǐn)?shù),并滿(mǎn)足對(duì)流輸運(yùn)方程：

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

連續(xù)性方程

雷諾時(shí)均方程

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

2 模型驗(yàn)證

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

圖1 不同時(shí)刻的潰壩波自由液面

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

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

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

圖2 模型三維尺寸圖

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

圖3 不同時(shí)刻的自由液面

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

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

圖4 梯級(jí)水庫(kù)壩頂?shù)膲簭?qiáng)和流速變化

4 結(jié)論

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