吳耕宇 白文 瞿麗霞

摘要：飛機設(shè)計過程中需要掌握各種舵面偏轉(zhuǎn)氣動特性。采用計算流體力學(xué)（CFD）技術(shù)計算舵面偏轉(zhuǎn)氣動力時，不同的舵面偏轉(zhuǎn)角度通常需要生成不同的網(wǎng)格、網(wǎng)格變形或重構(gòu)，這一過程需要耗費大量時間。一種簡化的計算方法是蒸騰邊界方法，通過在舵面的邊界條件中增加一個法向擾動速度以模擬舵面偏轉(zhuǎn)，從而可以在不進(jìn)行網(wǎng)格變形或重新生成網(wǎng)格的情況下計算不同舵面偏轉(zhuǎn)角度下的氣動力系數(shù)。本文采用跨聲速巡航標(biāo)模（TCR）鴨翼偏轉(zhuǎn)構(gòu)型驗證蒸騰邊界方法計算舵面偏轉(zhuǎn)氣動力的有效性，迎角范圍為-6°～10°，鴨翼偏轉(zhuǎn)角范圍為-15°～10°。通過對網(wǎng)格變形和蒸騰邊界兩種方法進(jìn)行氣動力計算對比，結(jié)果顯示，蒸騰邊界方法可以在網(wǎng)格保持不變的情況下獲得有效的氣動力計算數(shù)據(jù)，兩種方法多數(shù)工況的計算結(jié)果基本一致，只有在迎角和舵偏角產(chǎn)生疊加效應(yīng)使得舵面相對迎角較大的工況下，兩種方法的計算結(jié)果存在一定差異。

關(guān)鍵詞：計算流體力學(xué)；歐拉方程；舵面偏轉(zhuǎn)；蒸騰邊界條件方法；跨聲速巡航標(biāo)模

中圖分類號：V221.3文獻(xiàn)標(biāo)識碼：ADOI：10.19452/j.issn1007-5453.2021.06.003

在進(jìn)行飛機性能計算和操穩(wěn)特性分析時，需要大批量的舵面偏轉(zhuǎn)氣動力數(shù)據(jù)。在概念設(shè)計階段，通常采用DATCOM、渦格法、面元法等方法快速獲得批量數(shù)據(jù)。通過采用計算流體力學(xué)（CFD）計算，可以提高舵面偏轉(zhuǎn)氣動力計算結(jié)果的精準(zhǔn)度。進(jìn)行CFD計算首先需要生成計算網(wǎng)格，這也是最耗費人工的環(huán)節(jié)。飛機襟翼、升降舵、方向舵、鴨翼等各種舵面偏轉(zhuǎn)角及偏轉(zhuǎn)組合將產(chǎn)生上千個計算構(gòu)型，采用常規(guī)的CFD計算方法，每個計算構(gòu)型均需要重新生成計算網(wǎng)格、網(wǎng)格變形或重構(gòu)（以下統(tǒng)稱網(wǎng)格方法）；采用蒸騰邊界方法則無須重新生成計算網(wǎng)格，從而提高了獲得大批量舵面偏轉(zhuǎn)氣動力數(shù)據(jù)的效率。

CFD計算的控制方程主要有無黏流動的歐拉（Euler）方程和黏性流動的納維-斯托克斯（Navier-Stokes，N-S）方程，其中歐拉方程所需網(wǎng)格較小，計算速度較快，適用于氣動力的快速分析。傳統(tǒng)的歐拉方程CFD計算假定流體在物面的法向速度為0，即流體總是沿著物體表面流動，具有現(xiàn)實的物理意義。對于運動的物面邊界來說，流體微團(tuán)也相應(yīng)跟隨物面運動?；谶@一思路，Lighthill首先提出“蒸騰邊界”方法[1]，最初用于翼型厚度變化的模擬。國內(nèi)外多位專家學(xué)者對蒸騰邊界方法進(jìn)行了各種方面的應(yīng)用研究。N.L.Sanker等[2]使用具有后緣擺動襟翼NACA 64A006機翼（1°振幅）、F-5戰(zhàn)斗機機翼和大展弦比機翼（俯沖運動10°）；B.K.Bharadvaj等[3]采用大展弦比機翼算例進(jìn)行了6°副翼偏轉(zhuǎn)的計算；C.C.Fisher等[4]采用蒸騰邊界方法進(jìn)行氣動彈性分析；C.H. Stephens等[5]采用AGARD 445.6機翼和BACT機翼驗證蒸騰邊界方法在氣動彈性模擬中的有效性；J.C.Timothy等[6]將蒸騰邊界方法應(yīng)用于定常和非定常表面幾何擾動的模擬；R.M. Kolonay等[7]使用蒸騰邊界方法進(jìn)行幾何優(yōu)化使得構(gòu)型阻力減小；K.W. Huckriede等[8]通過蒸騰邊界方法模擬不穩(wěn)定氣流或水流；陸志良等[9]采用蒸騰邊界方法使用F-5機翼算例進(jìn)行舵面偏轉(zhuǎn)2°的模擬；宋萬強等[10]采用蒸騰邊界方法進(jìn)行了非定常動導(dǎo)數(shù)的計算；朱海濤等[11]曾經(jīng)使用TCR模型進(jìn)行非定常動導(dǎo)數(shù)的計算和分析，得到較為詳細(xì)的結(jié)果。

蒸騰邊界方法多用于氣動彈性變形研究領(lǐng)域，而在舵面偏轉(zhuǎn)氣動力計算領(lǐng)域研究較少，只有單一偏轉(zhuǎn)角度、單一工況的個別算例的研究[9]。本文采用TCR標(biāo)模研究多迎角、多舵面偏轉(zhuǎn)角工況的舵面偏轉(zhuǎn)氣動力計算的有效性及其適用范圍，為該方法應(yīng)用于飛機概念設(shè)計階段性能計算和操穩(wěn)特性分析提供參考依據(jù)。

1計算方法

實際計算中，根據(jù)舵面偏轉(zhuǎn)軸和偏轉(zhuǎn)角度，通過網(wǎng)格點原始坐標(biāo)和位移求出物面網(wǎng)格點法向ni，網(wǎng)格點法向速度值vni通過網(wǎng)格點運動速度得到。所有舵偏角工況計算過程相同。

2計算模型

跨聲速巡航標(biāo)模（transonic cruiser， TCR）是歐盟“飛機概念設(shè)計操控特性模擬”（SimSAC）項目[12]研制的用于評估舵面偏轉(zhuǎn)效應(yīng)和動導(dǎo)數(shù)計算的標(biāo)模。該模型在俄羅斯中央空氣流體動力學(xué)研究院（TsAGI）進(jìn)行了一系列風(fēng)洞試驗[13]，風(fēng)洞試驗來流速度為40m/s，雷諾數(shù)Re=7.9×105。風(fēng)洞試驗采用的模型的尺寸參數(shù)見表1。

在本文的研究中，首先進(jìn)行模型網(wǎng)格生成，分別生成鴨翼無偏轉(zhuǎn)的Euler網(wǎng)格和N-S方程網(wǎng)格，以及鴨翼偏轉(zhuǎn)-15°、-10°、-5°、5°和10°的Euler網(wǎng)格。各種計算網(wǎng)格的單元數(shù)量見表2。

TCR標(biāo)模鴨翼無偏轉(zhuǎn)和鴨翼偏轉(zhuǎn)-15°的表面網(wǎng)格分布如圖1所示。

計算所采用的來流條件與試驗相同，即馬赫數(shù)Ma= 0.1179，雷諾數(shù)Re=7.9×105。采用Euler方程進(jìn)行蒸騰邊界方法和網(wǎng)格方法CFD計算，計算迎角為-6°～10°（間隔為2°），鴨翼偏轉(zhuǎn)角為-15°～10°（間隔為5°）。

3計算結(jié)果分析

3.1定常氣動力計算結(jié)果

首先進(jìn)行定常流動計算以驗證計算網(wǎng)格和流場求解的正確性。定常流動法向力系數(shù)CN和俯仰力矩系數(shù)Cm計算結(jié)果及其與風(fēng)洞試驗的對比如圖2所示。

從圖中可以看出，法向力系數(shù)與風(fēng)洞試驗數(shù)據(jù)基本吻合。Euler方程計算的俯仰力矩系數(shù)在-6°～2°迎角工況下與風(fēng)洞試驗數(shù)據(jù)吻合較好；其他工況，Euler方程和RANS方程計算結(jié)果與試驗數(shù)據(jù)存在系統(tǒng)性偏差。鑒于本文主要是分析網(wǎng)格方法和蒸騰邊界方法計算舵面偏轉(zhuǎn)氣動力的差異，比對的是差量，計算結(jié)果滿足本文研究要求。

3.2舵面偏轉(zhuǎn)氣動力計算結(jié)果

采用風(fēng)洞試驗鴨翼偏轉(zhuǎn)氣動力數(shù)據(jù)進(jìn)行對比。不同的鴨翼偏轉(zhuǎn)角下，法向力系數(shù)的蒸騰邊界方法計算結(jié)果、網(wǎng)格方法計算結(jié)果及其與風(fēng)洞試驗數(shù)據(jù)對比結(jié)果如圖3所示。

可以看出，在不同的鴨翼偏轉(zhuǎn)角情況下，除鴨翼偏轉(zhuǎn)角較大且迎角較大時有微小差異以外，蒸騰邊界方法和網(wǎng)格方法計算得到的法向力系數(shù)均高度重合；蒸騰邊界方法計算得到的法向力系數(shù)與風(fēng)洞試驗數(shù)據(jù)也基本吻合。事實上，從圖3中也可以看出，法向力系數(shù)在不同的鴨翼偏轉(zhuǎn)角情況下的變化較小。

在飛機性能計算、操穩(wěn)分析等分析中，更關(guān)心舵面偏轉(zhuǎn)氣動力與相應(yīng)狀態(tài)下0°舵面偏轉(zhuǎn)角氣動力的差量。因此將計算所得到的鴨翼偏轉(zhuǎn)氣動力與相應(yīng)狀態(tài)下0°鴨翼偏轉(zhuǎn)角氣動力進(jìn)行差量的求取，進(jìn)行進(jìn)一步分析。在不同的鴨翼偏轉(zhuǎn)角下，法向力系數(shù)與相應(yīng)狀態(tài)下0°鴨翼偏轉(zhuǎn)角的法向力系數(shù)差量ΔCN的蒸騰邊界方法計算結(jié)果、網(wǎng)格方法計算結(jié)果及其與風(fēng)洞試驗數(shù)據(jù)對比結(jié)果如圖4所示。

以上結(jié)果顯示，在鴨翼偏轉(zhuǎn)角為負(fù)值時，網(wǎng)格方法和蒸騰邊界方法的法向力系數(shù)差量基本吻合；在鴨翼偏轉(zhuǎn)角為正值且迎角為正值時，隨著迎角增大，兩種方法計算的法向力系數(shù)存在明顯差異。原因是在這種情形下，迎角和舵偏角產(chǎn)生疊加效應(yīng)，使得舵面相對迎角較大，利用物面法向速度模擬物面變形對氣流造成的吸/吹影響，已難以真實反映氣流方向與物面法向角度的相對變化[10]。雖然兩種方法差量計算結(jié)果基本一致，但與風(fēng)洞試驗數(shù)據(jù)相比均存在一定的系統(tǒng)性偏差，考慮到法向力系數(shù)試驗和計算的精度一般為0.01量級，而且這種CFD與風(fēng)洞試驗的系統(tǒng)性偏差是一個專門的研究領(lǐng)域，本文不進(jìn)行進(jìn)一步分析。

不同迎角和不同鴨翼偏轉(zhuǎn)角下，俯仰力矩系數(shù)蒸騰邊界方法計算結(jié)果、網(wǎng)格方法計算結(jié)果及其與風(fēng)洞試驗數(shù)據(jù)的對比結(jié)果如圖5所示。

在不同的鴨翼偏轉(zhuǎn)角下，俯仰力矩系數(shù)與相應(yīng)狀態(tài)下0°鴨翼偏轉(zhuǎn)角的俯仰力矩系數(shù)差量ΔCm的蒸騰邊界方法計算結(jié)果、網(wǎng)格方法計算結(jié)果及其與風(fēng)洞試驗數(shù)據(jù)對比結(jié)果如圖6所示。

從圖5和圖6中可以看出，除鴨翼偏轉(zhuǎn)-15°且迎角為負(fù)的情況，以及鴨翼偏轉(zhuǎn)10°且迎角為10°的情況（這兩種情況下，鴨翼相對迎角分別為最大值-21°和20°），其他情況的俯仰力矩系數(shù)蒸騰邊界方法和網(wǎng)格方法均具有較高的吻合程度，計算數(shù)據(jù)與風(fēng)洞試驗數(shù)據(jù)吻合程度也較好。

4結(jié)論

采用蒸騰邊界方法計算了TCR標(biāo)模不同迎角下鴨翼偏轉(zhuǎn)法向力和俯仰力矩系數(shù)，與網(wǎng)格方法計算結(jié)果以及風(fēng)洞試驗數(shù)據(jù)進(jìn)行了對比分析，結(jié)論如下：

（1）蒸騰邊界方法可以在計算網(wǎng)格保持不變的情況下，在一定工況范圍內(nèi)獲得有效的舵面偏轉(zhuǎn)氣動力數(shù)據(jù)，大大提高采用CFD技術(shù)生成大批量數(shù)據(jù)的效率。

（2）除鴨翼相對迎角較大的情況，蒸騰邊界方法與網(wǎng)格方法計算的法向力系數(shù)和俯仰力矩系數(shù)絕對量和差量均高度吻合。

（3）本文主要研究工作采用Euler方程方法實施氣動力計算，法向力系數(shù)計算結(jié)果與風(fēng)洞試驗數(shù)據(jù)基本一致；除鴨翼相對迎角較大的情況，俯仰力矩系數(shù)計算結(jié)果隨飛機迎角增大與風(fēng)洞試驗數(shù)據(jù)存在系統(tǒng)偏差；兩種方法計算的法向力和俯仰力矩系數(shù)差量與風(fēng)洞試驗數(shù)據(jù)存在一定的系統(tǒng)偏差。

（4）在飛機概念設(shè)計階段需要快速迭代選型，因此需要快速生成大批量舵面偏轉(zhuǎn)氣動力數(shù)據(jù)，蒸騰邊界方法計算效率高，適用于Euler方程CFD計算方法。

參考文獻(xiàn)

[1]Lighthill M J. On displacement thickness[J]. Journal of Fluid Mechanics，1958，4（4）：383-392.

[2]Sankar N L. Ruo S Y，Malone J B. Application of surface transpiration in computational aerodynamics[C]//Proceedings ofAIAA24thAerospace Sciences Meeting，1986：1-10.

[3]Bharadvaj B K. Computation of steady and unsteady control surface loads in transonic flow[J]. AIAA Journal，1991，29（11）：1906-1911.

[4]Fisher C C，Arena A S. On the transpiration method for efficient aeroelastic analysis using an euler solver[C]//AIAA Meeting Papers on Disc，1996：623-631.

[5]Stephens C H，Arena A S. Application of the transpiration method for aeroservoelastic prediction using CFD[R]. AIAA Paper，1998：3092-3099.

[6]Timothy J C，Charles R.Transpiration boundary condition for computational fluid dynamics solutions in noninertial reference frames[J]. Journal ofAircraft，2003，41（5）：1252-1255.

[7]Kolonay R M，Thompson E D，Camberos J A. Active control of transpiration boundary conditions for drag minimization with an euler CFD solver[C]//Proceedings of 48th AIAA/ ASME/ASCE/AHS/ASC Structures，Structural Dynamics，and Materials Conference，Honolulu：AIAA，2007：1-14.

[8]Huckriede K W，Koop A H，Hospers J M. Finite-volume method with transpiration boundary conditions for flow about oscillating wings[C]//Proceedings of 48th AIAA/ASME/ASCE/ AHS/ASC Structures，Structural Dynamics，and Materials Conference，Orlando：AIAA，2010：1-20.

[9]陸志良.跨聲速機翼操縱面定常、非定常氣動力計算[J].南京航空航天大學(xué)學(xué)報， 1999， 31（1）： 92-96. Lu Zhiliang. Computation of steady and unsteady transonic flow for a wing with control surface deflections[J]. Journal of Nanjing University of Aeronautics & Astronautics， 1999， 31（1）： 92-96. （in Chinese）

[10]宋萬強，徐悅.一種基于歐拉方程的動導(dǎo)數(shù)簡化計算方法[J].航空科學(xué)技術(shù)， 2017， 28（2）： 39-42. Song Wanqiang， Xu Yue. A simplified method for computing dynamic derivatives based on Euler equations[J]. Aeronautical Science & Technology， 2017， 28（2）： 39-42. （in Chinese）

[11]朱海濤，白文.基于非定常CFD方法跨聲速民用飛機標(biāo)模俯仰氣動導(dǎo)數(shù)計算研究[J].航空科學(xué)技術(shù)， 2017， 28（7）： 19-24. Zhu Haitao， Bai Wen. Static and dynamic pitch derivatives of transonic full aircraft configuration with unsteady CFD Method[J]. Aeronautical Science & Technology， 2017， 28（7）： 19-24.（in Chinese）

[12]Rizzi A. Modeling and simulating aircraft stability and control. the SimSAC project[J]. Progress in Aerospace Sciences，2011，47：573-588.

[13]Khrabrov A，Kolinko K，Zhuk A，et al. Wind tunnel test

report[R]. SimSAC Report，2010.（責(zé)任編輯陳東曉）

An Efficient Approach for Computation of Control Surface Deflection Effect

Wu Gengyu，Bai Wen，Qu Lixia

Chinese Aeronautical Establishment，Beijing 100012，China

Abstract： Computing the aerodynamics of various control surface deflect configurations with CFD method generally needs to generate a new computational mesh， or deform， or reconstruct the mesh for each configuration， which is quite tedious and time consuming. An option method is to use transpiration boundary condition to simulate the movement of control surface deflection by adding a normal grid velocity on the boundary. This approach needs no extra mesh operations， which makes the large batch of computation tasks quite efficient. There is a lack of research on the applicability of transpiration boundary condition method for different angle of attack of aircraft and various control surface deflection angles. In this paper， Transonic Cruiser （TCR） model is selected to study the effectiveness of transpiration boundary condition method for the simulation of canard deflection effect. The angles of attack studied are from -6°to 10°， and canard deflect angles are from -15°to 10°. The computations are carried out using inviscid Euler flow equations. It is shown that the computed normal force and pitch moment coefficients using transpiration boundary condition method and mesh regeneration method match each other if the canard relative angle of attack is not so large. The computational results are also compared with wind tunnel test data， the agreement of normal force coefficients are acceptable， and the agreement of pitch moment coefficients are reasonable if the canard relative angle of attack is not so large. The comparison of aerodynamic coefficient delta value due to canard deflections further confirmed the effectiveness of the transpiration boundary condition method.

Key Words： computational fluid dynamics； Euler equations； control surface deflection； transpiration boundary condition method； transonic cruiser model