飛機不同攻角下近場飛行繞流場數值模擬

周鵬宇，柳文林，賀 俊，王 棟，謝名飛，李 洋

(1.海軍航空大學， 山東 煙臺 264001； 2. 92326部隊， 廣東 湛江 524005；3. 92236部隊，廣東 湛江 524005； 4. 94831部隊，福建 南平 354300； 5. 95172部隊，長沙 410000)

1 引言

飛機起飛和進場著陸階段是執行飛行任務的事故易發多發階段，據統計一半左右的空難發生在起飛和進場著陸階段[1]。高速軍用飛機因起降速度大，在特殊情況下中斷起飛和正常著陸時都需要放出減速傘以保障起降安全。據統計，高速戰機起降過程中沖出跑道事件多與減速傘異常有關。更進一步的分析表明，在放傘系統正常情況下，減速傘不能正常放出的事件多伴有放傘速度大、攻角小的現象，有資料提出這可能與飛機尾部的渦流、反壓區有關，但并無詳細的理論依據。隨著科學技術的進步，減速傘故障分析和放傘仿真研究也得到了進一步發展。沈燕良等[2]采用氣動控制技術，提出了減速傘放傘機構的氣壓和電控改進方法，研制了一種新型阻力傘機構，大大提高了放傘過程的可靠性，解決了無法拋傘的常見故障。朱武峰等[3]為了解決阻力傘氣動故障，構建了一套氣動仿真模型，分析了多狀態參數變化對放傘的影響作用，最后提出了一系列減少故障率的措施。

本文以計算流體力學為基礎，對不同速度和攻角下的飛機尾流場進行數值仿真，并綜合考慮地面效應和尾噴流的影響。通過仿真結果分析，驗證網格生成和仿真的合理性，為飛機安全起降、科研人員開展相關仿真研究提供理論參考。

2 幾何建模與網格劃分

2.1 幾何建模

飛機起落架細小精密部件較多，會導致網格質量降低，其結構對飛機整體氣動和流場影響很小[4]，為簡化研究，提高整體網格質量，網格劃分不考慮起落架結構。

原始模型包含了大量的幾何文件，實際網格劃分和數值仿真過程是在飛機表面殼體基礎上開展，因此需對體積模型進行表面提取，保證提取的面結構是一個連續的無縫面，便于后期網格劃分[5]。

對飛機模型前后襟翼做進一步的重構和修正。飛機在起飛和降落時通過放下襟翼增加機翼面積，改變翼型彎度，延緩機翼氣流分離，提高低速狀態下升力，有利于飛機快速起飛和緩慢著陸，保證飛行安全[6]。襟翼放下與否對飛機的外流場會產生一定影響，為提高仿真精度，模擬飛機真實地面滑行狀態，對機翼前后緣按實際比例分別分割出前緣襟翼和襟副翼，根據飛機起降規范，分別向下偏轉，偏轉角度為δ前襟=23°，δ襟副=18°。放襟翼前后的飛機模型如圖1。

圖1 飛機模型示意圖

2.2 網格劃分

網格是仿真計算域內一系列的離散點。計算流體力學通過離散控制方程，數值計算得到各個網格節點上的物理特性數據，得到的速度、壓力和溫度等數據即為需要的仿真的數值解[7]。網格劃分是數值計算的重要前處理工作，占據了整個仿真周期時間的80%左右，網格質量的高低是影響計算結果的關鍵因素[8]。該機型外形結構復雜，采用Delaunay非結構三角形網格生成技術，非結構網格對不規則復雜構型自適應能力較好。Delaunay三角形網格的顯著優點是能夠使每一個三角形網格最小角盡可能大，使之接近于等邊三角形單元，可顯著提高非結構網格生成效率[9]。

飛機表面劃分三角形非結構網格，網格尺寸為機身尺寸的0.05%左右。飛機尾椎部分是減速傘傘艙位置，是本文仿真關注的重點區域，需要單獨對該區域進行加密操作，保證后期流場仿真結果更加準確，尾椎部分網格尺寸值設為整機網格尺寸的一半。遠場部分的尺寸根據機身整體尺寸劃分。飛機長約21 m，翼展14 m，停機高度6 m。遠場尺寸取機身尺寸的20倍左右，設定為400 m×300 m×100 m的長方體。飛機正下方地面劃分60 m×60 m地效加密區，尾流場設置6個間距1 m的長方形縱切面加密區，每個長方形縱切面尺寸為4 m×2 m。圖2給出了機身和遠場各區域網格分布，表1給出了各區域面網格尺寸。

圖2 網格分布示意圖

表1 各區域面網格尺寸(m)Table 1 Mesh size of each area(m)

在使用求解器數值計算時，需要用湍流模型處理，這里所指的湍流通常是充分發展的湍流，多應用于高雷諾數流動中[10]。在近壁區，因黏性作用，雷諾數不大，湍流發展不充分，所以在近壁區需要采取特殊邊界層處理，在劃分網格時要完成此工作。本仿真網格邊界層第一層厚度取10-5m，最大邊界層網格數設為100，增長率1.3，邊界衰減系數0.9。最終生成的空間體網格類型包括四面體、金字塔和棱柱體等3大類。其中，金字塔和棱柱體網格集中在邊界層，四面體網格集中分布于整個遠場。進一步，以戰機模型及其表面網格為整體，攻角自5°起，以1°為單位依次抬升至10°攻角。每變化1°攻角，采用相同的體網格生成法重新劃分飛機邊界層網格和空間體網格，網格檢查無負體積，質量滿足計算要求。5°～10°攻角下Z=0 m縱切面體網格可視化如圖3所示，網格不同顏色代表網格體積大小。

圖3 Z=0 m縱切面體網格圖

5°～10°攻角下空間體網格類型及數量如表2所示。

表2 空間體網格統計(萬)

3 理論模型與數值方法

3.1 理論模型

三維雷諾平均Navier-Stokes(reynolds averaged navier-stokes,RANS)方程廣泛應用于各種計算流體力學工程實例中[11]。其方法可將控制方程的瞬時運動分為平均運動和脈動運動，根據不同的經驗和實驗方法假設雷諾應力項，進而求解RANS方程[12]。引入湍流模型的RANS方程，可較為準確地計算飛行器氣動力，是目前流場數值計算的主流方法。控制方程采用守恒形式，RANS方程如下：

(1)

(2)

(3)

p=ρRT

(4)

式(1)～(4)中：p為靜壓；τi j為應力張量；cp為比熱容；T為溫度；k為流體的傳熱系數；ST為粘性耗散項；R為普適氣體常數。

選擇能夠滿足對雷諾應力約束條件的可實現k-ε湍流模型，可保持與真實湍流一致。該湍流模型在實測數值計算中，性能表現更優于標準k-ε模型和RNGk-ε模型[13]，并加入了數學約束，Cμ不再是常量，改為隨應變變化率的函數，即：

(5)

3.2 數值方法

采用可較好計算可壓縮流動的密度基穩態求解器，其對高速流場的結構捕捉能力較強[14]。控制方程為守恒形式，采用耦合隱式算法和Roe離散格式。耦合求解器耦合了流動方程和能量方程，隱式算法可加速收斂，收斂速度快、效果好，相比于耦合顯式算法，需要的內存更大，對計算機的配置要求較高[15]。空間離散化方面，梯度項采用基于單元格的最小二乘法，流體對流項采用2階迎風格式，湍動能項和湍流耗散項均采用1階迎風格式。方程離散化后進行耦合求解，庫朗數取5，湍動能項和湍流耗散項的亞松弛因子均取0.8，湍流黏性項和固體壁面的亞松弛因子取1。實際仿真結果表明，以上參數設置效果良好，殘差收斂較快。

3.3 邊界條件

邊界條件方面，遠場尺寸約為飛機尺寸的20倍左右，可滿足仿真計算條件。戰機攻角自5°～10°變化，對應的滑跑速度自300～255 km/h變化。計算域遠場前面為遠場入口，設為壓強入口邊界條件，上面、左右面和后面為遠場出口，設為壓強出口邊界條件，遠場5個面均為自由空間；進氣道截面設為壓強出口邊界條件，尾噴口截面設為壓強入口邊界條件；機身表面和地面均設為無滑移固體壁面條件。環境氣壓為海平面標準大氣壓101 325 Pa，環境溫度為288 K。所有壓強入口條件的流動方向定義為垂直于邊界，壓強出口條件的流動方向定義為導自臨近單元，即緊鄰壓力出口的網格單元格流動方向。不同速度、攻角下各邊界條件具體參數見表3。

表3 邊界條件參數

執行仿真計算的工作站CPU為Intel Xeon，運行內存為192G，12核24線程，每度攻角計算耗時約1 d。求解結束后，將數據文件導入可視化軟件顯示云圖、描繪氣流流線和提取分析仿真數據。

4 仿真結果及對放傘影響分析

4.1 全局分析

在分析尾流場之前先對全機外流場進行仿真結果的可視化，便于把握整體流場分布情況。以7°攻角降落狀態為例，圖4為壓強全局視圖，可以看出:大于110 kPa的高壓區集中分布于尾噴口附近區域，此為尾噴口高壓氣流區；101～110 kPa的次高壓區主要集中在機頭、尾椎與尾椎下方地面、機身下腹部表面與其正下方地面、尾噴流噴射到地面的后方區域，其中飛機下方地面高壓區為地面效應作用的結果；小于101 kPa的低壓區主要集中在機翼上表面，機翼上下表面的壓差使飛機產生一定的升力。

圖4 壓強全局視圖

圖5為氣流流線圖。機身周圍可見明顯的機翼渦、邊條渦和翼尖渦。其中，邊條渦是該型戰機的一大特點，飛機前翼采用細長的三角邊條翼，邊條渦正是由邊條翼產生，是戰機飛行過程中前翼氣流分離所產生的渦流結構。該戰機飛行產生的邊條渦可以增大飛行攻角。翼尖渦的形成是由于飛機以一定的攻角滑行或飛行時會產生一定升力，飛機下翼面的壓強大于上翼面，下方的高壓氣流沿著翼尖位置向上滾卷流動到上方低壓區域，加上飛機本身具有前進速度，氣流向后流動，因此形成了較長的翼尖渦流。

圖5 氣流流線圖

4.2 軸向分析

軸向方面，給出尾噴口截面后X軸正方向4 m內特性，圖6給出了5°～10°攻角條件下Z=0 m縱切面壓強云圖。

圖6 Z=0 m縱切面壓強云圖

圖7給出了以兩尾噴口外截面Y對稱軸與尾椎表面交點為起點，沿X軸正向每隔0.15 m的表壓特性曲線。橫坐標為距離起點的X方向位置，縱坐標代表表壓大小。表壓也叫相對壓強，是以大氣壓為起點；絕對壓強是直接作用在物體表面的所有壓強，以絕對真空為起點。表壓和絕對壓強之間的關系為：

圖7 軸向表壓特性曲線

pg=pabs-p0

(6)

式中：pg為表壓；pabs為絕對壓強；p0為標準大氣壓。

圖7中，X=1.65 m之前為尾椎表面取點，X=1.65 m之后為X軸水平方向取點。從圖像和曲線分析可知，傘艙附近存在一定的反壓區(壓力與放傘方向相反，不利放傘即為反壓)，表壓范圍為-320～1 140 Pa(含尾椎傘艙表面)；不同攻角下，壓強均以尾椎末端點為反壓峰值點，沿軸向的兩側方向總體呈下降趨勢；在X<0.6 m的區域內為低于標準大氣壓的負壓區，且負壓值先下降后上升；攻角越大、速度越小，尾椎末端點反壓越小，不同攻角尾流場區域同位置點壓差值較小，范圍為-210～140 Pa，軸向反壓區大小與速度、攻角相關性不大。

4.3 橫向分析

橫向方面，給出以尾椎末端點為中心、Z軸正反兩方向各1.2 m范圍內特性，圖8給出了5°～10°攻角條件下過尾椎末端點橫切面壓強云圖。

圖8 橫切面壓強云圖

圖9給出了以尾椎末端點為中心、Z軸方向每隔0.15 m的表壓特性曲線。橫坐標為距離尾椎末端點的橫向位置，縱坐標為表壓大小。

圖9 橫向表壓特性曲線

從圖像和曲線分析可知，傘艙附近反壓區表壓范圍為360～1 140 Pa；不同攻角下，壓強均以尾椎末端點為反壓峰值點，沿Z軸兩側方向呈下降趨勢，兩側下降曲線較為對稱；攻角越大、速度越小，尾椎末端點反壓越小，不同攻角尾流場區域同位置點壓差值較小，范圍為-70～150 Pa，橫向反壓區大小與速度、攻角相關性不大。

4.4 垂向分析

垂向方面，給出以尾椎末端點為中心、Y軸上下兩方向各1.2 m范圍內特性，圖10給出了5°～10°攻角條件下過尾椎末端點垂切面壓強云圖。

圖10 垂切面壓強云圖

圖11給出了以尾椎末端點為中心、Y軸方向每隔 0.15 m的表壓特性曲線。橫坐標為距離尾椎末端點的垂向位置，縱坐標為表壓大小。

圖11 垂向表壓特性曲線

從圖像和曲線分析可知，傘艙附近反壓區表壓范圍為250～1 140 Pa；不同攻角下，壓強均以尾椎末端點為反壓峰值點，沿Y軸兩側方向呈下降趨勢，在尾椎末端點上方壓強持續下降直至接近外界大氣壓，在下方0.6 m后趨于平緩，受地效作用為穩定高壓區；攻角越大、速度越小，尾椎末端點反壓越小，不同攻角尾流場區域同位置點壓差值較小，范圍為-160～150 Pa，縱向反壓區大小與速度、攻角相關性不大。

4.5 氣流流線分布

圖12為5°～10°攻角條件下尾流場氣流流線。尾部流場區氣流平順，未見明顯分離渦流。

圖12 尾流場氣流流線示意圖

4.6 放傘影響分析

飛機的減速傘系統及其剎車制動系統可使戰機在著陸和中止起飛時縮短滑跑距離，保證人機安全。減速傘艙在戰機尾椎處，傘艙門面積約0.13 m2，放傘時艙門打開機構向上打開，如圖13所示。

圖13 傘艙門圖

傘艙下部裝有鎖閉機構，艙內固定減速傘鎖。放傘時，由飛行員操作座艙內放傘按鈕執行該動作。中斷起飛或降落時，要求飛機速度不大于300 km/h情形下，才允許放減速傘。放傘時先彈出1 m2引導傘，張滿的引導傘牽引拉出25 m2主傘，主傘張滿后使飛機減速，縮短滑跑距離。

需打開減速傘艙門時，通過按壓座艙內開關按鈕，電流進入艙門鎖電磁鐵線圈，扇形輪轉動，使搖臂脫開，傘艙門在艙門打開機構的作用下迅速向上打開。在艙門關閉時彈簧處于壓縮狀態，壓縮力為1 380 N。打開艙門鎖時，彈簧伸張，在彈力作用下打開艙門。

傘艙附近反壓區表壓為-320～1 140 Pa(含尾椎傘艙表面)，傘艙門面積約0.13 m2，取表壓最大值與傘艙門面積相乘得到約148 N的流場最大反壓力，而減速傘艙門打開機構結構的彈簧壓縮力為1380 N，尾流產生的最大反壓也遠小于艙門打開時的彈簧作用力，對放傘影響作用較小可忽略。

5 結論

1) 在戰機尾椎傘艙附近劃分的4 m×1.2 m×1.2 m尾流場區域內，存在一定的反壓區，表壓為-320～1 140 Pa(含尾椎表面)。

2) 在不同姿態的地面滑行條件下，壓強均以尾椎末端點為反壓峰值點，沿軸向、橫向和垂向的兩側方向呈下降趨勢。

3) 攻角越大、速度越小，尾椎末端點反壓越小，尾流場區域同位置點壓差值為-210～150 Pa，反壓區大小與速度、攻角相關性不大。

4) 飛機尾流產生的反壓遠小于艙門打開時的彈簧作用力，尾部流場反壓區對放傘影響作用較小可忽略。