飛行器翼面協(xié)同優(yōu)化中的氣動代理模型研究

席 睿，李 冰，賈宏光

（1.華北水利水電大學(xué) 機械學(xué)院，河南 鄭州 450011；2.中國科學(xué)院長春光學(xué)精密機械與物理研究所， 吉林 長春 130033；3.長光衛(wèi)星技術(shù)有限公司，吉林 長春 130033）

引言

在飛行器設(shè)計中，其尾翼的設(shè)計的水平對其性能有重要影響。對此類翼面結(jié)構(gòu)的設(shè)計涉及氣動、結(jié)構(gòu)、氣動彈性多個領(lǐng)域，傳統(tǒng)設(shè)計只在單一的氣動或結(jié)構(gòu)學(xué)科內(nèi)獨立優(yōu)化而忽略氣動彈性作用對性能參數(shù)的影響，或只考慮簡單的靜氣動彈性影響而不考慮其動態(tài)氣動彈性現(xiàn)象的作用。忽略這些彈性變形會影響翼面優(yōu)化設(shè)計結(jié)果的準(zhǔn)確性，使設(shè)計的翼面在實際飛行中很難達(dá)到預(yù)期的性能。因此，需要對其進(jìn)行多個學(xué)科的分析以研究其協(xié)同效應(yīng)并借助先進(jìn)的優(yōu)化策略和有效算法才能實現(xiàn)綜合性能的最優(yōu)化。近年來在飛行器的研發(fā)模式中，多學(xué)科優(yōu)化算法（Multidisciplinary Design Optimiza- tion, MDO）和快速迭代成為主要研究方向[1]。

目前，國內(nèi)外研究人員在航天領(lǐng)域的產(chǎn)品設(shè)計中，就多學(xué)科優(yōu)化理論的應(yīng)用進(jìn)行了廣泛的研究，多學(xué)科設(shè)計優(yōu)化（MDO）的研究起源于20世紀(jì)70年代，在20世紀(jì)80年代發(fā)展為一個單獨的研究領(lǐng)域。

到了20世紀(jì)90年代，MDO的研究已在世界范圍獲得了關(guān)注，不僅在理論研究上取得更多進(jìn)展，更是在工程上得到了實際應(yīng)用。比較典型的案例如Rockdyne公司的計劃用于X-33的塞式噴管發(fā)動機的設(shè)計，該設(shè)計以總升重比為目標(biāo)函數(shù)，采用耦合計算流體力學(xué)的模型、結(jié)構(gòu)有限元模型、計算熱力學(xué)模型等多學(xué)科的MDO模型，其結(jié)果不僅提高了設(shè)計質(zhì)量，還大大減少了設(shè)計時間。另一個計劃是高速民用飛機（HSCT）的設(shè)計，該研究利用分布式網(wǎng)絡(luò)和并行設(shè)計方法縮短了設(shè)計周期、獲得了更佳方案。國內(nèi)對這方面的研究是20世紀(jì)90年代中期開始的，通過向國外學(xué)習(xí)也取得了一些進(jìn)展，如文獻(xiàn)[2-4]研究了多學(xué)科設(shè)計算法在飛行器設(shè)計中的應(yīng)用。文獻(xiàn)[5]提出了混合多目標(biāo)粒子群算法在飛行器氣動布局設(shè)計中的應(yīng)用。文獻(xiàn)[6]對航天器性能樣機提出了一種氣動一體化多目標(biāo)優(yōu)化方法。

飛行器的氣動計算對整體方案設(shè)計影響較大，在傳統(tǒng)設(shè)計中，多采用基于有勢能、小擾動的工程算法，適用范圍小且需要后期修正，難以引入多學(xué)科優(yōu)化框架。現(xiàn)今對氣動特性計算主要依靠計算流體力學(xué)（CFD）仿真計算和風(fēng)洞實驗，這種高精度計算量很大，無法滿足設(shè)計中快速分析的需求。

目前在飛行器設(shè)計領(lǐng)域中針對多學(xué)科優(yōu)化和氣動外形優(yōu)化的問題多采用利用代理模型技術(shù)利用實驗設(shè)計法（DOE）獲得樣本并選擇近似方法來擬合氣動特性對外形參數(shù)的響應(yīng)，以達(dá)到避免復(fù)雜數(shù)值計算、縮短設(shè)計周期的目的。近年來代理模型在航空、汽車、船舶的設(shè)計和優(yōu)化中都得到了廣泛的應(yīng)用。

本文按照圖1的設(shè)計流程針對飛行器翼面的外形，對氣動特性參數(shù)進(jìn)行多種類型代理模型構(gòu)建，探索不同近似方法的精度和計算成本差別，選取最合適的近似算法參與翼面的氣動－結(jié)構(gòu)協(xié)同優(yōu)化。

圖1 代理模型構(gòu)建流程圖

1 設(shè)計空間

飛行器的飛行環(huán)境需要跨聲速飛行，應(yīng)采用大后掠翼，同時設(shè)計狀態(tài)屬于超聲速，應(yīng)采用小展弦比的薄彈翼，六角翼型，其結(jié)構(gòu)如圖2所示，其幾何參數(shù)含義見表1。

圖2 飛行器翼面幾何參數(shù)示意圖

表1 翼面幾何參數(shù)

翼面的建模幾何參數(shù)只適合對翼面進(jìn)行三維建模，對于描述一個固定的翼面形狀，往往需要將其換算為一系列無量綱參數(shù)即優(yōu)化所涉及的設(shè)計變量X={λ, χ,b1,b2,c’root, k}，其中λ代表展弦比、χ代表前緣后掠角、b1代表翼前緣特征、b2代表后緣特征、c’root代表翼根相對厚度和k代表翼梢收縮比。

參考飛行器設(shè)計工程經(jīng)驗[7]確定了設(shè)計變量的取值范圍，如表2所示。

表2 翼面設(shè)計變量區(qū)間

2 試驗設(shè)計

試驗設(shè)計技術(shù)（Design of Experiment, DOE）是多學(xué)科設(shè)計優(yōu)化中代理模型的取樣策略。試驗設(shè)計中，系統(tǒng)的輸入變量被稱為因素（Factor），其在樣本點處的值被稱為水平（Level），樣本點對應(yīng)對的輸出值被稱為響應(yīng)（Responce）。

2.1 全析因設(shè)計

全析因設(shè)計（Full factorial）是指在一次完全試驗中，系統(tǒng)的所有因素的所有水平可能的組合都要被研究到的一種試驗設(shè)計方法[8]。

假設(shè)系統(tǒng)輸入變量也即因素的個數(shù)為nv，每個因素對應(yīng)的水平數(shù)為ni（i =1,…,n），則對系統(tǒng)進(jìn)行全析因試驗所需的試驗次數(shù)為：

全析因試驗?zāi)軌蚍治鲆蛩貙ο到y(tǒng)影響的大小和分析因素間的交互作用，但當(dāng)系統(tǒng)的因素和水平比較多時，根據(jù)上式計算所得的試驗次數(shù)，即樣本點個數(shù)會將會是一個很大的數(shù)字，所以除了低維低水平的問題，全析因試驗很不適用。

2.2 中心復(fù)合設(shè)計

中心復(fù)合設(shè)計（Central Composite Design, CCD）[9]是一種針對二次多項式響應(yīng)面模型進(jìn)行分批實驗的一種試驗設(shè)計方法。該方法先挑出每個因素最大和最小的兩個水平值利用正交表Ln（2nv）安排n次試驗，試驗結(jié)束后在中心點做n0次重復(fù)試驗，最后在每個因素的坐標(biāo)軸上，取臂長為α的兩個對稱點作為樣本點，臂長的確定按下式：

若因素的個數(shù)為nv，三次試驗總樣本點數(shù)N=n+n0+2nv，N個點分布在以中心點為球心的兩個同心球上。

2.3 拉丁方設(shè)計

拉丁方（Latin Hypercube）[10]，它是一種分層抽樣法，將每個因素的設(shè)計空間均勻分成n份構(gòu)成矩陣，按水平數(shù)隨機組合成下標(biāo)，在設(shè)計空間矩陣上取樣本點，每個因素水平只可用一次。試驗點在設(shè)計空間內(nèi)均勻分布，且是隨機的，每次采樣的結(jié)果都不相同。這種設(shè)計均勻性比較好，適合因素數(shù)目較多的情況。

2.4 正交設(shè)計

正交試驗設(shè)計（Orthogonal design）也叫Taguchi設(shè)計，是應(yīng)用范圍最廣的試驗設(shè)計方法之一，由日本學(xué)者田口玄一于20世紀(jì)40年代后期首次應(yīng)用而得名。這種試驗設(shè)計方法按照一種己經(jīng)擬定好地滿足正交試驗條件的表格來安排試驗。表格稱為正交表（Orthogonal array），表示為LA（Pn），其中L代表正交表；下標(biāo)A表示表中橫行個數(shù)，即總共試驗次數(shù)；P為因素的水平數(shù)，n是因素的個數(shù)。當(dāng)遇到各因素水平數(shù)不等的試驗時，可使用不等水平正交表或采用擬水平法。正交表可從試驗設(shè)計參考書中獲得[11]。正交設(shè)計是多因子試驗中最重要的一種設(shè)計方法，設(shè)計變量最好不多于10個。

考慮到本文高精度模型的計算代價比較高，應(yīng)盡量減少分析的次數(shù)以提高設(shè)計效率，相比于全因子設(shè)計和中心復(fù)合設(shè)計，正交設(shè)計計算量較小，相比于均勻設(shè)計和拉丁方設(shè)計，正交設(shè)計精度較高[12-13]。鑒于此，本文選擇用正交設(shè)計技術(shù)來生成試驗點集合。對設(shè)計變量的水平數(shù)設(shè)置越多，樣本點數(shù)目也就跟著增加，這樣可以提高近似的精度。但考慮到相應(yīng)的計算量過大，對每個設(shè)計變量取5水平。

設(shè)計變量在設(shè)計空間上水平數(shù)的數(shù)據(jù)如表3所示。

表3 翼面設(shè)計變量水平數(shù)和水平值

取樣即為一組6因素5水平的試驗，可使用L25（56）正交表來安排試驗，即需要進(jìn)行25次不同輸入的氣動特性仿真，輸入安排如表1所示。

3 數(shù)值分析響應(yīng)值提取

本文對25組不同翼面形狀的飛行器進(jìn)行了基于CFD/ CSD緊耦合方式的高精度學(xué)科分析。對分析結(jié)果提取翼面優(yōu)化所需的參數(shù)即升力系數(shù)CL、阻力系數(shù)CD、結(jié)構(gòu)最大應(yīng)力σmax、翼梢變形δmax、彈翼重量Weight作為響應(yīng)值。

考慮到氣動彈性現(xiàn)象，與進(jìn)行了氣動彈性仿真后得到數(shù)值相比，后者更符合風(fēng)洞實驗值[14]。所以響應(yīng)值應(yīng)提取考慮了彈翼氣動彈性現(xiàn)象的仿真結(jié)果。

4 代理模型近似方案確定

近似技術(shù)也稱代理模型方法，是代理模型技術(shù)的核心，其本質(zhì)是利用已知樣本點以數(shù)學(xué)手段生成能夠反映設(shè)計變量與響應(yīng)值之間映射關(guān)系的數(shù)學(xué)模型。多項式響應(yīng)面（Polyno- mial Response Surface Method, PRSM）、徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)（Radial Basis Squares Method, RBF-ANN）和Kriging模型（Kriging Model, KM）是多學(xué)科設(shè)計優(yōu)化中常用的代理模型方法。

4.1 多項式響應(yīng)面

響應(yīng)面法是一種用簡單的代數(shù)函數(shù)來表示高精度模型分析信息的近似方法。工程中最常用的二階響應(yīng)面模型的數(shù)學(xué)表達(dá)式如下式所示：

式中n為設(shè)計變量的維數(shù)；xi,xj為設(shè)計變量；待定系數(shù)b0，bj，bjj，bij可由最小二乘法確定。為了更高的精度可建立階數(shù)更高的響應(yīng)面模型，但響應(yīng)面模型的階數(shù)越高，所需的樣本點即高精度模型的分析次數(shù)越多。如以二階響應(yīng)面模型為例，需要不低于（n+1）（n+2）/2個高精度模型的分析信息。

4.2 RBF神經(jīng)網(wǎng)絡(luò)

人工神經(jīng)網(wǎng)絡(luò)是一種模擬生物大腦結(jié)構(gòu)的信息處理系統(tǒng)，除了可以用來進(jìn)行函數(shù)逼近還可以進(jìn)行最近相鄰模型分類、概率密度估計等計算，在多個領(lǐng)域均有應(yīng)用。

神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)由輸入層、隱層和輸出層三層處理單元組成。信號輸入和輸出的單元層分別為輸入層和輸出層，輸入層和輸出層中間的單元層稱為中間層或隱層。

層中的數(shù)據(jù)處理單元稱神經(jīng)元。在神經(jīng)元中將輸入激勵轉(zhuǎn)化為輸出響應(yīng)的數(shù)學(xué)表達(dá)式稱為傳遞函數(shù)。

神經(jīng)網(wǎng)絡(luò)各神經(jīng)元通過權(quán)值連接，在各層類型決定后，通過已知樣本點的設(shè)計量及狀態(tài)量數(shù)據(jù)來調(diào)整確定權(quán)值即訓(xùn)練網(wǎng)絡(luò)。訓(xùn)練結(jié)束后，網(wǎng)絡(luò)的輸出即由輸入數(shù)據(jù)和各單元相連的各輸入量的權(quán)值來決定，此時即可用來模擬原有樣本和響應(yīng)值的映射關(guān)系。

BP神經(jīng)網(wǎng)絡(luò)和RBF神經(jīng)網(wǎng)絡(luò)是最常被用作代理模型的兩種模型，二者都是前饋網(wǎng)絡(luò)。相比之下后者的結(jié)構(gòu)更為簡單，網(wǎng)絡(luò)的訓(xùn)練過程更為快捷，并且在函數(shù)逼近和模式識別方面的表現(xiàn)也更為優(yōu)秀。所以對非線性函數(shù)的逼近，選用RBF神經(jīng)網(wǎng)絡(luò)較為合適。

圖3所示的是一個典型的r維輸入的單隱層RBF網(wǎng)絡(luò)。其中‖dist‖表示求輸入矢量p和權(quán)值矢量w的距離；b1,b2為閾值；n為隱層輸入：

圖3 r維輸入的單隱層RBF網(wǎng)絡(luò)

其中LW是隱層到輸出層的權(quán)值。

Jin等人[15]的研究表明RBF神經(jīng)網(wǎng)絡(luò)的近似精度和效率都不錯，適合在MDO中應(yīng)用。

4.3 Kriging方法

Kriging模型是Danie Krige于20世紀(jì)50年代提出的一種估計方差最小的無偏估計模型，由全局模型和局部偏差疊加而成。該模型早期主要應(yīng)用在地質(zhì)領(lǐng)域，到現(xiàn)在已經(jīng)成為MDO中比較有代表性的一種代理模型近似方法[16]，其具體模型為：

其中近似項g（x）是設(shè)計空間內(nèi)的全局模擬，g（x）可以看作一個常數(shù)β，β值可由已知響應(yīng)值進(jìn)行估計。估計項Z（x）期望為0，方差為σ2，協(xié)方差不為零的隨機函數(shù)，表示全局模擬的插值，是全局模擬基礎(chǔ)上的局部偏差。上式可變?yōu)椋?/p>

式中R是相關(guān)矩陣，矩陣為對角線元素為1的對稱矩陣；R是相關(guān)函數(shù)；ns為樣本點個數(shù)。R選擇高斯函數(shù)，可以表示為：

式中nv是設(shè)計變量個數(shù)，θk為未知相關(guān)參數(shù)向量，可取為常數(shù)以簡化運算。

根據(jù)Kriging理論，未知點x處的響應(yīng)值y的估計值y′可表示為：

式中y是樣本點響應(yīng)值組成的ns維列向量；g為長度為n的單位列向量；r（x）是未知向量x與樣本輸入數(shù)據(jù)之間的關(guān)系向量：

相關(guān)參數(shù)向量θk取為常數(shù)θ，可由極大似然估計變?yōu)橐痪S優(yōu)化問題：

4.4 代理模型綜合評估標(biāo)準(zhǔn)

為了判斷代理模型是否可以代替原有分析模型，需要一定的評估檢驗標(biāo)準(zhǔn)，包括精度評估、效率評估和實現(xiàn)難度評估。

考慮到Kriging和RBF這樣的近似方法具有插值特性。所以需要在構(gòu)造的樣本點以外重新選取另外一批樣本點作為測試樣本點來參與對精度的評價。本文采用均方根誤差（Root Mean Square Error, RMSE）和復(fù)相關(guān)系數(shù)（R Squared Coef- ficient, R2）兩個標(biāo)準(zhǔn)來檢驗該模型的預(yù)測值對真實值的代理精度：

其中y和yreg分別是設(shè)計空間內(nèi)每個樣本點的真實響應(yīng)值和模型預(yù)測值；為所有樣本點真實響應(yīng)值的平均值；Ngrid是樣本點的個數(shù)。RMSE值越接近于0則代理模型精度越高。復(fù)相關(guān)系數(shù)R2的取值越接近于1則近似精度越高。

效率評估考慮的是代理模型的計算效率，評價的指標(biāo)是構(gòu)造代理模型所需要的成本與使用代理模型預(yù)測新的設(shè)計點響應(yīng)值所需的成本。這個成本包括計算時間計算所需的PC機內(nèi)存等，一般用計算時間來評價。

代理模型的實現(xiàn)難度指的是在軟件上的功能實現(xiàn)。目前工程上最常用的實現(xiàn)手段是基于Matlab環(huán)境的各種函數(shù)工具箱來實現(xiàn)自行開發(fā)。也有一些優(yōu)化軟件自身帶有一套近似模型模塊，可以省去編程而方便工程人員使用。但是由于程序代碼固定，當(dāng)需要對程序進(jìn)行修改時往往比較困難。一般來說，近似方法的原理越復(fù)雜，軟件開發(fā)就越困難。但是為了獲得較高的精度，有時必須對模型反復(fù)修改使之復(fù)雜性增加。

本文以設(shè)計變量集合X={λ，χ，b1，b2，c’root，k}為輸入，分別以響應(yīng)值集合Y1={CL}，Y2={CD}，Y3={σmax}, Y4={δmax}，Y5={Weight}為輸出，采用多項式響應(yīng)面、RBF神經(jīng)網(wǎng)絡(luò)和Kriging方法為近似方法各自創(chuàng)建了近似模型。在設(shè)計空間內(nèi)隨機生成五個測試樣本點，以均方根誤差RMSE和復(fù)相關(guān)系數(shù)R2為代理模型精度評價標(biāo)準(zhǔn)。

用這三種模型的近似精度校驗結(jié)果和時間統(tǒng)計分別如表4和表5所示：

表4 近似精度校驗結(jié)果

表5 代理模型時間統(tǒng)計

5 結(jié)論

從表中數(shù)據(jù)可以看出，多項式響應(yīng)面得到得近似模型整體效果最差，所有響應(yīng)的精度都低于其他兩種方法，產(chǎn)生這種結(jié)果的原因是設(shè)計變量維數(shù)（即因素數(shù)目）有6個，要構(gòu)建一個6因素的二階多項式需要至少（6+1）（6+2）/2=28個樣本點，而樣本點數(shù)目只有25個，所以造成只能使用一階多項式來近似模型，這樣顯著降低了近似精度。不過由于許多研究以表明多項式響應(yīng)面法本身不適合高階非線性問題，所以這種方案對本文翼面優(yōu)化的代理模型來說是不可取的。

Kriging方法與RBF神經(jīng)網(wǎng)絡(luò)相比，對某些數(shù)據(jù)如CD和Weight的近似效果比較好，但是對CL和δmax的近似精度則較低，RSM的值都高于0.1，R2值也都低于0.7。

RBF神經(jīng)網(wǎng)絡(luò)對所有數(shù)據(jù)的RMSE值均低于0.01，R2值均大于0.8。所以其整體近似精度要更為優(yōu)秀，具有更好的魯棒性。

就構(gòu)造代理模型的時間而言，多項式響應(yīng)面所用的時間最少，Kriging方法需要的時間最多，這與這三者本身數(shù)學(xué)模型的復(fù)雜度相符合。在主頻2.6 GHz，單核的PC機上這三種方法所耗費的預(yù)測時間均在1 s以內(nèi)，而主頻2.6 GHz，8核的服務(wù)器上對任一個樣本點進(jìn)行高精度仿真的計算時間則在20 h以上。相比而言，使用代理模型大大提高了學(xué)科分析的效率。

通過上述研究可以看出，多項式響應(yīng)面雖然效率最高但精度過低，Kriging方法可以對測試問題取得較好的精度，但是魯棒性差，而且耗時較多，RBF神經(jīng)網(wǎng)絡(luò)各項性能較為平均，綜合性能最好。考慮到氣動優(yōu)化時所需要的迭代點很多，而且數(shù)據(jù)之間關(guān)系復(fù)雜，本文選用了RBF神經(jīng)網(wǎng)絡(luò)作為翼面優(yōu)化設(shè)計的代理模型近似方法。