基于數(shù)值模擬的小浪底龍?zhí)ь^式泄洪洞防洪安全分析

張春晉,孫西歡,李永業(yè),張學(xué)琴,張雪蘭,楊小妮, 5,李 飛

(1.黃河水利科學(xué)研究院水利部黃河泥沙重點(diǎn)實(shí)驗(yàn)室,河南 鄭州 450003;2.太原理工大學(xué)水利科學(xué)與工程學(xué)院,山西 太原 030024; 3.晉中學(xué)院,山西 晉中 030600; 4.章丘黃河河務(wù)局,山東 濟(jì)南 250200; 5.太原理工大學(xué)現(xiàn)代科技學(xué)院,山西 孝義 032300)

黃河小浪底水電站左岸山體布置3條龍?zhí)ь^式泄洪洞,在保障水利工程安全方面發(fā)揮著不可替代的作用,為此對(duì)該龍?zhí)ь^式泄洪洞的泄流能力與安全穩(wěn)定進(jìn)行分析具有重要意義。黃河小浪底水電站龍?zhí)ь^式泄洪洞屬于高壩泄洪洞,水頭高、落差大,形成的高速水流問題備受關(guān)注[1]。近年來,國外學(xué)者很少涉及高壩泄洪洞研究,而國內(nèi)學(xué)者卻對(duì)高壩泄洪洞的水力學(xué)問題開展了大量工作。由于模型試驗(yàn)測(cè)量時(shí)產(chǎn)生的負(fù)壓不符合重力相似,為此采用數(shù)值模擬來彌補(bǔ)模型試驗(yàn)的不足[2]。數(shù)值模擬憑借成本低、效率高、靈活性強(qiáng)及數(shù)據(jù)量大等優(yōu)勢(shì),越來越受到廣大研究人員的重視。郭軍等[3]和胡濤等[4]對(duì)龍?zhí)ь^式泄洪洞摻氣減蝕問題進(jìn)行分析。張春晉等[2]和徐國賓等[5]采用RNGk-ε湍流模型對(duì)高水頭泄洪洞紊流流場(chǎng)進(jìn)行仿真。陳瑞華等[6]和張宏偉等[7]對(duì)龍?zhí)ь^式泄洪洞反弧段流場(chǎng)特性進(jìn)行數(shù)值計(jì)算。范靈等[8]和蘇小麗等[9]采用特征線法對(duì)龍?zhí)ь^式泄洪洞閘門啟閉時(shí)的水擊問題進(jìn)行仿真。上述成果表明龍?zhí)ь^式泄洪洞水力特性的數(shù)值模擬已日臻成熟。

本文以黃河小浪底水電站左岸2號(hào)龍?zhí)ь^式泄洪洞為例,采用商用Fluent 12.0軟件對(duì)設(shè)計(jì)與校核洪水位條件下龍?zhí)ь^式泄洪洞泄流能力、斷面流速、洞頂余幅、測(cè)壓管水頭、空化數(shù)及挑流沖刷等特性進(jìn)行防洪安全分析。同時(shí)闡述了龍?zhí)ь^式泄洪洞存在的安全問題及解決措施。

1 物理模型試驗(yàn)

物理模型為正態(tài)模型,滿足重力相似準(zhǔn)則,比例尺為1∶40。進(jìn)口段、洞身段及挑流段采用有機(jī)玻璃制作。邊墻采用紅磚砌筑。上游水庫與沖刷段采用水泥砂漿面拉毛處理。結(jié)合伊茲巴申公式[10],根據(jù)設(shè)計(jì)挑流水舌落點(diǎn)處強(qiáng)風(fēng)化巖石抗沖流速,得到?jīng)_刷段選用的石料中值粒徑為12.5 mm。流量測(cè)量采用無側(cè)收縮矩形薄壁堰。進(jìn)口段分3層布置9個(gè)壓強(qiáng)測(cè)點(diǎn)。洞身段距底板1.0 cm處布置36個(gè)壓強(qiáng)測(cè)點(diǎn)。挑流段沿底板中軸線布置8個(gè)壓強(qiáng)測(cè)點(diǎn)。洞身段設(shè)置6個(gè)流速測(cè)試斷面。鋼直尺測(cè)量測(cè)壓管水頭。采用SW40型日記式水位計(jì)量測(cè)水面高程,LS300-A型便攜式流速流量?jī)x測(cè)定流速分布。

圖1 泄洪洞體型設(shè)計(jì)圖(單位：m)

2 數(shù)值模擬

2.1 體型設(shè)計(jì)

龍?zhí)ь^式泄洪洞位于小浪底水電站左岸山體,地質(zhì)條件復(fù)雜。泄洪洞由進(jìn)口段、洞身段及挑流段組成,如圖1所示。進(jìn)口段為深水式有壓進(jìn)口,進(jìn)口段斷面由13.6 m×9.83 m(寬×高)漸變?yōu)?.0 m×7.2 m,頂板與側(cè)墻邊界為1/4橢圓(橢圓方程為x2/36+y2/4=1),底板高程為490 m。洞身段由渥奇段、斜直段、反弧段及緩坡段組成,長(zhǎng)度為261 m。斷面為7.0 m×9.0 m的城門洞型,拱形頂角為130°。渥奇段底板遵循y=x2/260,斜直段坡度為1∶2,反弧段半徑為30 m,緩坡段坡度為1∶100。挑流段長(zhǎng)度為30 m,呈對(duì)稱布置。挑坎反弧半徑為40 m,底板遵循y=x2/240+x/100,挑射角為25°。挑坎頂部高程為468.39 m。下游河床為砂礫石動(dòng)床,長(zhǎng)度為150 m。設(shè)計(jì)洪水位為535.55 m,校核洪水位為538.35 m。

2.2 幾何模型

采用Pro/E 5.0軟件建立龍?zhí)ь^式泄洪洞的幾何模型。幾何模型由上游水庫、進(jìn)口段、洞身段、挑流段及沖刷段組成，長(zhǎng)度為540 m。上游水庫底板高程為490 m，沖刷段底板高程為450 m。上游水庫樁號(hào)為D0-037至D0-007。進(jìn)口段樁號(hào)為D0-007至D0+029。洞身段樁號(hào)為D0+029至D0+290。挑流段樁號(hào)為D0+290至D0+320。沖刷段樁號(hào)為D0+320至D0+470。

RNGk-ε湍流模型可更好地處理旋流、高應(yīng)變率流及流線彎曲程度較大的流動(dòng)問題[11-14]。為此,選擇RNGk-ε湍流模型進(jìn)行計(jì)算。采用PISO算法耦合求解[15]?？刂品匠滩捎萌[格式迭代求解,并設(shè)置合理的松弛因子[16]。壓力項(xiàng)采用PRESTO格式進(jìn)行離散,動(dòng)量方程、湍動(dòng)能方程和湍流耗散率方程采用二階迎風(fēng)格式離散,雷諾應(yīng)力項(xiàng)采用一階迎風(fēng)格式離散。時(shí)間步長(zhǎng)為1×10-5～1×10-3s。

龍?zhí)ь^式泄洪洞水面線隨時(shí)間發(fā)生變化,水面線的處理對(duì)于求解精度有較大影響。本文采用VOF法對(duì)水面線進(jìn)行計(jì)算[7]。為了解決空化問題,引入Singhal空化模型[17],該模型由氣泡動(dòng)力學(xué)方程與氣液兩相守恒方程聯(lián)合推導(dǎo)得出。

2.3 網(wǎng)格劃分

采用ICEM CFD 17.0軟件進(jìn)行網(wǎng)格劃分。為了提高網(wǎng)格劃分效率,將計(jì)算域劃分為3部分,即上游水庫、洞身段及沖刷段。由于洞身段體型不規(guī)則,采用較小的非結(jié)構(gòu)化四面體網(wǎng)格加密,而上游水庫和沖刷段體型規(guī)則,采用較大的結(jié)構(gòu)化六面體網(wǎng)格加密。通過網(wǎng)格無關(guān)性檢驗(yàn)可知,當(dāng)網(wǎng)格數(shù)高于384萬時(shí),網(wǎng)格數(shù)量對(duì)于斷面平均流速影響可忽略不計(jì)。考慮計(jì)算時(shí)間與仿真精度,認(rèn)為網(wǎng)格單元數(shù)為384萬時(shí)可滿足網(wǎng)格無關(guān)性檢驗(yàn)要求[18]。

2.4 邊界條件

上游水庫進(jìn)口采用壓力進(jìn)口邊界,壓力按照壩前特征水位計(jì)算。閘門室和上游水庫上邊界采用大氣壓進(jìn)口邊界。沖刷段出口采用自由出流邊界。挑流段和沖刷段上邊界采用大氣壓出口邊界。壁面采用無滑移邊界。近壁面區(qū)域采用標(biāo)準(zhǔn)壁面函數(shù)法[11]處理。上游水庫與進(jìn)口段和挑流段與沖刷段之間的連接斷面采用Interface邊界,可實(shí)現(xiàn)內(nèi)部不同流場(chǎng)數(shù)據(jù)的實(shí)時(shí)交換與調(diào)用。

3 數(shù)學(xué)模型驗(yàn)證與應(yīng)用

3.1 泄流能力

泄流能力是龍?zhí)ь^式泄洪洞體型設(shè)計(jì)的核心指標(biāo)。泄流能力采用閘孔出流流量系數(shù)衡量,流量系數(shù)m[19]可表達(dá)為

(1)

式中:Q為流量;A和B分別是有壓段出口的寬和高;Hs為工作閘門頂水頭;g為重力加速度。

圖2為流量系數(shù)模擬值與試驗(yàn)值對(duì)比。從圖2可見模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過0.34%。上游庫區(qū)達(dá)到校核水位時(shí),m=0.610 6,該值略高于設(shè)計(jì)流量系數(shù),表明泄洪洞泄流量已達(dá)到防洪要求。

圖2 流量系數(shù)模擬值與試驗(yàn)值對(duì)比

3.2 斷面流速與洞頂余幅

圖3為中軸線水面線模擬值與試驗(yàn)值對(duì)比。從圖3可知:模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過1.09%;流體在離心力的作用下,造成反弧段流速增大,水深減小;流線在渥奇段發(fā)生彎曲,引起水流高度湍動(dòng),使得流體質(zhì)點(diǎn)挾帶空氣進(jìn)入流體,形成摻氣水流,水深增大;挑流段采用連續(xù)擴(kuò)散結(jié)構(gòu)使得斷面增大,水深減小。

圖3 中軸線水面線模擬值與試驗(yàn)值對(duì)比

圖4為龍?zhí)ь^段與洞身水平段交接部位流態(tài)。從圖4可知:龍?zhí)ь^段與洞身水平段交接處流態(tài)平穩(wěn),未出現(xiàn)明滿流交替現(xiàn)象;反弧段上游引起流線彎曲,導(dǎo)致流體質(zhì)點(diǎn)將空氣挾帶進(jìn)入流體,形成摻氣水流,而洞身水平段流態(tài)穩(wěn)定。因此反弧段上游水深較大;反弧段下游受到離心力和重力作用,流速增加。

圖4 龍?zhí)ь^段與洞身水平段交接部位流態(tài)

圖5為斷面流速分布云圖，圖中不同樁號(hào)處斷面相同,邊墻高度為6.3 m,寬度為7.0 m。從圖5可知:渥奇段水位高于緩坡段水位。渥奇段流線發(fā)生彎曲,使得流體形成摻氣水流,水位升高。緩坡段流體湍流強(qiáng)度小,摻氣量較少,導(dǎo)致水位較低;由于反弧段末端曲率半徑趨于零,使得流體在離心力作用下對(duì)底板壓力增大,引起水位降低,流速升高。在緩

圖5 不同樁號(hào)斷面流速分布云圖

坡段,流體對(duì)于底板的壓強(qiáng)減小,使得該區(qū)域水位升高。因此,緩坡段水位隨樁號(hào)升高;斷面流速分布呈現(xiàn)近壁面流速低,而遠(yuǎn)離壁面流速均勻。由于流體具有黏滯性,使得近壁面形成了流速梯度變化明顯的邊界層。邊界層引起泄流能力降低。為此,體型優(yōu)化時(shí)需考慮邊界層影響;根據(jù)設(shè)計(jì)規(guī)范,洞頂余幅不低于洞深的25%,可避免高速水流在洞內(nèi)形成明滿流交替的水流流態(tài)[20]。經(jīng)計(jì)算,斷面洞頂余幅大于30%,說明洞頂余幅設(shè)計(jì)較為合理。

圖6 不同樁號(hào)中軸線流速模擬值與試驗(yàn)值對(duì)比

圖6為中軸線流速模擬值與試驗(yàn)值對(duì)比。從圖6可知:模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過5.84%;由于底板邊界層厚度較小,僅對(duì)邊界層以外流場(chǎng)的流速進(jìn)行測(cè)量;中軸線處流速分布分為2部分,即紊流邊界層和邊界層外勢(shì)流區(qū)。紊流邊界層最外側(cè)流速為0.99倍平均流速,通過計(jì)算,紊流邊界層厚度為0.2～0.3 m。邊界層外勢(shì)流區(qū)流速一致,而紊流邊界層流速梯度變化大且呈斜直線分布。因此,流速分布存在明顯的突變點(diǎn)或拐點(diǎn),且距泄洪洞底板約為0.2 m。

3.3 測(cè)壓管水頭與空化數(shù)

空化數(shù)是判斷洞內(nèi)是否發(fā)生空化空蝕的重要指標(biāo)[3]?？栈瘮?shù)[5]可表達(dá)為

(2)

式中:u0為斷面平均流速;γ為流體容重;Pb為相對(duì)壓強(qiáng);Pv為飽和蒸汽壓;H為當(dāng)?shù)馗叱獭Ｒ罁?jù)飽和蒸汽壓強(qiáng)表,水溫20℃時(shí),飽和蒸氣壓(水頭)為0.238 4 m。

圖7與圖8分別為測(cè)壓管水頭與空化數(shù)沿樁號(hào)模擬值與試驗(yàn)值對(duì)比。實(shí)際工程中,當(dāng)空化數(shù)小于0.3時(shí),洞內(nèi)必然發(fā)生空化空蝕。當(dāng)空化數(shù)在0.1～0.3之間,必須注意空化問題。

圖7 測(cè)壓管水頭模擬值與試驗(yàn)值的對(duì)比

圖8 空化數(shù)模擬值與試驗(yàn)值的對(duì)比

從圖7和圖8可知:模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過4.62%;由于渥奇段底板彎曲程度大,使得高速流體產(chǎn)生了邊界層分離,導(dǎo)致近壁面壓強(qiáng)降低。渥奇段沿樁號(hào)摻氣程度減小,使得該區(qū)域流速升高。兩者共同導(dǎo)致渥奇段空化數(shù)沿樁號(hào)降低。為此需在渥奇段進(jìn)行補(bǔ)氣,降低流體對(duì)壁面的汽蝕程度;挑流段反弧結(jié)構(gòu)使得流體在慣性力作用下對(duì)底板的壓強(qiáng)增大,使得該區(qū)域測(cè)壓管水頭升高;由于反弧段曲率半徑為零,使得流體在離心力作用下斷面流速和底板壓強(qiáng)增大。反弧段下游壓強(qiáng)較大,而緩坡段壓強(qiáng)較小,使得反弧段下游區(qū)域底板壓強(qiáng)急劇降低,加之流速變化小,引起反弧段下游附近流場(chǎng)發(fā)生空化空蝕,導(dǎo)致底板嚴(yán)重剝落與磨損。為此需對(duì)反弧段下游附近底板進(jìn)行加固;由于渥奇段的底板沿樁號(hào)彎曲程度較大,使得流體在渥奇段無動(dòng)能轉(zhuǎn)化為壓能,導(dǎo)致流體離開底板,并在底板近壁面形成空腔。隨后下游流體隨即填補(bǔ)空腔,形成旋渦損失,使得渥奇段底板測(cè)壓管水頭沿樁號(hào)降低。

3.4 挑流沖刷

表1為挑流射程模擬值與試驗(yàn)值對(duì)比。從表1可知:模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過0.73%;由于沖坑位置與挑坎距離較遠(yuǎn),使得挑流沖坑不會(huì)危害到壩體安全。

表1 挑流射程模擬值與試驗(yàn)值對(duì)比

根據(jù)泄洪洞下游砂卵石河床地質(zhì)情況,結(jié)合沖坑深度經(jīng)驗(yàn)公式可計(jì)算出沖坑深度模擬值,沖坑深度公式[2]可表達(dá)為

(3)

式中:ts為沖坑深度;ht為沖坑后的下游水深;.β為水舌的入射角；η為流速脈動(dòng)系數(shù),取1.5～2.0;vt為水舌進(jìn)入下游水面流速;ω為河床顆粒水力粗度;q為挑流段出口模擬計(jì)算的單寬流量。

圖9為水舌寬度沿樁號(hào)模擬值與試驗(yàn)值對(duì)比。從圖9可知模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過3.95%;挑流段結(jié)構(gòu)為進(jìn)口窄、出口寬,導(dǎo)致水舌寬度沿樁號(hào)增加;由于挑入空中的水舌失去了邊界約束,在空氣阻力作用下,發(fā)生摻氣和分散,引起水舌在水平方向動(dòng)能減小,導(dǎo)致水舌沿樁號(hào)向外側(cè)擴(kuò)散角度降低。

圖9 水舌寬度模擬值與試驗(yàn)值的對(duì)比

根據(jù)《水工(常規(guī))模型試驗(yàn)規(guī)程》,水工模型需連續(xù)沖刷2.0 h,方可測(cè)量沖坑深度。模型試驗(yàn)沖坑深度可通過鋼直尺直接測(cè)量。表2為沖坑深度模擬值與試驗(yàn)值對(duì)比。從表2可知:模擬值與試驗(yàn)值基本一致,相對(duì)誤差不超過1.81%;對(duì)沖坑深度分析得到?jīng)_坑對(duì)左岸坡腳有掏空作用,需對(duì)左岸山體采取錨桿加固等措施或改變挑坎形式使得流體遠(yuǎn)離山體。

表2 沖坑深度模擬值與試驗(yàn)值對(duì)比

4 安全分析

4.1 安全問題

龍?zhí)ь^式泄洪洞防洪存在以下安全問題:校核洪水位泄流能力不足;洞頂余幅不足25%,造成明滿流交替,引起洞頂混凝土剝落和洞內(nèi)補(bǔ)氣不足;高壩泄洪洞底板發(fā)生空化空蝕,造成洞內(nèi)混凝土表面剝蝕;反弧段流體在離心力作用下對(duì)底板壓力增大,引起該區(qū)域下游附近底板發(fā)生空化空蝕；反弧段下游底板剝落的混凝土將堵塞下游過流通道;底板粗糙度較大,引起流體向四周擴(kuò)散,增加了泄洪洞的不穩(wěn)定性;霧化現(xiàn)象導(dǎo)致高強(qiáng)度降雨和局部山體滑坡;挑流將對(duì)岸坡山體起到?jīng)_刷與掏空作用。

4.2 解決措施

可采用“龍落尾”布置方式,將發(fā)生空化空蝕的位置布置在下游區(qū)域;提高洞頂高度,避免洞身內(nèi)發(fā)生明滿流交替;對(duì)易發(fā)生空化空蝕區(qū)域采用抗沖耐磨混凝土材料,提高混凝土耐久性;對(duì)渥奇段采用摻氣減蝕措施;對(duì)沖刷山體岸坡錨桿加固,提高山體抗沖能力。

5 結(jié) 語

本文以黃河小浪底水電站2號(hào)龍?zhí)ь^式泄洪洞作為研究對(duì)象,采用數(shù)值模擬對(duì)龍?zhí)ь^式泄洪洞進(jìn)行了防洪安全分析，同時(shí)探討了龍?zhí)ь^式泄洪洞存在的防洪安全問題及解決措施。研究表明采用數(shù)值模擬研究龍?zhí)ь^式泄洪洞水力特性可行,可為今后龍?zhí)ь^式泄洪洞泄洪能力與安全穩(wěn)定性評(píng)估提供參考。