直升機機身氣動特性CFD計算研究進展

龍海斌，劉正勝，吳裕平

(1.中國直升機設計研究所 總體氣動研究室，江西 景德鎮 333001;2.32382部隊，北京 100072))

0 引言

在直升機研制過程中，機身氣動特性數據是飛行性能、品質和載荷等的設計輸入，因此準確地獲得機身氣動特性數據非常重要。目前主要通過風洞試驗和CFD計算兩種方法來獲得機身氣動特性數據。其中風洞試驗方法在航空航天等領域應用廣泛，具有很成熟的數據修正和準確性驗證方法，因此在工程應用領域有很高的可信度。但是風洞試驗之前需要設計和制造風洞試驗模型，風洞試驗過程中又受到風洞的檔期安排和測量系統的穩定性等的影響，因此成本比較高，獲得機身氣動特性數據的周期相對比較長。CFD計算是近幾十年來發展起來的數值模擬方法。得益于計算機技術和數值模擬方法的進步，機身氣動特性CFD計算的速度越來越快，可有效地加快直升機型號的研制進度。但是，目前尚未形成統一的機身氣動特性標準CFD計算方法，CFD計算結果的驗證與確認是公認的難題。針對近年來直升機機身氣動特性CFD計算方面的進展情況，本文首先介紹直升機機身計算模型的常見構型和氣動外形，之后分別從網格類型選擇、網格數量控制、求解方法選取和計算結果分析等方面進行整理與分析，最后對機身氣動特性CFD計算未來的發展進行展望。

1 機身模型簡介

在直升機機身氣動特性CFD計算中，計算模型的氣動外形通常與風洞試驗模型保持一致。用于理論研究和計算方法驗證的標準機身模型通常只包含機身和尾梁部分，如ROBIN機身模型、NUAA機身模型等。型號研制過程中的機身模型包含機身、主槳轂、起落架、尾梁、平尾、垂尾和尾槳轂等。其中主槳轂通常為靜止狀態，部分機身模型的主槳轂為旋轉狀態。部分機身模型的尾槳轂上帶有尾槳葉。武裝型直升機機身模型還包含短翼或掛梁、外掛航炮、導彈和火箭彈等。部分直升機機身模型還帶有天線、光電吊艙等外掛物。在部分機身氣動特性CFD計算模型中還包含風洞支撐機構，計算域的邊界為風洞試驗段的內表面。

目前已經在飛行和試飛過程中的直升機中，機身與尾梁的過渡情形大致可分為兩類：一是機身與尾梁平緩過渡，類似于ROBIN機身模型,如美國的UH-60“黑鷹”直升機、S-97共軸高速直升機、V-22傾轉旋翼機以及武裝直升機等；二是尾梁比較細長，或者機身尾部有艙門，如ROBIN mod7機身模型、小型無人機和運輸直升機等。這些機身的腹部到尾梁的過渡段有比較大的流動分離，因此壓差阻力比較大。機身模型如圖1-圖3。

圖1 ROBIN標準機身模型

圖2 ROBIN mod7標準機身模型(反裝)[4]

圖3 機身風洞試驗模型[5]

本文中統計的機身類型包括無人直升機、單旋翼運輸直升機、單旋翼武裝直升機、共軸式直升機、傾轉旋翼機等。按機身模型包含的部件來分，機身模型可分為簡單的標準機身模型和型號光機身模型，帶槳轂、起落架和平尾等部件的機身模型，部分機身模型的槳轂為旋轉狀態。

2 網格類型

對流體計算域進行網格劃分是機身氣動特性CFD計算的第一步，而在網格劃分之前首先要選取網格類型。目前常用的網格類型有結構網格、非結構網格和上述兩者組合的混合網格。結構網格在拓撲結構上相當于矩形區域內的均勻網格，每一層網格上的節點數都相等。針對外形簡單的流體域，采用結構網格可以生成比較好的貼體網格，同時求解過程中的數值耗散比較小，因此外形相對比較簡單的標準機身模型流場計算多用結構網格，或者進行理論研究需要對流場進行精細化計算時，也采用結構網格。工程應用中的直升機機身包含主槳轂、起落架、平尾和垂尾等部件，部分機型還包含短翼和外掛武器等，難以進行結構網格劃分，因此基本上都采用非結構網格。非結構網格包含四面體非結構網格、笛卡爾網格等。部分情況下在機身表面附近劃分邊界層網格，以進一步提高對機身表面附近復雜流動的模擬能力。但是劃分邊界層網格容易導致網格質量降低，同時網格數量增長很多。多塊網格綜合了結構網格和非結構網格的優點，但是各塊之間的數據傳遞存在一定的難度。在部分情況下，可以將非結構網格轉化為多面體網格，或者直接劃分得到多面體網格。對于相同的機身外形和計算域，多面體網格的數量比較少，可加快CFD計算速度。部分計算過程中需要考慮槳轂旋轉，采用嵌套網格方法來計算旋轉槳轂和機身周圍的流場。在統計的52個劃分網格的算例中，非結構網格、結構網格和混合網格的算例數量和所占百分比如表1所示。從表中可以看出，非結構網格占61.54%，說明非結構網格在機身氣動特性CFD計算中應用比較廣泛，其中采用非結構四面體的算例為22個。在網格劃分過程中ICEM軟件應用比較多。

表1 各類型網格數和所占比例

針對邊界層網格劃分，有13個算例劃分了邊界層網格，占總的算例數量的1/4，說明在直升機機身氣動特性CFD計算過程中邊界層網格劃分相對比較少。這是由于在直升機機身氣動特性計算過程中通常帶有主槳轂、尾槳轂和起落架等部件，有比較多的拉桿等細小零件，表面形狀比較復雜，生成邊界層網格的難度比較大；同時，在生成邊界層網格時，總的網格數量增長比較多，網格整體質量變差，后續求解計算時可能難以收斂。在劃分邊界層網格的算例中，第一層網格厚度多為0.05mm或0.005mm。結構網格、笛卡爾網格、非結構網格及其邊界層網格分別如圖4、圖5和圖6所示。

圖4 結構網格示意圖[6]

圖5 笛卡爾體網格示意圖[8]

圖6 非結構及邊界層網格示意圖[15]

3 網格數量

根據CFD計算方法的基本原理，網格數量越多，則對流體域的劃分越細，CFD計算模擬的流動與實際流動的情況越接近。因此在進行理論研究時，為了得到相關區域的流場細節、渦的變化情況等，通常將網格劃分得比較細密，網格數量可能達到千萬量級。在部分精細化計算時也不考慮網格數量的限制。但是在實際工程應用中，部分情況下機身氣動特性計算完成之后，還需要進行旋翼/機身干擾等更復雜的流場計算，加上項目研制周期、可用計算資源等的限制，對機身計算域劃分的網格相對比較粗糙。在部分情況下，劃分的網格數量多，由于相鄰區域的疏密控制等問題，在部分區域容易生成質量很差的網格，造成計算過程中收斂速度變慢甚至難以收斂。

在各類直升機機身氣動特性計算過程中，劃分的網格數量在幾十萬到幾千萬之間。在統計的38個公布的網格劃分算例中，超過1000萬的算例有12個，約占總數的31.58%，網格數量在100萬量級的占50%。各量級網格數和所占比例如表2所示。通常來講，計算模型的部件數量越多，則劃分的網格數量就相對比較多。

表2 各量級網格數和所占比例

部分研究人員針對同一機身模型，劃分不同的網格數量進行氣動特性CFD計算，以研究網格數量變化對CFD計算結果的影響。文獻[21]針對同一個計算模型，劃分了4套不同數量的網格。根據CFD計算結果，在劃分的網格數量范圍內，機身的升力系數和阻力系數均隨著網格數量的增加而不斷減小，見圖7。文獻[4]針對ROBIN mod7標準機身模型，劃分三套不同數量的結構網格，之后采用OVERFLOW求解器進行氣動特性計算，計算過程中選擇S-A湍流模型，得到的0°和-5°攻角時的阻力系數如表3所示。從表中的結果可以看出，隨著網格數量增加，ROBIN mod7標準機身模型的阻力系數略有減小。

圖7 網格數量對氣動特性影響示意圖[21]

表3 阻力網格數量相關性分析[4]

4 計算方法與湍流模型

在各類直升機機身氣動特性CFD計算過程中，常用的計算方法是求解N-S方程方法。早期部分研究人員常用面元法和求解Euler方程的方法。其中面元法是通過在機身表面網格上布置流動的奇點來求解氣動問題，計算結果的精度比較粗糙，目前應用比較少。求解Euler方程的方法沒有考慮空氣的粘性，因此求解速度相對比較快，對計算資源的需求也比較少。在工程應用中方便后續采用相同的方法對旋翼/機身干擾、旋翼/地面等復雜流場進行計算。但是由于沒有考慮空氣的粘性，因此不能計算空氣的摩擦阻力，對比較復雜的流動區域也難以模擬。求解N-S方程的方法主要有雷諾平均(RANS)方法、大渦模擬(LES)方法、直接數值模擬(DNS)方法等。近年來還發展了部分混合方法，比如RANS和LES混合方法等。此外還有格子玻爾茲曼方法(LBM)等可用于機身氣動特性CFD計算之中。受限于可利用的計算資源和求解周期限制等因素，目前雷諾平均(RANS)方法應用比較廣泛，而大渦模擬(LES)、直接數值模擬(DNS)和格子玻爾茲曼方法(LBM)等在工程領域很少使用。在采用雷諾平均(RANS)方法求解N-S方程的過程中，需要增加方程來使得方程組封閉，即引入湍流模型。目前常用的湍流模型有零方程B-L模型、一方程S-A模型、兩方程

k

-

ε

和

k

-

ω

模型等。其中一方程S-A和兩方程

k

-

ω

湍流模型應用比較多。在統計的50個采用雷諾平均(RANS)方法的算例中：采用S-A湍流模型的有22個，占44%；采用

k

-

ω

湍流模型的有20個，其中采用SST

k

-

ω

湍流模型的算例數量為13個，占26%。這是由于S-A湍流模型能給出比較準確的氣動特性結果，同時在求解過程中只增加了一個方程，求解的速度相對比較快。

表4 主要的湍流模型應用對比

文獻[22]分別采用面元法和CFD計算方法計算了ROBIN機身的繞流場，得到了機身表面壓力系數，并與風洞試驗結果進行了對比分析，部分結果如圖8所示。其中CFD計算方法采用結構網格，求解器為CFL3D，求解過程中采用B-L湍流模型。結果表明CFD計算方法在分離流動和粘性流動中的模擬能力更強。文獻[23]采用ICEM軟件對ANSAT-M光機身模型進行了以六面體為核心的四面體網格劃分，之中在求解N-S方程中分別采用Spalart-Allmaras(S-A), SST

k

-

ω

(SST)和改型

k

-

kl

-

ω

(TR)湍流模型，計算結果如表5所示。從表中與風洞試驗結果的對比情況來看，SST

k

-

ω

湍流模型得到的阻力系數和升力系數與風洞試驗結果比較接近。

圖8 0°攻角時結果對比圖[22]

表5 不同湍流模型計算結果與風洞試驗結果對比

5 與試驗結果對比

在直升機機身氣動特性CFD計算完成之后，需要對計算結果的準確性和有效性等進行驗證和確認。目前主要采用與風洞試驗結果進行對比分析的方法進行驗證與確認。在大部分情況下對比的是整個機身模型的力和力矩系數，部分重點關注位置的表面壓力系數等。在風洞試驗過程中可以采用PIV方法測量流場細節，同時在CFD計算結果也能比較方便地提取出流線和渦結構等流場細節信息。上述兩者可進行比較與分析。文獻[24]對某常規單旋翼直升機、某無人直升機和某共軸式直升機(圖9中分別用C、W和G表示)光機身氣動特性進行了CFD計算，并與風洞試驗結果進行了對比分析。其中攻角和側滑角的范圍均為-16°至16°。文獻[25]對某常規單旋翼直升機、某無人直升機和某共軸式直升機(圖10中分別用C、W和G表示)的機身大攻角和大側滑角氣動特性進行了CFD計算，并與風洞試驗結果進行了對比分析。文獻[26]對某型武裝直升機大側滑角機身氣動特性進行了CFD計算，并與風洞試驗結果進行了對比分析(圖11)。文獻[35]采用Helios求解器對ROBIN mod7機身表面壓力系數進行了CFD計算。計算過程中劃分了不同的網格，其中：網格1為棱柱網格，數量為870萬；網格2為四面體網格，網格數量為2400萬；網格5為棱柱網格，網格數量為1550萬。從目前CFD計算結果與風洞試驗值的對比情況來看(圖12)，當攻角或側滑角變化時，機身氣動特性CFD計算結果的變化趨勢與風洞試驗值基本一致，但在具體數值上有一定的誤差，而且部分情況下的誤差比較大。目前CFD計算方法可實現0°～360°范圍內攻角或側滑角變化時的機身氣動特性計算。機身表面外形變化比較平緩的區域，CFD計算得到的機身表面壓力系數與風洞試驗值比較接近；而在主減速器整流罩尾部、機身與尾梁的過渡段等外形變化比較劇烈的區域，機身表面壓力系數的CFD計算結果與風洞試驗值相差比較大。

圖9 光機身氣動特性CFD計算結果與風洞試驗對比圖[24]

圖10 機身大攻角氣動特性CFD計算結果與風洞試驗對比圖[25]

圖11 機身大側滑角氣動力CFD計算結果與風洞試驗對比圖[26]

圖12 機身表面壓力系數CFD計算結果與風洞試驗對比圖[35]

6 總結與討論

通過對多種類型直升機機身氣動特性的CFD計算過程進行整理與分析，包括網格類型選取、網格數量控制、求解方法選擇和計算結果分析等，可得出如下結論：

1)目前CFD計算方法已經應用到各種類型直升機機身氣動特性計算之中。計算狀態的攻角和側滑角等基本上覆蓋了直升機正常飛行中的大部分狀態。大部分機身氣動特性CFD計算的結果與風洞試驗結果變化趨勢一致，具有比較高的準確度和可信度。

2)在機身氣動特性CFD計算過程中，目前還沒有形成統一的標準計算方法。其中的網格類型選取、網格數量控制和求解方法的選擇等均沒有形成相對統一的標準方法。

3)提高機身氣動特性CFD計算結果的準確度需要進行綜合考慮，僅僅增加網格數量或選用新的湍流模型等單個因素的改進難以大幅度提高CFD計算結果的準確度。

未來直升機機身氣動特性CFD計算的發展趨勢如下：

1)計算速度進一步加快。隨著計算機技術、超級計算中心和數值模擬方法的發展，將來機身氣動特性CFD計算將達到一秒鐘計算一個狀態的速度量級。

2)CFD計算結果的驗證與確認工作進一步發展。機身風洞試驗過程中增加部分流場細節方面的測量工作，同時風洞試驗段的內表面尺寸和支撐機構等的外形數模將進一步公開，以便CFD計算結果有更多的參考與驗證。

3)形成統一的標準計算方法。國家或行業層面將形成相對統一的標準指導規范，促進機身氣動特性CFD計算標準化發展。