劉 鑫,邢玉明,杜 佶
(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院,北京100191)
飛機(jī)結(jié)冰是指飛機(jī)在飛行時(shí)其外表面上水分積聚凍結(jié)成冰的現(xiàn)象。它影響飛機(jī)的氣動外形、飛行品質(zhì)及飛行安全,給飛機(jī)飛行帶來了極大的危害。如:飛機(jī)風(fēng)擋結(jié)冰會影響駕駛員視野;測溫、測壓感頭結(jié)冰會導(dǎo)致儀表指示失真;機(jī)(尾)翼結(jié)冰將影響氣動外形,增大飛機(jī)阻力,減小升力,影響全機(jī)操縱性、穩(wěn)定性;螺旋槳、直升機(jī)旋翼結(jié)冰會降低升力效率;軸流式壓氣機(jī)進(jìn)氣部件結(jié)冰能使發(fā)動機(jī)熄火等。據(jù)統(tǒng)計(jì),大約有9%的飛行事故是由結(jié)冰造成的。長期以來,結(jié)冰機(jī)理和飛機(jī)的防/除冰系統(tǒng)一直是國內(nèi)外飛機(jī)系統(tǒng)設(shè)計(jì)的重要研究課題[1]。
結(jié)冰源于云層中的過冷水滴,即低于0℃的液態(tài)水滴,過冷水滴是極不穩(wěn)定的,受到擾動后就可能轉(zhuǎn)化為冰。當(dāng)過冷水滴撞擊到?jīng)]有采取防冰措施的飛行器表面后,它可能立即結(jié)成明冰或向后溢流再結(jié)成明冰。研究中,可以通過冰風(fēng)洞試驗(yàn)、飛行試驗(yàn)和數(shù)值計(jì)算等途徑獲取結(jié)冰量或結(jié)冰區(qū)域。但由于試驗(yàn)花費(fèi)昂貴、周期長且受限于設(shè)備條件,很難模擬結(jié)冰包線中一些特殊狀態(tài)。因此美法英等國都在積極研發(fā)飛機(jī)結(jié)冰的數(shù)值計(jì)算方法和軟件。
進(jìn)行結(jié)冰計(jì)算的第一步需要計(jì)算水滴在模型表面上的撞擊特性和收集系數(shù)。撞擊特性主要通過兩種方法來計(jì)算:一種是Lagrangian方法;另一種是Eulerian方法。多數(shù)結(jié)冰計(jì)算軟件采用Lagrangian方法,在所得流場計(jì)算結(jié)果的基礎(chǔ)上計(jì)算過冷水滴的運(yùn)動軌跡,用來確定水滴在模型表面上的撞擊點(diǎn),從而進(jìn)一步計(jì)算水滴收集系數(shù)。Eulerian法由于在計(jì)算主流場的同時(shí),也計(jì)算得到每一個(gè)網(wǎng)格內(nèi)過冷水滴的質(zhì)量系數(shù)(或體積系數(shù)),因此不需要計(jì)算過冷水滴在流場中的軌跡。Lagrangian方法對二維或幾何形狀簡單的表面應(yīng)用起來比較簡便,但對多段翼型或發(fā)動機(jī)進(jìn)口部件等三維復(fù)雜外形,由于粒子釋放位置較難確定,很難準(zhǔn)確確定撞擊極限,且計(jì)算時(shí)用時(shí)較長。而Eulerian方法則不涉及這些問題[2-4]。
本文應(yīng)用Eulerian方法計(jì)算水滴收集系數(shù),采用Fluent商業(yè)軟件計(jì)算空氣流場,即將水滴視為擬流體通過Eulerian方法實(shí)現(xiàn)水滴流場的求解。
在Eulerian方法中引入容積分?jǐn)?shù)α的概念。它表示單位體積的混合流體中某一相所占的體積,第q相的容積分?jǐn)?shù)可以表示為αq,q相的體積可以定義為
在結(jié)冰條件下,由于水滴的容積分?jǐn)?shù)較小,一般在10-6量級,可認(rèn)為空氣和水滴之間是單向耦合的,即空氣流場影響水滴,而水滴對空氣流場的影響可以忽略[5-6]。這使得空氣相和水滴相可以分離求解:先得到空氣流場,再使用空氣流場結(jié)果計(jì)算水滴相。
對Eulerian法中空氣和水滴的控制方程可以進(jìn)行以下假設(shè):
(a)過冷水滴在云層中均勻分布,且以球形存在,不凝聚、不分解、不變形、不破裂;(b)過冷水滴在云層中,不與空氣進(jìn)行熱、質(zhì)交換,即過冷水滴不蒸發(fā)、物性參數(shù)不變化;(c)忽略作用在水滴上的空氣阻力對空氣動力的影響;(d)作用在水滴上的作用力僅有阻力和重力;(e)作用在水滴上的阻力是穩(wěn)態(tài)的;(f)空氣的紊流脈動不影響水滴的運(yùn)動。
過冷水滴穩(wěn)態(tài)控制方程如下[5]:
其中,u是水滴的速度矢量,ua是空氣的速度矢量,F(xiàn)是作用于水滴上的除阻力之外的力。K為空氣-水滴交換系數(shù),由下式給出:
式中,μa為空氣的動力黏度,dd是水滴直徑,f為阻力函數(shù),對于球形物體的阻力函數(shù)可以用下式進(jìn)行計(jì)算:
式中,CD為阻力系數(shù),Re為相對雷諾數(shù),其具體計(jì)算公式如下:
式中,ρa(bǔ)為空氣密度。
本文借助Fluent軟件計(jì)算空氣流場,其中進(jìn)口采用速度進(jìn)口邊界條件,出口采用壓力出口邊界條件,無滑移壁面邊界條件,使用k-e湍流模型計(jì)算空氣流場。將水滴設(shè)為擬流體,通過UDS模塊求解水滴相控制方程。
水滴相的入口采用速度入口邊界條件,速度大小與空氣速度相等,出口采用壓力出口邊界條件。對于水滴相的壁面邊界條件需要特殊處理[7],固體壁面對于空氣來說只是普通的固壁,但是對于水滴相來說卻相當(dāng)于一個(gè)單向出口,這個(gè)出口只容許水滴流出計(jì)算域而不容許流入計(jì)算域。這是由于水滴在撞擊到壁面時(shí)會積聚在壁面上而不會從壁面上繞流經(jīng)過,從而再一次進(jìn)入計(jì)算域中,故而水滴在撞擊到壁面時(shí)只能流出計(jì)算域而不能再次進(jìn)入。水滴在壁面附近的運(yùn)動情況如圖1所示,當(dāng)水滴撞擊到壁面上時(shí),即水滴速度ud與壁面微元單位法向向量點(diǎn)積小于等于零時(shí),水滴的體積分?jǐn)?shù)與水滴速度保持不變;在非撞擊區(qū)域,壁面微元單位法向量與水滴速度ud點(diǎn)積大于零時(shí),水滴體積分?jǐn)?shù)取零,水滴速度取下游網(wǎng)格單元速度值。
圖1 水滴相壁面邊界條件
局部水滴收集系數(shù)β為當(dāng)?shù)乇砻嫖⒃獙?shí)際水滴收集量與最大可能收集量的比值,其計(jì)算公式為:
式中,Wβ為實(shí)際水滴收集量,Wβmax為最大可能水滴收集量,ud為當(dāng)?shù)厮嗡俣仁噶?,n為部件表面的法向方向矢量,αx為當(dāng)?shù)厮误w積分?jǐn)?shù),V0為水滴初始速度,ρ為水滴密度,α0為初始水滴體積分?jǐn)?shù),可由遠(yuǎn)場液態(tài)水含量計(jì)算得出。
本文以NACA0012機(jī)翼作為計(jì)算模型討論不同氣象條件及飛行速度對水滴撞擊特性的影響。計(jì)算模型中,翼型弦長1 000mm,在網(wǎng)格劃分過程中,在圓柱周圍劃分10倍遠(yuǎn)場網(wǎng)格[8],遠(yuǎn)場長度方向適當(dāng)延伸。計(jì)算網(wǎng)格如圖2-3所示。
圖2 流場計(jì)算網(wǎng)格
圖3 NACA0012計(jì)算網(wǎng)格
本文關(guān)注不同氣象條件及飛行速度對飛行器表面水滴撞擊特性的影響,為便于對比分析,并從分析中獲得一般規(guī)律,選取計(jì)算條件如表1所示。
表1 飛行速度及氣象條件表
利用上述計(jì)算方法對表1中所列計(jì)算條件下NACA0012機(jī)翼的水滴撞擊特性進(jìn)行了數(shù)值模擬計(jì)算。水滴在機(jī)翼空氣流場中運(yùn)動,其運(yùn)動軌跡受到空氣黏性阻力的影響。由于水滴具有慣性,因此其保持原有運(yùn)動形式的趨勢撞擊到機(jī)翼表面。水滴的慣性及受到的阻力與水滴的平均容積直徑、來流速度等有關(guān)。
選取狀態(tài)1,取距翼根50mm截面處的計(jì)算結(jié)果與文獻(xiàn)[9]中結(jié)果作對比,水滴收集系數(shù)分布曲線如圖4所示。由圖4可知計(jì)算結(jié)果與文獻(xiàn)結(jié)果走勢符合,因此證明本文計(jì)算方法是合理的。
選取飛行速度不同而其他條件相同的狀態(tài)4,6,7對水滴撞擊特性進(jìn)行分析,計(jì)算結(jié)果如圖5所示。由計(jì)算結(jié)果可知,當(dāng)其他飛行參數(shù)和氣象條件不變時(shí),飛行速度減小,水滴受到的慣性力減小,慣性對水滴的影響減小,使得阻力對水滴的運(yùn)動軌跡影響增大。因此,水滴撞擊到機(jī)翼表面的趨勢有所減弱,從而可知機(jī)翼水滴局部水收集系數(shù)逐漸減小。
圖4 狀態(tài)1條件下機(jī)翼水滴收集系數(shù)分布曲線
圖5 不同飛行速度下機(jī)翼水滴收集系數(shù)分布曲線
選取水滴平均容積直徑不同而其他條件相同的狀態(tài)1,2,3對水滴撞擊特性進(jìn)行分析。隨著水滴平均容積直徑的增大,水滴本身的慣性也相應(yīng)增大,慣性力對水滴的影響增加,黏性阻力對水滴軌跡的影響減小,直徑大的液滴由于慣性作用,在機(jī)翼壁面附近不能像空氣流場一樣很輕松地繞流,而是直接撞擊在壁面上,故而,有圖6中所示的計(jì)算結(jié)果,即水滴收集系數(shù)隨著水滴平均容積直徑的增大而增大,且撞擊范圍增大。值得注意的是,在水滴平均容積直徑超過一定數(shù)值時(shí)(一般認(rèn)為是50μm),水滴將會產(chǎn)生沿程參數(shù)變化,并在壁面附近產(chǎn)生飛濺、反彈等物理現(xiàn)象,使得撞擊特性發(fā)生變化[10]。
圖6 不同平均容積直徑下機(jī)翼水滴收集系數(shù)分布曲線
選取液態(tài)水含量不同而其他狀態(tài)相同的狀態(tài)3,4,5對水滴撞擊特性進(jìn)行分析,由于在控制方程中,液態(tài)水含量的值可以換算出相應(yīng)的體積分?jǐn)?shù)α值,當(dāng)水滴流場處于穩(wěn)態(tài)時(shí),即時(shí)間導(dǎo)數(shù)項(xiàng)為零,方程中體積分?jǐn)?shù)項(xiàng)可以消去,故而液態(tài)水含量對水滴撞擊特性影響不大,計(jì)算結(jié)果(見圖7)也證明了理論分析的結(jié)果,即當(dāng)液態(tài)水含量變化時(shí),機(jī)翼附近水滴收集系數(shù)基本不變。
圖7 不同液態(tài)水含量下機(jī)翼水滴收集系數(shù)分布曲線
本文基于歐拉方法計(jì)算了不同飛行速度和氣象條件下NACA0012機(jī)翼水滴撞擊特性,使用Fluent商用軟件計(jì)算空氣流場,并通過編制UDS程序求解水滴相控制方程,對水滴相撞擊邊界條件進(jìn)行了特殊處理使之更好地反映實(shí)際物理過程。從計(jì)算結(jié)果中可以得出:隨著飛行速度的增加,局部水滴收集系數(shù)明顯上升,水滴撞擊范圍加大;隨著水滴平均容積直徑的增加,水滴收集系數(shù)上升明顯,水滴撞擊范圍加大;而液態(tài)水含量對水滴撞擊特性影響不大。
由于水滴平均容積直徑對水滴收集系數(shù)影響很大,且在水滴直徑達(dá)到一定數(shù)值時(shí)(一般為50μm以上),水滴參數(shù)沿程會發(fā)生改變,并會在壁面上產(chǎn)生飛濺、反彈等現(xiàn)象。故而,在過冷大水滴條件下的撞擊特性還有待深入研究。
[1] 裘燮綱.韓鳳華.飛機(jī)防冰系統(tǒng)[M].北京:航空專業(yè)教材編審組出版,1985:44-58.
[2] Mark G Potapczuk.LEWICE/E:An Euler Based Ice Accretion Code[R].AIAA 92-0037,1992.
[3] Hedde T,Guffond D.Improvement of the ONERA 3D icing code comparison with 3Dexperimental shapes[R].AIAA 93-0169,1993.
[4] Morency F,Heloise Beaugendre,Wagdi G Habashi.FENSAP-ICE:a study of ice shapes on droplet impingment[R].AIAA 2003-1223,2003.
[5] Suikno Wirogo,Shashidhar Sriramblatla.An eulerian method to calculate the collection efficiency on two and three dimensional bodies[R].AIAA 2003-1073,2003.
[6] Chirag Bhargava,Eric Loth,Mark Potapczuk.Aerodynamic simulations of the NASA glenn icing research tunnel[R].AIAA 2003-566,2003.
[7] 楊勝華,林貴平,申曉斌.三維復(fù)雜表面水滴撞擊特性計(jì)算[J].航空動力學(xué)報(bào),2010,25(2):284-290.
[8] ANSYS Inc.ANSYS fluent12.0User’s manual[M].New Hampshire:ANSYS Inc,2009.
[9] Gray A Ruff,Brian M Berkowitz.User’s manual for the NASA lewis ice accretion prediction code(LEWICE)[R].NASA,1990:1-27.
[10] Brahimi M T,Tran P,Chocron D,et al.Effect of supercooled large droplets on ice accretion characteristics[R].AIAA 97-0306,1997.