陳興龍
(甘肅省水利水電勘測設(shè)計研究院有限責(zé)任公司,蘭州 730099)
自20世紀(jì)50年代初以來,國內(nèi)一直在使用物理模型研究溢洪道上的水流行為[1]。工程師可以使用一系列水力設(shè)計圖來設(shè)計任何給定設(shè)計洪水位的溢洪道剖面[2]。CFD技術(shù)在溢洪道水流分析中的應(yīng)用相對較新[3]。隨著計算機技術(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)查方法。隨著計算機技術(shù)的發(fā)展,數(shù)值方法得以使用,并且在降低成本和減少分析時間方面更具有優(yōu)勢。
在本次調(diào)查中,將分析各種洪水條件下的標(biāo)準(zhǔn)WES溢洪道剖面,研究二維和三維模型。為了驗證模型結(jié)果,將計算結(jié)果與公布的數(shù)據(jù)進(jìn)行比較,結(jié)果顯示一致性良好,為今后研究中應(yīng)用CFD提供了技術(shù)保障。
本研究使用的CFD程序為FLOW-3D,通過有限差分法求解Navier-Stokes方程。該算法是一種基于SOLA方法的擴展方法,該方法由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ù)流體體積給出新的流體配置。
考慮無墩的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個單元組成??v橫比盡可能保持統(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)。
穩(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)
計算假設(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 流量和平均流速的比較
現(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ù)壓。
鑒于計算結(jié)果與二維分析中公布的數(shù)據(jù)之間的良好一致性,使用相同的技術(shù)對混凝土重力壩中央溢洪道在洪水位增加的情況下進(jìn)行分析。首先建立一個物理模型,對溢洪道在不同洪水位下的特性進(jìn)行研究;繪制物理模型試驗的測量值,以進(jìn)行比較,沿溢洪道中心線計算的壩頂壓力分布見圖6。
圖5 顯示穩(wěn)定狀態(tài)下流過溢洪道的含水量和流速分布
圖6 不同洪水位的峰值壓力分布CFD分析和物理模型試驗結(jié)果的比較
由圖6可以觀察到,獲得了相對較好的一致性。
本文基于CFD模型,從二維和三維兩個方面對溢洪道洪水特性進(jìn)行了模擬分析,并將計算結(jié)果與模型試驗結(jié)果進(jìn)行了驗證,取得了良好的一致性。但計算結(jié)果高估了流速,低估了溢洪道沿線的壓力分布。