鄧凱文 陳海昕

(清華大學(xué)航天航空學(xué)院，北京100084)

基于差分進化和RBF響應(yīng)面的混合優(yōu)化算法1)

鄧凱文 陳海昕2)

(清華大學(xué)航天航空學(xué)院，北京100084)

針對氣動優(yōu)化等昂貴優(yōu)化問題，提出了一種基于差分進化和RBF響應(yīng)面的混合優(yōu)化算法HSADE，該方法結(jié)合了差分進化算法的強全局尋優(yōu)能力和RBF響應(yīng)面方法的快速局部搜索能力，能夠同時有效地提高算法的局部搜索效率和全局尋優(yōu)能力.對各子算法中的策略和邏輯進行了多項改進，提出和應(yīng)用了基于雙敗淘汰賽的競賽賽制和參數(shù)自適應(yīng)等改進策略.對HSADE使用多個典型算例進行了測試，并橫向?qū)Ρ攘薔SGA-II，MOPSO和多目標(biāo)差分進化算法.測試結(jié)果表明，在大多數(shù)問題中HSADE在以世代距離表征的局部搜索效率和以超體積比表征的全局尋優(yōu)能力兩項指標(biāo)上都優(yōu)于其他算法，證實了以上混合策略及算法改進的有效性.將該算法應(yīng)用于一個翼型優(yōu)化問題和一個二維超聲速噴管膨脹面優(yōu)化問題，并橫向?qū)Ρ任唇?jīng)改良的差分進化算法DE和另一種混合算法NARSGA，結(jié)果表明在接近1000次的函數(shù)評估下，HSADE能相對其他算法進一步對翼型減阻0.5 count，在噴管優(yōu)化中HSADE得到的結(jié)果也好于其他兩種算法，表明該方法具有較強工程應(yīng)用價值.

差分進化,響應(yīng)面方法,計算流體力學(xué),多目標(biāo)優(yōu)化,優(yōu)化設(shè)計

引言

隨著高性能計算技術(shù)的快速發(fā)展，高精度數(shù)值模擬在飛機氣動設(shè)計中已得到大量應(yīng)用.盡管基于高精度數(shù)值模擬的優(yōu)化設(shè)計已被廣泛應(yīng)用且取得了較好的效果[1-3]，但其進一步應(yīng)用仍受到較大的限制.尤其是面對“昂貴優(yōu)化問題”[4]，即目標(biāo)函數(shù)評估的計算耗費很大的時候.以飛行器氣動優(yōu)化為例，三維客機翼身組合體全模的雷諾平均方程(RANS)計算，需使用6000萬網(wǎng)格.利用64個核心的計算節(jié)點進行計算，一個工況耗時通常可達36h.在存在時限的工程設(shè)計中，以上復(fù)雜度的問題的總計算次數(shù)往往是嚴(yán)格受限的，這就要求優(yōu)化方法能夠在限制總函數(shù)評估次數(shù)下具有更好的優(yōu)化效果.

因此針對類似于氣動優(yōu)化的“昂貴優(yōu)化問題”，優(yōu)化算法應(yīng)有較高的優(yōu)化效率.局部搜索效率和全局尋優(yōu)能力是評價優(yōu)化算法效率最重要的兩大指標(biāo).通常而言，啟發(fā)式算法的局部優(yōu)化效率要大大低于梯度類算法，但具有更好的全局尋優(yōu)能力.為提高啟發(fā)類算法的效率，一種做法是與其他優(yōu)化方法組合得到“混合算法”，以提升其綜合性能.

近年來研究人員在混合優(yōu)化算法上開展了大量的改進和嘗試.提出了包括以啟發(fā)式算法結(jié)合基于梯度方法構(gòu)建的混合算法[5-6]、以代理模型方法為主[7-8]的混合算法、以及啟發(fā)式算法和代理模型方法組合的混合算法[9-11]等混合方法.在氣動優(yōu)化設(shè)計領(lǐng)域，研究者對遺傳算法[12-14]、代理模型方法[15-16[17-18]開展了大量測試和應(yīng)用.在此基礎(chǔ)上，混合算法也因其較好的綜合性能得到了該領(lǐng)域?qū)W者的關(guān)注，文獻[19]提出了一種基于粒子群算法和Kriging響應(yīng)面的應(yīng)用于氣動優(yōu)化問題的混合算法，文獻[20]提出了一種將文化基因算法(memetic algorithm)作為局部搜索算子以提升遺傳算法搜索效率的混合算法并應(yīng)用于氣動優(yōu)化，文獻[21]提出了一種基于NSGA-II[22]和Kriging響應(yīng)面的混合優(yōu)化算法用于多段翼型等昂貴問題的優(yōu)化.

目前廣為應(yīng)用的啟發(fā)式算法多具有精英保留機制[23]，即算法會保留精英個體而逐代淘汰較劣個體，典型的例子有SPEA2[24]，NSGA-II[22]等.因此在擁有精英保留機制的混合算法中，若引入輔助方法加強局部的深度搜索，其得到的個體往往趨近局部最優(yōu)因而性能指標(biāo)優(yōu)秀，它們在接下來的數(shù)代中會引導(dǎo)種群進化直至被后代超越；若引入輔助方法進行增強的廣度搜索，其最終得到的個體雖存在較好的多樣性，但由于性能指標(biāo)通常不佳而難以存活多代，在降低算法效率的同時并不能很好地維持種群的多樣性.基于優(yōu)勢互補的原則，一種理想的混合算法組合是一種具有強全局尋優(yōu)能力的主算法結(jié)合幫助其進行快速局部搜索的輔助方法.

差分進化 (di ff erential evolution,DE)是一種于1995年提出的啟發(fā)式算法[25]，聞名于其強魯棒性和對高維問題良好的適應(yīng)性，擁有很強的全局尋優(yōu)能力.徑向基函數(shù)(radial basis function,RBF)響應(yīng)面(response surface,RS)方法是一種使用徑向基函數(shù)進行插值的代理模型方法，也是現(xiàn)在人工神經(jīng)網(wǎng)絡(luò)最著名的模型之一[26]，在工程優(yōu)化中被廣為使用.這種方法相較于工程中常用的其他代理模型如Kriging響應(yīng)面，薄板樣條曲面等有著光順性好，數(shù)學(xué)描述簡單的特點，同時它的精度也比較高[27].

本文根據(jù)以上思路，提出了一種針對氣動優(yōu)化等昂貴優(yōu)化問題的基于差分進化和RBF響應(yīng)面的自適應(yīng)混合算法HSADE(hybrid self-adaptive di ff erential evolution)，通過在魯棒性較好的差分進化算法的主流程中嵌入高精度的RBF響應(yīng)面方法，利用代理模型擬合已計算數(shù)據(jù)，給出最佳個體預(yù)測，加入差分進化算法的種群，以提高算法整體的數(shù)據(jù)利用率；從而通過較小的響應(yīng)面計算量換來較大的搜索效率提升. HSADE以差分進化為主體，因此也是以代為單位進行迭代優(yōu)化的.

值得注意的是，在優(yōu)化過程中，響應(yīng)面插值并不取代CFD計算，它的功能是為主算法在設(shè)計變量空間內(nèi)搜索新的具有潛在最優(yōu)性能的待分析個體.因此整個優(yōu)化過程的仍基于高精度的CFD分析.而RBF響應(yīng)面會隨優(yōu)化進程中新的計算結(jié)果實時更新以提升最優(yōu)點預(yù)測能力.

此外，本文在提出兩種方法的混合策略同時也對它們各自進行了一些改進，以提高算法整體效能.

1 混合算法HSADE

1.1 算法整體結(jié)構(gòu)

HSADE的核心思想是利用多種優(yōu)化策略給出的蘊含不同特性的個體，確保種群多樣性，強化精英個體，以同時提高算法的局部和全局尋優(yōu)能力.

具體而言，HSADE的思路是以差分進化為主干，在每一代的進化結(jié)束階段向種群中依次嵌入由RBF響應(yīng)面給出的提升算法搜索效率的個體和由空間增量填充法(increment space filling ISF)給出的提升種群多樣性的個體，和原有差分進化的種群中通過擇優(yōu)機制保留的高度多樣的個體共同構(gòu)成子代個體進入下一代，以提升算法的整體效率.同時考慮到為提高種群目標(biāo)函數(shù)評估的計算效率，每代各個體的函數(shù)評估將同時并行執(zhí)行.

HSADE的優(yōu)化流程概括在圖1中，描述如下：

(1)初始化算法參數(shù)，置代數(shù)k=1，使用中心采樣方法和空間增量填充法為主的為實驗設(shè)計方法(design of experiments,DOE)產(chǎn)生初代種群Pk作為設(shè)計起點，并為待構(gòu)建響應(yīng)面提供一個分布均勻的樣本集；其中空間增量填充法旨在按順序找到一系列個體，其設(shè)計變量x?滿足如下特性

其中xi代表當(dāng)前種群中各個體的設(shè)計變量；

(2)將種群Pk按照差分進化的流程依次進行變異和交叉，產(chǎn)生試驗種群Uk；使用RBF響應(yīng)面擬合所有已計算點的數(shù)據(jù)，使用差分進化算法優(yōu)化響應(yīng)面參數(shù)以最小化擬合誤差，使用差分進化算法在響應(yīng)面上尋優(yōu)，得到潛在最優(yōu)個體集合Ok；利用空間增量填充法給出增加種群多樣性的個體集合Dk；

(3)并行評估Uk,Ok和Dk中各個體目標(biāo)函數(shù)；

(4)對Pk和Uk按照雙敗淘汰制執(zhí)行選擇操作并淘汰末位個體，得到獲勝種群Pk′；

(5)令Pk+1=Pk′∪Ok∪Dk，若未滿足結(jié)束條件，令k=k+1，重新回到步驟(2).

優(yōu)化流程中種群大小設(shè)為Np，總優(yōu)化代數(shù)為Ngene，其余各主要參數(shù)都經(jīng)自適應(yīng)改造以提高算法通用性.

1.2 子算法介紹

1.2.1 差分進化

圖1 HSADE優(yōu)化流程圖Fig.1 Optimization fl w chart of HSADE

差分進化算法是一種將實數(shù)編碼應(yīng)用于連續(xù)空間優(yōu)化的啟發(fā)式算法，和遺傳算法相似，它主要有三個過程：變異、交叉和選擇.為避免混淆，將優(yōu)化問題中一個設(shè)計變量和其對應(yīng)的目標(biāo)函數(shù)組合起來統(tǒng)稱為一個個體個體上標(biāo)表示代數(shù)，下標(biāo)表示其在種群中的位置.當(dāng)代的種群記作目標(biāo)函數(shù)對設(shè)計變量的映射記為按以上定義，差分進化的過程可被簡單歸納為以下步驟：

(1)讀取參數(shù)，置代數(shù)k=1，初始化種群Pk；

(2)令Pk為父代種群，對Pk中的每一個個體執(zhí)行(3)～(5)步；

(3)從當(dāng)前種群中挑選n個基矢量，n取決于變異方式

(6)令k=k+1，若滿足結(jié)束條件，輸出最終種群相關(guān)信息，否則回到步驟(2).

以上步驟(2)～(5)的具體形式均有多種變體，因此差分進化算法并不局限于單一的形式，它可以用統(tǒng)一的記號描述：DE/X/Y/Z[25].X表示變異過程中對基矢量的選擇方式，主要有/rand/，/best/和/targetto-best/三種；Y表示參與變異的矢量對數(shù)；Z表示交叉方式，詳細(xì)介紹見文獻[25].HSADE選用/target-tobest/1作為差分進化算的變異方式，選用二項式交叉作為其交叉方式.

1.2.2 徑向基函數(shù)響應(yīng)面

徑向基函數(shù)響應(yīng)面是一種基于徑向基函數(shù)的插值方法，假設(shè)已獲得一系列已知的采樣點，(xi,yi)(i=1,2,···,NS),dim(xi)=Nx,dim(yi)=Ny.對于一個待評估點(x,y)，在已知x的情況下由下式給出y的估計值

其中，λ為大小為NS×Ny的權(quán)重系數(shù)矩陣，φ(x)為大小為1×NS的行向量，它們分別按下式確定

φ(x)的元素為

Φ的元素為

F的元素為

式中，Φ和F分別為基函數(shù)矩陣和目標(biāo)函數(shù)矩陣，?為徑向基函數(shù)，其較為常用的類型被總結(jié)在表1中. HSADE采用的基函數(shù)為逆多項式函數(shù)(inverse multiquadric).

表1 徑向基函數(shù)類型Table 1 Di ff erent types of radial basis functions

2 子算法的改進

在引入基本算法形成混合算法的同時，本文也對基本算法分別提出了一些改進措施，以提升算法效率.

2.1 對于差分進化的改進

2.1.1 基于雙敗淘汰制的選擇算子

為進一步提高算法的整體效率，HSADE采用一種稱為雙敗淘汰賽的選擇算子，設(shè)每代的選擇過程在父代種群和以父代種群生成的試驗種群間進行，且父代種群中各個體已按非占優(yōu)排序方法[22]和最近鄰域距離進行了從最優(yōu)到最劣的排序，令Uk中的個體與Pk一一對應(yīng)，則選擇流程如下：

(1)在父代種群中抽取前Nb個個體，Nb為預(yù)定義的正整數(shù)滿足Nb＜Np，則子代種群的前Nb個個體由經(jīng)基于種群環(huán)境的擇優(yōu)準(zhǔn)則在2.1.2節(jié)詳述)競爭得出；若負(fù)于則將加入敗者組B，設(shè)B的最終大小為NB；

(3)將混合組C由非占優(yōu)排序方法和最近鄰域距離進行排序，選出最好的Np-Nb個個體，并將其保留，作為子代種群的后個個體

2.1.2 種群信息輔助的擇優(yōu)邏輯

在傳統(tǒng)差分進化中，擇優(yōu)在父代和子代種群的個體間一一進行，這種策略存在可能淘汰潛在的較優(yōu)個體的弊病.一種典型的情況是，當(dāng)使用差分進化算法求解一個雙目標(biāo)最小化優(yōu)化問題，假設(shè)此時種群內(nèi)正在進行個體A={xA,yA}和個體B={xB,yB}的競爭以選出勝者進入子代，個體A,B和當(dāng)代種群的Pareto前緣如圖2所示.

圖2 一種典型的擇優(yōu)情景Fig.2 A typical selection occasion

在這種情況中，因為個體A和B是Pareto互不占優(yōu)的，且個體B的擁擠程度優(yōu)于A(A距離上游點更近)，因此傳統(tǒng)的擇優(yōu)邏輯會傾向于保留B而淘汰A.但實際上個體A才是真正利于種群整體進化的，因為它位于當(dāng)前Pareto前緣下方，能將其向前推進.這是傳統(tǒng)差分進化中擇優(yōu)邏輯存在的漏洞，它忽略了競爭過程中種群環(huán)境的信息.在擇優(yōu)過程中這種邏輯對個體的優(yōu)劣判斷僅以其競爭對手作為唯一參照，這會導(dǎo)致潛在的判斷失準(zhǔn).對于這個問題，本文采用一種如圖3描述的擇優(yōu)邏輯，稱為種群信息輔助的擇優(yōu)邏輯，在擇優(yōu)過程中種群內(nèi)各個體信息會被引入擇優(yōu)決策過程以挑選出最利于種群進化的個體進入子代.

圖3 種群信息輔助的擇優(yōu)邏輯Fig.3 Population information enhanced selection logic

2.1.3 參數(shù)自適應(yīng)策略

在基本差分進化算法中，算法性能嚴(yán)重依賴于參數(shù)的選擇[28-29]. 為增強算法對問題的適應(yīng)性, HSADE中應(yīng)用了一種簡潔的參數(shù)自適應(yīng)策略，具體方法為：在算法開始階段，遍歷種群中所有個體，為其分配一個獨立的變異常數(shù)Fi和交叉概率因子CRi，其中i為個體編號.若在指定代數(shù)Nc內(nèi)位于種群相同位置的個體在種群中的排名沒有上升，則按下式更新其算法參數(shù)

其中，下標(biāo)best表示在指定代數(shù)Nc內(nèi)Pareto占優(yōu)等級上升最多的個體的參數(shù)，R為[0,1]內(nèi)的隨機數(shù)，以上參數(shù)自適應(yīng)邏輯是一個促使各個體參數(shù)逐漸向更適宜進化的參數(shù)靠近的過程.

2.2 對于徑向基函數(shù)響應(yīng)面的改進

2.2.1 響應(yīng)面參數(shù)優(yōu)化、誤差的度量與快速估計

在優(yōu)化進行的每一代中，HSADE會對響應(yīng)面的參數(shù)進行優(yōu)化以最小化插值誤差.由于RBF響應(yīng)面為完全插值，且在昂貴氣動優(yōu)化問題中測試點是稀缺的，本文采用文獻[30]的交叉驗證方法定義響應(yīng)面誤差以充分利用已計算數(shù)據(jù)，避免響應(yīng)面過擬合.具體而言，按式(7)定義誤差矩陣E

式中Si(xi)為響應(yīng)面在樣本集合中去掉第i個個體時,在設(shè)計變量為xi時給出的估計值.

如此定義的誤差矩陣E可按下式快速計算

式中λij為權(quán)重矩陣的分量，X=Φ-1，其余定義見式(2)～式(4).令nj為第j個目標(biāo)函數(shù)的歸一化因子，定義響應(yīng)面擬合誤差eRS為

2.2.2 RBF響應(yīng)面參數(shù)重分配策略

本文提出了一種對RBF響應(yīng)面參數(shù)進行重分配的策略.如上小節(jié)描述，經(jīng)典徑向基函數(shù)響應(yīng)面理論中，多目標(biāo)插值可以表示為

式中，每個目標(biāo)函數(shù)在插值時對應(yīng)每個采樣點的基函數(shù)?都是相同的，為一個統(tǒng)一預(yù)定義的徑向基函數(shù).由于采樣點在設(shè)計變量空間的分布存在疏密不同的區(qū)域，物理直觀上，這些區(qū)域中的采樣點對應(yīng)的基函數(shù)應(yīng)該具有差異化的函數(shù)參數(shù)以體現(xiàn)采樣點在空間的分布.

本文通過重分配響應(yīng)面方法中對應(yīng)于各樣本點的徑向基函數(shù)參數(shù)，以適應(yīng)不同區(qū)域樣本點的疏密度，達到提升插值精度的目的.以逆多項式函數(shù)為例，具體操作如下：

(1)對編號為i的采樣點，取ri為設(shè)計變量空間中其他采樣點與該采樣點的最小歐氏距離，用于描述其稀疏程度；

(2)定義一列差異化的徑向基函數(shù)

其中ci=kSri，kS為全局縮放因子，原徑向基函數(shù)響應(yīng)面的待優(yōu)化參數(shù)c被kS取代；

(3)重新定義插值公式和徑向基函數(shù)矩陣

(4)根據(jù)重新定義的徑向基函數(shù)矩陣由式(3)求新的權(quán)重矩陣.

響應(yīng)面參數(shù)經(jīng)重分配后，對應(yīng)每個采樣點的基函數(shù)參數(shù)都不相同，具體效果見于下節(jié)的測試部分.

3 算法測試

本節(jié)將使用一些典型測試函數(shù)驗證和評估RBF響應(yīng)面的擬合精度和HSADE的綜合優(yōu)化效能.

3.1 響應(yīng)面擬合精度驗證

本節(jié)將基本RBF響應(yīng)面，經(jīng)過2.2.2節(jié)參數(shù)重分配后的RBF響應(yīng)面(簡寫為SRBF)和Kriging響應(yīng)面在一些典型測試函數(shù)[31]上的擬合精度進行對比，對于每個測試函數(shù)，分別取NS=20,40,80,120個隨機采樣點，NT=1600個固定測試點.每種情況獨立運行80次.使用平均歸一化均方根誤差作為3種方法擬合誤差的評價指標(biāo)，定義為

每次測試中,每種待測試響應(yīng)面的參數(shù)均使用差分進化算法來最小化擬合誤差,以避免參數(shù)選取不當(dāng)導(dǎo)致的對比失真，對于SRBF和RBF兩種方法優(yōu)化過程的誤差度量采取2.2節(jié)的定義.

測試結(jié)果總結(jié)在表2中，可以看出：

(1)SRBF在所有測試問題中的擬合精度較RBF明顯提高；

(2)SRBF在部分測試問題中的擬合精度超過了以精度見長的Kriging響應(yīng)面，并且這種擬合精度的優(yōu)勢隨著采樣點數(shù)目的增加越發(fā)明顯.在采樣點個數(shù)為20個時，SRBF在9個測試問題中僅有4個優(yōu)于Kriging響應(yīng)面，而當(dāng)采樣點數(shù)目逐步增加至80時，SRBF在8個問題的測試結(jié)果都優(yōu)于Kriging響應(yīng)面.

可以看出本文提出的對RBF響應(yīng)面改進措施是可靠有效的.

表2 幾種響應(yīng)面插值精度驗證Table 2 Examination of interpolation error of di ff erent response surfaces

3.2 HSADE優(yōu)化效能驗證

3.2.1 評價指標(biāo)

(1)世代距離GD

世代距離表征了當(dāng)前最優(yōu)解集距理論最優(yōu)解集的距離，是算法局部搜索效率的有效度量，按下式定義[22]

式中，Npf為當(dāng)前優(yōu)化結(jié)果最優(yōu)解集中的解的數(shù)目，di為該最優(yōu)解集中第i個解在目標(biāo)函數(shù)空間距理論最優(yōu)解集的最小歐氏距離，世代距離越小越好.

(2)超體積比HVR

超體積比表征了當(dāng)前算法搜索到的結(jié)果占整個搜索空間的比例，是算法全局尋優(yōu)能力的一種度量，定義如下[32]

其中HV表示當(dāng)前算法搜索結(jié)果在目標(biāo)函數(shù)空間所能覆蓋的超體積，定義為

HVi為以其中第i個個體和該優(yōu)化問題各目標(biāo)函數(shù)的最劣值構(gòu)成的個體為對角的在目標(biāo)函數(shù)空間的超立方體的體積，HV?為理論最優(yōu)解的超體積，超體積比越大越好.

3.2.2 HSADE優(yōu)化效能驗證

本文選擇ZDT1-4，ZDT6，DTLZ1-4和DTLZ6問題[22-23]來驗證HSADE的優(yōu)化效能，作為對比參考的算法有MOPSO[34]，NSGA-II(二進制編碼)和基本DE.

考慮到對測試集不同特性的要求，ZDT1-3問題的設(shè)計變量分別取30維和80維的版本，考察算法對高維問題的優(yōu)化能力；ZDT4、ZDT6和DTLZ1-4、DTLZ6問題測試其10變量、2目標(biāo)的情況，考察算法對復(fù)雜問題的優(yōu)化能力，對于ZDT4問題，取其簡化版本，令xi∈[0,1],i=2,3,···,10，其余參數(shù)范圍不變.

參數(shù)設(shè)置方面，由于算法在不同問題中的表現(xiàn)依賴于參數(shù)選擇，因此對于每種算法，本文應(yīng)用多種參數(shù)組合以盡可能消除參數(shù)選擇不當(dāng)導(dǎo)致的性能衡量失準(zhǔn).對于DE，變異交叉模式采用/rand/1/Bi，其參數(shù)選擇兩組：①F=0.5,CR=0.2；②F=0.5,CR=0.7；對于NSGA-II，取其基因數(shù)恒定為8，變異率取4個不同值：0.01,0.02,0.05,0.10；交叉率恒定為50%，交叉方式為二進制交叉.對于MOPSO，令優(yōu)化起始點隨機生成，學(xué)習(xí)因子范圍定為[1.0,2.0]，慣性權(quán)重取3個不同值，分別為0.4，0.6和0.8.

對于 HSADE，變異交叉模式采用/rand/1/targetto-best，自適應(yīng)參數(shù)范圍為F∈ [0.1,0.9]，CR∈[0.1,0.5]，Nb取為16.每代中響應(yīng)面給出2個個體，空間增量填充給出1個個體.

為保證公平性以及減少隨機因素的影響，對于包括 HSADE在內(nèi)的每種算法采取隨機初始化種群，對每個問題獨立運行10次.每次優(yōu)化過程每種算法共迭代32代，種群大小皆為32，總共約1024次函數(shù)評價.對于每次運行結(jié)果計算其世代距離和超體積比，取10次運行結(jié)果的中位數(shù)進行比較.

表3和表4給出了各優(yōu)化問題的比較結(jié)果.

由表3可以看出，在世代距離表征的收斂性指標(biāo)中，除問題DTLZ4中HSADE劣于NSGA-II，其余問題中HSADE都顯著好于其他算法.

表3 幾種算法在測試問題上的收斂性指標(biāo)結(jié)果Table 3 Convergence criteria of competing algorithms upon tested problems

表4 幾種算法在測試問題上的分布性指標(biāo)結(jié)果Table 4 Distribution criteria of competing algorithms upon tested problems

由表4可以看出，在超體積比表征的分布性指標(biāo)中，HSADE除了在 80維 ZDT3問題中劣于MOPSO，其余都優(yōu)于其他算法.

從結(jié)果的比較可以看出,引入的 RBF響應(yīng)面方法可以有效地提高算法的局部尋優(yōu)能力，同時HSADE的混合結(jié)構(gòu)也讓其保留了較強的多樣性保持能力，其綜合優(yōu)化能力較未經(jīng)改良的DE有明顯提高.

4 應(yīng)用驗證

為進一步驗證HSADE在實際氣動優(yōu)化問題中的表現(xiàn)，本節(jié)展示應(yīng)用HSADE，NARSGA[21](一種基于NSGA-II和Kriging響應(yīng)面的混合算法)和傳統(tǒng)差分進化算法解決兩個典型氣動優(yōu)化問題的實際效果.

4.1 翼型優(yōu)化

翼型優(yōu)化是氣動優(yōu)化領(lǐng)域中較為基本和具有代表性的問題.本節(jié)展示應(yīng)用以上3種算法對一個10%厚度的二維翼型在Ma=0.72和Ma=0.75兩種工作狀態(tài)進行減阻的具體效果.

優(yōu)化問題中，翼型上下表面各設(shè)7個控制點，由CST(class shape function transformation)[35]曲線擬合生成.優(yōu)化問題的詳細(xì)描述總結(jié)于表5，其中Cd1和Cm1表示Ma=0.72工況下的翼型阻力系數(shù)和力矩系數(shù)，Cd2表示Ma=0.75工況下的翼型阻力系數(shù)，R為翼型前緣半徑.優(yōu)化起點選定10%厚度的RAE2822等一系列翼型，如圖4(a)所示，其中Ref1-2為從課題組翼型庫中選取的層流翼型.初始種群中，除手動輸入的優(yōu)化起點外，其余個體由1.1節(jié)所述的實驗設(shè)計方法進行采樣得到.

表5 翼型優(yōu)化問題的優(yōu)化變量、目標(biāo)和約束Table 5 Variables,targets and constraints of the airfoil optimization problem

圖4 優(yōu)化進程Fig.4 Optimization progress

CFD計算采取定升力系數(shù)0.30，起始攻角1.5°，雷諾數(shù)取為1.50×107，馬赫數(shù)取為0.72和0.75.湍流模式為兩方程SST，離散格式為Roe格式，網(wǎng)格數(shù)約為3萬，CFD計算使用NSAWET程序[36-37].

算法參數(shù)設(shè)置為：

(1)對于HSADE，CR定為0.2，F(xiàn)范圍為[0.1,0.9]，每代由DE給出個體29個，ISF給出1個，RBF響應(yīng)面給出2個，Nb為16；

(2)對于DE，CR定為0.2，F(xiàn)為0.5，其余設(shè)置等同于HSADE；

(3)對于NARSGA，取文獻中默認(rèn)參數(shù).

每種算法優(yōu)化總代數(shù)設(shè)為32代，種群大小皆為32，排除重復(fù)個體后各方法皆各進行了約1000次目標(biāo)函數(shù)評估.3種算法優(yōu)化過程中的世代距離收斂曲線，超體積比和響應(yīng)面誤差變化見圖4，計算世代距離的參考最優(yōu)解集取3種算法解集集合的非占優(yōu)解.

3種算法優(yōu)化得到的最優(yōu)翼型與優(yōu)化起點各翼型的幾何外形對比見圖5(b)～圖5(d)中，每種算法在2種工況下的最小阻力個體表面壓力分布對比見圖6，圖中c代表弦長.3種算法得到的Pareto前緣和優(yōu)化起點中的最優(yōu)性能翼型的對比見圖7，3種方法得出的兩種工況阻力最小值(圖表中標(biāo)注為minCd1和minCd2)和Pareto前緣中點個體(圖表中標(biāo)注為middle)的性能和對比優(yōu)化起點的優(yōu)化量見表6，從優(yōu)化結(jié)果可以看出：

(1)從圖4的優(yōu)化過程看，HSADE的世代距離收斂速度明顯強于未經(jīng)改良的DE，在末端雖然存在振蕩，但在優(yōu)化末期仍保持對于NARSGA的輕微優(yōu)勢；在超體積比發(fā)展上，HSADE在整個優(yōu)化過程中都明顯領(lǐng)先于NARSGA和DE，顯示出HSADE具有對目標(biāo)函數(shù)空間更強的探索能力.另外在優(yōu)化過程中HSADE的響應(yīng)面的精度經(jīng)校核處于較小的范圍內(nèi)，滿足使用標(biāo)準(zhǔn).

圖5 初始外形和3種算法得到的外形Fig.5 Shapes of initial airfoil and optimal airfoils obtained by optimizers

圖6 最優(yōu)翼型表面壓力分布形態(tài)Fig.6 Surface pressure distribution of optimal airfoils

(2)從圖5看，在約束了升力的情況下，優(yōu)化得出的最優(yōu)翼型中，HSADE的結(jié)果在翼型下表面存在更大的差異性，表明了HSADE能夠在優(yōu)化過程中很好地保持種群的多樣性，使優(yōu)化進程不易發(fā)生停滯或早熟收斂.

(3)從表6和圖7中看，HSADE得到的翼型在兩個工況的阻力指標(biāo)下都明顯優(yōu)于DE和NARSGA，在接近的CFD模擬次數(shù)下阻力優(yōu)化量比另兩種算法領(lǐng)先達0.5 count，HSADE得到的Pareto前緣明顯優(yōu)于另外兩種方法.

(4)從壓力分布上看，在兩個工況下，在翼型前緣，和初始翼型壓力分布相比3種算法得到的阻力最小個體的逆壓梯度區(qū)有所縮小或完全消除；DE的結(jié)果在Ma=0.72工況下在下表面出現(xiàn)了二次加速的現(xiàn)象；Ma=0.75工況下，在翼型前緣HSADE和NARSGA的結(jié)果出現(xiàn)了一道微弱的激波，但它們在最大厚度前的壓力分布較初始翼型更顯合理，較大地減小了形阻和總阻力；兩種工況下HSADE的結(jié)果壓力恢復(fù)最為平穩(wěn)，顯示出較小的分離可能；綜合來講，該問題中HSADE和NARSGA的結(jié)果在壓力分布上表現(xiàn)更好,相較于初始壓力分布也有明顯的改善.

圖7 3種算法的Pareto前緣對比Fig.7 Comparisons of Pareto fronts obtained by HSADE,NARSGA and DE

表6 3種算法得到的最優(yōu)翼型的性能參數(shù)和優(yōu)化量Table 6 Performance and optimization measurement of the optimal airfoils obtained by optimizers

4.2 二維超聲速噴管膨脹面優(yōu)化

本小節(jié)展示使用以上3種算法對一個用于超聲速推進系統(tǒng)的單斜面膨脹噴管(single expansion ramp nozzle,SERN)進行優(yōu)化并對比其效果.在這個問題中，初始種群的產(chǎn)生依靠算法自身的采樣方法進行.

圖8展示了3種算法在優(yōu)化過程中的世代距離收斂曲線，超體積比變化和HSADE的響應(yīng)面誤差變化，計算世代距離使用的參考解集仍取合并3種算法最優(yōu)解集中的非占優(yōu)部分.噴管構(gòu)型參考文獻[38].通過5個控制點使用NURBS(non-uniform rational B-spline)曲線擬合得到噴管上膨脹面的構(gòu)型，下?lián)醢?長度固定)的傾斜角可控，該構(gòu)型見圖9(a).噴管進口高度H0=0.5657m，出口高度H6=2.4m，長度Lt=3.294m，下?lián)醢彘L度Lp=0.7m.

圖8 優(yōu)化進程Fig.8 Optimization progress

設(shè)計變量為 5個控制點的歸一化高度Ci(i=1,2,3,4,5),下?lián)醢鍍A斜角α.Hi由Ci表示為Hi=H0+Ci(H6-H0)，各設(shè)計變量的范圍和起始個體A的參數(shù)由表7給出.

表7 噴管設(shè)計變量范圍Table 7 Variable ranges of the nozzle design

噴管來流馬赫數(shù)為 4.0，來流靜溫 2221.6K，來流靜壓 2516.6Pa，噴管入口馬赫數(shù) 1.011，噴管入口總溫和總壓分別為T0=1673.6K和p0=144645Pa，湍流模式為兩方程SST，離散格式為Roe格式，網(wǎng)格總量約為9萬.

優(yōu)化目標(biāo)分別為推力T和俯仰力矩M，均要求最大化.每種算法進化總代數(shù)皆為32代，種群大小皆為32，對于每種算法，均計算了約1000個有效個體，算法參數(shù)設(shè)置和4.1節(jié)完全相同.

圖9(b)～圖9(d)展示了3種算法優(yōu)化結(jié)果中3個比較有的代表性個體——推力最大個體(標(biāo)記為maxT)、力矩最大個體(標(biāo)記為maxM)和位于Pareto前緣的中點位置個體(標(biāo)記為middle)的膨脹面形狀.圖10為3種算法得到的Pareto前緣，表7總結(jié)了3種算法得出的代表性個體的相關(guān)性能參數(shù).

從優(yōu)化結(jié)果可以得到如下結(jié)論：

(1)從圖8看，在優(yōu)化過程中HSADE的世代距離和超體積比皆明顯領(lǐng)先于HSADE和DE，且HSADE的響應(yīng)面誤差盡管在優(yōu)化后期存在波動，仍然保持在可以接受的范圍內(nèi).

圖9 噴管構(gòu)型與3種算法獲得的代表性個體外形Fig.9 Nozzle concepts and representative optimal shapes obtained by HSADE,NARSGA and DE

圖10 3種算法的Pareto前緣Fig.10 Pareto fronts obtained by HSADE,NARSGA and DE

(2)從圖9看，從推力最大過渡到力矩最大，3種方法得出的噴管型面變化趨勢是一致的，同時HSADE和DE的結(jié)果中3種型面差異更大，可能表征了這兩種方法具有更強的多樣性保持能力.

(3)從圖10和表8看，HSADE得到的最優(yōu)解集最大推力和最大力矩均優(yōu)于其余兩種算法，且HSADE得到的Pareto前緣相較于另外兩種方法也更加靠前；NARSGA得到的Pareto前緣和DE互有交叉，但其前緣個體更充足，整體上略優(yōu)于DE.

表8 3種算法得到的代表個體的性能參數(shù)Table 8 Performance of representative optimal shapes obtained by HSADE,NARSGA and DE

5 結(jié)論

提出了一種混合優(yōu)化策略，基于差分進化和RBF響應(yīng)面形成了一種混合優(yōu)化算法.并引入雙敗淘汰選擇算子，參數(shù)自適應(yīng)等加以進一步改進.大量測試表明：差分進化與徑向基函數(shù)響應(yīng)面的組合給出了一種可靠且實用的優(yōu)化方法.能夠在保持全局尋優(yōu)能力的前提下在相近的總函數(shù)評估次數(shù)下達到更好的優(yōu)化結(jié)果，對昂貴優(yōu)化問題十分適用；與以往的混合優(yōu)化算法相比，本文的算法在全局尋優(yōu)能力和局部搜索能力均有明顯提高.

1 Cheung S,Aaronson P,Edwards T.CFD optimization of a theoretical minimum-drag body.Journal of Aircraft,2015,32(1):193-198

2 盧文書，王帥培，馬元春.基于CFD/CSD與Kriging插值模型的大展弦比復(fù)合材料機翼靜氣動彈性優(yōu)化設(shè)計.應(yīng)用力學(xué)學(xué)報，2015，32(4)：581-585(Lu Wenshu,Wang Shuaipei,Ma Yuanchun. Staticaeroelasticoptimizationofahigh-aspect-ratiocompositewing based on CFD/CSD and Kriging model.Chinese Journal of Applied Mechanics,2015,32(4):581-585(in Chinese))

3 Cazacu R,Grama L.Steel truss optimization using genetic algorithms and FEA.Procedia Technology,2014,12:339-346

4 Jones DR,Schonlau M,Welch WJ.Efficient global optimization of expensive black-box functions.Journal of Global Optimization, 1998,13(4):455-492

5 Su GS.Gaussian process assisted di ff erential evolution algorithm for computationally expensive optimization problems//Bilof R ed. Pacific-AsiWorkshop on Computational Intelligence and Industrial Application，PACIIA’08，Wuhan，2008.Los Alamitos：IEEE, 2008.272-276

6 Tabatabaei SME,Kadkhodaie-Ilkhchi A,Hosseini Z,et al.A hybrid stochastic-gradient optimization to estimating total organic carbon from petrophysical data:a case study from the Ahwaz Oilfield SW Iran.Journal of Petroleum Science&Engineering,2015,127(1): 35-43

7 孫美建,詹浩.Kriging模型在機翼氣動外形優(yōu)化中的應(yīng)用.空氣動力學(xué)學(xué)報,2011,29(6):759-764(Sun Meijian,Zhan Hao.Application of Kriging surrogate model for aerodynamic shape optimization of wing.Acta Aerodynamica Sinica,2011,29(6):759-764(in Chinese))

8 Huang D,Allen TT,Notz WI,et al.Global optimization of stochastic black-box systems via sequential kriging meta-models.Journalof Global Optimization,2006,34(3):441-466

9 Singh HK,Isaacs A,Ray T.A hybrid surrogate based algorithm (HSBA)to solve computationally expensive optimization problems.Evolutionary Computation,2014:1069-1075

10 Elsayed SM,Ray T,Sarker RA.A surrogate-assisted di ff erential evolution algorithm with dynamic parameters selection for solving expensive optimization problems.Evolutionary Computation,2014: 1062-1068

11 Liu B,Zhang QF,Gielen GGE.A Gaussian process surrogate model assisted evolutionary algorithm for medium scale expensive optimization problems.IEEE Transactions on Evolutionary Computation,2014,18(2):180-192

12 Antunes AP,Azevedo JLF.Studies in aerodynamic optimization basedongeneticalgorithms.JournalofAircraft,2014,51(3):1002-1012

13 Nam T,Chakraborty I,Gross J,et al.Multidisciplinary design optimization of a truss-braced wing concept//14th AIAA Aviation Technology,Integration,and Operations Conference，Atlanta，2014.Reston：AIAA,2014

14 Gibertini G.Aerodynamic shape optimisation of a proprotor and its validation by means of CFD and experiments.Aeronautical Journal, 2015,119(1120):1223-1251

15 Han ZH,Zimmerman R,Grtz S.Alternative cokriging method for variable-fidelitsurrogate modeling.AIAA Journal,2012,50(5): 1205-1210

17 Zingg DW,Nemec M,Pulliam TH.A comparative evaluation of genetic and gradient-based algorithms applied to aerodynamic optimization.European Journal of Computational Mechanics,2008, 17(1-2):103-126

18 Carrier G,Destarac D,Dumont A,et al.Gradient-based aerodynamic optimization with the elsA software//52nd Aerospace Sciences Meeting.2014,10:6.2014-0568

19 白俊強,王波,孫智偉等.基于松散式代理模型管理框架的亞音速機翼優(yōu)化設(shè)計方法研究.西北工業(yè)大學(xué)學(xué)報,2011,29(4):515-519(Bai Junqiang,Wang Bo,Sun Zhiwei,et al.Developing optimization design of subsonic wing with loose type of agent model.Journal of Northwestern Polytechnical University,2011,29(4): 515-519(in Chinese))

20 Kim HJ,Liou MS.Aerodynamic optimization using a hybrid moga-local search method//51st AIAA/ASME/ASCE/AHS/ASC Structures,Structural Dynamics,and Materials Conference 18th AIAA/ASME/AHS Adaptive Structures Conference 12th,2010: 2911

21倪昂修,張宇飛,陳海昕.NSGA-Ⅱ算法的改進及其在多段翼型縫道參數(shù)優(yōu)化中的應(yīng)用.空氣動力學(xué)學(xué)報,2014,32(2):252-257(Ni Angxiu,Zhang Yufei,Chen Haixin.An Improvement to NSGA-II algorithm and its application in optimization design of multi-element airfoil.Acta Aerodynamica Sinica,2014,32(2):252-257(in Chinese))

22 Deb K,Pratap A,Agarwal S,et al.A fast and elitist multiobjective genetic algorithm:NSGA-II.IEEE Transactions on Evolutionary Computation,2002,6(2):182-197

23 公茂果,焦李成,楊咚咚等.進化多目標(biāo)優(yōu)化算法研究.軟件學(xué)報, 2009,20(2):271-289(Gong Maoguo,Jiao Licheng,Yang Dongdong,et al.Research on evolutionary multi-objective optimization algorithms.Journal of Software,2009,20(20):271-289(in Chinese))

24 Zitzler E,Laumanns M,Thiele L.SPEA2:improving the strength pareto evolutionary algorithm.Eurogen,2001,3242(103):95-100

25 Storn R,Price K.Di ff erential evolution:a simple and efficient adaptive scheme for global optimization over continuous spaces.Journal of Global Optimization,1995,23(4):341-359

26 Park J,Sandberg IW.Universal approximation using radial-basisfunction networks.Neural Computation,1991,3(2):246-257

27 穆雪峰,姚衛(wèi)星,余雄慶等.多學(xué)科設(shè)計優(yōu)化中常用代理模型的研究.計算力學(xué)學(xué)報,2005,22(5):608-612(Mu Xuefeng,Yao Weixing,Yu Xiongqing,et al.A survey of surrogate models used in MDO.Chinese Journal of Computational Mechanics,2005,22(5): 608-612(in Chinese))

28 Wang L,Huang FZ.Parameter analysis based on stochastic model for di ff erential evolution algorithm.Applied Mathematics&Computation,2010,217(7):3263-3273

29 Brest J,Greiner S,Boskovic B,et al.Self-adapting control parameters in di ff erential evolution:a comparative study on numerical benchmark problems.IEEE Transactions on Evolutionary Computation,2007,10(6):646-657

30 Rippa S.An Algorithm for selecting a good value for the parameter c in radial basis function interpolation.Advances in Computational Mathematics,1999,11(2):193-210

31 吳亮紅.動態(tài)差分進化算法及其應(yīng)用.北京:科學(xué)出版社,2014: 190-220(Wu Lianghong.Dynamic Di ff erential Evolution and its Application.Beijing:Science Press,2014:190-220(in Chinese))

32 Nebro AJ,Luna F,Alba E,et al.Abyss:adapting scatter search to multiobjective optimization.IEEE Transactions on Evolutionary Computation,2008,12(4):439-457

33 Deb K,Thiele L,Laumanns M,et al.Scalable test problems for evolutionary multiobjective optimization//Abraham A,Jain L,Goldberg Reds.EvolutionaryMultiobjectiveOptimization.London:Springer London,2005.105-145

34 Coello CAC,Lechuga MS.MOPSO:a proposal for multiple objective particle swarm optimization//Proceedings of the 2002 Congress on Evolutionary Computation，CEC’02，Honolulu,Hawaii，2002. Los Alamitos：IEEE,2002.2:1051-1056

35 Kulfan BM,Bussoletti JE.Fundamental parametric geometry representations for aircraft component shapes//11th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference，Portsmouth, Virginia，2006.Reston：AIAA，2006.6948

36 Chen HX,Fu S,Li FW.Navier-Stokes simulations for transport aircraft wing/body high-lift configurationsJournal of Aircraft,2003, 40(5):883-890

37 Zhang YF,Chen HX,Fu S.Improvement to patched grid technique with high-order conservative remapping method.Journal of Aircraft,2012,48(3):884-893

38 陳兵,徐旭,蔡國飆.二維超燃沖壓發(fā)動機尾噴管優(yōu)化設(shè)計.推進技術(shù),2002,23(5):433-437(Chen Bing,Xu Xu,Cai Guobiao.Optimization design of two dimensional scramjet nozzle based on N-S equations.Journal of Propulsion Technology,2002,23(5):433-437(in Chinese))

HYBRID OPTIMIZATION ALGORITHM BASED ON DIFFERENTIAL EVOLUTION AND RBF RESPONSE SURFACE1)

Deng Kaiwen Chen Haixin2)
(School of Aerospace Engineering，Tsinghua University，Beijing100084，China)

A new hybrid optimization algorithm HSADE(hybrid self-adaptive di ff erential evolution)based on di ff erential evolutionandradialbasisfunctionresponsesurfacewasproposedaimingataerodynamicoptimizationproblems.Through combing the merits of response surface method’s fast local searching ability and di ff erential evolution’s powerful global searching ability,the overall local and global search efficiency of HSADE were simultaneously enhanced.Several improvements were made on certain logics and strategies embedded in the processes of each sub-algorithm by proposing and utilizing strategies such as selection strategy based on double elimination and self-adaptive parameters.Having applied HSADE and several other typical optimization algorithms—NSGA-II,MOPSO and multi-objective di ff erential evolution to several benchmark functions,the results indicated HSADE was superior to other algorithms in most of the cases regarding local search ability represented by generation distance and global search ability symbolled by hyper volume ratio,which validated the e ff ectiveness of above improvements.Applying HSADE along with basic DE and NARSGA to an airfoil optimization problem and a hypersonic nozzle expansion surface optimization problem,the results showed HSADE was able to obtain airfoils with extra 0.5 count drag reduction and nozzles with better performance than othertwo algorithms under approximately 1000 function evaluations,which indicated high engineering application potential of HSADE.

di ff erential evolution,response surfaces,computational flui dynamics,multiobjective optimization,optimization design

V211.3

A

10.6052/0459-1879-16-285

2016–10–17收稿，2017–01–16錄用,2017–01–20網(wǎng)絡(luò)版發(fā)表.

1)中航工業(yè)產(chǎn)學(xué)研專項(cxy2014QH14)和清華大學(xué)自主科研計劃(2015THZ0)資助項目.

2)陳海昕，教授，主要研究方向：空氣動力學(xué)，計算流體力學(xué).E-mail：chenhaixin@tsinghua.edu.cn

鄧凱文,陳海昕.基于差分進化和RBF響應(yīng)面的混合優(yōu)化算法.力學(xué)學(xué)報,2017,49(2):441-455

Deng Kaiwen,Chen Haixin.Hybrid optimization algorithm based on di ff erential evolution and RBF response surface.Chinese Journal of Theoretical and Applied Mechanics,2017,49(2):441-455