基于CFD模型的洪水條件下溢洪道水流特性模擬

陳興龍

(甘肅省水利水電勘測設(shè)計研究院有限責(zé)任公司，蘭州 730099)

0 引 言

自20世紀(jì)50年代初以來，國內(nèi)一直在使用物理模型研究溢洪道上的水流行為[1]。工程師可以使用一系列水力設(shè)計圖來設(shè)計任何給定設(shè)計洪水位的溢洪道剖面[2]。CFD技術(shù)在溢洪道水流分析中的應(yīng)用相對較新[3]。隨著計算機(jī)技術(shù)的進(jìn)步，更高效的CFD程序可以在三維中求解Navier-Stokes方程，并且也有許多湍流模型可供選擇。

在溢洪道設(shè)計中，剖面設(shè)計應(yīng)確保在最大洪水下水流過溢洪道結(jié)構(gòu)時，不會造成不利影響，如頂部和下游的氣穴。理想情況下，溢洪道表面應(yīng)剛好承受設(shè)計水頭下的大氣壓。當(dāng)水庫水位低于該洪水位時，溢洪道上方的壓力將高于大氣壓力。當(dāng)水庫水位遠(yuǎn)高于設(shè)計水頭時，沿溢洪道頂部將產(chǎn)生負(fù)壓，該壓力可能由于氣穴而損壞溢洪道的混凝土面，并對包括閘門結(jié)構(gòu)在內(nèi)的其他部件產(chǎn)生不利影響。

國內(nèi)的大多數(shù)大壩及其相關(guān)水工建筑物是在20世紀(jì)50年代和60年代早期設(shè)計的，以應(yīng)對當(dāng)時的洪水。此后，收集和處理了更可靠的長期水文數(shù)據(jù)。在許多情況下，修訂后的最大洪水等級有所增加。為了得出最佳補救設(shè)計，許多大壩需要考慮最具成本效益的方法，以調(diào)查最大洪水增加情況下溢洪道流量的行為。在過去的研究中，使用物理比例模型是唯一的調(diào)查方法。隨著計算機(jī)技術(shù)的發(fā)展，數(shù)值方法得以使用，并且在降低成本和減少分析時間方面更具有優(yōu)勢。

在本次調(diào)查中，將分析各種洪水條件下的標(biāo)準(zhǔn)WES溢洪道剖面，研究二維和三維模型。為了驗證模型結(jié)果，將計算結(jié)果與公布的數(shù)據(jù)進(jìn)行比較，結(jié)果顯示一致性良好，為今后研究中應(yīng)用CFD提供了技術(shù)保障。

1 數(shù)值計算方法

本研究使用的CFD程序為FLOW-3D，通過有限差分法求解Navier-Stokes方程。該算法是一種基于SOLA方法的擴(kuò)展方法，該方法由Hirt等在洛斯阿拉莫斯國家實驗室開發(fā)。流體體積(VOF)方法用于計算自由表面運動。

所有控制微分方程，如連續(xù)性和動量方程，均采用面積(2D)和體積(3D)孔隙度函數(shù)表示。該公式(分?jǐn)?shù)面積/體積障礙物表示法)用于模擬復(fù)雜的幾何區(qū)域。

任何復(fù)雜的障礙物幾何體都可以使用偏好技術(shù)來表示。每個單元(柵格)中障礙物占據(jù)的體積部分(或2D中的面積)在分析開始時定義。還計算了每個單元中的流體分?jǐn)?shù)，流體分?jǐn)?shù)的連續(xù)性方程、動量方程或輸運方程是使用偏好函數(shù)來表示的，使用有限差分近似對每個方程進(jìn)行離散。與某些有限元/體積或邊界擬合網(wǎng)格方法不同，這種網(wǎng)格技術(shù)不需要重新網(wǎng)格，并且在瞬態(tài)分析期間不會導(dǎo)致任何網(wǎng)格變形，因此可以很容易地應(yīng)用精確的求解算法。

以一個時間增量推進(jìn)解決方案的基本算法包括以下3個步驟：

步驟1：根據(jù)動量(Navier-Stokes)方程的顯式近似，使用所有平流、壓力和其他加速度的初始條件或以前的時間步長值計算速度。

步驟2：調(diào)整壓力以滿足連續(xù)性方程。

步驟3：更新流體自由表面或界面，根據(jù)流體體積給出新的流體配置。

2 二維溢洪道模型

2.1 模型網(wǎng)格及邊界條件

考慮無墩的WES溢洪道模型。之所以選擇這種方法進(jìn)行分析，是因為測量結(jié)果不受任何三維效應(yīng)的影響。因此，該模型接近于真實的二維流動問題。

溢洪道剖面的幾何形狀符合水力設(shè)計。它有一個垂直的上游面和一條曲線，由壩頂中心線前方的3個半徑(R=0.04Hd、R=0.20Hd和R=0.50Hd；Hd為設(shè)計水頭)定義。波峰中心線下游的剖面由以下方程式定義：

(1)

二維溢洪道頂部視圖見圖1。網(wǎng)格由X(水平)方向的95個單元和Z(垂直)方向的98個單元組成。縱橫比盡可能保持統(tǒng)一，特別是在敏感的區(qū)域，以提高求解精度和計算速度。對于該網(wǎng)格，X和Z方向的最大縱橫比分別為2.3和2.5。請注意，這里使用Z方向代替方程式(1)中定義的Y方向。

圖1 二維溢洪道頂部的網(wǎng)格特寫

本次調(diào)查的設(shè)計水頭取10 m。左邊界(上游)距離壩頂25 m，右邊界(下游)距離壩頂22 m。底部邊界位于壩頂下方18 m，頂部邊界位于壩頂上方14 m。假設(shè)以下邊界條件：

上游邊界：流速為零的靜水壓力；流體高度等于H。

下游邊界：流出邊界。

上游底部：無流動—被下方障礙物堵塞。

下游底部：流出邊界。

頂部邊界：對稱性—在這種情況下，由于重力，無影響。

設(shè)置初始條件，使水頭為H的流體體積位于溢洪道頂部。當(dāng)達(dá)到穩(wěn)定狀態(tài)時，在15 s的總時間內(nèi)進(jìn)行瞬態(tài)流分析，這是通過檢查系統(tǒng)的流速和動能等結(jié)果確定的。使用1 000 kg/m3的恒定水密度，這里假設(shè)水是不可壓縮的。在負(fù)Z方向施加9.81 m/s2的重力加速度。檢查3種不同的水頭工況(H/Hd=1.33、1.00和0.50)。

2.2 壓力分布

穩(wěn)態(tài)時沿波峰的表面壓力分布見圖2，這些壓力取自溢洪道附近的單元。圖2中還繪制了測量數(shù)據(jù)。可以觀察到，計算結(jié)果給出了略高的負(fù)壓，但總體趨勢和量級與測量數(shù)據(jù)非常一致；還可以看到一些壓力振蕩，它們可能歸因于局部網(wǎng)格效應(yīng)。

對于設(shè)計水頭工況(H/Hd=1.00)，水流沿溢洪道產(chǎn)生的壓力接近于零，即使在數(shù)值模擬中未引入曝氣。當(dāng)壓頭高于設(shè)計壓頭時，出現(xiàn)負(fù)壓；當(dāng)水頭低于設(shè)計水頭時，產(chǎn)生正壓。圖3顯示了穩(wěn)定狀態(tài)下溢洪道頂部H/Hd=1.33的壓力等值線，可以觀察到頂部上方的負(fù)壓區(qū)域。

圖2 不同水頭頂部壓力分布計算值與測量結(jié)果的比較

圖3 溢洪道頂部上方的負(fù)壓分布(單位：Pa)(H/Hd=1.33)

2.3 流速與流量分布

計算假設(shè)壁面完全光滑，且流動無黏性，湍流的影響將是未來研究的主題。根據(jù)試驗結(jié)果，尖頂堰/溢洪道上的流量可表示為：

Q=CLH1.5

(2)

式中：Q為流量，m3/s；C為流量系數(shù)；L為堰頂有效長度，m；H為壩頂上方的實測水頭，不包括流速水頭。

其中：流量系數(shù)近似為：

C=3.27+0.40(H/h)

(3)

式中：h為堰高。

穩(wěn)態(tài)下波峰上方的速度矢量見圖4，每種情況下的計算流量和平均速度(波峰中心線X方向的速度)見表1。為了便于比較，表1中還顯示了式(2)中的流量和相應(yīng)的平均流速。可以看出，計算值比經(jīng)驗結(jié)果高估了約10%～20%，這可能與分析中使用的無黏流條件有關(guān)。

圖4 溢洪道壩頂上方的速度矢量(單位：m/s)

表1 流量和平均流速的比較

3 三維溢洪道模型

現(xiàn)在考慮帶橋墩的WES溢洪道模型。橋墩的存在將產(chǎn)生三維效應(yīng)，此分析使用了相當(dāng)粗糙的網(wǎng)格。除了模型中包含橋墩和橋臺外，幾何結(jié)構(gòu)與二維情況類似。對稱條件應(yīng)用于兩個z-x邊界平面。上游z-y邊界面上保持恒定的靜水壓頭，允許水通過下游的z-y和x-y邊界平面流出。穩(wěn)態(tài)下的模型視圖見圖5。

對結(jié)果的初步檢驗表明，分析結(jié)果與公布的數(shù)據(jù)有較好的一致性，橋墩附近的負(fù)壓高于其余部分的負(fù)壓。

4 溢洪道結(jié)構(gòu)的應(yīng)用

鑒于計算結(jié)果與二維分析中公布的數(shù)據(jù)之間的良好一致性，使用相同的技術(shù)對混凝土重力壩中央溢洪道在洪水位增加的情況下進(jìn)行分析。首先建立一個物理模型，對溢洪道在不同洪水位下的特性進(jìn)行研究；繪制物理模型試驗的測量值，以進(jìn)行比較，沿溢洪道中心線計算的壩頂壓力分布見圖6。

圖5 顯示穩(wěn)定狀態(tài)下流過溢洪道的含水量和流速分布

圖6 不同洪水位的峰值壓力分布CFD分析和物理模型試驗結(jié)果的比較

由圖6可以觀察到，獲得了相對較好的一致性。

5 結(jié) 論

本文基于CFD模型，從二維和三維兩個方面對溢洪道洪水特性進(jìn)行了模擬分析，并將計算結(jié)果與模型試驗結(jié)果進(jìn)行了驗證，取得了良好的一致性。但計算結(jié)果高估了流速，低估了溢洪道沿線的壓力分布。