陳 功，唐志共，*，鄧 晨，王文正

(1. 國(guó)防科技大學(xué) 空天科學(xué)學(xué)院，長(zhǎng)沙 410073；2. 中國(guó)空氣動(dòng)力研究與發(fā)展中心，綿陽(yáng) 621000)

0 引 言

升力式再入飛行器由于具有可重復(fù)使用、短時(shí)多次空天往返、任務(wù)載荷靈活、可實(shí)現(xiàn)水平著陸等優(yōu)點(diǎn)，必將成為人類未來(lái)構(gòu)建空天一體能力的重要平臺(tái)。但與航空飛行器不同的是，升力式再入飛行器再入飛行速域大、飛行高度變化劇烈、飛行環(huán)境復(fù)雜，導(dǎo)致氣動(dòng)特性變化極大[1]。為滿足采用單一布局形式的再入飛行器在整個(gè)再入階段的飛行需求，必須對(duì)飛行器的氣動(dòng)布局進(jìn)行精細(xì)優(yōu)化，以滿足不同飛行階段的控制要求。實(shí)際中考慮到再入過(guò)程中的防熱、容積率和配平等約束，升力式再入飛行器通常采用大鈍頭體、厚翼型和小展弦比的升力體布局。這種布局形式極易導(dǎo)致側(cè)向氣動(dòng)力在跨聲速段出現(xiàn)非常強(qiáng)烈的非線性特性[2]，給飛行控制系統(tǒng)設(shè)計(jì)帶來(lái)困難，嚴(yán)重影響再入飛行器的水平起降能力和飛行安全[3]。為解決這一問(wèn)題，早在HL-20[4]研制過(guò)程中，氣動(dòng)工程師就采用了扭轉(zhuǎn)體襟翼來(lái)消除舵面偏轉(zhuǎn)對(duì)飛行器側(cè)向氣動(dòng)特性的不利干擾。Lee[5]和Chaudhary[6]等在研究X-33和X-37兩種飛行器的側(cè)向非線性動(dòng)力學(xué)特性的基礎(chǔ)上，采用Weissman判據(jù)分析了飛行器的側(cè)向穩(wěn)定性。日本Jeong等[7]通過(guò)研究HYFLEX[8]飛行器布局參數(shù)改變引起的飛行器穩(wěn)定性和操縱性的變化，分析了再入飛行器布局參數(shù)對(duì)控制性能的敏感性。國(guó)內(nèi)唐偉[9-10]、周軍[11]、王穎[12]等在分析高超聲速飛行器布局參數(shù)變化對(duì)飛行器側(cè)向氣動(dòng)特性影響的基礎(chǔ)上，提出了減小側(cè)向非線性的布局改進(jìn)和減小操縱耦合的建議。

本文以解決實(shí)際工程中飛行器布局設(shè)計(jì)關(guān)鍵問(wèn)題為目標(biāo)，發(fā)展了一種以減小跨聲速段側(cè)向氣動(dòng)非線性的EGO布局優(yōu)化設(shè)計(jì)方法。該方法以飛行器主要的氣動(dòng)敏感布局參數(shù)為優(yōu)化對(duì)象，通過(guò)建立Kriging代理模型來(lái)替代傳統(tǒng)的數(shù)值計(jì)算，實(shí)現(xiàn)對(duì)不同參數(shù)布局氣動(dòng)特性的快速預(yù)測(cè)，并可在優(yōu)化過(guò)程中對(duì)代理模型進(jìn)行自適應(yīng)修正，可以提高布局優(yōu)化效率。通過(guò)對(duì)一種升力式再入概念飛行器機(jī)翼布局的設(shè)計(jì)驗(yàn)證，優(yōu)化后布局明顯降低了飛行器在跨聲速段的側(cè)向氣動(dòng)非線性，滿足了飛行控制系統(tǒng)的設(shè)計(jì)要求。

1 概念飛行器及側(cè)向氣動(dòng)特性

1.1 概念飛行器基礎(chǔ)布局

如圖1所示，設(shè)計(jì)了一種常規(guī)布局形式的升力式再入飛行器為研究對(duì)象，進(jìn)行氣動(dòng)力計(jì)算與側(cè)向氣動(dòng)特性分析。該飛行器采用了鈍頭錐+細(xì)長(zhǎng)外形機(jī)體的“混合機(jī)體”氣動(dòng)布局，機(jī)翼為兼顧高低速的小展弦比切尖雙三角翼。為了提高升力，機(jī)腹設(shè)計(jì)接近于一塊平板，機(jī)身兩邊也向后以一定曲率延展成為豐滿橫向截面。機(jī)體后部設(shè)計(jì)有大面積體襟翼，用于縱向配平。機(jī)翼后緣布置有襟副翼用于俯仰和滾轉(zhuǎn)控制。機(jī)翼外側(cè)末端加裝了垂直尾翼，尾翼后緣有航向控制舵面。

圖1 再入飛行器氣動(dòng)布局Fig. 1 Aerodynamic configuration of a typical reentry vehicle

1.2 基礎(chǔ)布局氣動(dòng)特性計(jì)算

為獲得再入飛行器跨聲速段的氣動(dòng)特性，這里采用數(shù)值計(jì)算方法獲得了飛行器的側(cè)向氣動(dòng)數(shù)據(jù)。如圖2所示，計(jì)算采用的是非結(jié)構(gòu)混合網(wǎng)格，用三棱柱模擬附面層，四面體模擬空間流場(chǎng)同性區(qū)域，通過(guò)金字塔網(wǎng)格過(guò)渡三棱柱和四面體網(wǎng)格區(qū)域，總網(wǎng)格規(guī)模約為6百萬(wàn)。計(jì)算大氣模型為理想氣體模型，湍流模型取為SST湍流模型，解算器采用二階迎風(fēng)離散格式。

圖2 數(shù)值計(jì)算非結(jié)構(gòu)網(wǎng)格Fig. 2 Unstructured grid used for numerical simulations

橫航向靜穩(wěn)定系數(shù)導(dǎo)數(shù)Clβ、Cnβ是飛行器控制系統(tǒng)設(shè)計(jì)時(shí)主要關(guān)心的側(cè)向氣動(dòng)變化量。這里采用如下差分公式計(jì)算滾轉(zhuǎn)和偏航力矩系數(shù)導(dǎo)數(shù)：

式中Clβ,3°、Cnβ,3°分 別表示側(cè)滑 β= 3° 時(shí)計(jì)算得到的飛行器滾轉(zhuǎn)和偏航力矩系數(shù)。

圖3和圖4給出了基準(zhǔn)布局飛行器在跨聲速段、不同飛行速度時(shí)的橫航向靜穩(wěn)定系數(shù)導(dǎo)數(shù)的計(jì)算結(jié)果。由圖中的計(jì)算結(jié)果可知，在跨聲速區(qū)域再入飛行器基準(zhǔn)布局的側(cè)向靜穩(wěn)定系數(shù)導(dǎo)數(shù)隨著迎角的增大，顯示出明顯的非線性。特別是滾轉(zhuǎn)力矩系數(shù)導(dǎo)數(shù)在超聲速（Ma =1.2）時(shí)與高亞聲速（Ma =0.5）小迎角狀態(tài)的極性相反。為分析造成以上側(cè)向非線性的原因，圖5給出了飛行器在Ma =0.9不同飛行狀態(tài)時(shí)的壓力云圖：有側(cè)滑時(shí)，在機(jī)翼和翼尖垂尾的共同作用下，在靠近垂尾的機(jī)翼外側(cè)部分形成明顯的非對(duì)稱低壓區(qū)，且隨著迎角增大，低壓區(qū)明顯變化。從而導(dǎo)致飛行器側(cè)向氣動(dòng)特性出現(xiàn)較大的非線性。

圖3 偏航力矩導(dǎo)數(shù)變化曲線Fig. 3 Derivative of the yaw moment coefficient varies with the angle of attack

圖4 滾轉(zhuǎn)力矩導(dǎo)數(shù)變化曲線Fig. 4 Derivative of the roll moment coefficient varies with the angel of attack

圖5 基準(zhǔn)布局壓力云圖（Ma = 0.9，β = 3°）Fig. 5 Pressure contour of the baseline configuration(Ma = 0.9，β = 3°)

2 氣動(dòng)布局優(yōu)化方法

2.1 布局優(yōu)化變量

垂尾和機(jī)翼的干擾是導(dǎo)致再入飛行器基準(zhǔn)布局側(cè)向氣動(dòng)非線性的主要原因。因此，將小展弦比雙三角機(jī)翼的內(nèi)外側(cè)的兩個(gè)后掠角Φ1、Φ2以及垂直尾翼的后掠角Φ3作為布局優(yōu)化參數(shù)（圖6），以提高再入飛行器的側(cè)向氣動(dòng)線性度為目標(biāo)進(jìn)行布局優(yōu)化設(shè)計(jì)。

圖6 布局優(yōu)化參數(shù)Fig. 6 Schematics of the configuration optimization variables

2.2 布局優(yōu)化目標(biāo)函數(shù)

為在更大速度范圍內(nèi)提高飛行器的側(cè)向氣動(dòng)特性，選擇以下目標(biāo)函數(shù)：

式中ΔClβi和ΔCnβi（i=1,2,3）分別是在Ma=0.5、0.9、1.2時(shí)由式（3）計(jì)算得到的力矩系數(shù)線性偏差量：

其中Aα+B和Cα+D分別是由三個(gè)迎角α1、α2、α3狀態(tài)點(diǎn)處的Clβ和Cnβ采用最小二乘法擬合得到的線性函數(shù)。考慮到跨聲速段已經(jīng)接近再入飛行器的飛行末段，飛行器的實(shí)際飛行迎角較小，因此將三個(gè)計(jì)算的迎角樣本點(diǎn)設(shè)計(jì)為α1= 0°、α2= 5°、α3= 10°。為計(jì)算目標(biāo)函數(shù)，一種布局樣本就需要計(jì)算如表1所示的18個(gè)飛行狀態(tài)，這對(duì)于常規(guī)的Navier-Stokes計(jì)算代價(jià)很大。為減少計(jì)算資源消耗，提高設(shè)計(jì)效率，這里采用一種基于Kriging代理模型的自適應(yīng)優(yōu)化方法對(duì)飛行器布局進(jìn)行優(yōu)化。

表1 單樣本目標(biāo)函數(shù)計(jì)算狀態(tài)Table 1 Simulation parameters for the calculation of one sample’s objective function

2.3 布局優(yōu)化方法

多目標(biāo)優(yōu)化要求全局優(yōu)化算法具有目標(biāo)函數(shù)估值次數(shù)少的特點(diǎn)[13]。在飛行器氣動(dòng)布局設(shè)計(jì)時(shí)，由于飛行器典型特征點(diǎn)的氣動(dòng)數(shù)據(jù)計(jì)算時(shí)長(zhǎng)較多導(dǎo)致這個(gè)矛盾非常突出[14]?；贙riging代理模型的自適應(yīng)迭代優(yōu)化方法由于同時(shí)考慮模型預(yù)測(cè)值和預(yù)測(cè)方差來(lái)指導(dǎo)生成新的樣本，能夠大幅提高模型的最優(yōu)逼近精度，并可以大大減小目標(biāo)函數(shù)的估值次數(shù)。這里采用具有最大期望改善（Expected Improvement，EI）迭代準(zhǔn)則的EGO優(yōu)化方法完成對(duì)再入飛行器布局的優(yōu)化設(shè)計(jì)[15]。

EI是指在一個(gè)給定特征點(diǎn)的響應(yīng)值目標(biāo)改進(jìn)的概率。作為優(yōu)化搜索的質(zhì)量因素，來(lái)權(quán)衡局部搜索和全局搜索。設(shè)點(diǎn)x出的預(yù)測(cè)值為，標(biāo)準(zhǔn)差為s，預(yù)測(cè)結(jié)果滿足標(biāo)準(zhǔn)分布，則可以給出EI的公式為：

其中fmin=min(y1,···,yn)，φ和φ分別表示標(biāo)準(zhǔn)正態(tài)分布函數(shù)和正態(tài)分布密度函數(shù)。上式兼顧了搜索位置處的目標(biāo)預(yù)測(cè)值和預(yù)測(cè)值誤差：當(dāng)預(yù)測(cè)誤差很大時(shí)，s就很大，引導(dǎo)EGO算法向預(yù)測(cè)誤差大的位置搜索；當(dāng)目標(biāo)預(yù)測(cè)誤差小時(shí)，值就很大，引導(dǎo)EGO算法向極小值附近搜索。因此，EGO算法通過(guò)在EI值最大處增加迭代點(diǎn)，來(lái)兼顧全局性和收斂性[16]。對(duì)于優(yōu)化過(guò)程中的約束：

可以轉(zhuǎn)化為概率函數(shù)：

這里采用遺傳算法求解EGO模型的EI值，建立再入飛行器布局優(yōu)化的流程如圖7所示。

圖7 布局優(yōu)化方法流程Fig. 7 Flowchart of the aerodynamic configuration optimization process

3 氣動(dòng)布局優(yōu)化結(jié)果及分析

規(guī)定再入飛行器機(jī)翼的待優(yōu)化布局參數(shù)的變化范圍如下：

由拉丁方試驗(yàn)設(shè)計(jì)方法在以上三個(gè)待優(yōu)化參數(shù)構(gòu)成的三維參數(shù)空間內(nèi)確定12個(gè)初始布局樣本點(diǎn)。圖8給出了不同布局樣本在設(shè)計(jì)空間內(nèi)的分布，布局樣本較為均勻的分布在三維設(shè)計(jì)空間，并取中間樣本作為基準(zhǔn)布局。

圖8 優(yōu)化布局樣本分布Fig. 8 A scatter plot of samples for the optimization configuration

采用與上節(jié)相同的N-S解算器，分別計(jì)算所有樣本的側(cè)向氣動(dòng)特性，獲得每個(gè)布局樣本由式（2）計(jì)算的優(yōu)化目標(biāo)函數(shù)值，建立飛行器布局優(yōu)化參數(shù)和側(cè)向氣動(dòng)線性度偏差的預(yù)測(cè)代理模型。并考慮以下升阻比約束：

經(jīng)過(guò)兩輪迭代優(yōu)化，分別得到優(yōu)化布局Opt1和Opt2。圖9給出了上述兩種優(yōu)化布局由式（2）和式（3）計(jì)算得到的側(cè)向力矩系數(shù)線性偏差量。

圖9 優(yōu)化布局側(cè)向力矩系數(shù)線性偏差Fig. 9 Linear deviation of the lateral moment for optimized configurations

圖10和圖11分別給出了Opt1布局和Opt2布局的橫航向靜穩(wěn)定力矩系數(shù)導(dǎo)數(shù)與基準(zhǔn)布局的對(duì)比結(jié)果。表2分別給出了Ma =0.9時(shí)，基準(zhǔn)布局、Opt1布局和Opt2布局的機(jī)翼布局參數(shù)和升阻比結(jié)果。

圖10 布局滾轉(zhuǎn)力矩系數(shù)導(dǎo)數(shù)Fig. 10 Derivative of the roll moment coefficient

圖11 偏航力矩系數(shù)導(dǎo)數(shù)Fig. 11 Derivative of the yaw moment coefficient

由上述計(jì)算結(jié)果可知，Opt1布局與Opt2布局的側(cè)向氣動(dòng)線性度相較基準(zhǔn)布局明顯提升。由Opt1布局參數(shù)的優(yōu)化結(jié)果可知，增大內(nèi)側(cè)機(jī)翼后掠角、減小外段機(jī)翼后掠角和垂直尾翼后掠角可以明顯提高再入飛行器在跨聲速段的側(cè)向氣動(dòng)線性度，但飛行器的升阻比會(huì)略有下降。為此，Opt2布局減小內(nèi)側(cè)翼后掠角、增大外翼和垂尾后掠角，同時(shí)兼顧了側(cè)向氣動(dòng)線性度和升阻比要求。

圖12給出了Opt2布局Ma= 0.9、β= 3°時(shí)的壓力云圖。對(duì)比圖5和圖12，可知優(yōu)化后布局隨著迎角增大，在靠近垂尾的機(jī)翼外側(cè)仍會(huì)出現(xiàn)不對(duì)稱的低壓區(qū)。但相較基準(zhǔn)布局，左右機(jī)翼低壓區(qū)范圍的不對(duì)稱性明顯減弱，從而增加了飛行器側(cè)向氣動(dòng)特性的線性度，提高了飛行器的可控性。

表2 布局優(yōu)化結(jié)果Table 2 Aerodynamic configuration optimization results

圖12 Opt2布局壓力云圖（Ma = 0.9，β = 3°）Fig. 12 Pressure contour for the configuration Opt2(Ma = 0.9，β = 3°)

4 結(jié) 論

本文主要針對(duì)升力式再入飛行器跨聲速段側(cè)向氣動(dòng)非線性的問(wèn)題，發(fā)展了一種基于Kriging代理模型的自適應(yīng)布局優(yōu)化設(shè)計(jì)方法。通過(guò)對(duì)一種概念飛行器的機(jī)翼布局參數(shù)的優(yōu)化結(jié)果表明，優(yōu)化后布局可以明顯降低飛行器在跨聲速段的氣動(dòng)非線性，同時(shí)滿足升阻比約束，證明了所發(fā)展的方法是基本可行和有效的。為進(jìn)一步提高再入飛行器布局設(shè)計(jì)水平，未來(lái)可進(jìn)一步研究建立氣動(dòng)布局參數(shù)到控制性能要求的關(guān)聯(lián)模型和評(píng)估準(zhǔn)則，提高氣動(dòng)布局滿足控制性能要求的水平，實(shí)現(xiàn)氣動(dòng)與控制的一體化快速設(shè)計(jì)，實(shí)現(xiàn)飛行器綜合性能優(yōu)化。