柳海濤，徐建榮，孫雙科，彭 育

（1.流域水循環(huán)模擬與調(diào)控國家重點實驗室，中國水利水電科學(xué)研究院，北京 100038）（2.中國電建集團 華東勘測設(shè)計研究院有限公司，浙江 杭州 311122）

1 研究背景

大型水電站泄洪霧化對于工程安全與周邊環(huán)境有重大影響[1-4]，如何評估與防范霧化危害，是水利水電工程設(shè)計中的重要研究課題之一。泄洪霧化的本質(zhì)是水團、水滴以及水霧在大氣中的運動，其中水舌入水激濺占據(jù)主導(dǎo)地位，由此構(gòu)成了霧化降雨的主體。國內(nèi)目前主要采用隨機濺水?dāng)?shù)學(xué)模型對霧化降雨進行分析[5-8]，該類模型通過描述每個水滴在空氣中的運動過程，同時考慮空氣阻力、浮力、重力等條件，求得霧化濺水區(qū)內(nèi)下墊面降雨強度分布。上述模型在計算中主要考慮泄洪水力學(xué)條件與地形條件的影響，對于海拔高程與氣象條件的影響尚未考慮。

目前我國水電工程建設(shè)與運行已從低海拔地區(qū)發(fā)展到2000 m以上高海拔地區(qū)，泄洪霧化與氣象環(huán)境之間相互影響不容忽視。從定性的角度看，海拔高程對于霧化降雨的影響，主要體現(xiàn)在三個方面：（1）不同海拔地區(qū)大氣壓強不同，從而影響水流摻氣、擴散過程，引起雨霧總量發(fā)生改變；（2）不同海拔地區(qū)空氣密度與黏滯系數(shù)不同，使得霧化水滴在運動中所受阻力、浮力等作用發(fā)生改變，引起霧化降雨分布發(fā)生改變；（3）不同地域氣壓與空氣密度的差異，使得的泄洪雨霧輸運與沉降過程發(fā)生改變，最終影響外圍雨霧濃度分布規(guī)律。

針對上述問題，需通過不同海拔地區(qū)對比試驗、原型觀測、理論分析、數(shù)學(xué)模型等手段進行系統(tǒng)研究。鑒于問題的復(fù)雜性，本文僅針對海拔高程對于霧化降雨分布的影響，開展探索性研究。通過理論分析，在前人研究的基礎(chǔ)上，建立考慮當(dāng)?shù)睾０胃叱逃绊懙男购殪F化降雨隨機濺水?dāng)?shù)學(xué)模型，結(jié)合實際工程霧化原型觀測資料，進行數(shù)學(xué)模型計算驗證。然后，針對不同海拔高程條件下，泄洪霧化降雨分布的變化規(guī)律進行敏感性分析，為今后全面考慮氣象條件對于泄洪霧化過程的影響，奠定理論研究基礎(chǔ)。

2 隨機濺水?dāng)?shù)學(xué)模型的基本理論

將泄洪水舌入水噴射過程視為一種隨機噴射現(xiàn)象，采用拉格朗日方法，對噴射水滴在空氣中的運動進行跟蹤，然后通過統(tǒng)計方法求得濺水區(qū)內(nèi)下墊面上的降雨強度分布。

2.1 水滴運動的微分方程水滴在運動過程中，受到重力、浮力和空氣阻力的共同作用，由此建立水滴運動的力學(xué)微分方程：

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

2.2 水滴隨機噴射條件

（1）水滴直徑d采用概率密度函數(shù)分布[8]：

（2）水滴初始噴射速度u采用概率密度函數(shù)分布：

（3）水滴出射角θ采用概率密度函數(shù)分布：

（4）水滴出射偏轉(zhuǎn)角φ采用概率密度函數(shù)分布：

式中：μ為偏轉(zhuǎn)角眾值，依據(jù)水舌平面偏轉(zhuǎn)角度取值；σ為偏轉(zhuǎn)角的均方差，σ取值在20°～30°之間，本文取22.5°。

（5）水滴噴射速度與角度的眾值采用下式表示[10]：

（6）水滴噴射顆粒流量n：

水舌入水激濺主要發(fā)生于入水前緣，其噴射厚度h可以表示為[11]：

式中：η為系數(shù)，可取25；C為含水濃度，可根據(jù)水舌入水時斷面形態(tài)判斷，一般地，對于充分發(fā)展的摻氣水舌，一般在0.03左右；v為水的運動黏滯系數(shù)，m2/s；R為水力半徑，m；u*為摩阻流速，m/s，其表達(dá)式為其中τ為空氣阻力

水舌入水噴射的總流量可進一步表示為：

式中：ks為噴濺系數(shù)，ks=0.01～0.03，本文取0.01；l為水舌入水前緣總長度，m。

一旦噴射總流量確定，則對應(yīng)的水滴顆粒流量可以表示為：

式中dm為整個噴射過程的水滴平均粒徑，該粒徑值未知且不同于粒徑眾值。為此，本文模型首先計算粒徑眾值對應(yīng)的顆粒流量然后，由式（2）中概率密度函數(shù)，隨機生成個噴射水滴其對應(yīng)噴射總流量為最后，式（8）中總流量對應(yīng)的顆粒流量修正為

通過上式（2）—式（9），可以隨機生成一組噴射水滴，作為初始條件，運用式（1）計算每個水滴的運動過程。

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

（1）水滴運動阻力系數(shù) Cf采用下式表示[12]：

（2）空氣運動黏滯系數(shù) va采用下式表示[13]：

式中Tr為蘭氏溫度，oR，可表示為Tr=1.8Ta+491.67，Ta為空氣溫度，℃。

（3）空氣密度ρa與海拔高程、氣溫的關(guān)系可用下式表示：

式中E為海拔高程，m。

（4）當(dāng)?shù)貧鈮簆與海拔高程E、氣溫Ta的關(guān)系采用下式表示：

式中：p0為海平面標(biāo)準(zhǔn)氣壓，取101 325Pa；cp為空氣定壓比熱容，取1005 J/（kg·K）；M為干燥空氣的摩爾質(zhì)量，取0.02896kg/mol；R0為普適氣體常數(shù)，取8.315 J/（mol·K）。

將上式（10）—式（12）與式（1）聯(lián)立求解，可以考慮工程當(dāng)?shù)睾０闻c氣溫對于霧雨運動過程的影響。由于水舌風(fēng)與霧化降雨相互伴生，不可分割，以往的研究對于水舌風(fēng)場多采用正態(tài)分布假定或者均一化假定處理。對此，本文通過將三維水舌風(fēng)場轉(zhuǎn)換為結(jié)構(gòu)化數(shù)據(jù)，然后實時求解水滴瞬時鄰近風(fēng)速，用于閉合式（1）中水滴運動方程，從而較全面考慮河谷中復(fù)雜風(fēng)場對于雨霧運動過程的影響。

2.4 濺水模型的數(shù)值計算方法

（1）噴射時間的離散。在濺水過程中，水滴分布于整個空間，只有當(dāng)水滴到達(dá)地面時才能形成降水。為此，在本文模型中，首先，根據(jù)水滴自由拋射后在空中的最長停留時間t，確定濺水噴射總歷時T，一般應(yīng)滿足T≥2t；然后，將計算時間劃分為m個時間步長dt，假設(shè)在噴射歷時T內(nèi)的水滴總量為N=Tn，n為水滴噴射顆粒流量，則每個步長內(nèi)的水滴總量為Ni=N/m；第三，運用4階龍格-庫塔法求解每個時間步長內(nèi)水滴的運動方程，其中對于第i個步長內(nèi)的Ni個水滴，其最長飛行歷時為

（2）降雨強度的計算。降雨強度定義為一定時間內(nèi)穿過單位面積平面的總水量除以降雨歷時。為此，本文模型中，首先將噴射區(qū)域內(nèi)地面離散為小尺度的集雨網(wǎng)格，然后在每一個時間步長dt內(nèi)，判斷Ni個水滴的垂向位置，若其在n-1時刻位于地面高程以上，而n時刻位于該平面以下，則表明水滴已降落地面，此時根據(jù)其平面坐標(biāo)，將其水量計入相應(yīng)的集雨網(wǎng)格，即有Volume=Volume0+πd36，同時該水滴的飛行終止。由于采用隨機函數(shù)模擬水滴出射條件，模型需要重復(fù)M次噴射過程，當(dāng)所有M次噴射過程計算完成后，統(tǒng)計平面上每個網(wǎng)格內(nèi)積累的水滴總體積，再除以噴濺歷時T、網(wǎng)格垂直投影面積A與重復(fù)次數(shù)M，即可得到該網(wǎng)格內(nèi)時均降雨強度分布P=Volume/(M AT)。

（3）水滴運動過程中風(fēng)速與地形的影響。模型求解水滴運動方程中，需要獲取水滴鄰近風(fēng)速與地面高程，前者用于求解水滴所受風(fēng)場加速度，后者用于判斷水滴是否降落地面，從而反映兩者對于霧化降雨分布的影響。為此，本文模型中，首先采用相對成熟的Fluent軟件求解三維河谷風(fēng)場，然后通過相應(yīng)的轉(zhuǎn)換程序，將非結(jié)構(gòu)化的三維風(fēng)場與地形數(shù)據(jù)轉(zhuǎn)換為結(jié)構(gòu)化的風(fēng)場數(shù)據(jù)Wind(i , j,k)及地形數(shù)據(jù)[T o po(i , j)]，其中i、j、k為三維結(jié)構(gòu)化坐標(biāo)，為節(jié)省空間，上述數(shù)據(jù)以二進制文件存儲；隨機濺水模型運行初始，先讀取上述風(fēng)場與地形數(shù)據(jù)，并按照其結(jié)構(gòu)化坐標(biāo)存儲在相應(yīng)的數(shù)組單元中。對于空間任意水滴，若要獲取其鄰近風(fēng)速，可通過對其瞬時坐標(biāo)（x,y,z）數(shù)值取整，找到其相鄰的8個結(jié)構(gòu)化坐標(biāo)[Ii,Ji,Ki]，i=1～8。根據(jù)其結(jié)構(gòu)化坐標(biāo)，直接從對應(yīng)數(shù)組單元讀取相應(yīng)的風(fēng)速數(shù)據(jù)，然后插值得到水滴的鄰近風(fēng)速Ufki為空間插值函數(shù)。同樣地，采用二維空間插值方法可以求得水滴的鄰近地面高程。

（4）數(shù)值振蕩的抑制。從水滴運動微分方程得知，式（1）中右端項為水滴所受的運動加速度，其阻力項的量值與水滴粒徑成反比，由于采用隨機函數(shù)模擬，水滴噴射粒徑為一隨機變量，當(dāng)某一水滴粒徑趨近于0時，該水滴所受阻力加速度則趨于無窮大，水滴運動則產(chǎn)生數(shù)值振蕩。從理論上分析可知，單位時間步長內(nèi)，水滴所受的最大加速度可表示為UT為水滴終極速度，U0為水滴初始速度，dt為時間步長。在水平方向上，水滴終極速度UT=U0，即為當(dāng)?shù)仫L(fēng)速；在垂線方向上，水滴終極速度為即為重力、浮力、阻力達(dá)到平衡時速度。當(dāng)水滴初始時刻所受的運動加速度為a0≥amax，表明水滴運動速度在該時間步長中已經(jīng)達(dá)到穩(wěn)態(tài)，則在該時間步長中，令a0=amax，下一時間步長中，水滴初始速度U0=UT。

3 數(shù)學(xué)模型驗證分析

3.1 計算條件中電建昆明院科研所于2014年8月17日，針對小灣水電站泄洪洞全開泄洪過程，水舌入水形態(tài)與泄洪霧化進行了原型觀測。觀測過程中，壩上水位1236.32 m，實測泄洪流量3500 m3/s，當(dāng)?shù)睾０?000 m，氣溫20℃，觀測時間50 min。泄洪洞實地泄洪霧化情況見圖1，水舌入水形成的雨霧在泄洪風(fēng)場作用下沿兩岸爬升，高度接近300 m，通過沿岸布置霧化降雨測點，得到了下游降雨強度的實測點據(jù)。本次觀測工況中泄洪水舌入水條件如表1，水舌入水位置及下游河谷地形見圖2。

圖1 小灣電站泄洪洞全開霧化情況

3.2 水舌風(fēng)場計算結(jié)果泄洪水舌入水過程中，在河谷中形成局部風(fēng)場，對霧化降雨分布范圍產(chǎn)生影響。根據(jù)當(dāng)?shù)睾０闻c氣溫條件，由式（11）—式（13）求得當(dāng)?shù)卮髿怵禂?shù)1.6791×10-5m2/s，大氣密度1.0805 kg/m3，大氣壓強89 999.2 Pa。然后，基于Fluent軟件，求解泄洪洞水舌下游河谷風(fēng)場。計算中采用了RNGk-ε模型封閉連續(xù)與動量方程[14]，水舌入水?dāng)嗝娌捎盟俣冗吔?，周圍河谷地形按照?biāo)準(zhǔn)壁面函數(shù)進行處理，頂部為壓力出口邊界。圖3為計算得到的河谷地面附近風(fēng)速等值線圖。由圖可知，河谷泄洪風(fēng)場受到入水條件與地形的影響，在水舌入水附近風(fēng)速可達(dá)50 m/s量級，然后向外圍逐步衰減，縱向影響范圍可達(dá)下游1500 m，兩岸爬升高度超過300 m，其復(fù)雜形態(tài)對于霧化降雨分布將會產(chǎn)生影響。

3.3 降雨強度分布計算結(jié)果根據(jù)泄洪水舌入水條件、泄洪風(fēng)場以及當(dāng)?shù)貧庀髤?shù)，采用隨機濺水?dāng)?shù)學(xué)模型進一步計算河谷中霧化降雨分布形態(tài)。圖4為泄洪洞下游霧化降雨強度分布形態(tài)，由圖可知，霧化降雨分布受到風(fēng)場與地形的影響，霧化雨區(qū)在兩岸爬升高度約250 m，縱向分布范圍受到下游右岸河谷地形的阻擋，邊界位于泄洪洞下游880 m。圖5為泄洪洞下游霧化降雨強度等值線圖，為便于比較，將原型觀測點據(jù)一并匯出（圖中方框數(shù)據(jù)），結(jié)果表明兩者在分布規(guī)律上甚為吻合。

表1 泄洪洞水舌入水條件

圖2 小灣泄洪洞水舌入水位置與河谷地形等高線

圖3 小灣泄洪洞下游地面附近風(fēng)速分布

圖4 泄洪洞霧化降雨分布形態(tài)

圖5 泄洪洞霧化降雨強度分布

表2 不同海拔地區(qū)氣象因子計算結(jié)果

圖6 海拔高程50m條件下泄洪洞下游霧化降雨分布

3.4 海拔高程對于霧化降雨影響的敏感性分析小灣電站當(dāng)?shù)睾０胃叱碳s1000 m，為分析海拔高程的影響，本文在泄洪條件、河谷地形、當(dāng)?shù)貧鉁鼐蛔兊那疤嵯拢芯坎煌０胃叱?0 m、2000 m、3000 m情況下，泄洪霧化降雨分布的變化規(guī)律，上述海拔高程對應(yīng)的氣象參數(shù)見表2。圖6—圖8為不同海拔高程條件下，泄洪霧化降雨分布計算結(jié)果，結(jié)果分析表明：（1）隨著海拔高程的增加，在相同的泄洪條件下，泄洪洞下游霧化降雨區(qū)分布范圍會有所增加。當(dāng)海拔高程50 m時，兩岸雨區(qū)爬升高度約190～220 m，縱向邊界位于泄洪洞下游約850 m；而當(dāng)?shù)睾０胃叱?000 m時，兩岸雨區(qū)爬升高度可達(dá)230～280 m，縱向邊界位于泄洪洞下游950 m，若非受到右岸地形的阻擋，雨區(qū)縱向變化范圍會更大。（2）隨著海拔高程的增加，霧化降雨強度分布規(guī)律發(fā)生坦化。表3為霧化降雨強度沿河谷縱向的分布數(shù)據(jù)，由表可知，在大于400 mm/h的等值線區(qū)域，分布范圍有所減小，而在小于400 mm/h的區(qū)域，分布范圍有所增大。（3）海拔高程除了對霧化降雨分布產(chǎn)生影響以外，還會對泄洪雨霧總量以及雨霧輸運過程產(chǎn)生影響，此外氣溫條件也是重要的影響因素，對此需要進一步開展系統(tǒng)研究。

圖7 海拔高程2000m條件下泄洪洞下游霧化降雨分布

圖8 海拔高程3000m條件下泄洪洞下游霧化降雨分布

表3 泄洪霧化降雨強度沿河谷縱向的分布數(shù)據(jù)

4 結(jié)語

本文通過理論分析，建立了考慮海拔高程與氣溫條件影響的泄洪霧化隨機濺水?dāng)?shù)學(xué)模型，結(jié)合小灣電站泄洪洞霧化原型觀測數(shù)據(jù)進行了驗證，兩者數(shù)據(jù)吻合良好。在此基礎(chǔ)上，針對不同海拔高程條件下，泄洪洞下游霧化降雨分布變化規(guī)律進行敏感性分析。在泄洪條件、地形、氣溫不變前提下，隨著海拔高程的增加，泄洪洞下游霧化降雨區(qū)分布范圍會有所增加，同時霧化降雨強度分布規(guī)律發(fā)生坦化，證明海拔高程對于霧化降雨強度分布的確存在一定影響。本文研究僅針對泄洪霧化降雨過程，下一步擬對水舌激濺總源量以及雨霧輸運過程，繼續(xù)開展海拔高程與氣溫條件的影響分析。