彭燕祥，張 華，吳國(guó)華

（華北電力大學(xué) 水利與水電工程學(xué)院，北京 102206）

1 研究背景

隨著中國(guó)水電事業(yè)的不斷發(fā)展，水電站向高水頭、大流量的方向發(fā)展。泄洪霧化問題越發(fā)突出，尤其是采用挑流消能方式泄洪時(shí)，會(huì)在水電站下游局部區(qū)域形成大面積的雨霧。相關(guān)的研究表明：挑流泄洪霧化與水舌自身摻氣擴(kuò)散、水舌之間相互碰撞和水舌入水噴濺密切相關(guān)，要對(duì)挑流泄洪霧化進(jìn)行更深入的研究，就必須對(duì)挑流水舌的運(yùn)動(dòng)特性，以及水舌與空氣的耦合作用引起的水舌風(fēng)場(chǎng)進(jìn)行研究。

針對(duì)水舌的運(yùn)動(dòng)特性，許多學(xué)者建立不同的數(shù)學(xué)模型進(jìn)行了研究，劉宣烈等［1-3］通過實(shí)驗(yàn)對(duì)水舌運(yùn)動(dòng)過程的摻氣和擴(kuò)散過程進(jìn)行研究，得到了水舌運(yùn)動(dòng)的軌跡方程以及水舌橫、縱斷面的擴(kuò)散規(guī)律。梁在潮等［4］對(duì)水舌運(yùn)動(dòng)軌跡、水舌摻氣量以及霧源量等參數(shù)進(jìn)行了研究，并對(duì)這些參數(shù)進(jìn)行擬合。張華等［5］通過建立水舌運(yùn)動(dòng)微分方程，并用數(shù)值方法進(jìn)行求解，來研究空中水舌運(yùn)動(dòng)特性。劉士和等［6］運(yùn)用相似理論在考慮空氣阻力及卷吸摻氣的作用下，建立了挑流水舌的數(shù)學(xué)模型。劉繼廣等［7］依據(jù)摻氣濃度相似律準(zhǔn)則進(jìn)行氣水兩相流模擬研究。張貝貝等［8］提出一種在排放口近區(qū)考慮射流垂向效應(yīng)、在全計(jì)算域使用二維平面模型的計(jì)算方法。但是這些研究?jī)H僅是對(duì)水舌進(jìn)行了研究，由水舌引起的水舌風(fēng)場(chǎng)并沒有涉及。

為了研究水舌與空氣的耦合作用引起的水舌風(fēng)場(chǎng)。李敏等［9］采用非靜力閉合邊界層模式，模擬了水電站泄洪條件下水舌風(fēng)的發(fā)展過程。彭燕祥等［10］用CFD方法模擬了二灘水電站泄洪時(shí)水舌的運(yùn)動(dòng)及水舌風(fēng)場(chǎng)的空間分布。彭莘仔等［11］采用N-S方程模擬射流的水動(dòng)力過程。但是這些方法都是基于有限元網(wǎng)格的歐拉方法，這個(gè)方法對(duì)于流體運(yùn)動(dòng)的范圍有其局限性，網(wǎng)格的劃分密度以及范圍都對(duì)誤差的影響很大。

格子玻爾茲曼方法（Lattice Boltzmann Method，LBM）是一種介觀流體系統(tǒng)模擬方法，兼具宏觀流體連續(xù)性和微觀分子動(dòng)力學(xué)方法的優(yōu)點(diǎn)，而且計(jì)算效率較高，現(xiàn)已被廣泛應(yīng)用于多組分多相流、化學(xué)反應(yīng)流、氣固兩相流的研究［12-14］。在水汽兩相流兩方面，宋歌等［15］應(yīng)用LBM方法對(duì)離心式霧化噴嘴出口處的液膜厚度、霧化錐角以及其他霧化特性進(jìn)行了研究。葉永等［16］采用LBM對(duì)三維潰壩流進(jìn)行了數(shù)值模擬與驗(yàn)證。兩者都沒有對(duì)空氣流場(chǎng)進(jìn)行研究。張華等［17］應(yīng)用格子Boltzmann兩相流模型，對(duì)水舌在空氣中的運(yùn)動(dòng)進(jìn)行了數(shù)值模擬，但是模擬的水舌是二維的，且研究的水舌尺度較小。

為了分析三維挑流水舌的運(yùn)動(dòng)特性，以及兩股挑流水舌的碰撞特性，本文采用LBM方法對(duì)挑流水舌進(jìn)行三維數(shù)值模擬研究。

2 格子玻爾茲曼方法

粒子運(yùn)動(dòng)的速度矢量見圖1。

圖1 D3Q19網(wǎng)格模型與粒子速度向量

3 挑流水舌空中運(yùn)動(dòng)數(shù)值模擬

3.1 計(jì)算模型與參數(shù)設(shè)置 采用長(zhǎng)江水科院［20］對(duì)泄洪霧化概化模型，進(jìn)行數(shù)值模擬研究，如圖2所示，泄水建筑物長(zhǎng)2.3 m，高1.7 m，寬0.28 m，出口挑坎挑角20°。數(shù)值模擬工況采用了典型的實(shí)驗(yàn)工況，即水頭0.244 m，水頭落差2.05 m。

圖2 數(shù)值模擬幾何模型

3.2 模型參數(shù)設(shè)置 模型采用xflow兩相流計(jì)算模型，分別為液體相和空氣相，液體相密度為998.3 kg/m3，動(dòng)力黏度為0.001 Pa·s?？諝庀嗝芏葹?.225 kg/m3，動(dòng)力黏度為1.7894 10-5Pa·s。

由于時(shí)間步長(zhǎng)和解析度直接關(guān)系到計(jì)算結(jié)果的精度和計(jì)算效率。所以在保證計(jì)算結(jié)果精度的同時(shí)還必須考慮計(jì)算時(shí)間。經(jīng)多方面考慮以及多次試算之后，選擇時(shí)間步長(zhǎng)為0.001 s，解析度為0.03 m最為合適。

環(huán)境參考?jí)毫θ∫粋€(gè)標(biāo)準(zhǔn)大氣壓，邊界條件1、2截面采用壓力進(jìn)口條件，3截面采用壓力出口，相對(duì)壓力都設(shè)置為0。其余面均設(shè)置為自由滑移壁面邊界條件。

3.3 水舌風(fēng)場(chǎng)計(jì)算結(jié)果 圖3為t=2.98 s時(shí)的水舌風(fēng)場(chǎng)的空間分布計(jì)算圖，從圖中可以比較清楚的看出水舌與空氣運(yùn)動(dòng)交界面，挑流水舌隨著水舌的運(yùn)動(dòng)不斷的摻氣擴(kuò)大。水舌的運(yùn)動(dòng)帶動(dòng)周圍的空氣沿水舌方向運(yùn)動(dòng)，并且在水舌入水處，水舌表面空氣流速與水舌的流速相近風(fēng)速達(dá)到6 m/s，風(fēng)向沿水舌速度方向。

圖3 水舌風(fēng)場(chǎng)的空間分布

由于水舌下緣在臨近下游水面區(qū)域與下游水面形成一個(gè)相對(duì)閉塞的空間，特別時(shí)水舌入水附近，如圖3（c）所示，水舌下緣的風(fēng)場(chǎng)情況相較于水舌上緣相對(duì)更加的復(fù)雜，在水舌下沿出現(xiàn)較多的漩渦流場(chǎng)，如圖3（d）所示。從圖中可以看出水舌入水與下游水面碰撞，引起下游水面液滴飛濺，在下游水面形成一個(gè)與水舌入水面為邊界的水滴噴濺區(qū)，并且沿水舌外緣入水區(qū)域的噴濺水滴，由于水舌入水夾角的關(guān)系，沿x方向的速度分量比其他區(qū)域的大。這也是泄洪霧化的主要霧源。

3.4 水舌運(yùn)動(dòng)特性驗(yàn)證分析 數(shù)值模擬采用了典型工況：水頭0.244 m，水頭落差2.05 m。對(duì)水舌的內(nèi)緣和外緣軌跡和挑距，以及水舌的厚度沿水舌軌跡的變化等進(jìn)行分析。

采用2次多項(xiàng)式對(duì)水舌的內(nèi)緣、外緣以及中心內(nèi)線進(jìn)行擬合，分別得到：

內(nèi)緣曲線表達(dá)式為：

外緣曲線表達(dá)式為：

中心線曲線表達(dá)式：

挑流水舌的外緣和內(nèi)緣，分別與文獻(xiàn)［5］做了對(duì)比，如圖4所示。

圖4 挑流水舌外緣和內(nèi)緣的空間變化

從圖4中可以看出，兩者的水舌軌跡基本一致。在水舌的中間段，本文的外緣線比文獻(xiàn)［5］的結(jié)果偏高，中間段水舌的擴(kuò)散速度比文獻(xiàn)［5］更快。

圖5為水舌流速與下游水位截面位置的體積分布計(jì)算結(jié)果，水舌入水的速度在6 m/s左右。水舌內(nèi)緣挑距約為2.8 m，外緣挑距約為3.5 m，橫向?qū)挾燃s為0.75 m。與實(shí)測(cè)數(shù)據(jù)對(duì)比誤差分別為9.6%、7.9%、7.1%，如表1所示，可以看出挑距的模擬值相對(duì)實(shí)測(cè)數(shù)據(jù)偏低。水舌入水的形狀如圖5（b）所示，大致為一個(gè)圓角的矩形趨近于橢圓形，與文獻(xiàn)［21］中單孔泄洪實(shí)驗(yàn)的水舌落水形狀基本一致，這可以為水舌入水噴濺霧源的進(jìn)一步研究提供參考。

表1 水舌挑距實(shí)測(cè)值［17］與模擬值的對(duì)比

圖5 挑流水舌俯視圖與水舌入水截面圖

圖6為水舌模擬與理論計(jì)算的厚度與長(zhǎng)度關(guān)系變化圖，其中H為水舌厚度/m，H0為水舌初始斷面厚度/m，S為水舌長(zhǎng)度/m。通過挑流水舌厚度與水舌長(zhǎng)度的關(guān)系圖，可以將水舌的厚度變化可以分成三個(gè)階段，如圖的Ⅰ、Ⅱ、Ⅲ階段。

圖6 挑流水舌厚度與水舌長(zhǎng)度的關(guān)系

Ⅰ階段水舌擴(kuò)散緩慢，且水舌的厚度的增長(zhǎng)與水舌長(zhǎng)度呈線性關(guān)系，這是由于Ⅰ階段水舌尚有核心區(qū)存在，水舌摻氣量少，水體比較密實(shí)。

Ⅱ階段水舌的厚度擴(kuò)散速度先加快后慢慢變緩，這是因?yàn)橛捎谒嗨俣鹊脑黾?，空氣?duì)水舌的阻力增加，這樣水舌內(nèi)部的紊動(dòng)也就隨之增加，引起水舌的破碎和大量摻氣，使水舌厚度快速的增加，隨著水舌摻氣的增加逐漸趨于飽和，使得摻氣速度逐漸趨于變緩，水舌的厚度增加速度變緩。

Ⅲ階段為水舌的均勻擴(kuò)散段，水舌的厚度的增長(zhǎng)與水舌長(zhǎng)度呈現(xiàn)近似線性變化。這與文獻(xiàn)［1］中通過實(shí)驗(yàn)數(shù)據(jù)分析，將水舌厚度的變化分為三個(gè)不同變化區(qū)段，即：（1）出口緩擴(kuò)段；（2）水舌急劇擴(kuò)散段；（3）水舌均勻擴(kuò)散段；結(jié)論基本一致。

如圖7為水舌速度與水舌長(zhǎng)度的關(guān)系對(duì)比圖，其中的文獻(xiàn)［5］-v、文獻(xiàn)［5］-vx、文獻(xiàn)［5］-vy為采用文獻(xiàn)［5］的數(shù)學(xué)模型，忽略空氣阻力情況下的計(jì)算結(jié)果。從圖中可以看出本文-v和文獻(xiàn)［5］-v隨著S增加先減小后增加、本文-v和文獻(xiàn)［5］-v的上升階段，本文-v相對(duì)與文獻(xiàn)［5］-v的值偏小，這是因?yàn)椴捎梦墨I(xiàn)［5］數(shù)學(xué)模型計(jì)算水舌時(shí)，沒有考慮空氣阻力以及水舌風(fēng)場(chǎng)的影響。

圖7 水舌速度與水舌長(zhǎng)度的關(guān)系

本文-vx值隨著水舌挑S的增加而減小，而文獻(xiàn)［5］-vx水舌由于沒有考慮空氣阻力的影響，水舌速度文獻(xiàn)［5］-vx的值隨著水舌挑距S的增加幾乎不變。

本文-vy和文獻(xiàn)［5］-vy隨S增加先減小到0后增加，當(dāng)本文-vy和文獻(xiàn)［5］-vy等于0時(shí)，水舌上升到最高點(diǎn)。本文-vy與文獻(xiàn)［5］-vy隨S變化趨勢(shì)基本是一致，這是由于在y方向水舌雖然受重力和空氣阻力以及水舌風(fēng)場(chǎng)的共同作用，但是重力起主導(dǎo)作用。

3.5 水舌渦量分析 渦量為流體速度的旋度：

渦量是描寫旋渦運(yùn)動(dòng)的重要物理量之一，在xy平面取水舌的中心平面得水舌的中間截面的渦量云圖，即渦量在z軸的投影用Ωz表示，如圖8所示。從圖中可以看出水舌的上下表面在初始階段Ωz遠(yuǎn)遠(yuǎn)高于內(nèi)部水舌，這是因?yàn)樗喑跏茧A段摻氣并不充分，僅僅外部水舌與空氣作用，使得水舌的上下表面的Ωz遠(yuǎn)遠(yuǎn)高于內(nèi)部水舌。這也是水舌摻氣的重要影響因素，隨著水舌的摻氣量不斷變大，Ωz向水舌的內(nèi)部輸運(yùn)，在水舌的末端Ωz達(dá)到最大，Ωz達(dá)到139 s-1。

圖8 水舌中心截面渦量云圖

4 水舌碰撞的數(shù)值模擬

格子Boltzmann方法可以更加容易的實(shí)現(xiàn)兩股水舌的碰撞，以及由于水舌碰撞而形成的相對(duì)復(fù)雜的水汽交界面，如圖9所示。

圖9 兩股水舌碰撞水舌速度分布圖

圖9和圖10分別為水舌碰撞的數(shù)值模擬圖與文獻(xiàn)［22］中水舌碰撞的實(shí)驗(yàn)圖象。從圖中可以看出兩者的碰撞水舌形態(tài)基本一致，水舌經(jīng)過碰撞后，水花四濺，水舌在碰撞后形成一個(gè)錐形的射流形態(tài)，并且由于兩股水舌的碰撞形成的壓力波影響，使得上股水舌的下表面破碎形成一股反向微射流，微射流以液滴為主且運(yùn)動(dòng)狀態(tài)較為復(fù)雜，但其大多數(shù)液滴的運(yùn)動(dòng)方向與水舌的運(yùn)動(dòng)方向相反。這為進(jìn)一步研究多股水舌泄流，以及水舌入水噴濺等復(fù)雜情況提供參考。

圖10 兩股水舌碰撞的試驗(yàn)結(jié)果圖［22］

5 結(jié)論

本文采用格子Boltzmann方法針對(duì)三維挑流水舌以及兩股水舌的碰撞進(jìn)行數(shù)值模擬研究，得到如下結(jié)論：

（1）格子Boltzmann方法能很好的模擬泄洪挑流水舌和水舌風(fēng)場(chǎng)，水舌的運(yùn)動(dòng)界面清晰可見。并且不用對(duì)水舌的周圍的大片區(qū)域進(jìn)行網(wǎng)格劃分，計(jì)算更加便捷，效率更高。水舌內(nèi)緣挑距約為2.8 m，外緣挑距約為3.5 m，橫向?qū)挾燃s為0.75 m。與實(shí)測(cè)數(shù)據(jù)對(duì)比誤差分別為9.6%、7.9%、7.1%，挑距的模擬值相對(duì)實(shí)測(cè)數(shù)據(jù)偏低。為泄洪挑流水舌的流場(chǎng)研究提供新的方法。

（2）挑流水舌在初始階段由于水舌與空氣的相互作用，使得水舌的上下表面的Ωz遠(yuǎn)遠(yuǎn)高于內(nèi)部水舌，這也是水舌摻氣的重要影響因素，隨著水舌的摻氣量不斷變大，Ωz向水舌的內(nèi)部輸運(yùn)，在水舌的末端渦量達(dá)到最大，渦量達(dá)到139 s-1。

（3）格子Boltzmann方法可以對(duì)兩股水舌的碰撞進(jìn)行模擬，水舌碰撞的數(shù)值模擬與實(shí)驗(yàn)條件下水舌運(yùn)動(dòng)形態(tài)基本一致，水舌經(jīng)過碰撞后，水花四濺，水舌在碰撞后形成一個(gè)錐形的射流形態(tài)，并且由于兩股水舌的碰撞形成的壓力波影響，使得上股水舌的下表面破碎形成一股反向微射流。為多孔泄流的復(fù)雜情況，水舌的運(yùn)動(dòng)特性以及流場(chǎng)研究提供參考。