王海軍， 張 凡， 李會(huì)平

(1.天津大學(xué) 水利工程仿真與安全國(guó)家重點(diǎn)實(shí)驗(yàn)室， 天津 300350；2.天津大學(xué) 建筑工程學(xué)院， 天津 300350)

1 研究背景

隨著我國(guó)高壩建設(shè)的快速發(fā)展，泄水建筑物水頭越來(lái)越高。溢洪道是水利樞紐中重要的泄洪建筑物。隨著下泄流量增大，對(duì)溢洪道的過(guò)流能力及結(jié)構(gòu)的安全性要求越來(lái)越高。溢洪道的過(guò)流能力及其壁面壓強(qiáng)、臨底流速的分布情況等因素都對(duì)溢洪道泄水、消能有重要的影響[1]。臨底流速是底流消能工的一項(xiàng)重要水力學(xué)指標(biāo)，同時(shí)也是影響結(jié)構(gòu)安全的一項(xiàng)重要因素。傳統(tǒng)的溢洪道消力池水力特性分析主要是采用物理模型試驗(yàn)法[2]以及原型觀測(cè)法，尤其是前者，在已做的工程設(shè)計(jì)和水力靜動(dòng)力特性研究中起到了極為重要的作用。楊敏等[1]利用物理模型試驗(yàn)的方法對(duì)跌坎消力池的臨底流速進(jìn)行了探究，分析比較了跌坎消力池與傳統(tǒng)消力池的臨底流速分布規(guī)律；王立杰等[3]采用物理模型試驗(yàn)的方法針對(duì)新型射流消能的流態(tài)演變與極限臨底流速進(jìn)行了研究，同時(shí)也發(fā)現(xiàn)，在采用物理模型試驗(yàn)時(shí)，模型建立較為復(fù)雜、成本較高、測(cè)量精準(zhǔn)性難以保證，且存在比較大的簡(jiǎn)化，比尺效應(yīng)影響較大。隨著近年來(lái)計(jì)算機(jī)技術(shù)的發(fā)展及數(shù)值計(jì)算方法的完善，采用數(shù)值模擬方法分析研究復(fù)雜水流也成為經(jīng)濟(jì)有效的手段?；谟?jì)算流體動(dòng)力學(xué)的大型商業(yè)軟件 Fluent、Flow-3D等的開發(fā)利用，水利領(lǐng)域的數(shù)值仿真模擬也得到快速發(fā)展[4]。李玲等[5]研究了三維VOF模型及其在溢洪道水流計(jì)算中的應(yīng)用；薛宏程等[6]對(duì)溢洪道出口斜切型挑坎挑射水舌進(jìn)行了三維數(shù)值模擬；王青等[7]基于Flow-3D對(duì)陡坡彎道水流進(jìn)行了三維數(shù)值模擬；李樹寧等[8]對(duì)不同體型參數(shù)的跌坎型消力池進(jìn)行了數(shù)值模擬計(jì)算，為工程設(shè)計(jì)提供依據(jù)等。

本文采用 Fluent軟件，對(duì)某水電站溢洪道泄洪進(jìn)行精細(xì)數(shù)值模擬，結(jié)合該水電站原型觀測(cè)數(shù)據(jù)，驗(yàn)證所采取數(shù)值模擬方法的合理性。并在此基礎(chǔ)上探尋溢洪道水流臨底流速的分布規(guī)律及臨底流速梯度對(duì)溢洪道結(jié)構(gòu)安全的影響。

2 臨底流速模擬方法

臨底流速的數(shù)值模擬主要包括基本的流體控制方程、自由液面VOF模型、模型求解方法。

2.1 控制方程

將連續(xù)性方程和不可壓縮黏性流體運(yùn)動(dòng)的Navier-Stokes方程作為流體運(yùn)動(dòng)的控制方程，紊流模型采用的是RNGk-ε模型[9]。

(1) 連續(xù)性方程

(1)

式中：ui為i方向的速度分量，m/s；xi為坐標(biāo)方向；下標(biāo)i為變量的不同方向，取值1、2、3。

(2) 運(yùn)動(dòng)方程

(2)

式中：fi為作用于單位質(zhì)量水體的質(zhì)量力，m/s2；ρ為流體密度，kg/m3；p為壓強(qiáng)，N/m2；t為時(shí)間，s；uj為j方向的速度分量，m/s；xi、xj為坐標(biāo)方向(下標(biāo)i、j為變量的不同方向，取值1、2、3)。

(3)k方程

(3)

式中：Gk為剪切產(chǎn)生項(xiàng)，m2/s3；Gb為浮力產(chǎn)生項(xiàng)，m2/s3；k為單位質(zhì)量紊動(dòng)動(dòng)能，m2/s2；ε為紊動(dòng)動(dòng)能耗散率，m2/s3；ν為運(yùn)動(dòng)黏性系數(shù)，m/s2；νt為紊流運(yùn)動(dòng)黏性系數(shù)，m/s2。

(4)ε方程

(4)

式中：系數(shù)C1ε、C2ε、C3ε為常量；σk、σε為k方程和ε方程的Prandtl數(shù)；Rε為附加項(xiàng)，m2/s4，以適應(yīng)瞬變流和流線彎曲的影響[10]。

2.2 VOF模型

VOF模型(volume of fluid)基本原理是通過(guò)研究網(wǎng)格單元中流體與網(wǎng)格體積比函數(shù)來(lái)追蹤自由液面。相比于Mixture模型與歐拉模型更加適用于計(jì)算任意氣液分界面的分界面。VOF模型中不同的流體共用一套運(yùn)動(dòng)方程，通過(guò)引進(jìn)相體積分?jǐn)?shù)這一變量，實(shí)現(xiàn)對(duì)每一個(gè)計(jì)算單元相界面的追蹤。在每個(gè)控制容積內(nèi)，所有相體積分?jǐn)?shù)額總和為1，所有變量及其屬性在所控制容積內(nèi)各相共享，且代表了容積的平均值。在任何給定控制容積內(nèi)的變量和其屬性代表了一相或多相的混合。假設(shè)單元內(nèi)第m相流體的體積分?jǐn)?shù)為αm，則存在3種情況：

(1)αm=0，單元體積內(nèi)不存在第m相流體；

(2)0<αm<1，單元體積內(nèi)包含了第m相流體和其他一相或多相流體，且存在著不同相流體之間的交界面；

(3)αm=1，單元體積內(nèi)充滿第m相的流體。單元體積內(nèi)具有相同的速度和壓力場(chǎng)。

在VOF模型中，追蹤第m相流體自由水面的控制微分方程為：

(5)

通過(guò)對(duì)該方程的求解可以完成自由水面的追蹤[11-14]。

2.3 模型求解方法

流場(chǎng)計(jì)算采用的是有限體積法，其是一種分塊近似的計(jì)算方法。具體的將計(jì)算區(qū)域劃分為一系列不重復(fù)的控制體積，并使每個(gè)網(wǎng)格點(diǎn)周圍有一個(gè)控制體積，這樣計(jì)算區(qū)域被離散成許多小的體積單元，然后將待解的微分方程對(duì)每一個(gè)小單元體進(jìn)行積分，便得出一組離散方程解。其求解方法分為分離解法和耦合解法兩大類。具體的擴(kuò)展分為SIMPLE算法以及在其基礎(chǔ)上改進(jìn)的SIMPLEC、SIMPLER、PISO算法。由于PISO算法對(duì)瞬態(tài)問(wèn)題有明顯的優(yōu)勢(shì)[15]，故本模擬采用PISO算法。

3 臨底流速仿真模擬

3.1 工程簡(jiǎn)介

某工程采用了溢流表孔底流消能的泄洪消能方式。泄洪建筑物的泄水孔為開敞式表孔，總共分為5孔，單孔尺寸為13 m×20 m。邊墩、中墩厚度分別為4.0和5.0 m，溢流壩段總寬為93.0 m；孔口設(shè)有檢修平板門和弧形工作門，閘墩前緣呈尖圓流線形，溢流壩堰頂高程為1 398.0 m。堰面為WES曲線，曲線末端接1∶0.75直線段和半徑R=55 m的反弧段；溢流壩后為長(zhǎng)約200.0 m的泄槽和170.0 m的消力池，入池跌坎高6.0 m，池底高程1 276.5 m，消力池后設(shè)海漫與河床平順連接。溢洪道縱剖面示意圖如圖1所示。

圖1 溢洪道剖面示意圖

3.2 模擬步驟

模型計(jì)算求解步驟為: (1)在Auto CAD (前處理器) 中建立幾何模型。(2)將模型導(dǎo)入ICEM中，進(jìn)行網(wǎng)格的劃分。(3)將劃分完成的網(wǎng)格模型導(dǎo)入Fluent中，并設(shè)置相關(guān)參數(shù)。參數(shù)設(shè)置包括計(jì)算模型選取、材料屬性和邊界條件設(shè)定以及求解方法的選取和運(yùn)行環(huán)境的設(shè)置等。(4)對(duì)流場(chǎng)的初始化、求解。(5)進(jìn)行計(jì)算結(jié)果的后處理及結(jié)果的分析。

3.3 ICEM精細(xì)網(wǎng)格劃分

網(wǎng)格的劃分是流場(chǎng)數(shù)值模擬的重要環(huán)節(jié)。網(wǎng)格質(zhì)量的高低直接影響計(jì)算精度及收斂速度。溢洪道泄流的模擬范圍大、流道體型復(fù)雜，本文研究分析的

重點(diǎn)是臨底流速分布，因此在網(wǎng)格劃分時(shí)，近底處進(jìn)行了細(xì)化加密，上部網(wǎng)格適當(dāng)放大。計(jì)算模型約525×104個(gè)節(jié)點(diǎn)，約498×104個(gè)網(wǎng)格。

3.4 邊界條件

模型自由表面的處理采用VOF模型，邊界條件采用壓力進(jìn)出口邊界條件。模型進(jìn)流邊界為壓力進(jìn)口，給定進(jìn)水口水頭；出流邊界為壓力出口，給定出口水頭；固壁面采用無(wú)滑移壁面，法向和切向速度均為0。

4 模型驗(yàn)證

利用建立的溢洪道三維數(shù)學(xué)模型，對(duì)已進(jìn)行過(guò)原型觀測(cè)的工況進(jìn)行模擬計(jì)算。采用3#、4#孔全開工況(表1中工況1)的模擬結(jié)果進(jìn)行驗(yàn)證，此工況堰上水頭為117.15 m，下游水頭13.50 m，溢洪道流量3 336 m3/s。比較數(shù)值模擬計(jì)算值和原型觀測(cè)值，對(duì)所建立的數(shù)學(xué)模型進(jìn)行驗(yàn)證。

4.1 水流流態(tài)驗(yàn)證

溢洪道縱剖面水氣相分布如圖2所示。從圖2可以看出水流出閘室后，沿溢流面下泄，水流在溢流堰曲線段水流平順下泄，出閘室后的水流橫向擴(kuò)散。水流沿溢流堰面急速下泄至消力池內(nèi)，水流在消力池內(nèi)劇烈摻混，充分消能，出消力池后水流趨于平順泄向下游，消力池消能效果較好。整個(gè)過(guò)程水流在溢流堰反弧段與泄槽陡坡段水深均較小。模擬水流流態(tài)與原型觀測(cè)流態(tài)分布基本一致，原型觀測(cè)如圖3所示。

4.2 水面線驗(yàn)證

提取自由液面水深和實(shí)測(cè)水深對(duì)比，圖4為3#孔中心線處溢洪道水深模擬計(jì)算值與原型實(shí)測(cè)值比較圖(橫坐標(biāo)為分析位置距離壩趾的距離)，二者吻合度較好。溢流堰陡坡處水深下降迅速，至消力池后水深迅速增加，出消力池后水深開始減小。

通過(guò)流態(tài)和水面線驗(yàn)證，數(shù)值仿真模擬與原型觀測(cè)結(jié)果基本吻合，可見(jiàn)數(shù)值仿真模擬結(jié)果可作為原型觀測(cè)數(shù)據(jù)的補(bǔ)充。

圖2溢洪道縱剖面水氣相分布圖3原型觀測(cè)水流流態(tài)

圖4水深計(jì)算值與觀測(cè)值對(duì)比圖圖5縱向最大臨底流速變化趨勢(shì)圖

5 臨底流速分布規(guī)律

分別選取原型觀測(cè)中出現(xiàn)的最大流量、最小流量以及中間流量3種典型工況進(jìn)行數(shù)值模擬，從而分析臨底流速分布規(guī)律。具體工況如表1所示，分析3種模擬工況，提取3#孔中軸線的縱向最大臨底流速分布，如圖5所示。從數(shù)值仿真結(jié)果可知，縱向不同斷面臨底流速最大值順?biāo)鞣较虺尸F(xiàn)先增大后減小的趨勢(shì)。同時(shí)觀察到在1-1斷面(距離壩址102 m，見(jiàn)圖1)和2-2斷面(距離壩址200 m，見(jiàn)圖1)的臨底流速最大值在工況3時(shí)分別達(dá)到了42.24和44.95 m/s，如表2所示。

圖6～7為3種工況下，3#孔中軸線在1-1斷面與2-2斷面處距離底板不同位置臨底流速變化曲線圖，可以發(fā)現(xiàn)在這兩個(gè)斷面處近底側(cè)流速梯度均較大，最大值出現(xiàn)在工況3的斷面2處，達(dá)到了52.3 m/(s·m)，具體值見(jiàn)表2。

表1 3種模擬工況表

表2 不同斷面最大臨底流速和流速梯度

圖6 1-1斷面臨底流速曲線圖圖7 2-2斷面臨底流速曲線圖

3種工況下溢洪道中軸線處近底板處的湍流動(dòng)能如圖8所示，可見(jiàn)在1-1斷面和2-2斷面處湍流動(dòng)能值并不是最大值，其湍流動(dòng)能值相對(duì)于最大值處在一個(gè)較小的范圍內(nèi)。

結(jié)合原型觀測(cè)情況，發(fā)現(xiàn)在1-1斷面和2-2斷面處前后均發(fā)現(xiàn)了底板以及墻體的沖刷破壞，如圖9、10所示。對(duì)比發(fā)現(xiàn)與數(shù)值模擬發(fā)現(xiàn)的臨底流速變化規(guī)律相吻合。即原型底板發(fā)生破壞處臨底流速、流速梯度均較大。

圖8 中軸線處湍流動(dòng)能變化曲線

圖9工程反弧段破壞照片圖10 2-2斷面附近側(cè)壁破壞照片

由于臨底流速的增大，在近底處邊界層變薄，流速梯度增大，根據(jù)于洪銀等[16]的研究，過(guò)大的流速梯度是工程沖刷破壞的主要原因之一，可見(jiàn)該工程出現(xiàn)沖刷破壞的主要原因之一是局部臨底流速及其梯度過(guò)大。

6 結(jié) 論

(1)通過(guò)工程原型觀測(cè)數(shù)據(jù)與數(shù)值仿真結(jié)果對(duì)比表明，采用RNGk-ε和VOF模型仿真模擬溢洪道泄流是可行的。在此基礎(chǔ)上開展泄流臨底流速分布規(guī)律研究，可有效解決原型或模型臨底流速難以測(cè)量的問(wèn)題。

(2)工程實(shí)際沖刷破壞處與臨底流速數(shù)值仿真中最大臨底流速及近底側(cè)臨底流速最大梯度出現(xiàn)位置是吻合的，可見(jiàn)引起該工程沖刷破壞的主要原因之一是局部臨底流速及其梯度過(guò)大。