飛翼翼型高維目標(biāo)空間多學(xué)科綜合優(yōu)化設(shè)計(jì)

鄭傳宇, 黃江濤,*, 周 鑄, 劉 剛, 高正紅, 許 勇(. 中國(guó)空氣動(dòng)力研究與發(fā)展中心, 四川 綿陽(yáng) 6000; . 西北工業(yè)大學(xué) 翼型葉柵空氣動(dòng)力國(guó)防科技重點(diǎn)實(shí)驗(yàn)室, 陜西 西安 7007)

飛翼翼型高維目標(biāo)空間多學(xué)科綜合優(yōu)化設(shè)計(jì)

鄭傳宇1, 黃江濤1,*, 周 鑄1, 劉 剛1, 高正紅2, 許 勇1
(1. 中國(guó)空氣動(dòng)力研究與發(fā)展中心, 四川 綿陽(yáng) 621000; 2. 西北工業(yè)大學(xué) 翼型葉柵空氣動(dòng)力國(guó)防科技重點(diǎn)實(shí)驗(yàn)室, 陜西 西安 710072)

傳統(tǒng)的優(yōu)化方法在處理高維多目標(biāo)問(wèn)題上面臨多重困難，針對(duì)高維多目標(biāo)優(yōu)化問(wèn)題，進(jìn)行合理的目標(biāo)降維是一種實(shí)用高效的方法。本文首先討論、研究了PCA目標(biāo)降維算法在高維多目標(biāo)應(yīng)用中的思路，通過(guò)典型高維度目標(biāo)測(cè)試函數(shù)中的應(yīng)用，驗(yàn)證了降維方法的有效性，進(jìn)一步將該方法推廣應(yīng)用于翼型氣動(dòng)隱身多學(xué)科綜合設(shè)計(jì)中，對(duì)綜合設(shè)計(jì)高維度目標(biāo)空間進(jìn)行主成分分析。利用主成分分析提取主成分并辨識(shí)“冗余”或者不重要的目標(biāo)，將冗余目標(biāo)去除或者作為約束加入到重要目標(biāo)的優(yōu)化過(guò)程中。結(jié)果顯示，目標(biāo)空間降維以后的優(yōu)化設(shè)計(jì)結(jié)果滿足力矩、阻力發(fā)散、巡航升阻比、低速升力特性以及隱身特性等綜合設(shè)計(jì)的要求。進(jìn)一步探討、展望了該方法在飛行器多目標(biāo)、多學(xué)科優(yōu)化設(shè)計(jì)中的應(yīng)用前景。

主分量分析；冗余目標(biāo)；相關(guān)性；高維多目標(biāo)；綜合優(yōu)化設(shè)計(jì)

0 引 言

伴隨高性能計(jì)算機(jī)的發(fā)展，基于Pareto思想的數(shù)值進(jìn)化算法成為解決航空工程中復(fù)雜多目標(biāo)多約束問(wèn)題的最佳途徑之一。該類方法將當(dāng)前解集作為進(jìn)化群體，可以在指定的搜索范圍內(nèi)尋找最優(yōu)解集。具有適用范圍廣，易實(shí)現(xiàn)并行編程求解等優(yōu)勢(shì)。近年來(lái)，利用數(shù)值進(jìn)化算法來(lái)解決航空工程中的多目標(biāo)優(yōu)化問(wèn)題逐漸成為研究熱點(diǎn)。

該方法對(duì)于解決一般多目標(biāo)問(wèn)題效果顯著，但是對(duì)于諸如直升機(jī)旋翼翼型設(shè)計(jì)以及超臨界翼型設(shè)計(jì)等設(shè)計(jì)目標(biāo)數(shù)量較多的問(wèn)題便顯現(xiàn)出不足。總的來(lái)看，當(dāng)目標(biāo)空間維數(shù)達(dá)到4個(gè)及以上時(shí)，該方法便存在諸多局限。首先，Pareto最優(yōu)前沿點(diǎn)的數(shù)量會(huì)隨著目標(biāo)數(shù)量的增加而呈指數(shù)級(jí)增長(zhǎng)，這會(huì)大大增加算法的時(shí)間和空間復(fù)雜度。此外，目標(biāo)維數(shù)增加導(dǎo)致的非支配解數(shù)量劇增使得對(duì)優(yōu)化結(jié)果的辨別和篩選變得異常困難，在極多目標(biāo)優(yōu)化問(wèn)題中，甚至?xí)霈F(xiàn)種群中所有個(gè)體都是非支配解的情況。因此，傳統(tǒng)的基于Pareto思想的數(shù)值進(jìn)化算法在解決高維目標(biāo)空間中優(yōu)化設(shè)計(jì)問(wèn)題上存在搜索過(guò)程緩慢、難以收斂等問(wèn)題?？傮w來(lái)講，高維目標(biāo)空間的優(yōu)化設(shè)計(jì)問(wèn)題是當(dāng)前進(jìn)化多目標(biāo)優(yōu)化領(lǐng)域的研究難點(diǎn)[1-6]。

當(dāng)前解決該問(wèn)題的途徑主要包含以下三個(gè)方面：

1) 針對(duì)高維多目標(biāo)優(yōu)化問(wèn)題，改進(jìn)優(yōu)化算法。比如通過(guò)放寬算法中的Pareto占優(yōu)機(jī)制，或者調(diào)整對(duì)個(gè)體的挑選規(guī)則，來(lái)達(dá)到使算法快速收斂的目的。然而這些改進(jìn)在工程應(yīng)用中是否具有普適性還需要進(jìn)一步研究確定。此外，即便通過(guò)改進(jìn)算法能夠得出最優(yōu)解集，由于目標(biāo)數(shù)目過(guò)多帶來(lái)的計(jì)算量過(guò)大的問(wèn)題依然存在，而且高維目標(biāo)空間優(yōu)化結(jié)果的可視化難以實(shí)現(xiàn)，難以進(jìn)行進(jìn)一步?jīng)Q策。

2) 加權(quán)系數(shù)平均方法。利用靜態(tài)或動(dòng)態(tài)加權(quán)系數(shù)對(duì)設(shè)計(jì)目標(biāo)進(jìn)行綜合評(píng)估，將問(wèn)題轉(zhuǎn)化為單目標(biāo)優(yōu)化問(wèn)題。該方法對(duì)小于3個(gè)目標(biāo)的優(yōu)化設(shè)計(jì)有一定的可行性，但對(duì)于高維多目標(biāo)問(wèn)題來(lái)講可行性較差，主要原因在于眾多目標(biāo)之間存在錯(cuò)綜復(fù)雜的相關(guān)關(guān)系，造成權(quán)系數(shù)選擇的困難。

3) 引入數(shù)學(xué)分析中的降維思想，將高維多目標(biāo)優(yōu)化問(wèn)題進(jìn)行主要成分分析，提取決定問(wèn)題本質(zhì)的主要成分，將冗余目標(biāo)剔除、或者轉(zhuǎn)化為約束條件，在不失問(wèn)題主特征的前提下，將高維多目標(biāo)優(yōu)化轉(zhuǎn)化為低維優(yōu)化問(wèn)題。可以預(yù)見(jiàn)，該類方法對(duì)于實(shí)際工程復(fù)雜問(wèn)題具有重大的理論意義和工程應(yīng)用價(jià)值，該類方法的降維設(shè)計(jì)主要應(yīng)用于模式識(shí)別、信號(hào)和圖象處理、控制理論和其它領(lǐng)域中[7-9]，而在飛行器多目標(biāo)優(yōu)化設(shè)計(jì)中的應(yīng)用研究較少。

本文首先以典型高維多目標(biāo)測(cè)試函數(shù)DTLZ5(3,5)為例，采用結(jié)合小生境技術(shù)的多目標(biāo)粒子群算法，論述了基于主成分分析方法降維的有效性；進(jìn)一步以某飛翼布局飛行器翼型氣動(dòng)多目標(biāo)多約束優(yōu)化設(shè)計(jì)以及氣動(dòng)、隱身多目標(biāo)綜合設(shè)計(jì)為例，說(shuō)明了PCA方法在工程應(yīng)用中的合理性與有效性，進(jìn)一步展望了該類方法在飛行器氣動(dòng)多目標(biāo)優(yōu)化設(shè)計(jì)中的應(yīng)用前景。

1 主成分分析

主成分分析 (Principal Component Analysis， PCA) 是一種將主要因素從復(fù)雜事物中分辨出來(lái)的統(tǒng)計(jì)學(xué)方法，通過(guò)對(duì)復(fù)雜問(wèn)題的主成分分析可以降低問(wèn)題復(fù)雜度，幫助研究者掌握所研究問(wèn)題的主要矛盾。其主要操作方法為給出問(wèn)題中n個(gè)變量的m個(gè)觀察值，通過(guò)PCA規(guī)定的操作得到r(r

對(duì)于n個(gè)樣本，每個(gè)樣本具有m項(xiàng)指標(biāo)：

x1,x2,x3,…,xm

1) 對(duì)于m項(xiàng)指標(biāo)，為避免因量綱不同而可能產(chǎn)生的不合理影響，首先需要對(duì)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理:

2) 計(jì)算指標(biāo)的相關(guān)矩陣R，求解R矩陣的m個(gè)特征值：λ1≥λ2≥λ3≥…λm以及對(duì)應(yīng)的正交化特征向量：

νi=[νi1,νi2,νi3,…,νim],i=1,2,…,m

3) 設(shè)定截?cái)嚅撝礣C(文中取閾值為T(mén)C=0.97)。

yi=νi1x1+νi2x2+νi3x3+…νimxm,

5) 合理選擇各個(gè)主成分。

由上述分析步驟得到的第一主成分，可以被看作是通過(guò)高維目標(biāo)空間中心的一條直線，該直線與空間所有點(diǎn)的距離最小，第二主成分與第一主成分正交，同樣通過(guò)高維目標(biāo)空間中心。經(jīng)過(guò)PCA分析，可以得到數(shù)據(jù)集的相關(guān)矩陣，利用相關(guān)矩陣可以判斷兩個(gè)目標(biāo)之間的關(guān)系是相互一致或者相互沖突；進(jìn)一步基于特征向量的量值進(jìn)行冗余變量的分析、剔除，因此，基于PCA的目標(biāo)降維步驟可以被獲得[13]，如第二部分所示。

2 基于PCA方法的高維多目標(biāo)降維

綜合第一部分PCA算法的分析，基于PCA方法的高維多目標(biāo)降維分為以下幾步[14](圖1)：

1) 將迭代計(jì)數(shù)器初始化為I=0，設(shè)定初始目標(biāo)集M為空集，根據(jù)具體問(wèn)題要求設(shè)定一個(gè)閾值TC，本文中取TC=0.97。

2) 對(duì)于目標(biāo)集中的所有目標(biāo)，隨機(jī)初始化種群，對(duì)每一個(gè)個(gè)體計(jì)算其各個(gè)目標(biāo)對(duì)應(yīng)性能，得到一組數(shù)據(jù)集P。

3) 對(duì)數(shù)據(jù)集P進(jìn)行PCA處理。以上述設(shè)定的閾值為標(biāo)準(zhǔn)，對(duì)各目標(biāo)進(jìn)行取舍。具體操作方法如下：

① 對(duì)目標(biāo)矢量進(jìn)行歸一化處理，計(jì)算其相關(guān)矩陣及其各特征值對(duì)應(yīng)的特征矢量R(i,j)、v(i,j)。接下來(lái)即可根據(jù)得到的特征矢量v(i,j)來(lái)確定主要目標(biāo)。

② 對(duì)于第一特征矢量所包含的元素，選擇最正和最負(fù)的兩個(gè)元素所對(duì)應(yīng)的目標(biāo)進(jìn)入集合M；如果所有元素同號(hào)，選擇兩個(gè)絕對(duì)值最大元素所對(duì)應(yīng)目標(biāo)進(jìn)入集合M。

③ 接下來(lái)檢查閾值，如果閾值滿足預(yù)先設(shè)定的閾值要求，則結(jié)束分析過(guò)程。否則，檢查特征值大小:

若特征值小于0.1，選擇絕對(duì)值最大的元素所對(duì)應(yīng)的目標(biāo)進(jìn)入M。

否則

a、選擇最正和最負(fù)兩個(gè)元素所對(duì)應(yīng)的目標(biāo)進(jìn)入集合M。

b、如果所有元素同號(hào)，選擇兩個(gè)絕對(duì)值最大元素所對(duì)應(yīng)目標(biāo)進(jìn)入集合M。

4) 依次對(duì)各個(gè)特征矢量作分析，直至分析過(guò)的特征矢量所對(duì)應(yīng)權(quán)重滿足閾值停止，并輸出最優(yōu)解集。

對(duì)于高維多目標(biāo)問(wèn)題，經(jīng)過(guò)PCA分析，可以識(shí)別目標(biāo)之間的相互關(guān)系，從而提取出優(yōu)化目標(biāo)中的主要目標(biāo)，剔除冗余目標(biāo)，簡(jiǎn)化了原始設(shè)計(jì)問(wèn)題，為進(jìn)一步?jīng)Q策奠定了基礎(chǔ)。

圖1 PCA分析流程圖Fig.1 Flow chart of PCA

3 基于多目標(biāo)粒子群算法以及PCA方法的高維多目標(biāo)優(yōu)化

3.1 關(guān)于多目標(biāo)粒子群算法的說(shuō)明

本文研究工作采用的多目標(biāo)進(jìn)化策略為課題組開(kāi)發(fā)的基于小生境技術(shù)的粒子群優(yōu)化算法(Niche MOPSO)，在保證全局特性的前提下，為加速標(biāo)準(zhǔn)算法的收斂性，文中對(duì)粒子群學(xué)習(xí)模式引入了向群體中心學(xué)習(xí)模式、向非劣解集中心學(xué)習(xí)模式以及向自身歷史最優(yōu)中心學(xué)習(xí)模式。新型粒子群算法表達(dá)式如下[15]:

Vk+1= ωVk+c1×r1×(XpBest-Xk)+

c2×r2×(XgBest-Xk)+

c3×r3×(Xc-Xk)+

c4×r4×(Xpareto_c-Xk)+

c5×r5×(Xpbest_c-Xk)

其中XpBest_c表示自身歷史最優(yōu)中心，Xpareto_c表示非劣解集最優(yōu)中心，Xc表示群體中心，XgBest表示群體最優(yōu)位置，XpBest表示個(gè)體歷史最優(yōu)位置，Vk表示粒子運(yùn)動(dòng)變化的速度矢量，ω為慣性因子，c1、c2、c3、c4、c5為學(xué)習(xí)因子，r1、r2、r3、r4、r5為[0,1]之間均勻分布的隨機(jī)數(shù)。

對(duì)極小化問(wèn)題而言，Pareto解的定義[16-17]為：對(duì)于可行解X*，當(dāng)且僅當(dāng)不存在可行解X ，使

(a) fi(X)≤fi(X*), i∈{1,…,n};

(b) 至少存在一個(gè)j∈{1,…,n}，使fj(X)

兩個(gè)條件成立時(shí)，則可行解X*為一個(gè)Pareto最優(yōu)解。

需要指出的是為保持種群的多樣性，使Pareto前沿分布更均勻，本文的粒子群優(yōu)化算法中加入了小生境技術(shù)。在更新粒子最優(yōu)時(shí)，如果兩代粒子之間不互相支配，將由概率決定最佳位置是否更新，文中的概率函數(shù)取值為：

P=0.9(I-Gen/Genmax)2

其中，Gen為種群進(jìn)化代數(shù)，Genmax為總進(jìn)化代數(shù)，即選取概率在初始階段較大，加速Paerto前沿的推進(jìn), 在進(jìn)化中逐級(jí)減小，有利于優(yōu)化過(guò)程最終收斂。

3.2 PCA高維目標(biāo)降維有效性測(cè)試

DTLZ測(cè)試函數(shù)集是由Deb特殊設(shè)計(jì)的，其具有已知的Pareto最優(yōu)解，用來(lái)測(cè)試優(yōu)化算法的尋優(yōu)性能[18]。本文采用DTLZ5(3,5)測(cè)試函數(shù)對(duì)基于粒子群算法及PCA主成分分析方法的優(yōu)化效果進(jìn)行有效性測(cè)試。DTLZ5測(cè)試函數(shù)包含若干冗余目標(biāo)，用來(lái)測(cè)試算法對(duì)存在冗余目標(biāo)情況下的尋優(yōu)能力。DTLZ5(3,5)表示該函數(shù)有5個(gè)目標(biāo)，其中3個(gè)為非冗余目標(biāo)，即函數(shù)的真實(shí)Pareto前沿的維度數(shù)。首先基于這5個(gè)目標(biāo)進(jìn)行多目標(biāo)優(yōu)化，得到優(yōu)化結(jié)果對(duì)應(yīng)的三個(gè)非冗余目標(biāo)上的顯示如圖2(a)所示。顯然由于目標(biāo)數(shù)量較多，優(yōu)化結(jié)果難以收斂到真實(shí)最優(yōu)解。使用PCA主成分分析方法處理優(yōu)化結(jié)果，按照第2節(jié)中描述的過(guò)程提取主成分，去除冗余成分。最后，以得到的主成分為目標(biāo)進(jìn)行優(yōu)化，優(yōu)化結(jié)果如圖2(b)所示。

(a) Result of direct optimization

(b) Result of optimization after PCA

DTLZ5(3,5)測(cè)試函數(shù)的真實(shí)Pareto前沿為球面，測(cè)試結(jié)果表明，經(jīng)過(guò)PCA主成分分析提取主成分并去除冗余成分后，算法更易收斂到全局最優(yōu)解，表明PCA主成分分析能夠有效辨識(shí)主成分，達(dá)到對(duì)高維目標(biāo)有效降維的目的。

4 高維多目標(biāo)氣動(dòng)優(yōu)化

文中的研究工作是基于課題組開(kāi)發(fā)的飛行器多學(xué)科、多目標(biāo)集成化設(shè)計(jì)平臺(tái)(AMDEsign)開(kāi)展的，該設(shè)計(jì)平臺(tái)中包含了氣動(dòng)多目標(biāo)模塊(AMP)、多學(xué)科模塊(MDP)、伴隨優(yōu)化模塊(Adjoint)，嵌入氣動(dòng)分析、隱身計(jì)算以及結(jié)構(gòu)重量估算等高、低精度學(xué)科分析模塊，具備氣動(dòng)、隱身、結(jié)構(gòu)等多學(xué)科優(yōu)化設(shè)計(jì)能力。對(duì)于氣動(dòng)設(shè)計(jì)中高維度多目標(biāo)優(yōu)化問(wèn)題，文中基于氣動(dòng)多目標(biāo)模塊(AMP)開(kāi)展研究。

4.1 參數(shù)化方法

在優(yōu)化設(shè)計(jì)中，對(duì)設(shè)計(jì)對(duì)象進(jìn)行參數(shù)化描述是設(shè)計(jì)的基礎(chǔ)。本文選用CST(class function/shape function transformation, CST)參數(shù)化方法作為翼型描述方法。CST參數(shù)化方法中各參數(shù)有明確的幾何意義，并且具有控制參數(shù)少、適應(yīng)性強(qiáng)、精度高等優(yōu)點(diǎn)。翼型CST參數(shù)化方法可分為直接CST方法、擾動(dòng)CST方法以及中弧線疊加厚度CST方法。其中，擾動(dòng)CST方法是將CST方程作為幾何擾動(dòng)疊加到初始翼型上得到新的翼型，由于本文是在初始翼型的基礎(chǔ)上，在給定設(shè)計(jì)范圍內(nèi)進(jìn)行變形以尋找翼型的最優(yōu)幾何外形，所以擾動(dòng)CST方法更適用于本文研究。

CST方法的一般表達(dá)式為[19]：

其中，C(x)為類函數(shù)，S(x)為型函數(shù)，yTE為翼型后緣的y坐標(biāo)。類函數(shù)C(x)的定義如下：

其中，對(duì)于一般的翼型N1和N2分別取0.5和1.0。型函數(shù)S(x)的定義如下：

其中，Si(x)為Bernstein多項(xiàng)式，N為多項(xiàng)式的階數(shù)，Ai為待定系數(shù)。

擾動(dòng)CST方法的一般表達(dá)式為：

式中y0(x)為原來(lái)基準(zhǔn)翼型的坐標(biāo)值，即將CST表達(dá)式得到的值作為擾動(dòng)量疊加到初始翼型得到新外形。

4.2 流場(chǎng)數(shù)值模擬方法及代理模型

對(duì)翼型流場(chǎng)的高可信度模擬是對(duì)翼型氣動(dòng)性能進(jìn)行優(yōu)化的基礎(chǔ)和前提。綜合考慮本文研究的飛翼翼型氣動(dòng)優(yōu)化目標(biāo)，本文選用基于雷諾平均的Navier-Stocks(Reynolds Averaged Navier-Stocks, RANS )方程，并采用Roe格式為空間離散格式，湍流模型采用剪切應(yīng)力輸運(yùn)SST模型。流場(chǎng)計(jì)算網(wǎng)格由程序自動(dòng)生成，并使用多重網(wǎng)格技術(shù)提高收斂速度。此外，為減少優(yōu)化過(guò)程中的計(jì)算量，本文采用Kriging代理模型用于對(duì)優(yōu)化迭代過(guò)程中種群個(gè)體氣動(dòng)性能的預(yù)測(cè)。

4.3 翼型高維多目標(biāo)氣動(dòng)優(yōu)化

對(duì)于高亞聲速巡航的飛行器設(shè)計(jì)來(lái)說(shuō)，一副綜合性能優(yōu)異的機(jī)翼是整個(gè)氣動(dòng)設(shè)計(jì)環(huán)節(jié)中最為重要的部分，而用來(lái)進(jìn)行三維機(jī)翼壓力分布形態(tài)控制的典型剖面翼型設(shè)計(jì)是機(jī)翼設(shè)計(jì)的關(guān)鍵。本文以某飛翼布局飛行器翼型作為基礎(chǔ)翼型進(jìn)行多目標(biāo)氣動(dòng)優(yōu)化設(shè)計(jì)研究，進(jìn)一步說(shuō)明基于PCA降維方法的有效性與可行性。該飛翼布局飛行器翼型具有高升阻比、高隱身特性、良好自配平特性、良好失速特性，翼型形狀如圖3所示。本文在初始翼型較為優(yōu)良的綜合性能基礎(chǔ)上進(jìn)行優(yōu)化得到進(jìn)一步的性能提升。

圖3 初始翼型Fig.3 Initial airfoil

在該飛翼翼型優(yōu)化設(shè)計(jì)中主要關(guān)注的性能指標(biāo)有：巡航狀態(tài)下的阻力系數(shù)CD,M=0.695、阻力發(fā)散馬赫數(shù)下的阻力系數(shù)CD,M=0.715、俯仰力矩系數(shù)Cm以及低速失速特性等。其中低速失速特性主要體現(xiàn)在起飛狀態(tài)失速迎角相比初始翼型要有所提高，在本文的研究中，若優(yōu)化得到的優(yōu)化翼型在失速點(diǎn)附近迎角對(duì)應(yīng)的升力系數(shù)高于初始翼型在相同迎角下對(duì)應(yīng)的升力系數(shù)，則可認(rèn)為該優(yōu)化翼型的低速失速特性有所提高。對(duì)于本文研究的翼型，其巡航馬赫數(shù)為0.695，阻力發(fā)散馬赫數(shù)為0.715，低速狀態(tài)下失速迎角為11.5°。故本文的研究中選定的優(yōu)化目標(biāo)為巡航阻力系數(shù)、阻力發(fā)散阻力系數(shù)、俯仰力矩系數(shù)以及低速狀態(tài)11°迎角和12°迎角下升力系數(shù)共計(jì)五個(gè)優(yōu)化目標(biāo)。此外，需要強(qiáng)調(diào)的是，本節(jié)針對(duì)翼型進(jìn)行氣動(dòng)優(yōu)化是以翼型的隱身性能相對(duì)于初始翼型不變差為前提進(jìn)行的，故在本節(jié)的優(yōu)化過(guò)程中，需要將隱身性能作為約束加入到優(yōu)化程序中。

對(duì)確定的五個(gè)優(yōu)化目標(biāo)進(jìn)行降維處理。首先，利用拉丁超立方方法在指定設(shè)計(jì)空間中隨機(jī)取樣，利用RANS方法對(duì)采集到的樣本點(diǎn)對(duì)應(yīng)翼型的氣動(dòng)性能進(jìn)行綜合評(píng)估，得到一個(gè)包含各樣本點(diǎn)對(duì)應(yīng)上述五個(gè)優(yōu)化目標(biāo)性能的數(shù)據(jù)集。本文中選定樣本點(diǎn)共計(jì)570個(gè)，進(jìn)行PCA主成分分析，得到各目標(biāo)對(duì)應(yīng)各主成分的元素如表1所示。其中百分?jǐn)?shù)即為各主成分對(duì)應(yīng)的權(quán)重，每一列中的數(shù)字即為該主成分中各目標(biāo)對(duì)應(yīng)的元素。

表1 主成分元素對(duì)應(yīng)情況Table 1 Principal component’s elements

表中百分?jǐn)?shù)表示每一個(gè)主成分對(duì)應(yīng)的權(quán)重，按照上文所述基于PCA方法的高維多目標(biāo)降維的規(guī)則可篩選得到三個(gè)主要目標(biāo)，分別為巡航阻力系數(shù)CD,M=0.695、巡航俯仰力矩系數(shù)Cm以及低速狀態(tài)12°迎角下升力系數(shù)CL,α=12°。由表中數(shù)據(jù)可以看出，兩個(gè)狀態(tài)下的阻力系數(shù)以及低速狀態(tài)不同迎角下的升力系數(shù)對(duì)應(yīng)的主成分元素之間具有相同的符號(hào)并且絕對(duì)值大小接近，表明這兩對(duì)目標(biāo)內(nèi)部具有很強(qiáng)的相關(guān)性，所以其中存在冗余目標(biāo)。此外，巡航狀態(tài)的阻力系數(shù)與低速狀態(tài)的失速性能指標(biāo)在第一主成分中具有相反的符號(hào)，即表明巡航性能與低速性能往往相互矛盾，這一點(diǎn)也是符合工程經(jīng)驗(yàn)的。在實(shí)際的優(yōu)化設(shè)計(jì)中，通過(guò)PCA主成分分析辨別出冗余目標(biāo)及各目標(biāo)之間的關(guān)聯(lián)關(guān)系后，可以采用去除冗余目標(biāo)或者將冗余目標(biāo)轉(zhuǎn)化為約束條件的方法達(dá)到降維的目的。此外，對(duì)于設(shè)計(jì)中不主要關(guān)注的性能指標(biāo)，也可以將其設(shè)定為約束條件，這樣既可保證優(yōu)化結(jié)果至少滿足該指標(biāo)最低要求，也可達(dá)到降維的目的。故本文對(duì)該飛翼翼型的優(yōu)化設(shè)計(jì)中，選定巡航狀態(tài)下阻力系數(shù)CD,M=0.695、低速12°迎角下升力系數(shù)CL,a=12為優(yōu)化目標(biāo)，將其他指標(biāo)作為約束條件，即優(yōu)化翼型相對(duì)初始翼型對(duì)應(yīng)性能不變差。優(yōu)化得到一組Pareto最優(yōu)解，從優(yōu)化結(jié)果中選取性能各具特點(diǎn)的三副翼型Opt 1、Opt 2以及Opt 3，如圖4所示，其中翼型Opt 1具有較好的低速升力特性，翼型Opt 2在巡航狀態(tài)下減阻明顯，翼型Opt 3兼顧巡航阻力特性及低速升力特性。優(yōu)化翼型與初始翼型綜合性能對(duì)比如圖5所示。

由優(yōu)化結(jié)果可看出，優(yōu)化翼型相對(duì)初始翼型具有

圖4 優(yōu)化得到的Pareto解分布Fig.4 Pareto solutions of optimized result

(a) Polar

(c) Drag divergence

(d) Pitching moment

(e) Low speed lift coefficient

更好的巡航升阻力特性、低速升力特性并保持了俯仰力矩特性，這表明通過(guò)PCA主成分分析去除冗余目標(biāo)后，對(duì)主要目標(biāo)進(jìn)行優(yōu)化所得出的結(jié)果在主要目標(biāo)及冗余目標(biāo)上均有提升。表2中列出了從優(yōu)化結(jié)果中挑選出的三副翼型相對(duì)于初始翼型各主要目標(biāo)參數(shù)的變化情況。從優(yōu)化翼型可見(jiàn)，通過(guò)PCA主成分分析確定的與主要目標(biāo)呈現(xiàn)較強(qiáng)正相關(guān)關(guān)系的冗余目標(biāo)獲得了與該主要目標(biāo)相同方向且幅度相當(dāng)?shù)膬?yōu)化效果，例如作為冗余目標(biāo)的阻力發(fā)散馬赫數(shù)下的阻力系數(shù)與作為主要目標(biāo)的巡航狀態(tài)阻力系數(shù)獲得了相近的優(yōu)化效果。而PCA主成分分析中確定的存在相互矛盾的兩個(gè)目標(biāo)，如巡航狀態(tài)阻力系數(shù)與低速狀態(tài)升力性能，在優(yōu)化結(jié)果中其矛盾性也得到體現(xiàn)。例如，翼型Opt 3在所有優(yōu)化翼型中具有最佳的巡航阻力性能，但是其低速升力性能卻是優(yōu)化翼型中較差的；而翼型Opt 1具有所有優(yōu)化翼型中最佳的低速升力性能，但是其巡航阻力性能是優(yōu)化翼型中表現(xiàn)較差的。

表2 優(yōu)化翼型相對(duì)初始翼型性能變化Table 2 Performance improvement of several optimized airfoil

圖6為優(yōu)化前后翼型幾何外形的變化示意圖，可以看出從優(yōu)化結(jié)果中挑選的三副翼型都朝著相似的方向變化，即前緣彎度變大，這有利于提高翼型失速性能。圖7給出了翼型優(yōu)化前后壓力分布云圖，圖8

(a) Optimized airfoil Opt 1

(b) Optimized airfoil Opt 2

(c) Optimized airfoil Opt 3

(a) Initial airfoil

(b) Optimized airfoil Opt 1

(c) Optimized airfoil Opt 2

(d) Optimized airfoil Opt 3

所示為優(yōu)化前后翼型巡航狀態(tài)下壓力分布的對(duì)比。從圖7及圖8的結(jié)果可見(jiàn)，優(yōu)化后翼型產(chǎn)生的激波發(fā)生后移，且激波強(qiáng)度明顯減弱，這種壓力分布的改變對(duì)巡航減阻有所幫助。

圖8 壓力分布對(duì)比Fig.8 Comparison of pressure coefficient

5 氣動(dòng)、隱身高維度多目標(biāo)綜合優(yōu)化

以該飛翼布局飛行器翼型氣動(dòng)、隱身多學(xué)科高維多目標(biāo)綜合設(shè)計(jì)為例，進(jìn)一步研究PCA方法在多學(xué)科綜合設(shè)計(jì)中的有效性與合理性。

5.1 翼型隱身性能計(jì)算方法

5.1.1 隱身性能工程估算法

結(jié)合上一部分的研究結(jié)果，進(jìn)一步考慮隱身性能對(duì)該型翼型進(jìn)行優(yōu)化設(shè)計(jì)。為減少計(jì)算量，本文對(duì)樣本點(diǎn)對(duì)應(yīng)的隱身性能計(jì)算采用工程估算方法，即高頻近似物理光學(xué)法以及蘇聯(lián)學(xué)者尤費(fèi)賽夫提出的物理繞射理論[20]。為了避免軸向焦散，使用等效電磁流(MEC)方法，使得不在Keller錐上的散射方向也能計(jì)算繞射場(chǎng)，計(jì)算中物體的楔脊被等效絲狀的線電流和線磁流代替[21]。

物理光學(xué)近似下面元的鏡面散射電磁場(chǎng)表示為：

除鏡面散射電磁場(chǎng)對(duì)RCS的貢獻(xiàn)外，邊緣的繞射電磁場(chǎng)也是一個(gè)重要的散射源，由等效絲狀線電流和線磁流代替的等效電磁流分別為：

邊緣的繞射電場(chǎng)表示為：

最后對(duì)所有的照明面元和照明楔，合成散射場(chǎng)和繞射場(chǎng)得到的物體的雷達(dá)散射面積：

5.1.2 基于Maxwell方程組的電磁場(chǎng)求解方法

為保證隱身性能工程估算與優(yōu)化設(shè)計(jì)相結(jié)合的方法對(duì)隱身性能的提升行之有效，本文最后采用了高精度求解麥克斯韋方程組方法驗(yàn)證優(yōu)化翼型隱身性能優(yōu)化效果。任何電磁問(wèn)題的電磁場(chǎng)解都滿足如下時(shí)變麥克斯韋方程組：

式中，B是磁感應(yīng)強(qiáng)度矢量，H是磁場(chǎng)強(qiáng)度，D是電位移矢量，E是電場(chǎng)強(qiáng)度矢量，J是自由電流密度，ρ是自由電荷密度。無(wú)源情況下如在自由空間中，J=0，ρ=0。

在直角坐標(biāo)系下，無(wú)源條件下的麥克斯韋方程組可寫(xiě)為守恒形式：

其中：

本文使用的課題組研發(fā)的電磁場(chǎng)求解程序采用了多塊并行的時(shí)域有限體積法(FVTD)，并使用雙時(shí)間步推進(jìn)。

5.2 翼型氣動(dòng)、隱身高維度多目標(biāo)綜合優(yōu)化

依然采用上文采樣得到的570個(gè)樣本點(diǎn)，使用工程估算計(jì)算樣本點(diǎn)對(duì)應(yīng)的隱身性能，得到包含氣動(dòng)性能及隱身性能在內(nèi)的共六個(gè)目標(biāo)的數(shù)據(jù)集，對(duì)其進(jìn)行PCA主成分分析，得到各目標(biāo)對(duì)應(yīng)各主成分的元素如表3所示。

表3 主成分元素對(duì)應(yīng)情況Table 3 Principal component’s elements

根據(jù)PCA分析結(jié)果，可篩選得到四個(gè)主要目標(biāo)，分別為巡航阻力系數(shù)CD,M=0.695、巡航俯仰力矩系數(shù)Cm、低速狀態(tài)12°迎角下升力系數(shù)CL,a=12以及雷達(dá)散射面積RCS。由于該飛翼布局飛行器設(shè)計(jì)目標(biāo)主要為增加巡航半徑，提高隱身能力，故本文中選定巡航狀態(tài)阻力系數(shù)及雷達(dá)散射面積作為優(yōu)化目標(biāo)，其他目標(biāo)轉(zhuǎn)化為約束，優(yōu)化得到Pareto最優(yōu)解集分布如圖9所示。從Pareto最優(yōu)解集中挑選出三副綜合性能優(yōu)良的翼型進(jìn)行氣動(dòng)隱身綜合性能驗(yàn)證，優(yōu)化翼型各目標(biāo)性能相對(duì)初始翼型提升幅度由表4給出。

圖9 優(yōu)化得到的Pareto解分布Fig.9 Pareto solutions of optimized result

ParameterOpt1Opt2Opt3CD，M=0．695-6．7%-7．1%-8．7%CD，M=0．715-9．8%-10．1%-12．7%abs(Cm)-22．5%-19．7%-9．9%lowspeedCL，a=113．1%2．6%1．8%lowspeedCL，a=126．3%5．9%3．2%RCS-12．3%-9．8%-8．9%

當(dāng)把隱身性能作為優(yōu)化目標(biāo)之一時(shí)，氣動(dòng)性能的提升幅度相對(duì)第4節(jié)的研究有一定程度減小，尤其是翼型的低速升力特性優(yōu)化幅度減少最為明顯，這表明翼型的氣動(dòng)性能與隱身性能在一定程度上存在相互矛盾，為獲得隱身性能上的提升，在氣動(dòng)性能上需要作出一定妥協(xié)。

對(duì)于隱身性能的驗(yàn)證，采用了基于高精度求解麥克斯韋方程組的方法，對(duì)初始翼型及優(yōu)化翼型前向單站雷達(dá)-30°～30°角度范圍照射下雷達(dá)散射面積RCS值進(jìn)行了計(jì)算，圖10所示為計(jì)算所得初始翼型及優(yōu)化翼型散射磁場(chǎng)Dsz分布情況，圖11所示為優(yōu)化翼型與初始翼型在-30°～30°之間單站RCS空間分布的對(duì)比?？梢钥闯觯瑑?yōu)化翼型在單站雷達(dá)負(fù)角度入射時(shí)雷達(dá)散射面積相對(duì)初始翼型有明顯降低，即對(duì)地隱身性能有明顯改善。由于本文優(yōu)化目標(biāo)即為對(duì)地隱身性能，故該性能提升的同時(shí)，相應(yīng)的付出了雷達(dá)正角度入射時(shí)散射面積有所增加的代價(jià)。綜合來(lái)看，使用基于PCA主成分分析的多目標(biāo)優(yōu)化方法對(duì)解決綜合考慮氣動(dòng)及隱身極多目標(biāo)優(yōu)化問(wèn)題具有良好效果，能夠獲得兼顧氣動(dòng)及隱身性能的綜合性能良好的優(yōu)化結(jié)果。

(a) Initial airfoil

(b) Optimized airfoil Opt_1

(c) Optimized airfoil Opt_2

(d) Optimized airfoil Opt_3

圖11 單站RCS空間分布對(duì)比Fig.11 Comparison of RCS space distribution in single azimuth

6 結(jié) 論

文中研究了PCA方法在高維多目標(biāo)優(yōu)化降維研究中的應(yīng)用，并且對(duì)綜合考慮氣動(dòng)、隱身性能的翼型多目標(biāo)設(shè)計(jì)開(kāi)展研究：

1) DTLZ5典型測(cè)試函數(shù)表明，PCA主成分分析方法能夠有效降維，有效識(shí)別出冗余目標(biāo)。

2) 翼型氣動(dòng)、隱身多目標(biāo)設(shè)計(jì)優(yōu)化研究表明，PCA分析方法對(duì)冗余目標(biāo)的識(shí)別，與相關(guān)性分析及流動(dòng)物理機(jī)制保持兼容一致。

3) 所設(shè)計(jì)的氣動(dòng)外形滿足多個(gè)設(shè)計(jì)狀態(tài)下阻力發(fā)散特性、升阻比特性、力矩特性以及隱身特性等多個(gè)需求。

4) PCA主成分分析方法對(duì)于高維氣動(dòng)多目標(biāo)優(yōu)化設(shè)計(jì)問(wèn)題是有效的，在不失問(wèn)題主特征前提下進(jìn)行有效降維，降低問(wèn)題的復(fù)雜程度，提高非劣解集的可視化水平，更加有利于設(shè)計(jì)人員對(duì)多目標(biāo)問(wèn)題做出合理有效的決策。

未來(lái)飛行器外形設(shè)計(jì)趨勢(shì)必然是多學(xué)科、多目標(biāo)、多約束綜合優(yōu)化?？梢灶A(yù)見(jiàn)，PCA主成分分析方法將發(fā)揮重要作用?；诂F(xiàn)有的飛行器多目標(biāo)集成化設(shè)計(jì)平臺(tái)，開(kāi)展氣動(dòng)、隱身、結(jié)構(gòu)等學(xué)科的多學(xué)科、高維度多目標(biāo)優(yōu)化，是該方法具體應(yīng)用的進(jìn)一步研究方向。

[1]Knowles J, Corne D. Quantifying the effects of objective space dimension in evolutionary multi-objective optimization[M]. Fourth International Conference on Evolutionary Multi-Criterion Optimization (EMO-2007), 2007, 757-771

[2]Purshouse R C. Evolutionary many-objective optimisation: an exploratory analysis[C]//Proceedings of the 2003 IEEE Congress on Evolutionary Computation, IEEE Press, Piscataway, 2003: 2066-2073

[3]Schutze O, Lara A, Coello C. On the influence of the number of objectives on the hardness of a multiobjective optimization problem[J]. IEEE Transactions on Evolutionary Computation, 2011, 15(4): 444-455

[4]Huang L F. Research on evolutionary multiobjective algorithms[D]. Hefei: University of Science and Technology of China, 2009. (in Chinese)黃林峰. 多目標(biāo)優(yōu)化算法研究[D]. 合肥: 中國(guó)科學(xué)技術(shù)大學(xué), 2009

[5]Lei Y Y, Jiang W Z, Liu L J, et al. Many-objective optimization based on sub-objective evolutionary algorithm[J]. Journal of Beijing University of Aeronautics and Astronautics, 2015, 41(10): 1910-1917. (in Chinese)雷宇曜, 姜文志, 劉立佳, 等. 基于子目標(biāo)進(jìn)化的高維多目標(biāo)優(yōu)化算法[J]. 北京航空航天大學(xué)學(xué)報(bào), 2015, 41(10): 1910-1917

[6]Chen X H, Li X, Wang N. Objective reduction with sparse feature selection for many objective optimization problem[J]. Acta Electronica Sinica, 2015, 43(7): 1300-1307. (in Chinese)陳小紅, 李霞, 王娜. 高維多目標(biāo)優(yōu)化中基于稀疏特征選擇的目標(biāo)降維方法[J]. 電子學(xué)報(bào), 2015, 43(7): 1300-1307

[7]Turk M, Pentland A. Eigenfaces for recognition[J]. Journal of Cognitive Neuroscience, 1991, 3(1): 71-86

[8]Zhang D, Zhou Z H. (2D)2PCA: 2-directional 2-dimensional PCA for efficient face representation and recognition[J]. Neurocomputing, 2005, 69: 224-231

[9]Kim K I, Park S H, Kim H J. Kernel principal component analysis for texture classification[J]. IEEE Signal Processing Letters, 2001, 8(2): 39-41

[10]Abdi H, Williams L J. Principal component analysis[J]. Wiley Interdisciplinary Reviews: Computational Statistics, 2010, 2(4): 433-459

[11]Zhu C M. Research on integrated learning algorithm for developmental robots[D]. Harbin: Harbin Engineering University, 2008. (in Chinese)朱長(zhǎng)明. 發(fā)育機(jī)器人集成學(xué)習(xí)算法研究[D]. 哈爾濱: 哈爾濱工程大學(xué), 2008

[12]周草臣. 基于PCA的高維多目標(biāo)優(yōu)化算法研究[D]. 重慶: 重慶大學(xué), 2014. Zhou C C. The Study on many-objective optimization algorithms based on PCA[D]. Chongqing: Chongqing University, 2014. (in Chinese)[13]Deb K, Saxena D K. Searching for Pareto-optimal solutions through dimensionality reduction for certain large-dimensional multi-objective optimization problems[C]//IEEE Congress on Evolutionary Computation, 2006, 3353-3360

[14]Zhao K. Complex aerodynamic optimization and robust design method based on computation fluid dynamics[D]. Xian: Northwestern Polytechnical University, 2014. (in Chinese)趙軻. 基于CFD的復(fù)雜氣動(dòng)優(yōu)化與穩(wěn)健設(shè)計(jì)方法研究[D]. 西安: 西北工業(yè)大學(xué), 2014

[15]Kennedy J, Eberhart R. Particle swarm optimization[C]//Proceedings of IEEE International Conference on Neural Networks Perth, 1995: 1942-1948

[16]Jackson I R H. Convergence properties of radial basis function[J]. Control Approximation, 1988, 4(2): 243-264

[17]Wang R W, Gao Z H. The structural topology optimization based on ant colony optimization[J]. Chinese Journal of Applied Mechanics, 2011, 28(3): 232-237. (in Chinese)王榮偉, 高正紅. 基于改進(jìn)粒子群算法的翼型多目標(biāo)優(yōu)化研究[J]. 應(yīng)用力學(xué)學(xué)報(bào), 2011, 28(3): 232-237

[18]Khare V, Yao X, Deb K. Performance scaling of multi-objective evolutionary algorithms[M]. Evolutionary Multi-Criterion Optimization. Springer Berlin Heidelberg, 2003: 376-390

[19]Kulfan B M. A universal parametric geometry rep-resentation method-‘CST’[R]. AIAA Paper 2007-0062, 1-36

[20]Deans S R. Hough transform from the Radon trans-form. IEEE Trans[J]. Pattern Analysis and Machine Intelligence, 1981, 3(2): 185-188

[21]Altintas A, Russer P. Time-domain equivalent edge currents for transient scattering[J]. IEEE Trans, Antennas Propag, 2001, 49(4): 602-606.

Multidisciplinary optimization design of high dimensional target space for flying wing airfoil

ZHENG Chuanyu1, HUANG Jiangtao1,*, ZHOU Zhu1, LIU Gang1, GAO Zhenghong2, XU Yong1
(1. China Aerodynamics Research and Development Center, Mianyang 621000, China; 2. National Key Laboratory of Aerodynamic Design and Research, Northwestern Polytechnical University, Xi’an 710072, China)

There are many difficulties when traditional optimization methods are used to deal with high-dimensional multi-objective problems. For the high-dimensional multi-objective optimization problem, it is a practical and effective method to reduce dimensionality. In this paper, we first discuss the idea of principal component analysis (PCA) target dimension reduction method in high-dimensional multi-objective applications. The effectiveness of the dimension reduction method is verified by the application of the test function of typical high dimensional objects. After that, the method is applied to the aerodynamics and stealth multidisciplinary integrated design of airfoils, and principal component analysis is carried out on the integrated design of the high dimensional target space. The PCA is used to extract principal components and identify redundant or unimportant targets. Redundant targets are removed or added as constraints to the optimization process of important targets. The results show that the optimal design results of the target space after dimension reduction satisfy the requirements of torque, drag divergence, cruise lift-drag ratio, and low-speed lift characteristics as well as stealth characteristics. The application prospect of this method in the multi-objective and multi-disciplinary optimization design of aircrafts is further discussed.

principal component analysis; redundant target; relevance; high-dimensional multi-objective; integrated optimization design

0258-1825(2017)04-0587-11

2017-04-02;

2017-06-29

國(guó)家自然科學(xué)基金( 11402288)，國(guó)家重點(diǎn)研發(fā)計(jì)劃"數(shù)值飛行器原型系統(tǒng)"重點(diǎn)專項(xiàng)(基金號(hào)2016YFB0200704)，裝備預(yù)研基金重點(diǎn)項(xiàng)目資助( 9140A13021015KG29038)

鄭傳宇(1993-)，男，碩士研究生，主要研究方向：飛行器氣動(dòng)設(shè)計(jì)與計(jì)算空氣動(dòng)力學(xué). E-mail: zhengcy11@163.com

黃江濤(1982-)，男，副研究員，研究方向：飛行器氣動(dòng)外形多學(xué)科優(yōu)化與計(jì)算空氣動(dòng)力學(xué). E-mail：hjtcyf@163.com

鄭傳宇, 黃江濤, 周鑄, 等. 飛翼翼型高維目標(biāo)空間多學(xué)科綜合優(yōu)化設(shè)計(jì)[J]. 空氣動(dòng)力學(xué)學(xué)報(bào), 2017, 35(4): 587-597.

10.7638/kqdlxxb-2017.0079 ZHENG C Y, HUANG J T, ZHOU Z, et al. Multidisciplinary optimization design of high dimensional target space for flying wing airfoil[J]. Acta Aerodynamica Sinica, 2017, 35(4): 587-597.

V211.3

A doi: 10.7638/kqdlxxb-2017.0079