非結(jié)冰氣象條件下機(jī)翼熱氣防冰系統(tǒng)數(shù)值模擬

郁嘉，卜雪琴，2，*，林貴平，2，申曉斌，2，馬文濤，2

(1.北京航空航天大學(xué)航空科學(xué)與工程學(xué)院，北京100191; 2.北京航空航天大學(xué)人機(jī)工效與環(huán)境控制重點(diǎn)學(xué)科實驗室，北京100191)

非結(jié)冰氣象條件下機(jī)翼熱氣防冰系統(tǒng)數(shù)值模擬

郁嘉1，卜雪琴1，2，*，林貴平1，2，申曉斌1，2，馬文濤1，2

(1.北京航空航天大學(xué)航空科學(xué)與工程學(xué)院，北京100191; 2.北京航空航天大學(xué)人機(jī)工效與環(huán)境控制重點(diǎn)學(xué)科實驗室，北京100191)

為了在自然結(jié)冰飛行前預(yù)估飛機(jī)防冰系統(tǒng)的性能，減少飛行風(fēng)險，開展了非結(jié)冰氣象條件下飛機(jī)機(jī)翼熱氣防冰系統(tǒng)的數(shù)值模擬研究。數(shù)值模擬采用了計算流體力學(xué)方法;外部對流傳熱系數(shù)的計算采用附面層積分方法。機(jī)翼熱氣防冰系統(tǒng)表面溫度的計算是將蒙皮外部散熱熱流、內(nèi)部熱氣加熱熱流以及蒙皮導(dǎo)熱三者進(jìn)行熱耦合，并進(jìn)行了改進(jìn)，改進(jìn)的方法不要求防冰腔內(nèi)、外兩套網(wǎng)格在重合面網(wǎng)格處一致，是通過雙向線性插值將一方面網(wǎng)格信息插值傳遞到另一方表面網(wǎng)格上，提高了計算效率。

熱氣防冰系統(tǒng);機(jī)翼;數(shù)值模擬;非結(jié)冰氣象條件;表面溫度

0 引言

為了驗證飛機(jī)防冰系統(tǒng)的性能以及各部件的有效性，防冰系統(tǒng)性能的理論分析和試驗在飛機(jī)適航認(rèn)證中缺一不可［1］。試驗一方面驗證理論計算，另一方面直接有力地檢驗飛機(jī)防除冰系統(tǒng)及其部件的有效性，而理論計算可用于任意條件下的防冰系統(tǒng)性能分析。自然結(jié)冰條件下的飛行試驗是適航中必不可少的一個步驟，但風(fēng)險性很高。因此在自然結(jié)冰飛行試驗之前一般都要進(jìn)行理論分析和計算，還要在冰風(fēng)洞內(nèi)或非結(jié)冰氣象條件下開展防除冰系統(tǒng)性能試驗，以初步地檢驗防除冰系統(tǒng)的性能。非結(jié)冰氣象條件下防冰系統(tǒng)的理論計算為其在自然結(jié)冰條件下的理論計算提供依據(jù)。國內(nèi)較少開展飛機(jī)防冰系統(tǒng)在非結(jié)冰氣象條件下的數(shù)值模擬［2-3］，也未曾得到試驗驗證，因此開展這方面的研究具有重要的意義。

本文以某型飛機(jī)為背景，開展了非結(jié)冰氣象條件下機(jī)翼熱氣防冰系統(tǒng)的數(shù)值模擬研究，為此方法擴(kuò)展到結(jié)冰氣象條件下防冰系統(tǒng)的模擬提供依據(jù)，有助于提高飛機(jī)在自然結(jié)冰條件下的飛行安全。

1 數(shù)值模擬方法

歐美國家由于有冰風(fēng)洞試驗數(shù)據(jù)的支持，已經(jīng)開發(fā)了比較成熟的結(jié)/防冰計算軟件，例如:美國的LEWICE［4］，加拿大的FENSAP-ICE［5］和CANICE［6］，英國的TRAJICE2［7］，法國的ONERA［8］，意大利的CIRAMMIL［9］等。本文采用內(nèi)外傳熱耦合的方法進(jìn)行防冰系統(tǒng)表面溫度計算，實現(xiàn)了非結(jié)冰氣象條件下防冰系統(tǒng)性能的預(yù)估。

1.1 計算模型

外部流場計算選取三維典型機(jī)翼的一段，并選取防冰腔某一段進(jìn)行腔內(nèi)熱氣流動與換熱計算及防冰表面溫度的計算，如圖1所示。選取的機(jī)翼段展長1.76m;防冰腔段展向長0.1m，蒙皮厚度1.6mm。

圖1 機(jī)翼及防冰腔計算模型位置Fig.1 Positions of the model of the wing and anti-icing cavity for simulation

飛機(jī)機(jī)翼熱氣防冰系統(tǒng)比較復(fù)雜，本文研究的熱氣防冰系統(tǒng)示意圖如圖2所示，是雙排陣列笛形管噴口、雙蒙皮熱氣通道結(jié)構(gòu)。發(fā)動機(jī)引氣出來的熱氣經(jīng)過管路流到前緣縫翼內(nèi)的笛形管內(nèi)，熱氣從笛形管上的雙排交錯小孔射流出到前腔，同時加熱蒙皮表面，熱氣經(jīng)過雙蒙皮通道后流到后腔，并從后腔的排氣孔排出到環(huán)境中。

圖2 熱氣防冰腔示意圖Fig.2 Schematic of hot air anti-icing cavity

1.2 數(shù)值模擬總體思想

機(jī)翼熱氣防冰系統(tǒng)性能評價的重要指標(biāo)之一是表面溫度分布。本文數(shù)值模擬的最終目標(biāo)是得到防冰系統(tǒng)工作時防冰表面溫度結(jié)果。

熱氣防冰系統(tǒng)表面溫度的計算需要考慮:①防冰腔外部即飛機(jī)機(jī)翼外部氣流對蒙皮產(chǎn)生的換熱;②防冰腔內(nèi)熱氣對蒙皮產(chǎn)生的換熱;③蒙皮本身的導(dǎo)熱。簡單來說，非結(jié)冰氣象條件下防冰腔可認(rèn)為是無相變氣-氣熱交換器;結(jié)冰氣象條件下防冰腔則是有相變(撞擊水蒸發(fā))氣-氣熱交換器。復(fù)雜之處在于防冰腔內(nèi)、外部空氣的流動比較復(fù)雜，要想準(zhǔn)確地預(yù)測傳熱現(xiàn)象比較困難。

防冰表面溫度計算步驟如下:

步驟1計算翼身組合體外部流場。

步驟2設(shè)置表面溫度初始值，計算防冰腔內(nèi)部流動傳熱。

步驟3利用附面層積分方法計算外部對流傳熱系數(shù)，并分析外部散熱熱流。

步驟4將外部熱流作為防冰腔外表面熱流邊界條件，進(jìn)一步計算蒙皮導(dǎo)熱以及防冰腔內(nèi)部流動與換熱，得到新的表面溫度。

步驟5重復(fù)步驟3和步驟4直到所有控制容積相鄰兩步的表面溫度結(jié)果的最大差值小于小量ε。

整個計算基于計算流體力學(xué)(CFD)的方法，最終得到非結(jié)冰氣象條件下防冰腔表面溫度結(jié)果。

1.3 外部流場計算

利用CFD方法計算機(jī)翼外部流場，采用Fluent計算工具。圖3給出了外部三維流場計算網(wǎng)格，劃分工具采用Gambit軟件。網(wǎng)格總數(shù)為147900個網(wǎng)格，為結(jié)構(gòu)化的六面體網(wǎng)格。

圖3 外部流場計算網(wǎng)格Fig.3 Grid for external flow field of wing-body combination

機(jī)翼外部流場計算采用壓力遠(yuǎn)場邊界條件，無黏流模型。壓力修正采用分離式求解方法，即Simple算法。方程組的離散采用二階迎風(fēng)格式。外部流場計算收斂后可近似得到機(jī)翼附面層外邊界處速度，為附面層積分方法計算外部對流傳熱系數(shù)做準(zhǔn)備。

監(jiān)控速度及能量的殘差值來檢驗計算是否收斂。計算收斂后，將空氣速度及溫度結(jié)果導(dǎo)出，作為計算外部對流傳熱系數(shù)的輸入條件。

1.4 外部對流傳熱系數(shù)

計算外部對流傳熱系數(shù)時，根據(jù)機(jī)翼表面的空氣流動特點(diǎn)，機(jī)翼可分成駐點(diǎn)、上表面和下表面，其中:上下表面均有層流區(qū)、湍流區(qū)和過渡區(qū)。本文總結(jié)了前人通過理論推導(dǎo)和大量試驗得到的各區(qū)域的對流傳熱系數(shù)計算公式，采用的是附面層積分方法。

在駐點(diǎn)位置，根據(jù)Smith和Spalding經(jīng)過大量試驗后得出的結(jié)論，駐點(diǎn)區(qū)的對流傳熱系數(shù)計算式為［10］

式中:Nustag為駐點(diǎn)努塞爾數(shù);hstag為駐點(diǎn)對流傳熱系數(shù);c為機(jī)翼弦長;kair為空氣導(dǎo)熱系數(shù);Re∞=V∞c/vair為基于弦長的遠(yuǎn)場雷諾數(shù);V∞為遠(yuǎn)場空氣速度;vair為空氣運(yùn)動黏度;ue為附面層外邊界處切向方向速度; s為弧長。

層流區(qū)域的對流傳熱計算式為［11］

式(2)考慮了壓力梯度及非等溫表面情況。式中:Res=ues/vair為基于弧長的當(dāng)?shù)乩字Z數(shù);Ve為附面層外界處速度;ΔT=Ts－Te為當(dāng)?shù)乇砻婧透矫鎸油膺吔缣帨夭睢?/p>

湍流區(qū)域的對流傳熱系數(shù)計算式為［11］

式中:St為斯坦頓數(shù);ρa(bǔ)ir為空氣密度;cp，air為空氣定壓比熱;ReΔ2，trub=ueΔ2/vair為基于焓厚度的當(dāng)?shù)乩字Z數(shù);Pr為空氣普朗特數(shù)。ReΔ2，trub可近似為如下積分公式［11］

式中:ReΔ2，tr為轉(zhuǎn)捩起始位置的焓厚度雷諾數(shù)，計算表達(dá)式為

過渡區(qū)域的對流傳熱系數(shù)計算引入間斷因子γ(s)，代表在弧長s處流動呈現(xiàn)湍流的概率，其表達(dá)式為［12］

式中:str為過渡區(qū)起始位置距離駐點(diǎn)位置的弧長; sE為過渡區(qū)終止位置對應(yīng)的弧長。那么:

通過以上積分公式的計算，可得到各區(qū)域的對流傳熱系數(shù)。

1.5 內(nèi)部流場及蒙皮導(dǎo)熱初始計算

防冰腔內(nèi)部熱氣流動與換熱比較復(fù)雜，將內(nèi)部熱氣流動換熱及蒙皮導(dǎo)熱問題一起考慮，本文利用計算流體力學(xué)軟件Fluent來開展計算。首先需要將防冰腔內(nèi)的空間及蒙皮離散化即劃分網(wǎng)格。防冰腔曲面數(shù)量較多，結(jié)構(gòu)復(fù)雜，不利于網(wǎng)格結(jié)構(gòu)化劃分，因此對模型進(jìn)行簡化如下:

1)由于笛形管噴口噴射方向基本垂直于縫翼前緣，因此選擇防冰腔段時垂直于前緣縫翼截取。

2)根據(jù)雙蒙皮通道位置，對雙蒙皮曲面進(jìn)行切割，構(gòu)建通道輪廓，將部分小曲面改為平面，例如通道末端銑刀形狀導(dǎo)致的小曲面。

3)忽略雙蒙皮通道內(nèi)蒙皮的厚度。

4)根據(jù)笛形管噴口位置，對笛形管圓柱曲面進(jìn)行切割，構(gòu)建與圓形噴口等面積的方形熱氣噴口，并忽略笛形管管壁厚度。

完成以上簡化工作后，對模型進(jìn)行網(wǎng)格劃分，采用Gambit網(wǎng)格劃分工具。圖4為防冰腔網(wǎng)格，包括蒙皮固體區(qū)域網(wǎng)格和流體流動區(qū)域網(wǎng)格。由于前腔流動變化劇烈，劃分結(jié)構(gòu)化的網(wǎng)格及蒙皮內(nèi)表面附面層網(wǎng)格有助于提高計算的精度，為了保證前腔結(jié)構(gòu)化網(wǎng)格的生成，前腔進(jìn)行分區(qū)處理。固壁區(qū)域的網(wǎng)格數(shù)為46104，腔內(nèi)氣體流動區(qū)域的網(wǎng)格數(shù)為393924。

圖4 防冰腔內(nèi)部流場計算網(wǎng)格Fig.4 Grid for internal flow field of anti-icing cavity

防冰腔內(nèi)熱氣流動換熱以及蒙皮的導(dǎo)熱計算需要聯(lián)立求解熱氣的N-S方程、能量方程以及固體區(qū)域的導(dǎo)熱方程。另外，湍流模型采用Spallart-Allmaras一方程模型，比較適合沖擊射流曲面的流動計算。采用壓力入口、壓力出口邊界條件，蒙皮外表面采取等溫邊界條件。

防冰腔內(nèi)流動變化比較劇烈，特別是在壁面附近，嘗試計算后，發(fā)現(xiàn)壓力的修正采用Coupled算法在計算收斂方面更具優(yōu)勢。方程組的離散采用二階迎風(fēng)格式。

1.6 內(nèi)外傳熱耦合計算表面溫度

為了得到蒙皮表面溫度，本文將外部熱流、內(nèi)部熱流以及蒙皮自身導(dǎo)熱三者一起強(qiáng)固耦合計算，有別于結(jié)/防冰軟件FENSAP-ICE所采用的松散耦合方法［13-14］。本文在計算導(dǎo)熱與防冰腔內(nèi)部流動與換熱時，將外部熱載荷作為邊界條件加載到防冰腔蒙皮外表面，每個迭代步不斷地交替更新表面溫度和外部熱載荷，通過監(jiān)控相鄰迭代步的表面溫度變化來判斷迭代計算的收斂，最終得到表面平衡溫度。

這里將外部流場計算網(wǎng)格簡稱為外網(wǎng)格;內(nèi)部流場計算網(wǎng)格簡稱為內(nèi)網(wǎng)格;外網(wǎng)格和內(nèi)網(wǎng)格的重合面稱為界面。內(nèi)外傳熱耦合計算表面溫度存在的難點(diǎn)之一是:如何在界面處給內(nèi)、外網(wǎng)格相互傳遞變量值的信息，達(dá)到傳熱耦合求解的目的。

在文獻(xiàn)［3］的研究中，外部流場計算區(qū)域在機(jī)翼展向方向的范圍較小，僅為所研究防冰腔的展向長度，外網(wǎng)格和內(nèi)網(wǎng)格在界面處的面網(wǎng)格一致，因此界面網(wǎng)格處的變量值傳遞時只需要尋找相同位置的網(wǎng)格即可，如圖5(a)所示，找到相同位置的界面網(wǎng)格后，就可以方便地將計算所得的外網(wǎng)格界面處的熱載荷傳遞給內(nèi)網(wǎng)格界面，作為計算新的表面溫度的輸入;將內(nèi)網(wǎng)格計算所得的表面溫度傳遞給外網(wǎng)格界面，作為計算新的防冰熱載荷的輸入。

由于本文外流場計算針對整個機(jī)翼進(jìn)行，計算區(qū)域相對較大，若采用上述方法來實現(xiàn)界面的數(shù)據(jù)傳遞，外網(wǎng)格的界面網(wǎng)格必然會較密，網(wǎng)格數(shù)量的增多將降低計算速度。本文對界面間數(shù)據(jù)傳遞方法進(jìn)行改進(jìn)，采用了線性插值的方法進(jìn)行雙向插值，實現(xiàn)不一致界面間的數(shù)據(jù)傳遞，如圖5(b)所示。這樣外網(wǎng)格的界面網(wǎng)格可以疏于內(nèi)網(wǎng)格的界面網(wǎng)格，大大提高了計算速度。

圖5 內(nèi)、外網(wǎng)格界面數(shù)據(jù)傳遞Fig.5 Data transfer between the internalexternal grid interfaces

根據(jù)外網(wǎng)格界面網(wǎng)格節(jié)點(diǎn)a的坐標(biāo)找到內(nèi)網(wǎng)格界面上的4個最近的網(wǎng)格，插值計算得到a'的表面溫度，將a'的表面溫度結(jié)果傳遞給外網(wǎng)格節(jié)點(diǎn)a。根據(jù)內(nèi)網(wǎng)格界面網(wǎng)格節(jié)點(diǎn)b'的坐標(biāo)找到外網(wǎng)格界面上的4個最近的網(wǎng)格，插值計算得到b的防冰熱載荷，將b的熱載荷結(jié)果傳遞給內(nèi)網(wǎng)格節(jié)點(diǎn)b'。其他界面網(wǎng)格如法炮制，即可完成內(nèi)、外網(wǎng)格界面間的數(shù)據(jù)傳遞。插值計算時根據(jù)展向z坐標(biāo)以及弧長s進(jìn)行二維線性插值，定義見圖5(b)。

編寫子函數(shù)來尋找每個界面網(wǎng)格插值計算所需的4個網(wǎng)格，并在表面溫度迭代計算之前運(yùn)行，將網(wǎng)格序號存貯，為每次迭代計算中插值及傳遞數(shù)據(jù)使用。

利用上述內(nèi)外傳熱耦合思想以及改進(jìn)后的內(nèi)、外網(wǎng)格界面數(shù)據(jù)傳遞方法，對防冰腔在非結(jié)冰氣象條件下防冰表面溫度開展了計算。

2 數(shù)值模擬結(jié)果與分析

在計算表面對流傳熱系數(shù)時，則根據(jù)雷諾數(shù)大小來劃分層流區(qū)、湍流區(qū)和過渡區(qū)。當(dāng)雷諾數(shù)小于5× 105，認(rèn)為層流區(qū)［15］;雷諾數(shù)大于2×106，認(rèn)為湍流區(qū)［15］;雷諾數(shù)在以上兩者之間為過渡區(qū)。

計算狀態(tài)如表1所示，表中p為壓力，α為迎角，T∞為環(huán)境溫度，Ma為馬赫數(shù)。

表1中狀態(tài)1和狀態(tài)2的防冰系統(tǒng)引氣狀態(tài)分別為:壓力0.26 MPa，溫度200℃;壓力0.25 MPa，溫度200℃。蒙皮固體區(qū)域?qū)嵯禂?shù)為121W/ (m·K)。

表1 計算狀態(tài)Table 1 Calculation cases

圖6和圖7給出了狀態(tài)1的數(shù)值模擬三維云圖結(jié)果。從防冰腔三維結(jié)果中抽取中間截面的溫度與對流傳熱系數(shù)結(jié)果，如圖8所示。

圖6 防冰腔外表面對流傳熱系數(shù)三維云圖Fig.6 Contours plot of the convective heat transfer coefficient of the anti-icing cavity external surface

圖7 防冰表面溫度三維云圖(單位:K)Fig.7 Surface temperature distribution of the anti-icing cavity(unit:K)

圖8 二維截面對流傳熱系數(shù)及表面溫度模擬結(jié)果Fig.8 2D surface heat transfer coefficient and tem perature simulation results

由圖6可以看到，機(jī)翼實際駐點(diǎn)處的對流傳熱系數(shù)達(dá)到一個峰值;氣流在駐點(diǎn)處分流，氣流往上表面流動時，對流傳熱系數(shù)沿著氣流流動方向首先在層流區(qū)逐漸下降，然后在過渡區(qū)逐漸升高，完全轉(zhuǎn)變成湍流后對流傳熱系數(shù)又開始下降，這一點(diǎn)從圖8中二維截面的對流傳熱系數(shù)曲線圖也可看出;氣流往下表面流動時，由于流動速度較低，下表面一直處于層流狀態(tài)，對流換熱系數(shù)一直下降。氣流往上表面流動時對流傳熱系數(shù)變化并不是一直都平滑，某些地方出現(xiàn)拐點(diǎn)，如圖8中二維截面s=0.1附近，分析發(fā)現(xiàn)s=0.1附近也是溫度變化的拐點(diǎn)位置，對流傳熱系數(shù)與表面溫度有關(guān)系，因為表面溫度影響著附近空氣的物性參數(shù)，從而影響對流傳熱系數(shù)。

由圖7可知，射流駐點(diǎn)附近的表面溫度明顯高于其他部分的溫度，呈現(xiàn)出笛形孔交錯排列的形狀，上下表面溫度逐漸降低。從圖8可知狀態(tài)2的表面溫度整體低于狀態(tài)1的，這主要是因為狀態(tài)2的環(huán)境溫度低于狀態(tài)1，另外還因為飛行速度高于狀態(tài)1，這樣導(dǎo)致對流散熱量會高于狀態(tài)1。對比圖8中狀態(tài)1和狀態(tài)2的對流換熱系數(shù)可看出狀態(tài)2駐點(diǎn)處hs要高于狀態(tài)1的，且上表面層流向紊流的轉(zhuǎn)變提前，導(dǎo)致上表面后部狀態(tài)2的hs要高于狀態(tài)1。

從本文計算結(jié)果來看，在非結(jié)冰氣象條件下飛行時，此段防冰系統(tǒng)工作時能夠使表面溫度最高達(dá)到380K(107℃)左右，上表面最低達(dá)到310 K(37℃)左右，下表面最低達(dá)到335 K(62℃)左右。在結(jié)冰氣象條件下，由于環(huán)境溫度低加上撞擊水的蒸發(fā)效果，表面整體溫度會低于非結(jié)冰氣象條件下的表面溫度。但從此處非結(jié)冰氣象條件下較高的表面溫度結(jié)果來看，即使在結(jié)冰氣象條件下，駐點(diǎn)附近以及下表面的溫度能夠滿足防冰要求。因為在結(jié)冰氣象條件下，上表面結(jié)冰相對于下表面結(jié)冰對飛機(jī)性能的影響要小，且一般來說機(jī)翼上表面水滴撞擊范圍要小于下表面的水滴撞擊范圍，上表面溢流水極限也會遠(yuǎn)小于上表面防護(hù)范圍，因此在上表面后部即使溫度低(在結(jié)冰氣象條件下會遠(yuǎn)低于此處計算結(jié)果37℃)，由于沒有水的存在不會發(fā)生結(jié)冰。本文非結(jié)冰氣象條件下機(jī)翼熱氣防冰系統(tǒng)性能的預(yù)測為飛機(jī)防冰系統(tǒng)在非結(jié)冰氣象/結(jié)冰氣象條件下的性能飛行試驗提供了參考依據(jù)。

3 結(jié)論

本文提出并詳細(xì)介紹了內(nèi)外傳熱耦合計算熱氣防冰系統(tǒng)表面溫度的方法。該方法基于Fluent軟件進(jìn)行二次開發(fā)，可擴(kuò)充性較好。

通過改進(jìn)內(nèi)外網(wǎng)格數(shù)據(jù)傳遞的方法，采用雙向線性插值，實現(xiàn)了不一致界面間的數(shù)據(jù)傳遞。該方法可減小外流場計算網(wǎng)格，相對以往內(nèi)外網(wǎng)格在界面網(wǎng)格處須一致的情況，有效地提高了計算速度。

采用該方法對某型飛機(jī)某段防冰系統(tǒng)在非結(jié)冰氣象條件下的表面溫度開展計算，結(jié)果表明整個防冰表面溫度較高。射流駐點(diǎn)附近表面溫度明顯高于其他部分的溫度，分布呈笛型孔交錯排列形狀。上表面表面溫度低于下表面表面溫度。

［1］Civil Aviation Administration of China.Civil aviation regulations of China，Part 25:airworthiness standards for transport category airplanes.4th ed.［S］.Beijing:Civil Aviation Administration of China，2001:144.(in Chinese)中國民航總局.中國民用航空規(guī)章第25部:運(yùn)輸類飛機(jī)適航標(biāo)準(zhǔn)第四次修訂版［S］.北京:中國民用航空總局，2011:144.

［2］Bu Xueqin，Lin Guiping.Predictions of collection efficiency and anti-icing surface temperature［J］.Journal of Beijing University of Aeronautics and Astronautics，2007，33(10):1182-1185.(in Chinese)卜雪琴，林貴平.基于CFD的水收集系數(shù)及防冰表面溫度預(yù)測［J］.北京航空航天大學(xué)，2007，33(10):1182-1185.

［3］Bu Xueqin，Lin Guiping，Yu Jia.Three dimensional conjugate heat transfer simulation for the surface temperature of wing hot-air antiicing system［J］.Journal of Aerospace Power，2009，24(11): 2495-2500.(in Chinese)卜雪琴，林貴平，郁嘉.三維內(nèi)外熱耦合計算熱氣防冰系統(tǒng)表面溫度［J］.航空動力學(xué)報，2009，24(11):2495-2500.

［4］Wright W B.User manual for the NASA Glenn ice accretion code LEWICE version 2.2.2［R］.NASA CR-2002-211793，2002.

［5］Morency F，Beaugendre H，Baruzzi G S，et al.FENSAP-ICE:a comprehensive 3D simulation system for in-flight icing.AIAA 2001-2566［R］.Reston:AIAA，2001.

［6］Saeed F，Paraschivoiu I.Optimization of a hot-air anti-icing system.AIAA 2003-0733［R］.Reston:AIAA，2003.

［7］Gent R W.TRAJICE2，A combined water droplet and ice accretion prediction program for aerofoil［R］.DRA Technical Report No.TR90054，1990.

［8］Hedde T，Guffond D.ONERA three-dimensional icing model［J］.AIAA Journal.1995，33(6):1038-1045.

［9］Fortin G，Perron J，Mingione G，et al.CIRAAMIL ice accretion code improvement.AIAA 2009-3968［R］.Reston:AIAA，2009.

［10］Smith A G，Spalding D B.Heat transfer in a laminar boundary layer with constant fluid properties and constant wall temperature［J］.Journal of the Royal Aeronautical Society，1958，62(1):60-64.

［11］Ambrok G S.Approximate solution of equations for the thermal boundary layer with variations in boundary layer structure［J］.Soviet Physics-Technical Physics，1957(2):1979-1986.

［12］Abu-Ghannam B，Shaw R.Natural transition of boundary layers——The effects of turbulence，pressure gradient and flow history［J］.Journal of Mechanical Engineering Science，1980，22 (5):213-228.

［13］Croce G，Habashi W G，Guevremont G，et al.3D thermal analysis of an anti-icing device using［R］.FENSAP-ICE.AIAA-1998-0193，1998.

［14］Wang H，Tran P，Habashi W G，et al.Anti-icing simulation in wet air of a piccolo system using FENSAP-ICE［R］.SAE paper 2007-01-3357，2007.

［15］Qiu Xiegang，Han Fenghua.Aircraft anti-icing system［M］.Beijing:Compilation and Examination Group of Aero Specialized Teaching Materials，1985:166.(in Chinese)裘燮綱，韓鳳華.飛機(jī)防冰系統(tǒng)［M］.北京:航空專業(yè)教材編審組，1985:166.

Numerical simulation of a wing hot air anti-icing system in dry air conditions

Yu Jia1，Bu Xueqin1，2，*，Lin Guiping1，2，Shen Xiaobin1，2，Ma Wentao1，2
(1.School of Aeronautic Science and Engineering，Beihang University，Beijing 100191，China;2.Fundamental Science on Ergonomics and Environment Control Laboratory，Beihang University，Beijing 100191，China)

In order to evaluate the effectiveness of an ice protection system in natural atmospheric icing conditions and reduce the flight test risks in icing conditions.Numerical analysis and validation for a wing hot-air anti-icing system in dry air conditions are investigated.The computational fluid dynamics (CFD)method was used for the simulation of the wing anti-icing system.The thermal boundary layer integral method was applied to acquire the external convective heat transfer coefficient.Surface equilibrium temperatures were obtained by coupling the external convective heat loss，the internal heat gain and the thermal conductivity through the skin.An improved method was proposed to exchange the values of the surface temperature and heat loads between the internal grid interface and the external grid interface.The method reduced the total cell number of the external flow field and accelerated the calculation speed.

hot air anti-icing system;wing;numerical simulation;dry air condition;surface temperature

V244.1+5

A

10.7638/kqdlxxb-2012.0212

0258-1825(2016)05-0562-06

2012-12-20;

2013-03-28

國家自然科學(xué)基金(51206008)

郁嘉(1979-)，男，上海人，博士，研究方向:飛機(jī)防除冰，飛行力學(xué)與飛行安全.E-mail:Yujia@buaa.edu.cn

卜雪琴*(1982-)，女，江西萍鄉(xiāng)人，講師，研究方向:飛機(jī)防除冰.E-mail:buxueqin@buaa.edu.cn

郁嘉，卜雪琴，林貴平，等.非結(jié)冰氣象條件下機(jī)翼熱氣防冰系統(tǒng)數(shù)值模擬［J］.空氣動力學(xué)學(xué)報，2016，34(5):562-567.

10.7638/kqdlxxb-2012.0212 Yu J，Bu X Q，Lin G P，et al.Numerical simulation of a wing hot air anti-icing system in dry air conditions［J］.Acta Aerodynamica Sinica，2016，34(5):562-567.