龔 淼，黃 文

（中國民航大學(xué)航空工程學(xué)院，天津 300300）

0 引言

航空發(fā)動(dòng)機(jī)在慢車狀態(tài)時(shí)，尾流是影響地面作業(yè)安全的重要因素，尤其是對(duì)于民航局正在大力推進(jìn)的飛機(jī)慢車除冰作業(yè)，合理的慢車安全區(qū)規(guī)劃是影響作業(yè)安全和除冰效率的關(guān)鍵。此外，開展對(duì)航空發(fā)動(dòng)機(jī)尾部流場的研究，準(zhǔn)確分析尾部熱流場的特性對(duì)除冰作業(yè)軌跡規(guī)劃和航空發(fā)動(dòng)機(jī)設(shè)計(jì)都具有重要的工程意義。

在基于航空發(fā)動(dòng)機(jī)二維模型的尾流數(shù)值分析方面，目前研究方向有飛機(jī)內(nèi)墻壁和空氣之間傳熱的緊密關(guān)系、發(fā)動(dòng)機(jī)高速時(shí)的燃?xì)馍淞髁鲌龇植家?guī)律、不同飛行參數(shù)下的離散相模型發(fā)動(dòng)機(jī)尾流場流場分布［1－3］。在基于航空發(fā)動(dòng)機(jī)三維模型尾流數(shù)值分析方面，研究方向包括發(fā)動(dòng)機(jī)排氣后的紅外損傷位置和程度、不同形狀的噴管尾流的沖擊特性、橢圓噴管對(duì)自由雙噴流混合和湍流特性的影響、飛機(jī)推進(jìn)系統(tǒng)和發(fā)動(dòng)機(jī)尾流結(jié)構(gòu)對(duì)各性能參數(shù)變化、收斂噴管和錐形噴管拋物線方程對(duì)流場影響、噴氣式飛機(jī)的流動(dòng)RANS模型優(yōu)化［4－10］。

國內(nèi)外對(duì)高馬赫數(shù)和一定飛行高度下的尾流數(shù)值模擬和紅外輻射計(jì)算較多，對(duì)慢車時(shí)熱流場研究較少。大部分學(xué)者對(duì)比二維和三維模型模擬結(jié)果，發(fā)現(xiàn)三維模型更加符合實(shí)際情況。目前多數(shù)文獻(xiàn)主要研究二維和三維合成噴流上的流動(dòng)行為和特性［11］。

由于大型發(fā)動(dòng)機(jī)建模和計(jì)算量龐大，相應(yīng)數(shù)據(jù)還在獲取中，而簡化模型又存在較大的誤差，因此本文采用SW120B發(fā)動(dòng)機(jī)為研究對(duì)象，在慢車狀態(tài)下分析尾部流場的溫度和流速特性，獲得地面作業(yè)的安全范圍。以真實(shí)小型發(fā)動(dòng)機(jī)為基礎(chǔ)模型的數(shù)值分析有助提升后期對(duì)大型發(fā)動(dòng)機(jī)建模分析的準(zhǔn)確率。

1 噴管幾何及數(shù)值模擬

1.1 噴管幾何模型及計(jì)算域

以玄云SW120B發(fā)動(dòng)機(jī)數(shù)據(jù)為理論研究對(duì)象，尾噴管模型的實(shí)際模型如圖1（a）所示。簡化尾噴管模型，建立尾噴管三維幾何模型，如圖1（b）所示［12］。建立三維內(nèi)外流場計(jì)算模型，外部計(jì)算域?yàn)殚L方體，沿z軸方向尺寸為噴管長度的35倍，沿y軸和沿x軸方向尺寸均為噴管入口直徑的10倍，如圖2（a）所示［11］。計(jì)算域視圖如圖2所示，區(qū)域1為尾噴管，區(qū)域2為流體待求解區(qū)，a-b為入口，e-f為流場出口，c-e和d-f為流場外邊界。

圖1 尾噴管模型

圖2 計(jì)算域簡化幾何結(jié)構(gòu)

1.2 網(wǎng)格劃分及邊界條件

在尾噴管出口附近流域細(xì)化網(wǎng)格，采用自由四面體，在流場邊界處設(shè)置邊界層網(wǎng)格。數(shù)值計(jì)算網(wǎng)格數(shù)量約為667萬，模型分別使用了標(biāo)準(zhǔn)偏度、增長率和條件數(shù)來評(píng)估網(wǎng)格質(zhì)量，3種評(píng)估的平均網(wǎng)格質(zhì)量均在0.7以上，滿足流體計(jì)算質(zhì)量要求［12］，網(wǎng)格劃分和單元質(zhì)量分布如圖3所示。圖4選取了噴口中心軸線上網(wǎng)格數(shù)分別為6.18×106、6.67×106、7.72×106時(shí)的尾流溫度，3種不同網(wǎng)格數(shù)的計(jì)算結(jié)果差異小于0.8%。經(jīng)過綜合考慮，最終選擇網(wǎng)格數(shù)為6.67×106的模型進(jìn)行數(shù)值模擬分析。

圖3 尾噴管計(jì)算域網(wǎng)格劃分和質(zhì)量分布

圖4 不同網(wǎng)格數(shù)的溫度

邊界條件：噴管入口給定進(jìn)氣流總溫、流速，流動(dòng)方向與噴管軸向一致；計(jì)算域出口指定靜壓；噴管壁面采用絕熱無滑移固體壁面邊界條件；外邊界給定壓力遠(yuǎn)場條件［26］，尾噴管入口為弱可壓縮湍流［14］。

發(fā)動(dòng)機(jī)尾部氣流組成成分復(fù)雜，涉及噴流和自由流的混合，本研究簡化地面和流場的影響［15］：

1）假設(shè)發(fā)動(dòng)機(jī)燃?xì)鉃槔硐霘怏w，在做功過程中完全燃燒，各組成分化學(xué)性質(zhì)和氣體定壓比不發(fā)生改變；

2）假設(shè)尾噴管內(nèi)部氣體為純氣相流動(dòng)，等熵膨脹；

3）忽略噴管壁的傳熱作用和空氣重力影響，不考慮燃?xì)饧t外輻射作用。

1.3 控制方程

采用雷諾平均Navier-Stokes（RANS）方程，忽略黏性耗散導(dǎo)致的發(fā)熱，渦流黏度可以增大動(dòng)量的黏性傳遞，從而計(jì)算出無法求解的小尺度渦流傳遞的動(dòng)量，在合理的計(jì)算成本下，恰好彌補(bǔ)了標(biāo)準(zhǔn)方程計(jì)算成本高的問題［16－18］。

不可壓縮流動(dòng)的RANS方程為：

式中：ρ為流體的密度；Cρ為流體的比熱容；μ為動(dòng)力黏度；uˉ為平均速度向量；μT為湍流黏度；為平均重力體積力和其他體積力。

對(duì)于湍流黏度μT，需要引入湍流動(dòng)能k和耗散率ε來增加湍流變量。k-ε標(biāo)準(zhǔn)模型具有較好的穩(wěn)定性，在許多工業(yè)應(yīng)用中都頗受歡迎［19－20］。本數(shù)值模擬選擇k-ε標(biāo)準(zhǔn)模型，壁面處理采用標(biāo)準(zhǔn)壁面函數(shù)法，保證在結(jié)果收斂的情況下提高算法精度［21－22］。

k-ε標(biāo)準(zhǔn)模型控制方程為：

式中：ep為湍流耗散率；μT為湍流黏度；Cμ為隨應(yīng)變率變化的量，簡化計(jì)算Cμ＝0.09；Pk為生成相，Cε1＝1.44，Cε2＝1.92，σε＝1.3，σk＝1.0。

數(shù)值模擬是流體傳熱研究，主要求解非定常的物理項(xiàng)，采用2階離散線性單元對(duì)動(dòng)量守恒方程里的對(duì)流項(xiàng)和擴(kuò)散項(xiàng)處理，對(duì)流體傳熱的控制方程為：

式中：ρ為流體的密度；Cρ為流體的比熱容；u為速度向量；T為熱力學(xué)溫度；k為熱傳導(dǎo)系數(shù)；Qted為熱彈性阻尼熱源。

流體計(jì)算域中，入口處控制其法向流入速度，出口處控制其壓力，并抑制回流：

流入為：

流出為：

式中：n為法向矢量；q為熱傳導(dǎo)通量；k為湍流動(dòng)能；ε為耗散率。

2 計(jì)算結(jié)果

分析發(fā)動(dòng)機(jī)慢車狀態(tài)下，尾流的溫度場和速度場的數(shù)值模擬結(jié)果。這些結(jié)果從不同方向和截面進(jìn)行比較，用三維云圖的形式展現(xiàn)。由于SW120B怠速時(shí)轉(zhuǎn)速：38 000 r／min，推力：12 kg（環(huán)境溫度15℃）。在COMSOL軟件進(jìn)行模擬時(shí)，取飛行高度H＝0，外流馬赫數(shù)Ma＝0，尾噴管進(jìn)口總溫Tin＝973 K，流速Vin＝57.5 m／s，遠(yuǎn)場壓力條件為標(biāo)準(zhǔn)大氣壓，環(huán)境溫度為室溫293 K。

2.1 尾流溫度場結(jié)果分析

圖5（a）為整個(gè)計(jì)算域的溫度場分布以及尾噴管噴口處的局部放大云圖，可以看出：溫度場是沿軸向方向不均勻擴(kuò)展的穩(wěn)態(tài)場，氣流經(jīng)流出尾噴管后，一定距離內(nèi)沿軸向是旋轉(zhuǎn)對(duì)稱，并且中心溫度最高。隨著軸向距離增大，溫度逐漸降低，溫度變化較大的區(qū)域比較集中，溫度范圍為315.4～976.0 K。尾噴口附近流場溫度最高可達(dá)976.0 K，尾流的核心高溫溫度范圍為569.1～976.0 K。

圖5 發(fā)動(dòng)機(jī)尾流溫度分布云圖

圖5（b）顯示了尾流溫度場各切面的分布云圖，可以看出：在0.3＜z＜0.9 m、y＜0.2 m、x＜0.2 m時(shí)，溫度范圍為329.2～573.1 K；溫度在z＞1 m、y＞0.3 m、x＞0.3 m時(shí)，溫度等值線較少，在y軸和x軸方向的影響范圍較小。不同切面的尾流溫度區(qū)域在向z軸方向不斷變大，各切面的溫度呈正態(tài)分布［23］，中心最高溫度點(diǎn)為圓心。在一定范圍內(nèi)，溫度往外降低的速率一致，越往外越低，其射流方向與z軸有一定的向下偏角。

數(shù)值仿真過程中采用簡化后的湍流模型并忽略環(huán)境重力影響。湍流在高溫快速流動(dòng)時(shí)，會(huì)產(chǎn)生不同方向的渦流，尾流中心線的方向會(huì)發(fā)生變化［24］。每個(gè)切面的中心溫度在逐漸向y軸負(fù)方向和z軸正方向靠近，發(fā)動(dòng)機(jī)尾流屬于不穩(wěn)定合成射流，傳熱系數(shù)會(huì)受到一定的影響［25］。這說明尾流噴射出來后，形成的溫度射流場并不是完全的二維軸對(duì)稱分布，采用三維模型進(jìn)行數(shù)值模擬更加貼近實(shí)際情況。由于圖5（b）中，每個(gè)切面的溫度擴(kuò)散程度和具體變化并不明顯，為了更加清晰地分析尾流溫度場特性，將各切面的溫度分布云圖以三維映射的形式分開表示。

圖6為截取的yz平面三維溫度分布圖?？梢钥闯?，尾流的分布符合高斯分布，中心溫度最高，逐漸往四周擴(kuò)散，溫度梯度在z＜0.5 m內(nèi)變化明顯，快速下降；在z＞2.5 m的區(qū)域，溫度變化緩慢，屬于低溫區(qū)域，最低可達(dá)325.0 K；尾流溫度在y軸和x軸方向溫度擴(kuò)散范圍比較小，高溫區(qū)域集中在以y軸和x軸±0.2 m以內(nèi)。

圖6 發(fā)動(dòng)機(jī)尾流yz截面三維溫度分布

在圖7中展示了距離噴口不同位置的溫度分布，分別選取了z軸方向距離噴口0.3、0.9、1.5和2.1 m的截面。由圖可知：z＝0.3 m時(shí)，高溫危險(xiǎn)區(qū)域分布以核心高溫為圓心，越靠近噴口的位置，云圖尖端越突出，溫度變化越劇烈，往四周降低的速率均勻且快速，作用區(qū)域較小。相反，z＞0.9 m以后，越遠(yuǎn)離噴口，溫度變化梯度越小，作用范圍變大，中心溫度不均勻地向四周降低擴(kuò)散。

圖7 發(fā)動(dòng)機(jī)尾流xy截面距噴口不同位置的溫度分布云圖

2.2 尾流流速場結(jié)果分析

圖8為數(shù)值模擬計(jì)算所得的尾流流速分布云圖，其中圖8（a）和（b）分別為整個(gè)計(jì)算域的流場，尾噴管出口局部放大和不同切面的流速分布云圖?？梢钥闯?，流場等溫線在噴口處分布密集，流速變化區(qū)域集中，符合在低雷諾數(shù)的混合流中自由紊動(dòng)射流的流動(dòng)特征［39］。Weiss等［40］研究表明了混合流的數(shù)值模擬更加符合實(shí)際情況。流速作用在y軸和x軸的范圍較小且相似，z軸范圍較大。流速沿著噴口逐漸降低，z＜0.2 m，y＜0.1m，x＜0.1m內(nèi)為高速區(qū)域，流速范圍為77～110 m／s。從局部放大圖可以看出，在尾噴管內(nèi)部流速由小變大，尾噴管為收斂型，在流量保持不變的情況下，出口截面面積減少，導(dǎo)致流速增加；低流速區(qū)域0.5＜z＜2.1 m，流速最高不超過20 m／s。

圖8 發(fā)動(dòng)機(jī)尾流流速分布云圖

從圖8（b）可以看出，尾流的流速分布區(qū)域?yàn)椴痪鶆虻恼龖B(tài)分布，各切面的云圖形狀變化較大，射流的整體是向y軸負(fù)方向和z軸正方向有一定的偏轉(zhuǎn)角度。流速場射流z＞0.9 m后，其作用的氣體流速基本和環(huán)境混合，無明顯速度變化。Sharma等［28］研究發(fā)現(xiàn)，合成射流在中心線速度下降速率方面和穩(wěn)定射流有一定的偏差，射流擴(kuò)散呈非線性增長。這表明流體在經(jīng)過高溫高速作用后，其內(nèi)部結(jié)構(gòu)會(huì)發(fā)生變化，導(dǎo)致其中心線發(fā)生偏轉(zhuǎn)，這與本文中的數(shù)值模擬計(jì)算結(jié)果相吻合。部分學(xué)者數(shù)值模擬了尾流的特性，得到結(jié)論：流體內(nèi)部形狀相互作用導(dǎo)致漩渦偏移，將渦流的軸向速度轉(zhuǎn)變?yōu)閲娚淞?，渦旋向后擴(kuò)散增加［29］。實(shí)際的流速場射流仍然是動(dòng)態(tài)復(fù)雜的，需要不斷修正優(yōu)化數(shù)值模型，并且與二維模型的研究結(jié)果有較大的變化，三維模擬仿真更加符合實(shí)際情況。

圖9為發(fā)動(dòng)機(jī)尾流yz截面流速分布云圖，可以清楚地看出，氣流流速在很短距離內(nèi)驟降，從110 m／s降低到32 m／s，此后的溫度梯度變化緩慢，基本保持無差別，沿｀軸方向的高流速作用半徑范圍為0.1 m。結(jié)合圖8（a）和圖9發(fā)現(xiàn)，氣體流速的整體作用影響范圍不大。為了更直觀地解釋和分析流速特性，截取靠近和遠(yuǎn)離噴口的yx兩個(gè)切面進(jìn)行分析。圖10（a）和（b）分別為距離噴口0.2 m和2.1 m的流速分布云圖，可知：靠近噴口的尾流流速徑向作用范圍較小，為細(xì)小的高速射流，以高速流為圓心向四周快速降低；遠(yuǎn)端流速云圖分布較為平緩，流速在每個(gè)方向擴(kuò)展均勻，流速大小無明顯差別。

圖9 發(fā)動(dòng)機(jī)尾流yz截面流速分布云圖

圖10 發(fā)動(dòng)機(jī)尾流xy截面距噴口不同位置的流速分布云圖

結(jié)合以上對(duì)尾流溫度場和流速場的分析，從不同方向的截面和不同角度的視角出發(fā)，參考《MH5001-2013民用機(jī)場飛行區(qū)技術(shù)標(biāo)準(zhǔn)》中關(guān)于除冰安全區(qū)規(guī)劃標(biāo)準(zhǔn)和美國相關(guān)尾流分離技術(shù)標(biāo)準(zhǔn)。其中，人體耐熱溫度為低于323 K，特種設(shè)備作業(yè)安全壞境溫度不超過313 K，最大承受風(fēng)速不超過29 m／s。得出小型發(fā)動(dòng)機(jī)的熱流場分布和人體作業(yè)安全區(qū)域示意圖，如圖11所示［30］。分析可知，尾流的溫度和流速均在y方向和x方向上作用距離相似，y方向作用距離比x方向大，選取y方向結(jié)果代替x方向結(jié)果。圖11為yz平面的尾流熱流場區(qū)域分布圖，由于尾流在y軸方向的作用范圍較小，以與z軸夾角為30°的扇形劃分。

圖11 發(fā)動(dòng)機(jī)尾流熱流場區(qū)域示意圖

由圖11可以看出，溫度區(qū)域劃分為z＜0.7 m，y＜0.35 m，為高溫區(qū)域，溫度范圍為569～976 K；0.9＜z＜2.5 m，y＜0.5 m，為較高溫區(qū)域（人體灼傷危險(xiǎn)區(qū)域），溫度范圍為323～569 K；z＜2.5 m，y＜0.7 m，為低溫區(qū)域（人體安全作業(yè)區(qū)域），溫度小于313 K。流速區(qū)域劃分為z＜0.2 m，為高流速區(qū)域，流速范圍為77～110 m／s；0.2＜z＜2.1 m，y＜0.3 m，為低流速區(qū)域，流速小于20 m／s；z＜2.1 m，沒有流速作用，為環(huán)境流速區(qū)域。人體安全可作業(yè)區(qū)域?yàn)閦＞2.1 m，y＞0.5 m。高溫高流速危險(xiǎn)區(qū)域?yàn)閥＜0.2 m，z＜0.2 m。

3 結(jié)論

1）尾流出口熱流場呈逐漸向外擴(kuò)散的錐體的分布趨勢，xy截面場的分布符合正態(tài)分布，中心為最高點(diǎn)，逐漸向四周降低。溫度場的作用范圍相對(duì)流速場較為明顯；流體的流速經(jīng)尾噴管后，射流主要作用范圍在噴口附近，隨后降低梯度較大。

2）尾流溫度場的高溫核心范圍為z＜0.7 m，y＜0.35 m，x＜0.2 m，核心溫度為569.1～976.0 K；人體安全作業(yè)范圍為z＞2.5 m，y＞0.5 m，x＞0.5 m，溫度不超過323 K。

3）尾流流場的流速核心高速區(qū)域?yàn)閦＜0.2 m，y＜0.3 m，x＜0.2 m，流速范圍為77～100 m／s；人體安全作業(yè)區(qū)域?yàn)?.5＜z＜2.1 m，y＞0.3 m，流速最高不超過20 m／s；z＞2.1 m的流體與環(huán)境混合，無明顯影響。

本研究方法可為民航大型客機(jī)發(fā)動(dòng)機(jī)尾流場研究提供參考。接下來，根據(jù)波音737－800及A320等主流機(jī)型發(fā)動(dòng)機(jī)參數(shù)，開展對(duì)主流民航發(fā)動(dòng)機(jī)的建模和尾流特性分析，深入探究民航客機(jī)在慢車除冰狀態(tài)下尾流場安全區(qū)的動(dòng)態(tài)分布。