沖擊式水輪機噴射機構泥沙磨損特性的數值模擬

葛新峰，孫 潔，李 陽，吳 丹，張 雷，化洪昌

（1.河海大學能源與電氣學院，江蘇南京 210098；2.江蘇理工學院，江蘇常州 213001；3.黃河水利科學研究院，河南鄭州 450003）

1 研究背景

在沖擊式水輪機運行的過程中，水流會攜帶一定量的泥沙通過機組。對于一些泥沙含量較高的河流而言，噴射機構和轉輪都會受到較為嚴重的磨損，從而造成射流質量下降、功率下降，運行維護成本增加，嚴重時還可能威脅到機組的安全穩定運行［1］。

對于沖擊式水輪機磨損的研究，主要有實驗研究和數值計算兩種。但是由于實驗研究的一些局限性，現在數值計算的方法更加普遍［2］。大多數學者將注意力放在了轉輪斗葉內表面磨損的研究上，關注噴射機構磨損的學者較少［3-5］；然而，噴射機構的高速射流質量，會直接影響到轉輪的動特性，因此，噴射機構也應是磨損研究應該關注的對象。Thapa［6］研究了泥沙對水力機械的沖蝕機理，認為泥沙造成磨損的根本原因是其攜帶的大量動能，泥沙在流速最高或是水流加速區域的磨損情況最為嚴重。曹永等［7］經過研究發現，泥沙顆粒對壁面的沖蝕磨損率與顆粒直徑呈一定的函數關系。Messa等［2］研究發現，噴嘴和噴針是噴射機構中最容易受到磨損的部件，且噴嘴出口處的磨損更大。曾崇濟等［1］對噴射機構自由射流的水固氣三相流進行了三維非定常模擬，得到了泥沙直徑與隨流性、磨損位置的關系。Benzon等［8］和Jo等［9］研究了噴針頂角，噴嘴收縮角對噴射機構效率以及性能的影響。

在數模計算磨損的過程中，計算結果很大程度上依賴于選用的經驗或半經驗的磨損公式。本文的計算求解是基于Fluent19.2平臺的，在Fluent中，有多種磨損模型，如Generic，Finnie，McLaury和Oka模型。其中，Generic模型為廣義磨損模型，可以允許用戶自己定制磨損模型；Finnie模型最早于1960年提出［10］，應用較廣，但往往會出現過度預測磨損的情況［11］；McLaury模型是一種預測水中砂沖蝕率的模型，該模型主要用于模擬泥漿侵蝕過程中的侵蝕速率［12］；Oka模型由Oka等在總結大量實驗數據的基礎上提出的［13-14］，之后也得到一些學者的應用。蘇佳慧等［15］對比了4種模型對彎管沖蝕的適用性研究，發現McLaury 模型也會出現過度預測磨損的現象。相比于Oka 模型，Generic 模型留給用戶的自由度更大，更便于調整參數。為此，本文磨損模型選用Generic模型。

為探究泥沙顆粒對噴射機構的影響，本文基于Fluent19.2平臺，對型號為CJC601-L-45/2X3.5的模型沖擊式水輪機單噴嘴進行三維建模，研究對象為泥沙濃度較低、顆粒之間較少發生相互碰撞、且顆粒大小分布較為均勻的含沙水流，采用VOF多相流模型模擬水氣兩相流，在兩相流非定常計算穩定之后，加入離散相，進行離散相模型和磨損模型的計算，研究泥沙的直徑、濃度對噴射機構特性的影響。

2 數學模型

2.1 VOF多相流模型VOF模型是一種固定在歐拉網格下的表面跟蹤方法，適用于跟蹤多種互不相融的流體交界面。在VOF模型中，不同的流體組分共用一套動量方程，計算時記錄各流體組分所占有的體積率。本文將水和空氣視為連續相，用VOF模型可以很好的跟蹤到射流的自由液面。內部的控制方式如式（1）—（7）［16］：

式中：ρq為第q相的物理密度；為第q相的速度；aq為第q相的體積分數；ρ為密度。

2.2 RNG k-ε 模型RNGk-ε模型在標準k-ε模型的基礎上做出改進，其計算功能更強。它的湍動能和耗散率方程如下［17］：

2.3 離散相模型Fluent 中的離散相模型適合于水流中稀疏顆粒的計算，在計算的過程中，忽略粒子間的碰撞，且顆粒體積分數必須低于10%，才能保證計算的準確性，此外，此模型也不適用于無限期懸浮的顆粒流問題。基于以上特點，Fluent中的離散相模型與本文的研究對象是相適應的。在計算過程中，離散相的慣性、曳力、重力都會在拉式公式中考慮到。由于在計算過程中，不考慮泥沙顆粒形狀對磨損的影響，將其簡化為圓球形，因此曳力模型的公式選擇球形顆粒。顆粒的平衡方程在笛卡爾坐標系下的形式為［16］：

式中：u為流體相速度，m/s；up為顆粒速度，m/s；μ為流體動力黏度，N·s/m2；ρp為顆粒密度，kg/m3；dp為顆粒直徑，m；CD為曳力系數，對于球形顆粒可取常數；Re為顆粒雷諾數，定義如下［16］：

顆粒在運動時，其周圍流體由于加速作用，會在顆粒上附加作用力，其表達式為［16］：

流場中存在的流體壓力梯度引發的附加作用力為［16］：

在計算時，考慮到湍流對顆粒隨機性的影響，還應同時開啟隨機游走模型，使顆粒與流體的離散渦相互作用。此時，流體隨機脈動所引起的瞬時速度如下［16］：

2.4 磨損模型含沙水流流過噴射機構時，會對其表面造成切削磨損。應用磨損模型可以監測到泥沙顆粒在所有壁面的磨損和沉積情況。磨損通常與壁面材料、泥沙硬度、水流速度、沖擊角度等因素有關［18-19］。本文計算選用Generic磨損模型，粒子在幾何壁面上的反射取決于顆粒與固體表面的性質，此次計算選用沙粒對碳鋼的表面反射值。沖擊角函數采用線性分段函數；壁面恢復系數采用多項式函數；顆粒的直徑函數和速度指數函數均設置為常量。這些量的具體設置值都可以在參考文獻中找到［20］。沖蝕速率可由下式定義［16］：

式中：Rerosion為泥沙的沖蝕速率，kg/m2·s，表示壁面材料在單位時間單位面積上損失的質量；為顆粒的徑粒函數；a為顆粒對壁面的沖擊角；為沖擊角函數；v為顆粒相對于壁面的速度（m/s）；代表顆粒相對速度函數；Aface為壁面面積（m2）。

3 物理模型及網格劃分

3.1 物理模型及計算設置對CJC601-L-45/2X3.5 的模型沖擊式水輪機單噴嘴進行三維建模。該模型水輪機參數如下：設計水頭為32 m，設計流量為0.038 m3/s，額定出力為10 kW。為了使數值模擬的結果更加明顯，采用小開度進口。進口直徑為125 mm，進口面積為0.012 m2，為了深入研究自由射流的流動特性，將噴嘴出口的射流域以圓柱代替。

計算設置：水和空氣視為連續相，水氣兩相流進行三維瞬態計算，時間步長設置為0.005 s。采用RNGk-ε紊流模型，開啟VOF模型捕捉水氣邊界，使用隱式算法，通過隱式體積力（implicit body force）來平衡壓力梯度和動量方程中的體積力，以提高解的收斂。求解方法采用SIMPLIC算法，離散格式選用一階迎風格式。在設置松弛因子時，將Momentum 改為0.2，以提高殘差中連續性一相的收斂性，其余值保持默認。由于非定常流動的Fr和We都比較大，重力和表面張力對主流特性的影響小［21］，因此本文在計算過程中忽略這些力的影響。

當水氣兩相流完全收斂，流場穩定之后，增加泥沙顆粒，泥沙顆粒形狀默認為球形顆粒，選擇面入射，入射方向為法向。泥沙密度為2650 kg/m3，體積分數為1%，對應的泥沙入口質量為0.263 kg/s。對水氣沙三相流進行非穩態計算。

歐拉-拉格朗日法的計算過程：水相和氣相視作連續相，在歐拉場中用統一的運輸方程計算，在拉格朗日坐標系下對流場中的每一個顆粒軌跡進行追蹤。每一步，離散相求解器都計算顆粒從當前狀態起在積分時間內的運動軌跡及動量、質量和能量損益，并更新顆粒狀態。每一個連續相時間步對顆粒進行一次更新計算，連續相迭代與離散相計算交替進行，顆粒不斷向前推進［16］。在求解過程中，連續相與離散相是單項耦合的，即連續相在當前時間步的計算結果會作為離散相下一個時間步計算時的輸入條件，從而影響離散相的分布和流場，而連續相的流動則不受泥沙顆粒存在的影響。

3.2 邊界條件水氣兩相流邊界條件及泥沙顆粒入射的邊界條件設置如圖1所示。

（1）入口邊界。水流入口邊界設置為速度入口，流速大小根據模型實驗，設置為0.857 m/s。空氣入口邊界設置為壓力進口，相對壓力為0。泥沙入口邊界亦設為速度入口，速度大小與清水相同。

（2）出口邊界。出口邊界采用壓力出口邊界，相對壓力為0。

（3）壁面。所有壁面均為無滑移固體邊界，并使用標準壁面函數法來模擬近壁面區域的流動。

圖1 模型邊界條件設置

3.3 網格劃分及網格無關性分析在ICEM中對模型采用六面體結構化網格劃分，如圖2所示，對噴嘴出口位置進行局部加密。在計算磨損的過程中，磨損量會受到網格數量的影響，因此需要進行網格無關性驗證。表1給出了模型驗證的幾種不同網格數量。

圖2 模型的六面體結構化網格

表1 網格劃分

對表1中的5個模型先進行兩相流計算，流場穩定后加入離散相繼續計算5 s，可得到如圖3所示的網格單元數與磨損率的關系。隨著網格單元數的增加，入射面的網格數也隨之增加，發射的粒子數量也相應增多，從而影響噴射機構的磨損率。隨著網格數量的增加，磨損率呈下降趨勢；當網格數增加到某一定值后，磨損率則趨于穩定。

除了對磨損率進行網格無關性驗證外，還需要Roache［22］引入提出的網格收斂指數（Grid Convergence Index，GCI）對網格的獨立性做進一步檢查。對model 1，model 2，model 3的網格進行驗證，具體驗證過程參考文獻中的方法［23-25］，驗證結果如表2。

其中，下標1 ～3分別表示model1，model2，model3三套由疏到密的網格。rk，k+1為網格細化比，為網格控制單元數；δk，k+1為不同網格之間的相對誤差，；fk取不同網格下質量流量的收斂解（kg/s）；p為收斂精度，可用不動點迭代法求解；在計算GCI時，Fs為安全因子，基于三套網格取值，Fs取1.25。

可以看出，model 2和model 3之間的網格收斂指數為0.63%，說明model 3已經滿足精度要求。結合對磨損率進行的無關性驗證，加之考慮到計算資源的問題，在本文之后的計算中，都選用model 4，即2 586 116的網格單元數進行計算，此時，進口面的網格共發射粒子數為66 000。

圖3 網格無關性驗證

表2 網格獨立性驗證（GCI）

4 結果及討論

4.1 流動數值結果分析圖4為非穩態氣液兩相流的數值模擬結果。在圖4（a）中可以看出，水流在流動到接近噴嘴出口位置時，流速開始增加；流速在噴嘴出口處達到最大值，為36 m/s；之后，水流與射流機構碰撞，流入大氣。隨著水流與空氣不斷接觸，部分空氣卷吸、混摻入高速射流中，并獲得動能，隨水流一起向前流動，使得流速降低［26］。圖4（b）為水進入射流域后的體積分數云圖，不考慮重力作用，水流呈軸對稱分布。這與Zeng等［27］的研究結果較為一致。水流的對稱性與射流質量關系密切，偏心流會使射流發散，降低機組效率［28］。圖4（c）中，進口處壓力為566 kPa，隨著水流接近噴嘴出口，壓力迅速減小。這表明在噴嘴出口位置，大部分壓能轉換成了動能。

圖4 兩相流的速度、水氣體積分數及壓力分布云圖

4.2 顆粒直徑對磨損的影響水流在流經噴嘴收縮段時會加速，高速射流中所攜帶的泥沙顆粒會撞擊噴嘴和噴針表面，造成材料剝蝕。噴針的磨損會增加水流在一定開度下的過流面積，當噴針磨損比較嚴重的時候，即使機組是關閉的，也會有水流泄出，從而對整個機組的控制系統造成負面影響。泥沙顆粒直徑與噴射機構表面的撞擊位置和磨損量有著直接的關系，因此本文選擇了幾種不同的泥沙直徑作為研究對象，探討其對噴射機構磨損的影響。

如圖5所示，圖5（a）為清水的流場流線，圖5（b）—（f）為不同泥沙顆粒直徑的軌跡線。在噴嘴進口至過渡處之前的水流段，流場流線與噴針軸線平行；加入泥沙之后，軌跡線與流線并不重合，且隨著泥沙顆粒直徑的增加，軌跡線的不平行數越來越多，說明泥沙顆粒直徑越大，隨流性越差。在過流面積減小的過渡段，由于噴嘴突然收縮，而顆粒軌跡又與流線不重合，此時，泥沙顆粒會與噴針和噴嘴表面發生劇烈碰撞，導致噴嘴收縮段更易受到磨損。在射流流出噴嘴之后，泥沙顆粒會受到慣性力和離心力的影響，隨著泥沙直徑的增加，軌跡線向著噴針軸線的彎曲傾斜程度也隨之增加。

圖5 流場流線及不同泥沙顆粒直徑的軌跡線

經過以上分析，可知泥沙顆粒與噴嘴的碰撞磨損主要發生在噴嘴的過渡段。圖6為不同泥沙徑粒下的噴嘴表面磨損率云圖。從圖6可以看出，當泥沙顆粒直徑為0.01 mm和0.05 mm時，噴嘴的磨損范圍較小；當泥沙顆粒直徑增大到0.4 mm和1 mm時，噴嘴的磨損范圍顯著增加。但是噴嘴的磨損位置都是相似的。這表明泥沙顆粒的大小雖然會影響噴嘴磨損的范圍，但對噴嘴表面的磨損位置影響較小。且噴嘴的磨損形態主要為點狀和片狀磨損。

圖7（a）為某電站噴嘴磨損實物圖［6］，可以看出，噴嘴接近出口處的位置更易受到磨損，數值模擬的結果和現場磨損的結果較為一致，驗證了數模的準確性。

圖6 不同泥沙直徑下的噴嘴磨損率云圖

圖7 某電站噴嘴及碰嘴泥沙磨損圖［6］

圖8為噴針受到不同徑粒泥沙磨損的磨損率分布云圖。當泥沙徑粒為0.01 mm時，由于細小的顆粒隨流性較強，在噴針直桿段、過渡段和噴針頭部的針尖位置均有磨損；當泥沙直徑開始增加，其隨流性越來越差，泥沙顆粒在流動過程中，會不斷撞擊噴針桿部和噴針頸部，而在噴針頭部的針尖位置，由于大顆粒受到的離心力更大，將顆粒甩離針尖表面，因此針尖部位反而不易受到磨損；隨著顆粒直徑的增加，噴針直桿段的磨損會由于顆粒的離心力使顆粒偏離而減弱，而噴針頸部過渡段的磨損則會由于占主導因素的慣性力作用而增加。因此，噴針頭部的針尖位置更容易受到細小顆粒的磨損，而噴針過渡段的磨損與顆粒直徑大小呈正相關。對于模型沖擊式水輪機而言，水頭和流速都遠小于真機，此時的數值模擬，顯示的噴針位置磨損形態主要為點狀。對比圖7（b）（c）（d），為某一高水頭電站噴射機構磨損的實物圖［6］。由于高水頭、高流速，針尖位置磨損更為嚴重，且磨損形態主要為波紋狀。

4.3 顆粒濃度對磨損的影響河流中的攜沙量常受到徑流量、干旱指數和人為活動等因素的影響，因此，研究泥沙顆粒濃度對射流機構磨損的影響有著重要的意義。

表2 不同泥沙濃度下噴射機構的磨損率

在計算穩定后的連續相中分別加入不同濃度的泥沙顆粒（泥沙的體積分數分別為1%、2%、3%，4%和5%），繼續計算0.45 s 后，可以得到如圖10 所示的噴射機構磨損率變化曲線。從圖10 可以發現，隨著泥沙濃度的增加，噴嘴和噴針受到的最大磨損率幾乎呈線性增加；結合表2，可以看出，在同一泥沙濃度下，噴嘴受到的磨損要大于噴針受到的磨損；結合圖6（c）、圖8（b）和圖9可知，泥沙的濃度增加，對泥沙顆粒的運動軌跡以及噴針表面磨損的位置影響不大。

圖8 不同泥沙直徑下的噴針磨損率云圖

圖9 D=0.05mm，泥沙濃度為5%時的泥沙顆粒軌跡線及噴針磨損云圖

圖10 噴嘴和噴針最大磨損率隨泥沙濃度的變化

5 總結

本文對沖擊式水輪機的噴射機構進行了水氣沙三相非定常計算，研究了泥沙顆粒的直徑和濃度，對噴嘴和噴針磨損位置、磨損率的影響。數值模擬的結果和現場磨損的結果較為一致，驗證了數模的準確性。本文研究為機組的檢修維護和安全穩定運行提供了參考和理論依據。主要結論如下：（1）不考慮重力，水流呈軸對稱分布。噴嘴出口處流速達到了最大值，同時壓力出現最小值，這表明噴嘴出口處的壓能大部分轉化為了動能。（2）泥沙顆粒直徑越大，隨流性越差。泥沙顆粒的大小會影響噴嘴的磨損范圍，但噴嘴表面的磨損位置是相似的，噴針頭部的針尖位置更容易受到細小顆粒的磨損，噴針過渡段的磨損程度與顆粒直徑大小呈正相關。（3）噴嘴和噴針受到的最大磨損率會隨著泥沙濃度的增加而增大；在同一泥沙濃度下，噴嘴受到的磨損要大于噴針；泥沙濃度的增加，會影響磨損率，但對泥沙顆粒的運動軌跡及噴針磨損的位置影響不大。