呂志斌 萬(wàn)志強(qiáng) / LYU Zhibin WAN Zhiqiang

(北京航空航天大學(xué)，北京 100191)

0 引言

飛機(jī)是高度綜合的現(xiàn)代科學(xué)技術(shù)的體現(xiàn)，飛機(jī)結(jié)構(gòu)設(shè)計(jì)又是飛機(jī)設(shè)計(jì)的主要階段[1]。但隨著現(xiàn)代飛機(jī)尤其是民用客機(jī)與運(yùn)輸機(jī)追求高性能和低結(jié)構(gòu)重量等方面的要求，飛機(jī)結(jié)構(gòu)柔性較大，導(dǎo)致飛機(jī)結(jié)構(gòu)變形增大。加之近代工業(yè)在復(fù)合材料、智能材料等新型材料研究上的突破和應(yīng)用，使得氣動(dòng)彈性問(wèn)題在飛機(jī)設(shè)計(jì)中越來(lái)越凸顯出來(lái)[2]。傳統(tǒng)的“設(shè)計(jì)-校核-修改”的結(jié)構(gòu)設(shè)計(jì)模式越來(lái)越難以滿(mǎn)足現(xiàn)代飛機(jī)對(duì)于高性能和低重量的雙重要求。

基于此，隨著氣動(dòng)彈性學(xué)科的發(fā)展，以克服氣動(dòng)彈性效應(yīng)所帶來(lái)的不利影響的同時(shí)盡可能地降低飛機(jī)結(jié)構(gòu)重量為目標(biāo)的氣動(dòng)彈性結(jié)構(gòu)優(yōu)化技術(shù)已經(jīng)逐漸成為一種主動(dòng)設(shè)計(jì)的手段，應(yīng)用于飛機(jī)設(shè)計(jì)的各個(gè)階段[3]。

然而目前工程上在進(jìn)行氣動(dòng)彈性分析求解時(shí)氣動(dòng)力數(shù)據(jù)來(lái)源仍然是采用低階面元法，該方法計(jì)算效率高但精度低；而CFD(Computational Fluid Dynamics，簡(jiǎn)稱(chēng)CFD)方法雖然計(jì)算精度高，但效率太低，而且難以實(shí)現(xiàn)優(yōu)化過(guò)程的迭代計(jì)算。當(dāng)前的氣動(dòng)彈性?xún)?yōu)化算法則主要集中于單一的優(yōu)化算法，如可行方向法、遺傳算法等。這些算法均有其自身的局限性，如敏度信息不易獲得、搜索精度低、收斂速度慢等問(wèn)題[4]。同時(shí)考慮到氣動(dòng)彈性結(jié)構(gòu)優(yōu)化本身往往具有計(jì)算量大、效率低、時(shí)間長(zhǎng)、設(shè)備要求高等特點(diǎn)，難以滿(mǎn)足方案設(shè)計(jì)初期對(duì)時(shí)間的要求。

針對(duì)上述問(wèn)題，本文提出了采用基于高階面元法的靜氣動(dòng)彈性分析方法來(lái)求解約束響應(yīng)，使用遺傳/混合算法作為尋優(yōu)算法，并引入Kriging(克立格)代理模型來(lái)減少計(jì)算量和優(yōu)化時(shí)間。最后通過(guò)一個(gè)成熟的大展弦比機(jī)翼為算例，驗(yàn)證本文提出的靜氣動(dòng)彈性結(jié)構(gòu)優(yōu)化方法的可行性與可靠性，從而為飛機(jī)方案設(shè)計(jì)階段提供一種新的優(yōu)化方法。

1 基于高階面元法的靜氣動(dòng)彈性分析

1.1 高階面元法

高階面元法與低階面元法一樣都是線性方法，都是通過(guò)在物面網(wǎng)格上進(jìn)行布置偶極子或者源來(lái)模擬氣流流經(jīng)物體表面的流場(chǎng)，都是對(duì)線化后的小擾動(dòng)位勢(shì)流方程進(jìn)行求解[5]：

(1)

(2)

其中，M∞為來(lái)流馬赫數(shù)，?為擾動(dòng)位函數(shù)。

而與低階面元法不同的是，低階面元法在平板氣動(dòng)力網(wǎng)格上的壓力點(diǎn)需布置在每個(gè)氣動(dòng)網(wǎng)格的四分之一弦線上，高階面元法則沒(méi)有類(lèi)似限制。除可以按照幾何外形建立三維的氣動(dòng)力計(jì)算模型外，高階面元法在氣動(dòng)力計(jì)算模型的表面網(wǎng)格上布置連續(xù)的點(diǎn)源，在結(jié)點(diǎn)上布置連續(xù)的偶極子，線性分布奇點(diǎn)強(qiáng)度，如圖1所示[5]。

圖1 網(wǎng)格元素奇點(diǎn)分布示意圖

1.2 靜氣動(dòng)彈性分析框架

與常用的基于低階面元法的Nastran軟件不同，本文使用的靜氣動(dòng)彈性分析方法的氣動(dòng)力數(shù)據(jù)來(lái)源于基于高階面元法的Panair開(kāi)源程序，結(jié)構(gòu)計(jì)算基于柔度法的Nastran靜力分析模塊，剛性配平為求解基本的飛行力學(xué)平衡方程，載荷導(dǎo)數(shù)計(jì)算采用差分法。圖2所示即為基于高階面元法的靜氣動(dòng)彈性分析方法框架。

圖2 基于高階面元法的靜氣動(dòng)彈性分析方法框架

2 代理模型方法

代理模型方法主要包含兩方面的內(nèi)容：一是構(gòu)建模型的樣本點(diǎn)的取樣方法，又稱(chēng)實(shí)驗(yàn)設(shè)計(jì)方法；二是進(jìn)行數(shù)據(jù)擬合和預(yù)測(cè)的模型構(gòu)建方法，又稱(chēng)近似方法[6]。

2.1 拉丁超立方抽樣

拉丁超立方取樣方法于1979年由McKay提出，該方法最大的優(yōu)點(diǎn)是不需考慮樣本點(diǎn)的維度和取樣數(shù)目，通過(guò)控制取樣點(diǎn)的位置可以盡可能避免取樣點(diǎn)在小領(lǐng)域內(nèi)重合，從而實(shí)現(xiàn)取樣點(diǎn)在設(shè)計(jì)空間的均勻分布[7]。

假設(shè)需要在m維向量空間中抽取n個(gè)樣本，其具體步驟為：

1)將m維向量空間中的每一維均按策略(一般采用均勻分布的策略)分為互不重疊的n個(gè)區(qū)間，使得每個(gè)區(qū)間有相同的被選擇概率；

2)從向量空間每一維的每個(gè)區(qū)間中隨機(jī)抽取一個(gè)點(diǎn)作為該維該空間的值，按照劃分的區(qū)間共抽取m×n個(gè)值；

3)再?gòu)拿恳痪S中隨機(jī)抽取2種選擇的點(diǎn)組成m維向量，共抽取n個(gè)，至此樣本點(diǎn)抽取完成。

如圖3所示為用Matlab仿真出來(lái)的二維向量空間的拉丁超立方取樣的樣本點(diǎn)分布圖，取樣點(diǎn)數(shù)目為8個(gè)。

圖3 二維空間的拉丁超立方取樣

2.2 Kriging代理模型方法

Kriging模型是一種改進(jìn)的線型回歸模型，其函數(shù)表達(dá)式中既包含了線性回歸部分，又包含了非參數(shù)部分[8]：

y=F(x)+z(x)

(3)

式中，F(xiàn)(x)為線性回歸部分，具有全局近似的特點(diǎn)；z(x)為非參數(shù)隨機(jī)部分，具有局部偏差的特點(diǎn)。

線性回歸部分可以表示為：

F(x)=f1(x)β1+…+fp(x)βp

(4)

按照f(shuō)(x)形式的不同，可分為常數(shù)型、線性型和二次型。本文選擇如下所示的二次型，式中β可通過(guò)廣義最小二乘估計(jì)得到：

f1(x)=1

f2(x)=x1,…,fn+1(x)=xn

(5)

非參數(shù)部分具有如下的統(tǒng)計(jì)學(xué)特性：

E[z(x)]=0

Var(z(x)]=σ2

E[z(xi),z(xj)]=σ2G(xj,xi)

(6)

式中，G(θ,xj,xi)為i、j兩點(diǎn)的空間相關(guān)方程，該空間相關(guān)方程也有很多不同種類(lèi)，對(duì)模型構(gòu)造影響很大。本文取Guass函數(shù)作為空間相關(guān)函數(shù)：

(7)

式中，ndv為樣本點(diǎn)的維度，θk是空間相關(guān)函數(shù)在樣本點(diǎn)第j個(gè)方向的參數(shù)，為待求解的量。

對(duì)于一組給定樣本點(diǎn)X=[x1…xn]T及其響應(yīng)值Y=[y1…yn]T，由樣本點(diǎn)得到待測(cè)點(diǎn)的響應(yīng)值為：

(8)

式中c是待求權(quán)系數(shù)向量。根據(jù)無(wú)偏條件及估計(jì)方差最小條件，結(jié)合拉格朗日乘子，可求得c：

(9)

具體推導(dǎo)過(guò)程見(jiàn)參考文獻(xiàn)[8]。

3 優(yōu)化算法

與其他優(yōu)化問(wèn)題一樣，氣動(dòng)彈性結(jié)構(gòu)優(yōu)化問(wèn)題也可以用如下的公式進(jìn)行描述，即在設(shè)計(jì)空間中尋找一組變量使目標(biāo)函數(shù)(一般為結(jié)構(gòu)重量)最小化[9]：

Min.F(v)

(10)

S.T.gj(v)≤0j=1,2,…,ncon

(11)

(12)

其中，式(10)為目標(biāo)函數(shù)，本文為結(jié)構(gòu)重量最?。皇?11)為約束條件，包括靜氣彈約束(如發(fā)散速度約束、變形約束等)、動(dòng)氣彈約束(如顫振速度約束等)、強(qiáng)度約束、操穩(wěn)約束以及氣動(dòng)約束(如阻力約束等)等；式(12)為設(shè)計(jì)變量的上下限，指定設(shè)計(jì)空間的值域，又叫邊界約束。

3.1 遺傳算法

遺傳算法是一種借鑒自然選擇和遺傳進(jìn)化機(jī)理而發(fā)明出來(lái)的自適應(yīng)概率搜索算法，具有全局尋優(yōu)的能力[10]。作為進(jìn)化算法的一種，遺傳算法從初代個(gè)體開(kāi)始，以適應(yīng)度函數(shù)為評(píng)價(jià)標(biāo)準(zhǔn)，通過(guò)遺傳算子實(shí)現(xiàn)對(duì)個(gè)體的自然選擇和遺傳進(jìn)化，從而不斷改良個(gè)體直至滿(mǎn)足條件。

圖4所示為遺傳算法流程圖[9]，其基本步驟為：

圖4 遺傳算法流程圖

1)確定優(yōu)化策略，定義優(yōu)化參數(shù)，包括群體規(guī)模、交叉概率、變異概率、程序終止條件等遺傳參數(shù)，選擇編碼方式、選擇方式、交叉方式、變異方式等；

2)隨機(jī)產(chǎn)生指定規(guī)模的初代群體(即第一代父群體)；

3)計(jì)算群體中每個(gè)個(gè)體的目標(biāo)函數(shù)值以及約束響應(yīng)值，從而計(jì)算個(gè)體的適應(yīng)度值，從而對(duì)個(gè)體進(jìn)行評(píng)估；

4)判斷是否滿(mǎn)足程序終止條件，如果滿(mǎn)足則計(jì)算結(jié)束，否則繼續(xù)進(jìn)行計(jì)算；

5)通過(guò)基本遺傳操作(選擇、交叉和變異)產(chǎn)生新一代子群體，用子群體代替父群體形成新一代父群體，重復(fù)步驟3和步驟4。

3.2 分形算法

分形算法是一種借鑒樹(shù)的生長(zhǎng)來(lái)實(shí)現(xiàn)搜索尋優(yōu)的啟發(fā)式算法，具有較強(qiáng)的局部搜索能力。圖5所示為分形算法流程圖，其中選點(diǎn)、分枝、剪枝為分形算法的三個(gè)主要操作運(yùn)算[11]。分形算法的原理很簡(jiǎn)單，通過(guò)選點(diǎn)來(lái)確定當(dāng)前區(qū)域的最優(yōu)解，為分枝和剪枝提供參考依據(jù)；通過(guò)分枝來(lái)產(chǎn)生新的樣本點(diǎn)，計(jì)算樣本點(diǎn)的適應(yīng)度(為了統(tǒng)一個(gè)體評(píng)價(jià)標(biāo)準(zhǔn)，此處依舊引入遺傳算法中的適應(yīng)度概念代替目標(biāo)函數(shù)值)以判斷是否更新當(dāng)前全局和局部最優(yōu)解；通過(guò)剪枝來(lái)縮小可行域范圍，使可行域逐漸逼近最優(yōu)解。在分形優(yōu)化算法中，分枝提高了算法的精度，剪枝加快了算法的收斂速度。

3.3 遺傳/分形算法

遺傳/分形混合算法結(jié)合了兩種算法的優(yōu)點(diǎn)，遺傳算法用于全局搜索以避免優(yōu)化算法陷入局部最優(yōu)解，分形算法用于最優(yōu)個(gè)體周?chē)男》秶植繉?yōu)。圖6所示為遺傳/分形混合算法的流程圖，若該流程中不適用分形算法，則圖6為遺傳算法的流程圖。

確定遺傳策略和優(yōu)化參數(shù)后，首先使用遺傳算法對(duì)優(yōu)化問(wèn)題進(jìn)行一次基本遺傳操作。選擇遺傳操作后的最優(yōu)個(gè)體，以該個(gè)體為中心確定分形算法的初始可行域以及初始最優(yōu)解(即為該最優(yōu)個(gè)體適應(yīng)度)。然后對(duì)上述可行域使用分形算法進(jìn)行尋優(yōu)，并用分形算法得到的最優(yōu)個(gè)體代替上一次遺傳算法得到的最優(yōu)個(gè)體。完成上述操作后重新評(píng)估該代群體，判斷是否滿(mǎn)足算法的終止條件。

圖6 遺傳/分形混合算法流程圖

3.4 靜氣動(dòng)彈性結(jié)構(gòu)優(yōu)化框架

基于前文提到的基于高階面元法的靜氣動(dòng)彈性分析方法、代理模型以及遺傳/分形算法，本文建立了一種高效高精度的靜氣動(dòng)彈性結(jié)構(gòu)優(yōu)化框架，如圖7所示。在該優(yōu)化框架中，遺傳/分形算法為優(yōu)化算法(優(yōu)化時(shí)不考慮分形算法即為遺傳算法)，代理模型代替基于高階面元法的靜氣動(dòng)彈性分析方法作為個(gè)體目標(biāo)函數(shù)和約束響應(yīng)的求解器。在實(shí)施本文提出的優(yōu)化方法時(shí)，首先應(yīng)該完成代理模型的構(gòu)建和精度校驗(yàn)，如圖8所示。

圖7 靜氣動(dòng)彈性結(jié)構(gòu)優(yōu)化框架

圖8 代理模型構(gòu)建與精度校驗(yàn)流程圖

4 算例分析

4.1 算例描述

本文所用算例的氣動(dòng)和結(jié)構(gòu)模型如圖9和圖10所示，機(jī)翼半展長(zhǎng)16.27 m，翼尖弦長(zhǎng)1.62 m。優(yōu)化時(shí)，本文僅對(duì)機(jī)翼進(jìn)行優(yōu)化，而不改變機(jī)身的結(jié)構(gòu)重量。

圖9 算例氣動(dòng)模型

圖10 算例結(jié)構(gòu)模型

表1給出了算例的優(yōu)化方案，優(yōu)化狀態(tài)為飛行高度在11 200 m、飛行馬赫數(shù)0.8時(shí)的定直平飛狀態(tài)，副翼和方向舵偏角均為0，依靠升降舵和攻角進(jìn)行配平。

表1 算例的優(yōu)化方案

為了機(jī)翼的蒙皮滿(mǎn)足從翼根到翼尖、從后緣到前緣逐漸遞減的設(shè)計(jì)準(zhǔn)則，優(yōu)化時(shí)分別使用不同的線性函數(shù)來(lái)表示蒙皮的厚度變化。圖11和圖12為初始模型機(jī)翼上下表面蒙皮的厚度分布。

圖11 機(jī)翼上表面蒙皮厚度分度

圖12 機(jī)翼下表面蒙皮厚度分度

4.2 代理模型構(gòu)建與精度校驗(yàn)

在優(yōu)化前，首先需要完成代理模型的構(gòu)建和精度校核。參照?qǐng)D8的流程首先在機(jī)翼的設(shè)計(jì)變量空間中使用拉丁超立方取樣法選取4 000個(gè)樣本點(diǎn)作為樣本輸入用于模型構(gòu)建，使用隨機(jī)方法選取500個(gè)樣本點(diǎn)作為測(cè)試點(diǎn)用于精度校核，樣本點(diǎn)和測(cè)試點(diǎn)均采用基于高階面元法的靜氣動(dòng)彈性分析方法進(jìn)行目標(biāo)函數(shù)和約束響應(yīng)計(jì)算。

(13)

(14)

復(fù)相關(guān)系數(shù)R的值介于0和1之間，其值越大，說(shuō)明預(yù)測(cè)值和真實(shí)值相關(guān)性越大即越接近。一般而言，0.5～1.0即為強(qiáng)相關(guān)。相對(duì)均方根誤差RMSE與響應(yīng)值無(wú)關(guān)，值介于0和1之間，其值越小，說(shuō)明預(yù)測(cè)精度越高[12]。

圖13～圖16為校驗(yàn)點(diǎn)的Kriging模型預(yù)測(cè)結(jié)果與高階面元法的靜氣動(dòng)彈性程序計(jì)算結(jié)果的對(duì)比，表 2說(shuō)明了Kriging模型的預(yù)測(cè)表現(xiàn)。從圖表的結(jié)果來(lái)看，Kriging模型的預(yù)測(cè)精度很高，可以完全滿(mǎn)足工程應(yīng)用和科學(xué)研究的精度要求。

圖13 結(jié)構(gòu)重量預(yù)測(cè)誤差百分比

圖14 翼尖位移預(yù)測(cè)誤差百分比

圖15 翼尖扭角預(yù)測(cè)誤差百分比

圖16 副翼效率預(yù)測(cè)誤差百分比

表2 Kriging模型預(yù)測(cè)表現(xiàn)

4.3 優(yōu)化結(jié)果對(duì)比與分析

本文使用了遺傳算法和遺傳/分形算法，結(jié)合Kriging代理模型對(duì)機(jī)翼進(jìn)行了優(yōu)化，優(yōu)化共計(jì)進(jìn)行了200代，圖17顯示了兩種算法的目標(biāo)函數(shù)隨迭代代數(shù)變化的曲線，表3為算例的優(yōu)化結(jié)果。由于算例使用的機(jī)翼是一個(gè)成熟的機(jī)翼，因此可優(yōu)化的程度并不大。

圖17 優(yōu)化目標(biāo)函數(shù)曲線

從目標(biāo)函數(shù)曲線來(lái)看，遺傳/分形混合優(yōu)化算法相比于單一的遺傳算法具有更快收斂速度和更好的全局尋優(yōu)能力，更容易收斂到最優(yōu)解，且最優(yōu)解也更接近全局最優(yōu)解。從算例的優(yōu)化結(jié)果來(lái)看，相比于初始模型，結(jié)構(gòu)重量減小，翼尖變形和扭角有一定程度的提高，副翼效率略微減小，但都在約束條件范圍內(nèi)(由于存在約束條件閾值，使得在優(yōu)化時(shí)約束響應(yīng)能略微超過(guò)約束條件)。

表3 優(yōu)化結(jié)果

5 結(jié)論

本文建立了一種基于高階面元法和遺傳/分形算法的靜氣動(dòng)彈性結(jié)構(gòu)優(yōu)化算法，并引入Kriging代理模型方法來(lái)加快優(yōu)化時(shí)間，降低優(yōu)化成本。其中，基于高階面元法的靜氣動(dòng)彈性?xún)?yōu)化方法用于約束響應(yīng)的求解，遺傳/分形算法用于目標(biāo)函數(shù)尋優(yōu)。最后以一個(gè)大展弦比客機(jī)機(jī)翼為例，通過(guò)與遺傳算法的對(duì)比，說(shuō)明了本文提出方法的實(shí)用性和更強(qiáng)的尋優(yōu)效率。

值得一提的是，本文在對(duì)算例進(jìn)行優(yōu)化時(shí)，大部分的優(yōu)化時(shí)間花費(fèi)在了樣本數(shù)據(jù)的計(jì)算上，而代理模型的構(gòu)建和優(yōu)化計(jì)算總耗時(shí)不超過(guò)2h，也顯示了本文的方法在時(shí)間上的高效性。