徐建榮，柳海濤，彭 育，孫雙科，余凱文，3

（1.中國電建集團 華東勘測設(shè)計研究院有限公司，浙江 杭州 311122；2.中國水利水電科學研究院，北京 100038；3.南京水利科學研究院，江蘇 南京 210029）

1 研究背景

高壩工程泄洪霧化伴隨大壩全生命周期，是一種非自然的災(zāi)害性氣象事件[1]。近30年來我國大批高壩工程相繼開始建設(shè)與運行，并且開始向高海拔地區(qū)發(fā)展，泄洪環(huán)境更加復(fù)雜，泄洪霧化引發(fā)的安全與環(huán)境問題更加凸顯[2～5]。

高壩泄洪霧化的影響因素眾多，既包括樞紐布置、結(jié)構(gòu)型式、調(diào)度方式等主觀條件[6～8]，也包括風場、地形、氣壓、氣溫等客觀條件[9～10]。泄洪建筑物布置與調(diào)度方式，主要影響泄洪水舌擴散與運動過程，從而改變霧化源噴射條件。在挑射水流運動模擬方面，物理模型方法仍然是目前較為可靠的研究手段[11～12]，可以正確反映復(fù)雜水舌入水參數(shù)及分布形態(tài)，據(jù)此換算得到霧化源噴射條件[13]。風場與河谷地形，主要影響霧化雨滴的運動過程與降落位置，物理模型中水滴粒徑、表面張力、大氣黏滯性以及空氣阻力等因素無法與原型相似[14～16]，故應(yīng)以數(shù)學模型作為主要研究手段，由于噴射水相以水團、水滴形式運動，宜采用拉格朗日方法，通過追蹤水滴在風場中的運動軌跡直至下墊面，從而統(tǒng)計得到地面降雨量分布[17～18]。氣象條件與泄洪霧化過程相互影響，一方面，在高海拔的低溫、低氣壓地區(qū)，水滴運動阻力、浮力以及自身裂散、碰并、蒸發(fā)過程均發(fā)生改變，對此宜采用理論分析、原型驗證、模型試驗相結(jié)合的方法[19]，建立和完善理論模型；另一方面，泄洪霧化也會引起當?shù)仫L場、濕度、氣溫變化[20～21]，此時水相以霧滴與水汽形式運動，可采用歐拉方法，建立水氣兩相流擴散模型[22～23]，也可借鑒氣象學方面的相關(guān)知識開展研究[24～25]。上述問題的研究中，尚存一些理論與技術(shù)難題亟待解決：首先，水舌風伴隨泄洪全過程，受泄洪方式與河谷地形影響，對霧化雨滴運動及降落具有至關(guān)重要的作用，因此在復(fù)雜風場及下墊面地形條件下粒子運動模型的耦合求解方法，成為霧化數(shù)學模型的關(guān)鍵技術(shù)之一。最初是將水舌風場及地形概化為經(jīng)驗分布函數(shù)，直接求解水滴運動方程，但其得到的霧化降雨分布誤差甚大，后又提出采用人工神經(jīng)網(wǎng)絡(luò)方法，先對復(fù)雜風場與地形進行學習形成預(yù)報模塊，在求解水滴運動過程中，調(diào)用預(yù)報模塊實時計算水滴鄰近風速與下墊面高程[18]。實際計算表明，對于復(fù)雜三維風場，神經(jīng)網(wǎng)絡(luò)模塊學習難度甚大，同時若在每個時間步長內(nèi)調(diào)用神經(jīng)網(wǎng)絡(luò)模塊進行風場與地形計算，程序運算量巨大而至今無法應(yīng)用。第二是關(guān)于氣象因素對于霧化降雨分布的影響，由于其細觀機理復(fù)雜，是當前霧化研究的熱點問題之一，國內(nèi)學者針對高海拔低氣壓的影響開展了試驗室減壓模擬研究[19]，成果離實際應(yīng)用尚有差距。第三是泄洪運行調(diào)度也會改變霧化源與降雨分布形態(tài)，以往研究重點關(guān)注泄洪流量、入水流速與角度等宏觀影響因素，對于采用不同孔口泄洪組合時，引起的水舌落點、入水形態(tài)、水舌風等細觀差異考慮不足，無法分辨其對霧化降雨分布的影響。由于實際工程泄洪時，海拔高程與氣溫、樞紐布置與泄洪調(diào)度、下游風場與地形條件等影響因素往往同時存在，對霧化分析方法的精度與適應(yīng)性，提出了更高的要求。

針對上述問題，本文作者開展了探索研究，通過將復(fù)雜霧化源離散為一系列噴射點源，同步運算，解決了各種泄洪組合條件下入水條件差異問題[26]；研究提出了一種復(fù)雜風場與地形數(shù)據(jù)的結(jié)構(gòu)化存取方法，實現(xiàn)了其與巨量水滴運動的快速耦合計算；借鑒氣象學方面的知識，提出了氣壓、氣溫對于霧化降雨影響的分析公式[27]；基于上述成果，提出了改進的霧化濺水數(shù)學模型。本文在此基礎(chǔ)上，運用該模型對白鶴灘水電站在中小洪水條件下，壩身孔口不同泄洪組合，進行霧化影響敏感性分析，驗證其模擬精度與適應(yīng)性，為今后實際工程霧化影響的分析與預(yù)測，提供可靠技術(shù)手段。

2 隨機濺水數(shù)學模型的基本方法

2.1 水滴運動的基本方程霧化雨滴在運動過程中，受到重力、浮力和空氣阻力的共同作用，由此建立水滴運動的力學微分方程[18]：

式中，u、v、w為水滴在空間坐標[x、y、z]處的運動速度，m/s，uf、vf、wf分別為水滴鄰近風速，m/s；Cf為阻力系數(shù)；d為水滴粒徑，m；ρa為空氣密度，kg/m3；ρw為水的密度，kg/m3；g為重力加速度，m/s2。求解上述方程組需給定水滴隨機噴射條件[27]，并采用4階Runge-Kutta法進行數(shù)值求解。

2.2 海拔高程與氣溫對于水滴運動的影響海拔高程與氣溫的變化，會引起空氣密度、氣壓、黏滯系數(shù)、阻力系數(shù)的變化，最終引起泄洪霧化降雨分布發(fā)生改變。

2.3 數(shù)值模型計算流程水滴隨機濺水模型的數(shù)值計算流程見圖1。

圖1 隨機濺水數(shù)學模型計算流程

首先，根據(jù)海拔高程與氣溫條件，由式（2）—式（5）求得大氣密度、運動黏滯系數(shù)、阻力系數(shù)、當?shù)貧鈮旱葪l件，用于求解水舌風場作用下水滴運動方程；根據(jù)物理模型試驗結(jié)果，求得水舌入水條件，包括入水位置、速度與角度、前緣長度、含水量等變量，復(fù)雜水舌入水時具有多個離散點源；水舌風場與河谷地形，需要通過空間插值方法，轉(zhuǎn)化為結(jié)構(gòu)化數(shù)據(jù)格式[27]。在模型運行初始，先讀取上述風場與地形數(shù)據(jù)，按照其結(jié)構(gòu)化坐標存儲在相應(yīng)的數(shù)組單元中。

第二，根據(jù)水舌入射條件，給定隨機噴射眾值條件，包括每個噴射源的粒徑、噴射速度、垂直角、擴散角、顆粒流量等變量。將噴射歷時T劃分為m個時間步長，每個步長初始，根據(jù)噴射眾值條件，隨機生成一組噴射水滴[19]，水滴總量為Ni=nT/m，n為噴射源顆粒流量，時間步長序號i=1，m。然后，由式（1）計算此Ni組噴射水滴的運動過程，最大飛行時長為ti=(m -i+1)Tm，對應(yīng)計算時間步長序號j=1，m-i+1。

第三，在求解水滴運動過程中，對于任意水滴，需要判別其鄰近風速與地面高程。在時間步長初始，根據(jù)水滴空間位置先找到其相鄰的8個結(jié)構(gòu)化坐標，直接從對應(yīng)數(shù)組單元讀取相應(yīng)的風速數(shù)據(jù)[uf，vf，wf]，然后插值得到鄰近風速；同樣地，在時間步長末，根據(jù)水滴空間位置，插值得到其下墊面高程若飛行高度小于下墊面高程，則終止飛行，水滴體積累積到相應(yīng)的地面網(wǎng)格。

第四，整個噴射過程計算完成后，根據(jù)地面網(wǎng)格累積水量、投影面積與持續(xù)時間，計算得到降雨強度分布。上述隨機噴射過程需重復(fù)計算M次，為保證計算精度一般應(yīng)不小于50次。

3 白鶴灘電站壩身泄洪方式霧化影響分析

3.1 數(shù)學模型計算驗證中電建昆明院科研所于2014年8月17日，針對小灣水電站泄洪洞全開泄洪過程，水舌入水形態(tài)與泄洪霧化進行了原型觀測。觀測過程中，壩上水位1236.32 m，實測泄洪流量3500 m3/s，當?shù)睾０?000 m，氣溫20℃，作者運用改進后的隨機濺水數(shù)學模型，對于泄洪洞下游風場與霧化降雨分布進行了驗證計算[27]。根據(jù)當?shù)睾０闻c氣溫條件，求得當?shù)卮髿怵禂?shù)1.6791×10-5m2/s，大氣密度1.0805 kg/m3，大氣壓強89 999.2Pa?；谏鲜鰠?shù)，運用Fluent軟件求解泄洪洞水舌下游河谷風場。然后，根據(jù)泄洪水舌入水條件、泄洪風場以及當?shù)貧庀髤?shù)，計算得到下游河谷中霧化降雨分布形態(tài)。圖2為泄洪洞下游霧化降雨強度分布與實測點據(jù)對比分析結(jié)果，降雨強度計算值與實測值的相對誤差在10%以內(nèi)。

圖2 小灣電站泄洪洞下游霧化降雨強度計算結(jié)果

3.2 白鶴灘水電站壩身泄洪調(diào)度方案白鶴灘水電站壩身泄洪建筑物由6個表孔、7個深孔組成，壩下設(shè)長約400 m水墊塘、二道壩消能。電站壩頂海拔高程834 m，最大壩高289 m，下游河谷狹窄，泄洪流量大，霧化降雨分布范圍大。本文主要針對電站中小洪水下壩身孔口泄洪調(diào)度方案，開展霧化影響分析研究。表1為正常蓄水位下壩身孔口開啟方式，由此可知，不同開啟方式泄洪流量與水位落差基本相同、僅泄洪水舌入水條件略有不同。各工況對應(yīng)的水舌入水位置及形態(tài)見圖3，其中陰影部分為表孔入水形態(tài)。計算中，當?shù)睾０胃叱倘?00 m，當?shù)貧鉁厝?0℃，大氣密度1.1047 kg/m3，空氣黏滯系數(shù)1.6421×10-5m2/s，當?shù)貧鈮?2 190 Pa，作為泄洪風場計算與霧化降雨分析的基本環(huán)境參數(shù)。

表1 白鶴灘電站正常蓄水位下壩身孔口開啟方式與泄流條件

圖3 白鶴灘電站壩身孔口泄洪水舌入水位

3.3 典型工況泄洪霧化風場與降雨分布結(jié)果分析圖4與圖5分別為電站壩身1#與4#表孔、3#與6#表孔，2種泄洪組合的水舌風場與霧化降雨分布。當1#、4#表孔組合泄洪時，壩下水舌風場明顯偏向左岸，同時霧化暴雨區(qū)也偏向左岸；當3#、6#表孔組合泄洪時，壩下水舌風場與霧化暴雨區(qū)基本位于河道中央，此泄洪組合為相對較優(yōu)。究其原因在于，4#表孔水舌靠近水墊塘中部，呈對稱擴散，滿足水墊塘消能防沖要求，但1#表孔水舌落點靠近水墊塘左側(cè)邊緣，兩者聯(lián)合泄洪時，水舌風場與霧化降雨集中于左岸一側(cè)。

圖4 1#與4#表孔泄洪組合壩下水舌風場與霧化降雨分布

圖5 3#與6#表孔泄洪組合壩下水舌風場與霧化降雨分布

圖6與圖7分別為電站壩身1#與7#深孔、3#與5#深孔，2種泄洪組合的水舌風場與霧化降雨分布形態(tài)。當1#、7#深孔組合泄洪時，壩下風場與降雨分布明顯偏向右岸；當3#、5#深孔組合泄洪時，壩下風場與霧化雨區(qū)位于河道中央。究其原因在于，1#與7#深孔泄洪中心線呈向心分布，使得壩下風場與霧化降雨相對集中，橫向上1#深孔風場占優(yōu)，使得霧雨區(qū)偏向右岸一側(cè)。

圖6 1#與7#深孔泄洪組合壩下水舌風場與霧化降雨分布

圖7 3#與5#深孔泄洪組合壩下水舌風場與霧化降雨分布

本文所采用的數(shù)學模型，可以精確分辨泄洪孔口組合、水舌落點、分布形態(tài)及兩岸岸坡，對于泄洪霧化風場與降雨分布的影響，可用于今后泄洪調(diào)度方案的深化研究。

3.4 各種工況泄洪霧化降雨強度等值線分布結(jié)果分析據(jù)表深孔泄洪組合下霧化降雨點據(jù)，繪制出降雨強度等值線圖，具體分布見圖8與圖9。計算結(jié)果表明，表孔泄洪組合霧化降雨縱向范圍明顯小于深孔泄洪組合；表孔泄洪以3#與6#孔組合最優(yōu)，2#與5#孔組合次之，1#與4#孔組合較差；深孔泄洪以3#與5#孔組合最優(yōu)，2#與6#孔組合次之，1#與7#孔組合較差；當壩身孔口采用單獨開啟方式時，以3#表孔為最優(yōu)，4#深孔次之。

圖8 表孔泄洪組合壩下霧化降雨強度等值線分布

圖9 深孔泄洪組合壩下霧化降雨強度等值線分布

4 結(jié)論

白鶴灘水電站壩身布置6個表孔與7個深孔，在常遇洪水下頻繁采用局部開啟方式泄洪。由于各孔間布置體型與位置存在差異，在下游形成的霧化降雨分布明顯不同，為優(yōu)化壩身孔口泄洪調(diào)度方案，保證電站長期運行安全，本文針對其霧化影響進行了敏感性分析。首先，通過對現(xiàn)有數(shù)學模型的數(shù)值分析方法與計算流程進行改進，開發(fā)了可以綜合考慮水舌入水形態(tài)、當?shù)貧庀髼l件、復(fù)雜河谷地形的泄洪霧化數(shù)值分析模型。然后，根據(jù)小灣泄洪洞下游霧化原型觀測數(shù)據(jù)進行計算驗證，結(jié)果表明，數(shù)學模型在水舌風場、降雨強度分布方面與原型觀測結(jié)果較為吻合，在0～800 mm/h降雨強度范圍內(nèi)，相對誤差小于10%。在此基礎(chǔ)上，針對白鶴灘電站壩身泄洪8種孔口開啟方式，定量分析不同孔口水舌落點與入水形態(tài)差異，所引起的水舌風場與霧化降雨變化，通過對比降雨強度分布規(guī)律的差異，進一步判斷各種孔口開啟方式的優(yōu)劣。對比分析表明，在中小洪水泄洪條件下，表孔泄洪組合明顯優(yōu)于深孔泄洪組合；表孔泄洪以3#與6#孔組合最優(yōu)，2#與5#孔組合次之，1#與4#孔組合較差；深孔泄洪以3#與5#孔組合最優(yōu)，2#與6#孔組合次之，1#與7#孔組合較差；當壩身孔口單獨開啟時，以3#表孔為最優(yōu)，4#深孔次之。上述研究為今后實際工程泄洪安全問題分析，優(yōu)化泄洪建筑物設(shè)計與運行調(diào)度方案，提供了借鑒。