渦輪葉柵通道內(nèi)顆粒物沉積過程的數(shù)值模擬

楊曉軍， 祝佳雄

1.中國民航大學(xué) 航空工程學(xué)院， 天津 300300 2.中國民航大學(xué) 中歐航空工程師學(xué)院， 天津 300300

渦輪葉柵通道內(nèi)顆粒物沉積過程的數(shù)值模擬

楊曉軍1,2，*， 祝佳雄2

1.中國民航大學(xué) 航空工程學(xué)院， 天津 300300 2.中國民航大學(xué) 中歐航空工程師學(xué)院， 天津 300300

為了更加準(zhǔn)確獲得顆粒物在渦輪中的沉積分布，以某渦輪葉片為模型，選用最接近航空發(fā)動機(jī)內(nèi)部顆粒組成的Jim Bridger Power Station (JBPS)顆粒為污染物，同時，利用C++編寫合適的User Defined Function (UDF)經(jīng)過調(diào)試來分析顆粒沉積后葉片邊界的復(fù)雜變形和邊界網(wǎng)格依賴于時間變化的重構(gòu)生成，在考慮每個時間步長顆粒沉積在葉片上從而改變?nèi)~片幾何特性和換熱特性的情況下，來深入研究顆粒物沉積在葉片的整個過程，最終得出了沉積的分布情況，并且通過數(shù)值研究結(jié)果與實驗結(jié)果的對比，驗證了網(wǎng)格重構(gòu)與融合程序的合理性、準(zhǔn)確性。根據(jù)葉片變形情況預(yù)測腐蝕的發(fā)生情況。隨后，改變主流溫度、顆粒直徑來研究顆粒沉積特性。結(jié)果表明：顆粒主要沉積在葉片壓力面中部，但會使得葉片前緣和壓力面中部均產(chǎn)生明顯變形；葉片前緣由于顆粒沉積使得粗糙度增加形成鋸齒形，最先遭受腐蝕；顆粒直徑影響顆粒沉積的分布與沉積生長速度；只影響沉積速度，并不改變沉積分布。

沉積； 變形； UDF； 腐蝕； 渦輪葉片； 粒徑； 數(shù)值模擬

隨著燃?xì)廨啓C(jī)在能源供應(yīng)上的提高以及天然氣資源陸續(xù)被耗盡，尋找替代燃料開始變得非常重要。然而，替代燃料含有微量的煤灰和其他雜質(zhì)。這些雜質(zhì)會沉積在葉片和渦輪表面從而影響它們的熱交換特性，并且會降低氣膜冷卻的效率。同時伴隨著世界工業(yè)化的進(jìn)步，相當(dāng)多的地區(qū)包含了更多雜質(zhì)顆粒，從而引發(fā)嚴(yán)重的安全以及性能問題。因此，顆粒物在航空發(fā)動機(jī)內(nèi)部的沉積以及由此帶來的發(fā)動機(jī)性能的損失相關(guān)問題和所引發(fā)的環(huán)境問題[1]已經(jīng)成為國際上高度重視的一個研究課題[2-5]。

為了提高渦輪的工作效率，預(yù)測渦輪葉片的壽命，也為了優(yōu)化葉片的設(shè)計，必須深入研究葉片表面的換熱特性，在眾多影響葉片換熱特性的因素當(dāng)中，顆粒的沉積是唯一隨著葉片工作時間的增長而對葉片換熱影響愈加重要的關(guān)鍵因素，這是因為顆粒的沉積主要影響葉片表面的粗糙度，甚至?xí)斐蓪θ~片的熱腐蝕，因此研究顆粒沉積過程對于葉片設(shè)計與維護(hù)都有重要的意義[6]。

顆粒沉積與其影響在國外研究得非常早。在沉積機(jī)理的研究方面，David和Karen[7]進(jìn)行實驗研究發(fā)現(xiàn)，熔融程度不同的粒子沉積對于冷卻效果的影響可以忽略。Wright等[8]發(fā)現(xiàn)顆粒很有可能通過相互碰撞或者撞擊葉片而重新分配并且沉積在葉片上。同時，Wright通過實驗證明了處于端壁處的前緣氣膜冷卻區(qū)域更易吸引顆粒沉積。對于顆粒黏附的主要影響因素的實驗研究方面，David和Kaven[7]主要研究了冷卻孔的結(jié)構(gòu)對于粒子沉積帶來影響的抵御作用。根據(jù)實驗結(jié)果，槽式孔可以有效地降低由于顆粒沉積而損失的冷卻效率。同時，他們通過對前緣氣膜冷卻與端壁冷卻進(jìn)行實驗，得出橫向槽可以更有效地改善端壁冷卻的效率。Ahluwalia等[9]研究了煤灰在渦輪葉片上面的沉積和黏附現(xiàn)象。他們發(fā)現(xiàn)壓力面的黏著系數(shù)在0.000 3～0.11范圍內(nèi)變化并且與主燃?xì)夂捅诿鏈囟扔兄軓?qiáng)的聯(lián)系，但是沖擊角的變化對于黏著系數(shù)的影響是可以忽略的。對于沉積與腐蝕的關(guān)系以及影響沉積因素的研究方面，Crosby等[10]利用渦輪加速沉積設(shè)備(Turbine Accelerated Deposition Facility，TADF)研究顆粒的大小、氣體的溫度和沉積面的溫度如何影響沉積過程。他們發(fā)現(xiàn)當(dāng)氣體溫度跨過1 230 K 這個門檻，沉積增長會越來越快。同時發(fā)現(xiàn)，保持主流氣體溫度不變，隨著沉積面溫度的降低，沉積的粗糙度會下降(即沉積更加均勻)。

國內(nèi)方面航空發(fā)動機(jī)制造起步較晚，渦輪進(jìn)口溫度不夠高，因此研究顆粒物沉積起步較晚，關(guān)于沉積物的數(shù)值與實驗研究很少?；鸺虥_發(fā)動機(jī)研究者研究了顆粒組分對于沉積情況的影響，主要針對液體燃料。王德全等[11]研究了固沖發(fā)動機(jī)中顆粒沉積的主要分布位置，并無其他細(xì)致研究。關(guān)于航空發(fā)動機(jī)方面，對顆粒沉積生長研究甚少，國內(nèi)大多數(shù)學(xué)者研究與顆粒沉積相似的機(jī)翼結(jié)冰生長過程，以及結(jié)冰對于機(jī)翼的影響等，其對于顆粒沉積的研究有著借鑒作用。2004年，張大林和陳維建[12]應(yīng)用Navier-Stokes方程和紊流模型計算結(jié)冰表面外的空氣流場。在研究過程中發(fā)現(xiàn)數(shù)值模擬結(jié)果與實驗結(jié)果總是有相對較大的誤差，經(jīng)過仔細(xì)研究之后，發(fā)現(xiàn)是由于網(wǎng)格的重構(gòu)與融合方面沒有定義。隨后2005年，陳維建和張大林[13]提出了一種可考慮粗糙度影響的瘤狀冰結(jié)冰過程數(shù)值模擬方法。并且與國外文獻(xiàn)實驗數(shù)據(jù)對比，結(jié)果吻合較好。2011年，孫志國等[14]專門對于機(jī)翼結(jié)冰過程中二維以及三維的網(wǎng)格重構(gòu)進(jìn)行了研究，提出了網(wǎng)格重構(gòu)與融合的思想，并且將網(wǎng)格模塊NI-GRID集成于結(jié)冰計算軟件NUAA-ICE3D中，并對二維翼型和三維機(jī)翼的網(wǎng)格生成與重構(gòu)進(jìn)行了驗證。

這些研究為顆粒物沉積生長方面的探索提供了重要的理論依據(jù)和實驗數(shù)據(jù)，但是這些研究對于全面探索顆粒物沉積的過程機(jī)理還有一定的不足，有必要進(jìn)一步展開數(shù)值模擬研究來彌補(bǔ)實驗研究過程中所不能覆蓋的部分。本文嘗試?yán)糜袔缀螖?shù)據(jù)以及沉積實驗數(shù)據(jù)的葉片來進(jìn)行顆粒物在渦輪葉片方面的沉積數(shù)值研究，與實驗結(jié)果相結(jié)合檢驗數(shù)值研究的可靠性，并分析影響顆粒物沉積的因素與顆粒物沉積對于渦輪葉片的影響。

1 顆粒沉積數(shù)值模擬方法

1.1 控制方程

關(guān)于流體方面的三大方程，根據(jù)流體力學(xué)中的無量綱化結(jié)果，可以得知，在航空發(fā)動機(jī)工作的情況下，可以得到

(1)

由于渦輪內(nèi)部流動的復(fù)雜性和高速度性，在整個渦輪葉片表面溫度的梯度的分布，在微米以及亞微米級別的顆粒，其主要受4個力的作用: 阻力，Saffman升力、熱泳力(熱遷移力)以及布朗力。

這其中，在一定時間內(nèi)起主導(dǎo)作用的是阻力，隨著時間延長，Saffman升力占的比重會越來越大，而熱泳力只在溫度梯度大的地方稍微對固體的運(yùn)動有些影響，就葉柵通道全局來看，影響不大。布朗力是分子隨機(jī)運(yùn)動而產(chǎn)生的力，只會對整體的運(yùn)動有輕微的波動性影響，并不會影響固體顆粒宏觀的運(yùn)動趨勢。進(jìn)一步簡化之后，對固體顆粒列三大方程最終無量綱化的結(jié)果為

(2)

(3)

(4)

式中：St為Stokes數(shù)，其表示顆粒松弛時間和流體特征時間的比；下標(biāo)p表示顆粒，下標(biāo)conv表示對流，下標(biāo)k表示動量。

由式(2)～式(4)可以看出，固體顆粒的運(yùn)動受到Stokes數(shù)、主流溫度和主流速度的影響。其中，Stokes數(shù)主要是由顆粒的直徑來決定的，主流溫度與主流速度都是流場特性，通過改變進(jìn)口溫度與進(jìn)出口壓強(qiáng)比實現(xiàn)。

1.2 沉積模型

Tafti和Sreedharan[15]基于顆粒的黏性發(fā)展了一種沉積模型。顆粒的黏性隨著溫度改變并且可以根據(jù)煤灰顆粒的性質(zhì)進(jìn)行預(yù)測。煤灰顆粒的軟化溫度被認(rèn)為是臨界的沉積溫度，記為Ts。高于臨界溫度的顆粒100%會沉積，低于臨界溫度的顆粒有根據(jù)以下函數(shù)的沉積概率

(5)

式中：μcrit為顆粒在臨界溫度時候的黏度；μTp為熔融顆粒在目前溫度的黏度。

Senior和Srinivasachar[16]提出了一種根據(jù)煤灰顆粒的組分來計算煤灰顆粒黏性的方法。隨后，N’Dala等[17]顯示了溫度和黏性之間的關(guān)系式

(6)

式中：A和B是依賴于化學(xué)組分的常數(shù)，Tp為顆粒溫度。Senior和Srinivasachar通過實驗做出了A和B的變化曲線。更詳細(xì)的情況以及機(jī)理在文獻(xiàn)[15]中。

由于研究對象是第一級渦輪，渦輪進(jìn)口溫度較高，并且渦輪溫度場分布較速度場梯度變化更加劇烈，因此，數(shù)值模擬基于此沉積模型進(jìn)行模擬研究，并且與真實沉積結(jié)果進(jìn)行對比，來驗證在渦輪工作條件下，溫度場對于沉積的影響是主要的。

2 網(wǎng)格的重構(gòu)與融合

2.1 網(wǎng)格生成思想

首先，利用橢圓型偏微分方程的方法生成網(wǎng)格。因為橢圓形偏微分方程生成的網(wǎng)格有優(yōu)秀的光滑性，并且在橢圓偏微分方程的右側(cè)可以加入源項來控制調(diào)節(jié)生成網(wǎng)格的正交性。同時，橢圓型偏微分方程生成的網(wǎng)格質(zhì)量更加容易調(diào)節(jié)，利用Hilgenstock法源項可以有效地提高網(wǎng)格的質(zhì)量，從而保證迭代的穩(wěn)定性。

ξi=Pi

(7)

另一種方法是利用求解擴(kuò)散方程來進(jìn)行網(wǎng)格邊界的調(diào)整，擴(kuò)散光順調(diào)整即通過求解擴(kuò)散房產(chǎn)成對于網(wǎng)格的重構(gòu)進(jìn)行光順處理，從而優(yōu)化網(wǎng)格的質(zhì)量，減少產(chǎn)生負(fù)體積的危險，其擴(kuò)散方程為

(8)

式中：

Υ=1/Vα

(9)

其中：α為擴(kuò)散參數(shù)；V為正則體積，其分別代表了取決于壁面距離的擴(kuò)散光順與取決于單元體積的擴(kuò)散光順。

2.2 網(wǎng)格的重構(gòu)與融合思想

網(wǎng)格重構(gòu)和沉積模型是構(gòu)成User Defined Function (UDF)的主要部分，其編寫和選取方法對結(jié)果的影響非常大。目前發(fā)現(xiàn)主流的重構(gòu)方法主要包括以下2種：

1) 法線重構(gòu)法

將顆粒的體積分別投影到3個正方向上面，分別計算3個正方向上面的位移量，然后再將其合成為總的變形量，以二維來舉例，則可利用以下關(guān)系式進(jìn)行網(wǎng)格重構(gòu)與融合：

(10)

(11)

(12)

(13)

三維的情況只是多了z方向，此種方法依然需要對網(wǎng)格進(jìn)行光順處理，否則依然會很快出現(xiàn)負(fù)體積。

2) 順序重構(gòu)法[13]

首先，假設(shè)不存在顆粒沉積之后由于形狀而產(chǎn)生的滾轉(zhuǎn)以及其余位移形態(tài)選取貼近葉片表面的兩層網(wǎng)格為例，粗線表示葉片邊緣。虛線表示該時間步長之內(nèi)顆粒沉積的總量。對于計算網(wǎng)格邊界而言，只需要確定陰影邊界即圖1中虛線與法線部分的交點坐標(biāo)即可，這是由于顆粒沉積沿著貼面網(wǎng)格的法線方向，而結(jié)構(gòu)網(wǎng)格即六面體網(wǎng)格的邊恰好是貼面的法線方向，因此極大簡化了網(wǎng)格的重構(gòu)過程。具體的計算過程如下：

① 根據(jù)顆粒的沉積模型來設(shè)定顆粒沉積的條件。

② 計算顆粒單位時間步長內(nèi)的沉積數(shù)量和沉積質(zhì)量，隨后根據(jù)顆粒的密度求解變形量。

③ 尋找最近的沒有顆粒沉積的控制體，如圖1(b)中的控制體8。

④ 從控制體7開始計算各控制體內(nèi)顆粒沉積的量與形狀。具體方法為：根據(jù)控制體7內(nèi)的沉積量計算出顆粒在控制體7內(nèi)沉積的體積，然后從a點引射線到控制體6與7的交界處，即d點，使得Δacd與顆粒沉積的體積相等，即可求解出d點坐標(biāo)；之后的方法與計算d點相同，在控制體6中以d作射線與線段ef的交點記為g，將其類比為梯形面計算，則可計算出梯形面積，從來可以得出控制體6內(nèi)的沉積體積；

⑤ 控制體1～5內(nèi)的顆粒沉積體積與控制體6中的計算方法相同。

圖1 等體積順序重構(gòu)法 Fig.1 Equal-volume sequence method

2.3 重構(gòu)方法的比較與選擇

順序重構(gòu)方法的優(yōu)點在于網(wǎng)格的重構(gòu)有著順序性，因此在編輯UDF的工作量上面要小很多，邏輯光面的思考量比較小，同時，這種順序方法重構(gòu)的網(wǎng)格對于負(fù)體積的接納能力要優(yōu)于法線重構(gòu)方法。但是它有著一個致命的缺陷，經(jīng)過計算發(fā)現(xiàn)，此種方法的計算速度慢，并且對處理器和內(nèi)存的要求高。

同順序重構(gòu)方法相比，法線重構(gòu)方法屬于非結(jié)構(gòu)方法，F(xiàn)luent本身就是非結(jié)構(gòu)計算軟件，不適用規(guī)定順序的結(jié)構(gòu)化計算，這種方法首先不需要尋找出現(xiàn)沉積量較少的網(wǎng)格，更不需要再去尋找臨近網(wǎng)格進(jìn)行計算，從而節(jié)省了給沉積網(wǎng)格排序的時間，在網(wǎng)格數(shù)量達(dá)到一定程度的情況下，這2種方法的計算結(jié)果是一致的，并且法線重構(gòu)方法比較節(jié)約時間和計算內(nèi)存。它也有自身的缺陷，即法線重構(gòu)法對于負(fù)體積的接納能力要差于順序重構(gòu)方法，不適合長時間的或者濃度很大的顆粒沉積的數(shù)值模擬。

綜合起來，根據(jù)航空發(fā)動機(jī)的工作狀況，顆粒直徑和顆粒濃度都非常小，而且每經(jīng)過1 000飛行小時基本上就會進(jìn)行洗發(fā)工作，因此適合使用法線重構(gòu)方法來模擬顆粒的沉積過程。

3 沉積過程的模擬與驗證

3.1 模型與網(wǎng)格

首先，建立三維葉片模型，如圖2所示。進(jìn)行邊界條件設(shè)置，隨后設(shè)置并且編譯UDF進(jìn)入Fluent，進(jìn)行模擬并且與實驗結(jié)果對比從而驗證數(shù)值模擬和UDF的正確性。葉片弦長為120 mm，葉高為30 mm，葉片材料使用鋁合金，其密度為 2 719 kg/m3，比熱容為871 J/(kg·K)，其物理屬性見文獻(xiàn)[19]。計算域包括葉柵通道，前后延長大于一倍弦長。

圖2 葉片三維模型圖 Fig.2 3D modeling schema of blade

由于葉片表面是個曲面，因此用近似六面體的網(wǎng)格，并不完全是六面體，在壁面邊界(左右壁面)第1層使用非結(jié)構(gòu)化網(wǎng)格，其余網(wǎng)格均使用六面體網(wǎng)格，主要保證葉片表面的網(wǎng)格正交性比較好就可以了，計算網(wǎng)格如圖3所示。

圖3 葉片計算網(wǎng)格 Fig.3 Computing grid for blade

網(wǎng)格的單元總數(shù)為170 880，為了延長負(fù)體積的時間，可以進(jìn)一步進(jìn)行優(yōu)化，通過mesh/repair-improve/improve-quality指令進(jìn)行進(jìn)一步的優(yōu)化，此指令在優(yōu)化正交性的同時會降低扭曲度。

采用歐拉-拉格朗日方法進(jìn)行描述，選取離散相模型(Discrete Phase Model，DPM)，因為顆粒的Stokes數(shù)變化范圍比較大，歐拉模型和混合物模型都只適用于Stokes數(shù)一定范圍內(nèi)的多相模擬。關(guān)于網(wǎng)格的重構(gòu)與融合部分，利用法線重構(gòu)法，在UDF里面對于變形進(jìn)行第一次光順，完成UDF輪詢之后，在時間步長結(jié)尾對網(wǎng)格進(jìn)行第2次光順，此次光順通過求解擴(kuò)散方程來進(jìn)行光順。

顆粒選擇為Jim Bridger Power Station(JBPS)顆粒，其組成成分如表1所示，其顆粒的特性如表2所示。

表1 JBPS顆粒組成Table 1 JBPS composition

表2 JBPS顆粒特性Table 2 JBPS properties

3.2 軌跡模擬

在Fluent中，可以畫出顆粒的軌跡，在同樣進(jìn)出口溫度和進(jìn)出口壓比的條件下，分別取射入直徑d為1，10和100 μm的顆粒，觀察它們的曲線，如圖4所示。

圖4中，紅色軌跡代表100 μm顆粒，綠色代表10 μm顆粒，藍(lán)色代表1 μm顆粒，由圖中可以看出，隨著顆粒直徑的增大，在遇到葉片壓力面的情況下，主流體沿著壓力面表面流走，而顆粒卻由于慣性使得運(yùn)動狀態(tài)不能及時改變，因此撞擊在壓力面上面。同時，由于主流體與顆粒之間有了相對速度，這種相對速度相對于顆粒來說是向下的，根據(jù)動量定理，此時顆粒會受到向下的作用力，因此，這也在一定程度上促使顆粒撞擊在葉片表面上。比如紅色顆粒由于直徑過大，因此撞擊區(qū)域在葉片前緣，進(jìn)行一次反彈后再次撞擊到相鄰葉片的前緣和中部，再次進(jìn)行反射，由于葉片幾何構(gòu)型的原因，才使得其順利流走，同時，10 μm的綠色顆粒同樣撞擊在葉片尾緣區(qū)域，而1 μm的藍(lán)色顆?；旧暇透S著主流運(yùn)動，不與葉片進(jìn)行接觸。隨著顆粒直徑的減小，慣性力的作用就會越來越小，顆粒傾向于與主流體一起順著壓力面流走，因此，小直徑的顆粒更不容易撞擊在葉片上，連葉片都接觸不到，沉積下來的小直徑顆粒自然更低。

圖4 顆粒軌跡示意圖 Fig.4 Schematic of particle tracking

3.3 網(wǎng)格變形模擬

進(jìn)口壓強(qiáng)約選為1.8 MPa，出口為1 MPa，進(jìn)口溫度選為1 650 K，顆粒的入射速度為81 m/s，沉積模型選取臨界黏度模型，A與B如表2所示，模擬沉積生長情況從而得到葉片的變形情況，如圖5所示。

從圖5中可以看出，編譯過UDF之后的Fluent可以模擬在沉積情況下的葉片變形情況，前緣的沉積大大增加了葉片的粗糙度，使得葉片形成鋸齒狀，而壓力面中部均勻沉積，因此粗糙度改變并不大，只是增加了厚度。這種現(xiàn)象與葉片的幾何構(gòu)型和溫場分布有關(guān)，在前緣溫場梯度較密，因此沉積量不均勻，而到了葉片壓力面中部，溫度場分布較為均勻，葉片開始變得接近平板，因此更接近均勻沉積。

圖5 葉片變形情況 Fig.5 Deformation of blade

3.4 葉片沉積生長模擬與驗證

進(jìn)口壓強(qiáng)約選為1.8 MPa，出口為1 MPa，進(jìn)口溫度選為1 650 K，顆粒的入射速度為81 m/s，沉積模型選取臨界黏度模型，A與B如表2所示，圖6和圖7表示壓力面顆粒的沉積變形量與沉積質(zhì)量。

由圖6與圖7可以看出，葉片發(fā)生了明顯變形，沉積變形量分布基本上與沉積質(zhì)量保持一致，對于葉片壓力面來說，沉積使得葉片呈現(xiàn)橫向條紋式的變形，在壓力面中部，也就是直徑10～100 μm 顆粒的主要沉積區(qū)域，沉積使得縱向呈條紋式的變形，并且其隨著時間的推進(jìn)程度呈現(xiàn)拋物線形式從葉根與葉冠向著葉片中部發(fā)展。

圖7中，對于直徑在10～100 μm之間的顆粒，其沉積主要從葉片壓力面的中部進(jìn)行生長，逐漸向葉片前緣和尾緣生長，并且由于顆粒直徑的原因，前緣和尾緣很難造成大面積的沉積，這是因為在尾緣區(qū)域葉片壓力面法線與顆粒速度方向角度變小，使得大直徑的顆粒更加容易撞擊進(jìn)而反彈，前緣法線幾乎與顆粒入射速度角度重合，因此更加容易造成顆粒的反彈，這是由于顆粒的直徑變大，慣性力和動量也會增大，其黏性不足以平衡反作用力，從而導(dǎo)致顆粒無法沉積，只能反射。

圖6 葉片變形量 Fig.6 Blade deformation

圖7 沉積質(zhì)量 Fig.7 Deposition mass

圖8 實驗結(jié)果與數(shù)值模擬對比 Fig.8 Comparison of experiment and CFD simulation

圖8中的沉積實驗是在TADF裝置中進(jìn)行的快速沉積實驗，實驗葉片為CFM56系列高壓渦輪葉片，其入射顆粒的流量是200 μg/s，基本上接近于北京地面的PM2.5的濃度，入射顆粒為JBPS顆粒，入射粒徑估算為1～30 μm，進(jìn)口壓強(qiáng)為1.8 MPa，出口壓強(qiáng)為1 MPa，入射溫度為1 500 K，作用3 h。

實驗與數(shù)值研究所選用的沉積物顆粒組分完全相同，但是其質(zhì)量流量增大了萬倍；另一方面，由于真實的入射顆粒直徑無法控制，因此只能通過統(tǒng)計學(xué)理論來估算平均顆粒直徑為10 μm，其真實入射顆粒直徑為正態(tài)分布，期望值為10 μm。

至于沉積量的測量方面，在實驗中測量沉積量和沉積速度是非常難辦的，往往采用估算的方法，誤差非常大，因此主要對比實驗與數(shù)值研究的沉積分布，只要分布一致，即以相同條件下數(shù)值研究量化分析來彌補(bǔ)實驗不足。

從圖8可以看出，數(shù)值模擬和實驗結(jié)果[20]的沉積分布量非常一致，均分布在葉片壓力面中部(數(shù)值研究沉積程度為：藍(lán)<綠<紅)，其次是前緣，尾緣沉積非常少。并且實驗與數(shù)值研究壓力面中部的沉積變形均為縱向擴(kuò)展，而前緣部分變形為橫向鋸齒。在前緣變形的體現(xiàn)上面，實驗與數(shù)值計算均不能明顯體現(xiàn)出變形，這說明前緣沉積量至少與壓力面中部的沉積量積累差一個等級。

數(shù)值計算相對于實驗的優(yōu)勢即在于可以完美控制顆粒直徑，并不需要利用統(tǒng)計學(xué)估計出平均直徑；另外，數(shù)值研究可以無休止地計算下去，而不僅僅是實驗所能進(jìn)行的幾小時的沉積實驗，這避免了實驗單方面增大太多濃度而造成的誤差。

圖8驗證了在氣體剛剛從燃燒室進(jìn)入渦輪中，其氣動速度比較低，并沒有經(jīng)過多級渦輪加速，并且同時發(fā)動機(jī)內(nèi)部的雜質(zhì)顆粒直徑比較小，根據(jù)Tafti[15]的研究，在較小顆粒直徑的沉積模擬方面，臨界速度模型與臨界黏度模型的沉積效率相差不多，因此為研究燃燒室沉積而利用的臨界黏度模型完全可以勝任前幾級渦輪的沉積生長與分布數(shù)值研究，并且與真實的沉積情況有著驚人的一致性。

3.5 腐蝕預(yù)測

首先，由圖9可以看出，前緣的沉積呈現(xiàn)鋸齒狀，這樣會大大增大葉片前緣的粗糙度，從而增大顆粒的碰撞概率，最終加劇顆粒的磨蝕作用，使得葉片前緣最先遭受腐蝕。

在圖8實驗的基礎(chǔ)上，利用同樣的葉片進(jìn)行長時間的沉積實驗，將顆粒入射流量提高到200 mg/s，

溫度提高至1 700 K，其條件仍然與圖8相同，提高顆粒流量與主流溫度可以使得顆粒盡快沉積，造成腐蝕的時間縮短，實驗6 h之后的結(jié)果如圖10 所示。

圖9 葉片前緣變形 Fig.9 Deformation of blade leading edge

圖10 沉積分布與腐蝕實驗結(jié)果 Fig.10 Experimenal results of deposition and errosion

下面將葉片沉積的模擬與葉片腐蝕的實驗結(jié)果進(jìn)行對比：

從圖9中可以看出，X和Y分別為基準(zhǔn)坐標(biāo)軸，沉積會使得葉片表面尤其是前緣粗糙度增加，而腐蝕最嚴(yán)重的區(qū)域恰恰是沉積變形量分布最不均勻的區(qū)域，即粗糙度最大的區(qū)域葉片前緣，這是與腐蝕的形成機(jī)理有著重要關(guān)系的。正是因為粗糙度的增加，使得腐蝕速度大大加快，因此前緣和尾緣角落處最先受到腐蝕。腐蝕形成的原因是多方面的，其中比較重要的兩個方面是熔融顆粒物沉積并且滲透到葉片涂層中和顆粒的磨蝕撞擊，這2個腐蝕的主要形成原因造成了這種腐蝕與多步長法計算的沉積分布有著驚人相似這種結(jié)果。

從沉積變形圖與腐蝕實驗圖的對比來看，沉積不均勻的區(qū)域預(yù)示著腐蝕嚴(yán)重的區(qū)域，因此，沉積變形可以預(yù)測腐蝕的區(qū)域。

4 顆粒沉積特性研究

4.1 粒徑對于沉積的影響

顆粒直徑影響粒子的Stokes數(shù)，粒子的直徑越大，其對應(yīng)的Stokes數(shù)也越大，而Stokes數(shù)是影響沉積分布與沉積速度的重要因素。其余邊界條件同3.4節(jié)相同，改變粒徑，從而研究顆粒直徑對于沉積分布與沉積生長速率的影響。對于各直徑區(qū)間內(nèi)部顆粒在葉柵通道內(nèi)的流動特性3.2節(jié)已經(jīng)進(jìn)行了分析，選用直徑1～10 μm與直徑10～100 μm的顆粒進(jìn)行數(shù)值模擬。

圖11表示進(jìn)行360 s沉積之后，直徑為10 μm 的顆粒與直徑為100 μm的顆粒沉積情況的對比，橫坐標(biāo)表示從壓力面前緣至尾緣(0～1)，從圖11中可以看出，小直徑的沉積分布更傾向于葉片尾緣，而大直徑顆粒的沉積分布有向葉片前緣移動的趨勢，在直徑增大10倍的情況下已經(jīng)移動到了葉片中部。如果使用更大直徑的顆粒進(jìn)行沉積生長的數(shù)值與實驗研究會發(fā)現(xiàn)，隨著顆粒直徑的增大，沉積區(qū)域會從葉片尾緣向葉片前緣移動，這正是顆粒的運(yùn)動機(jī)理造成的結(jié)果。

圖12還表明了顆粒直徑對于沉積速度也是有很大影響的，同樣進(jìn)行360 s的沉積，10 μm的顆粒沉積速度為5×10-10kg/s，而100 μm顆粒沉積速度為2×10-6kg/s，并且100 μm顆粒沉積密度遠(yuǎn)大于10 μm顆粒。從1～100 μm選取多組不同粒徑進(jìn)行模擬。

以Stokes數(shù)為100的顆粒沉積速度為標(biāo)準(zhǔn)進(jìn)行無量綱化，即可得到一個無量綱數(shù)，表征為Stokes數(shù)小于100的顆粒沉積效率，發(fā)動機(jī)的內(nèi)的顆粒Stokes一般不超過10。從圖12可以看出，直徑的增大會使得沉積速度增大，從而導(dǎo)致沉積效率增加，當(dāng)顆粒直徑增大到一定程度，由于慣性力的作用，顆粒會100%沉積下來，達(dá)到最大沉積速度，之后不再改變。

圖11 粒徑對沉積的影響 Fig.11 Effects of diameter on deposition

圖12 Stocks數(shù)對沉積效率的影響 Fig.12 Effects of Stocks number on deposition efficiency

4.2 主流溫度對沉積的影響

保持其余邊界條件不變，與3.4節(jié)相同，改變主流溫度，從而研究主流溫度對于沉積的影響，將主流溫度1 200 K提高1 600 K，每隔100 K進(jìn)行一次研究，將各個溫度點對應(yīng)的沉積速度均除以1 800 K的沉積速度得到無量綱沉積速度，得到如下結(jié)果：

主流的改變不影響沉積的分布，顆粒沉積的分布如圖13所示，其變化的只是沉積速度的大小，換而言之，主流溫度只影響沉積的生長速度，主流溫度在一定范圍內(nèi)升高，沉積速度加快；根據(jù)顆粒的性質(zhì)以及入射速度的不同，存在一個臨界主流溫度，超過臨界溫度之后主流溫度的提高對沉積無影響。

圖13 主流溫度對沉積的影響 Fig.13 Effects of mainflow temperature on deposition

5 結(jié) 論

1) 基于臨界黏度模型與法線網(wǎng)格重構(gòu)法編寫了顆粒沉積計算與葉片變形方面的程序模塊，將其用UDF嵌入Fluent數(shù)值模擬，并且與國外已有實驗結(jié)果進(jìn)行對比，驗證本文計算方法的合理性。

2) 由沉積使得葉片變形的結(jié)果來預(yù)測了葉片腐蝕的情況，同時進(jìn)行關(guān)聯(lián)原因分析，并且與國外已有腐蝕方面的實驗進(jìn)行對比，從而驗證了預(yù)測方案的可行性。

3) 顆粒直徑既影響顆粒分布區(qū)域也影響顆粒沉積速度，這是由于顆粒直徑影響其受力狀態(tài)，慣性力在顆粒沉積方面有著很大影響。

4) 主流溫度對顆粒沉積分布無影響，與沉積速度正相關(guān)，其通過改變顆粒的物性從而對沉積速度與效率進(jìn)行影響。

[1] DONKELAAR A, MARTIN R, BRAUER M, et al. Global estimates of ambient fine particulate matter concentrations from satellite-based aerosol optical depth: Development and application[J]. Environmental Health Perspectives， 2010, 118(6): 847-855.

[2] BRUN K, NORED M， KURZ R. Particle transport analysis of sand ingestion in gas turbine engines[J]. Journal of Engineering for Gas Turbines and Power, 2012, 134(1): 012402.

[3] DUNN M G, PADOVA J, MOLLER J E， et al. Performance deterioration of a turbofan and a turbojet engine upon exposure to a dust environment[J]. Journal of Engineering for Gas Turbines and Power, 1987, 109(3): 336-343.

[4] KIM J, DUNN M G， BARAN A J， et al. Deposition of volcanic materials in the hot sections of two gas turbine engines[J]. Journal of Engineering for Gas Turbines and Power, 1993, 115(3): 641-651.

[5] SCHNEIDER O, DOHMEN H J, BENRA F K， et al. Investigations of dust separation in the internal cooling air system of gas turbines： GT-2003-38293[R]. New York: ASME, 2003.

[6] EKKAD S, HAN J C. Detailed heat transfer distributions on a cylindrical model with simulated TBC spallation: AIAA-1997-1595[R]. Reston: AIAA, 1997.

[7] DAVID G, KAREN A. Simulating particle deposition and mitigating deposition degradation effects in film cooled turbine sections: DE-FC21-92MC29061[R]. Washington, D.C.: UTSR, 2010.

[8] WRIGHT I G, LEVENS C, PINT B A. An analysis of the potential for deposition, erosion, or corrosion in gas turbines fueled by the products of biomass gasification or combustion: GT-2000-0019[R]. New York: ASME, 2000.

[9] AHLUWALIA R K, IM K H, WENGLARZ R A. Flyash adhesion in simulated coal-fired gas turbine environment[J]. Journal of Engineering for Gas Turbines and Power, 1989, 111(4): 672-678.

[10] CROSBY J M, LEWIS S, BONS J P, et al. Effects of temperature and particle size on deposition in land based turbines[J]. Journal of Engineering for Gas Turbines and Power, 2008, 130(5): 819-825.

[11] 王德全, 夏智勛, 胡建新, 等. 固沖發(fā)動機(jī)沉積數(shù)值模擬與試驗研究[J]. 固體火箭技術(shù)， 2009， 32(1)： 38-42.

WANG D Q, XIA Z X, HU J X, et al. Numerical simulation and experimental study of deposition for ducted rockets[J]. Journal of Solid Rocket Technology, 2009, 32(1): 38-42 (in Chinese).

[12] 張大林, 陳維建. 飛機(jī)機(jī)翼表面霜狀冰結(jié)冰過程的數(shù)值模擬[J]. 航空動力學(xué)報, 2004, 19(1): 138-141.

ZHANG D L, CHEN W J. Numerical simulation of the rime ice accretion process on airfoil[J]. Journal of Aerospace Power, 2004, 19(1): 138-141 (in Chinese).

[13] 陳維建, 張大林. 飛機(jī)機(jī)翼結(jié)冰過程的數(shù)值模擬[J]. 航空動力學(xué)報, 2005, 20(6): 1010-1017.

CHEN W J, ZHANG D L. Numerical simulation of ice accretion on airfoils[J]. Journal of Aerospace Power, 2005, 20(6): 1010-1017 (in Chinese).

[14] 孫志國, 朱程香, 付斌, 等. 結(jié)冰計算中翼面網(wǎng)格重構(gòu)研究[J]. 航空學(xué)報, 2011, 32(2): 231-241.

SUN Z G, ZHU C X, FU B,et al. Research of grid reconstruction on wing surface for icing calculation[J]. Acta Aeronautica et Astronautica Sinica, 2011, 32(2): 231-241 (in Chinese).

[15] TAFTI D K, SREEDHARAN S S. Composition dependent model for the prediction of syngas ash deposition with the application to a leading edge turbine vane, GT-20110-23655[C]//ASME Turbo Expo 2010: Power for Land, Sea, and Air. New York: ASME, 2010: 615-626.

[16] SENIOR C J, SRINIVASACHAR S. Viscosity of ash particles in combustion systems for prediction of particle sticking[J]. Energy and Fuels, 1995, 9(2): 277-283.

[17] N’DALA I, CAMBIER F, ANSEAU M R, et al. Viscosity of liquid feldspars. Part I: Viscosity measurements[J]. Transactions & Journal of the British Ceramic Society, 1984, 83(4)： 105-107.

[18] SHIH T I P, BAILEY R T. GRID2D/3—A computer program for generating grid systems in complex-shaped two- and three-dimensional spatial domains: NASA-TM-102453[R]. Washington, D.C.: NASA,1990.

[19] Allegheny Ludlum Corporation Company. Stainless steel types 309 and types 310[S]. Pittsburgh, PA: Allegheny Lulum Corporation Company, 2002.

[20] WEBB J B, CASSADAY B. Coal ash deposition on nozzle guide vanes: Part II—Computational modeling[J]. Journal of Turbomachinery, 2011, 135(1): 1757-1767.

(責(zé)任編輯： 張晗)

URL：www.cnki.net/kcms/detail/11.1929.V.20161221.1520.004.html

Numericalsimulationofparticledepositionprocessinsideturbinecascade

YANGXiaojun1,2,*,ZHUJiaxiong2

1.CollegeofAeronauticalEngineering,CivilAviationUniversityofChina,Tianjin300300,China2.Sino-EuropeanInstituteofAviationEngineering,CivilAviationUniversityofChina,Tianjin300300,China

AnumericalsimulationisconductedtoobtainmoreaccuratedistributionofparticledepositiononturbinebladewithJimBridgerPowerStation(JBPS).C++isusedtocodesuitableUserDefinedFunction(UDF),whichreconstructsgridstoaccomdatecomplexboundarydeformationineverytimestep.Insuchawaythatdepositioncanchangegeometryfeaturesandheattransfercharacteristics,andthewholedepositionprocessonbladecanbeanalyzedandtheexactdepositiondistrbutioncanbeobtained.Afaircomparisonofnumericalresultswithexistingexperimentalresultsshowstherationalityandaccuracyofgridreconstructioncode.Onthebasisofdepositiondistributionandgeometrydeformationofturbineblade,theerosionareaonthebladeispredicted.Thereafter,theeffectofmainstreamtemperatureandparticlediameterisanalyzed.Theresearchresultsshowthatthemaindepositionoccursonthecentralareaofbladepressureside,whereasthereisobviousdeformationinleadingedgeandcentralarea.Theobviousdeformationrevealsthatroughnessisincreasedanderosionismoresevereonleadingedge.Thediameterofparticlehasimportantinfluenceonparticledepositionanddepositionrate,whereastemperatureofmainstreamonlyhasinfluenceondepositionrate.

deposition;deformation;UDF;erosion;turbineblade;particlediameter;numericalsimulation

2016-06-14;Revised2016-07-18;Accepted2016-08-15;Publishedonline2016-12-211520

s：JointFundsoftheNationalNaturalScienceFoundationofChinaandCivilAviationAdministrationofChina(U1633113);TianjinResearchProgramofApplicationandAdvancedTechnology(14JCQNJC06800);theFundamentalResearchFundsfortheCentralUniversities(ZXH2012H004)

.E-mailxiaojunyoung@hotmail.com

2016-06-14；退修日期2016-07-18；錄用日期2016-08-15； < class="emphasis_bold">網(wǎng)絡(luò)出版時間

時間：2016-12-211520

www.cnki.net/kcms/detail/11.1929.V.20161221.1520.004.html

國家自然科學(xué)基金委員會與中國民用航空局聯(lián)合資助項目 (U1633113)； 天津市應(yīng)用基礎(chǔ)與前沿技術(shù)研究計劃 (14JCQNJC06800)； 中央高?；究蒲袠I(yè)務(wù)費(fèi) (ZXH2012H004)

.E-mailxiaojunyoung@hotmail.com

楊曉軍， 祝佳雄. 渦輪葉柵通道內(nèi)顆粒物沉積過程的數(shù)值模擬J. 航空學(xué)報,2017,38(5):120530.YANGXJ,ZHUJX.NumericalsimulationofparticledepositionprocessinsideturbinecascadeJ.ActaAeronauticaetAstronauticaSinica,2017,38(5):120530.

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2016.0237

V231.1

A

1000-6893(2017)05-120530-12