基于代理模型的大型民機機翼氣動優(yōu)化設(shè)計

韓忠華，張瑜，許晨舟，王凱，吳猛猛，朱震，宋文萍

西北工業(yè)大學(xué) 航空學(xué)院，西安 710072

大型民機的研究是關(guān)系中國國民經(jīng)濟發(fā)展和科學(xué)技術(shù)進(jìn)步的一項重大課題。未來大型民機設(shè)計將朝著更高效、更安全、更經(jīng)濟、更舒適、更環(huán)保的方向發(fā)展，這對其氣動性能提出了更高的要求。近年來，隨著計算流體力學(xué)(CFD)技術(shù)的飛速發(fā)展，基于CFD的氣動分析與優(yōu)化設(shè)計，在提高大型民機氣動與綜合性能方面正發(fā)揮越來越重要的作用[1-3]。在當(dāng)今美國航空航天領(lǐng)域，CFD約占?xì)鈩釉O(shè)計工作量的70%[4]。開展基于CFD的氣動優(yōu)化設(shè)計方法研究，對于提升大型民機氣動與綜合性能具有重要意義。國內(nèi)外在研發(fā)C919、B787、A350和A380等客機過程中，均大量采用了基于CFD的分析與設(shè)計技術(shù)。

增升減阻是大型飛機氣動設(shè)計的核心問題。研究表明，大型民機阻力系數(shù)每下降 1 count (1 count=0.000 1)，載重將可能增加 7.56%[4]；文獻(xiàn)[5]指出，對于某典型雙發(fā)動機的大型民機，起飛升阻比提高1%，意味著可增加14個乘客；著陸最大升力系數(shù)增加1%，則意味著可增加22個乘客。機翼作為提供升力的主要部件，其阻力占飛機總阻力的一半以上。因此，開展機翼氣動優(yōu)化設(shè)計，對于提升大型民機整體性能有著至關(guān)重要的作用。

在氣動設(shè)計領(lǐng)域，傳統(tǒng)的依賴經(jīng)驗、以試湊法(Cut-and-Try)為主的氣動設(shè)計方法，計算成本高、設(shè)計周期長，且無法保證獲得最優(yōu)的氣動外形。因此，國內(nèi)外學(xué)者對基于CFD的氣動優(yōu)化設(shè)計方法開展了大量研究?，F(xiàn)有的優(yōu)化設(shè)計方法主要分為3類：① 遺傳算法、粒子群算法等啟發(fā)式優(yōu)化算法[6]。這類方法具有較好的全局性，但優(yōu)化過程可能需要成千上萬次(甚至更多)地調(diào)用CFD分析，計算量巨大[7]。更嚴(yán)重的是，針對機翼、翼身組合體等復(fù)雜幾何外形，隨著設(shè)計變量數(shù)(維數(shù))的增加，導(dǎo)致計算量劇增，出現(xiàn)了所謂的“維數(shù)災(zāi)難”[8]，從而大大限制了其在復(fù)雜外形氣動優(yōu)化設(shè)計中的應(yīng)用。② 梯度優(yōu)化算法[9]。結(jié)合Jameson教授發(fā)展的Adjoint方法[10-11]，能夠顯著提高針對大規(guī)模變量優(yōu)化問題的效率，但該方法容易陷入局部最優(yōu)，且設(shè)計結(jié)果與初始外形密切相關(guān)。研究表明，即使采用“多起點”的梯度優(yōu)化策略，其優(yōu)化效果也可能難以與全局優(yōu)化算法媲美[12]。③ 代理優(yōu)化算法[13]。這是近年來興起的一種基于代理模型的高效優(yōu)化算法。所謂代理優(yōu)化算法，是指通過構(gòu)造目標(biāo)和約束函數(shù)的代理模型，再求解由“優(yōu)化加點準(zhǔn)則”[14]定義的子優(yōu)化問題，指導(dǎo)新的樣本點自適應(yīng)循環(huán)加入，不斷更新代理模型，直到樣本點序列收斂于全局最優(yōu)解的一種優(yōu)化算法。代理優(yōu)化算法具有全局優(yōu)化能力，同時優(yōu)化效率比遺傳算法等高1～2個量級，是一種有望解決航空航天復(fù)雜工程設(shè)計領(lǐng)域計算量過大和計算時間過長問題的高效全局優(yōu)化方法[15-22]。

本文研究將代理優(yōu)化算法應(yīng)用于大型民機超臨界機翼氣動優(yōu)化設(shè)計，通過梳理設(shè)計經(jīng)驗，探討了大型民機機翼氣動優(yōu)化設(shè)計思想及原則，并提出了一套面向工程應(yīng)用的多輪次氣動優(yōu)化設(shè)計方法，旨在為中國大型民機自主研發(fā)提供一定的技術(shù)支持。

1 大型民機機翼氣動優(yōu)化設(shè)計的基本準(zhǔn)則和要點

目前，國際上主流的先進(jìn)大型民機均采用超臨界機翼。超臨界機翼的基本設(shè)計思想是通過減小前緣負(fù)壓峰，削弱前緣氣流加速，從而降低最大當(dāng)?shù)伛R赫數(shù)，同時盡可能維持上表面的超聲速區(qū)和平緩的壓力分布，來獲得較小的激波阻力和較高的阻力發(fā)散馬赫數(shù)。此外，通過下表面后加載來彌補前半段彎度小的不足，以保證足夠的升力。

壓力分布形態(tài)是評判氣動性能的重要依據(jù)。對于無限翼展后掠機翼，可以通過后掠角轉(zhuǎn)換為二維翼型設(shè)計來獲得理想的三維機翼壓力分布。然而對于后掠角較大的有限翼展超臨界機翼設(shè)計，其二維/三維之間的轉(zhuǎn)換關(guān)系不是很明確，翼型的氣動特性很難直接反映在機翼上，需要直接開展機翼三維幾何外形的氣動設(shè)計。而且對于超臨界機翼，其壓力分布吸力峰值、負(fù)壓平臺的高度長度及壓力恢復(fù)速率等因素均相互影響，與幾何外形變化呈現(xiàn)很強的非線性關(guān)系，規(guī)律不易掌握[23]。因此，基于CFD的氣動優(yōu)化設(shè)計方法仍然面臨很大的困難與挑戰(zhàn)。對于基于數(shù)值優(yōu)化算法的超臨界機翼設(shè)計而言，要點如下：

1) 最小激波阻力設(shè)計。在高馬赫數(shù)條件下，盡可能減小激波阻力是優(yōu)化設(shè)計的重要目標(biāo)。在機翼平面外形固定的情況下，主要通過剖面設(shè)計來實現(xiàn)較小的激波阻力。由于超臨界機翼氣動特性對幾何外形十分敏感，各剖面氣動特性又相互影響，必須開展以機翼整體氣動性能為目標(biāo)的多剖面氣動優(yōu)化設(shè)計。此外，單點優(yōu)化設(shè)計的機翼，在偏離設(shè)計點時激波阻力可能會急劇增大，需要進(jìn)行多點設(shè)計。

2) 減小誘導(dǎo)阻力設(shè)計。誘導(dǎo)阻力也是機翼阻力的重要組成部分。優(yōu)化過程中需要在兼顧機翼結(jié)構(gòu)特性的前提下，通過優(yōu)化幾何扭轉(zhuǎn)和氣動扭轉(zhuǎn)等參數(shù)，改變展向載荷分布，使誘導(dǎo)阻力盡可能小。

3) 等壓線平行原則。研究發(fā)現(xiàn)，當(dāng)機翼表面等壓線形態(tài)整齊平行時，機翼會具有較好的氣動特性。在設(shè)計過程中，可以通過合理配置與設(shè)計各個站位處的翼型來實現(xiàn)等壓線平行。

4) 幾何兼容性原則。出于機翼加工的角度考慮，機翼剖面沿展向應(yīng)當(dāng)具有較好的幾何兼容性，這同樣也有利于實現(xiàn)各剖面壓力分布的形態(tài)相似，實現(xiàn)等壓線平行。

5) 一體化優(yōu)化設(shè)計。針對大后掠角和幾何扭轉(zhuǎn)變化較大的機翼，三維效應(yīng)顯著，不能采用傳統(tǒng)的先設(shè)計剖面翼型、再配置到機翼進(jìn)行三維設(shè)計的方法，而是需要直接開展三維流動環(huán)境下的機翼多剖面一體化設(shè)計。另外，對于大型民機，機身、發(fā)動機短艙及掛架對機翼的氣動性能有著十分明顯的影響。如果僅對單獨機翼進(jìn)行設(shè)計，獲得氣動特性提升，在考慮機身和發(fā)動機影響后，氣動性能可能顯著下降。因此，必須開展考慮機身、發(fā)動機短艙和掛架影響的一體化機翼氣動優(yōu)化設(shè)計。

6) 多點多目標(biāo)優(yōu)化設(shè)計。除了保證設(shè)計外形在巡航狀態(tài)下均具有較優(yōu)的氣動性能，還需要保證起飛和著陸等其他設(shè)計點的氣動優(yōu)化性能，以及阻力發(fā)散馬赫數(shù)、最大馬赫數(shù)及抖振邊界等綜合氣動性能。在巡航馬赫數(shù)和巡航升力系數(shù)附近，要求氣動性能應(yīng)具有較好的穩(wěn)定性，也需要通過多點多目標(biāo)優(yōu)化設(shè)計來實現(xiàn)。

大型民機機翼的這些設(shè)計準(zhǔn)則對氣動優(yōu)化設(shè)計方法提出了很高的要求。所采用的氣動優(yōu)化設(shè)計方法，必須能夠高效處理具有多個設(shè)計目標(biāo)、多個約束、大量設(shè)計變量和高度非線性的優(yōu)化問題。由于代理優(yōu)化算法的高效性和全局優(yōu)化特性，本文采用該方法開展大型民機機翼氣動優(yōu)化設(shè)計的研究。

2 基于代理模型的氣動優(yōu)化設(shè)計方法

針對大型民機超臨界機翼氣動優(yōu)化設(shè)計，本文基于代理優(yōu)化算法，提出了一種面向工程應(yīng)用的多輪次氣動優(yōu)化設(shè)計方法。將優(yōu)化設(shè)計與人工修型相結(jié)合，引入設(shè)計者的經(jīng)驗，提高了實用性。

2.1 代理優(yōu)化算法

代理優(yōu)化的基本原理如圖1所示，其過程可描述為：① 對設(shè)計空間進(jìn)行試驗設(shè)計，抽取樣本點并運行精確數(shù)值模擬分析獲得響應(yīng)值，建立初始代理模型；② 基于代理模型，按照一定的優(yōu)化加點準(zhǔn)則(如MSP、EI、PI、MSE、LCB[13]等)，采用傳統(tǒng)優(yōu)化算法求解相應(yīng)的子優(yōu)化問題，以很小的計算代價，對最優(yōu)解進(jìn)行預(yù)測；③ 對這些預(yù)測的最優(yōu)解再次運行數(shù)值分析，并將結(jié)果作為新的樣本數(shù)據(jù)添加到現(xiàn)有數(shù)據(jù)集中，不斷更新代理模型，直到所產(chǎn)生的樣本點序列收斂于局部或全局最優(yōu)解。

圖1 代理優(yōu)化算法框架Fig.1 Framework of surrogate-based optimization algorithm

2.2 面向工程應(yīng)用的多輪次氣動優(yōu)化設(shè)計方法

代理優(yōu)化算法應(yīng)用于工程實踐中，仍然有一個關(guān)鍵問題需要解決——如何確定設(shè)計空間。對于一個最優(yōu)解未知的優(yōu)化問題，理論上講，設(shè)計空間應(yīng)當(dāng)盡可能大，以包含優(yōu)化問題的最優(yōu)解。然而，過大的設(shè)計空間會導(dǎo)致優(yōu)化變得十分困難。一方面，在一個較大的設(shè)計空間內(nèi)建立足夠精度的代理模型需要大量的樣本點，因此計算量更大；另一方面，在設(shè)計空間中抽樣時，一組設(shè)計變量可能對應(yīng)非常“奇異”甚至完全不可行的氣動外形，導(dǎo)致無法獲得其響應(yīng)值，影響代理模型的建立。同時，隨著優(yōu)化的進(jìn)行，設(shè)計目標(biāo)和約束需要不斷地調(diào)整。因此，開展多輪優(yōu)化具有很強的工程實用性。

本文基于代理優(yōu)化算法，提出了一種面向工程應(yīng)用的多輪次氣動優(yōu)化設(shè)計方法，如圖2所示。該方法的基本原理是在代理優(yōu)化算法基礎(chǔ)上，執(zhí)行多輪次優(yōu)化：① 在下一輪次優(yōu)化時，采用上一輪次優(yōu)化的最優(yōu)結(jié)果作為基準(zhǔn)外形，重新設(shè)定設(shè)計空間及設(shè)計目標(biāo)等參數(shù)；② 每一輪優(yōu)化時，重新采用試驗設(shè)計進(jìn)行抽樣，并將上一輪的最優(yōu)結(jié)果作為初始樣本點之一；③ 為了加速優(yōu)化收斂過程，節(jié)省計算量，將上一輪次樣本點中靠近當(dāng)前最優(yōu)解的部分樣本點也作為初始樣本點加入到下一輪優(yōu)化中。

圖2 多輪次氣動優(yōu)化設(shè)計方法示意圖Fig.2 Schematic of multi-round aerodynamic optimization design method

多輪次氣動優(yōu)化設(shè)計方法，不僅能夠有效避免設(shè)計空間過大使得優(yōu)化效率降低的問題，也可以盡可能避免優(yōu)化過程中產(chǎn)生“奇異”的幾何外形導(dǎo)致CFD計算難以收斂的問題。此外，該方法還有2個方面的作用：① 避免矩陣剛性的問題。優(yōu)化設(shè)計在一個較大的設(shè)計空間內(nèi)進(jìn)行，如果樣本分布不均，會導(dǎo)致代理模型相關(guān)矩陣條件數(shù)大、代理模型近似精度降低的問題。代理模型保持了大范圍的分辨率，在所關(guān)心的最優(yōu)解附近的分辨率會顯著不足。② 避免樣本多樣性不足的問題。隨著優(yōu)化進(jìn)行，樣本點會逐漸聚集于當(dāng)前最優(yōu)解的附近而出現(xiàn)“早熟”。通過在每一輪次優(yōu)化中加入新的試驗設(shè)計樣本點，能夠增加設(shè)計空間中樣本點的多樣性，從而有利于建立更精確的代理模型，尋找全局最優(yōu)解。

2.3 數(shù)值優(yōu)化與人工修型相結(jié)合

大型民機機翼設(shè)計需要滿足等壓線平行和幾何兼容性準(zhǔn)則。盡管可以通過幾何約束來對氣動外形進(jìn)行控制，但對于大型民機機翼需要添加成千上萬的幾何約束，僅僅依靠氣動優(yōu)化設(shè)計方法難以獲得滿足要求的氣動外形。而僅依賴人工經(jīng)驗的“試湊法”也難以獲得最優(yōu)的設(shè)計方案。因此，有必要開展數(shù)值優(yōu)化與人工修型相結(jié)合的設(shè)計方法。一方面，氣動優(yōu)化設(shè)計的結(jié)果雖然可能不完全滿足設(shè)計準(zhǔn)則，但可以為設(shè)計者提供指導(dǎo)，給出正確的設(shè)計方向。另一方面，人工修型能夠在多輪次優(yōu)化設(shè)計過程中，不斷調(diào)整設(shè)計目標(biāo)與約束，并通過研究對比剖面壓力分布為優(yōu)化設(shè)計提供合理的設(shè)計空間。總之，數(shù)值優(yōu)化與人工修型相結(jié)合，可大大增強優(yōu)化設(shè)計的實用性。

3 氣動優(yōu)化設(shè)計方法驗證

3.1 Branin-Hoo函數(shù)優(yōu)化算例

Branin-Hoo函數(shù)是帶約束的全局優(yōu)化算法的經(jīng)典測試算例之一。優(yōu)化問題數(shù)學(xué)模型為

(1)

該優(yōu)化問題具有3個局部最優(yōu)解，如圖3所示，虛線以上為滿足約束的可行域，“×”代表滿足約束的全局最優(yōu)解。分別采用基于Kriging的代理優(yōu)化算法和SQP(Sequential Quadratic Programming)梯度優(yōu)化算法開展優(yōu)化，并進(jìn)行比較,優(yōu)化收斂曲線如圖4所示。代理優(yōu)化算法分別采用10組不同初始樣本點執(zhí)行優(yōu)化，均收斂于全局最優(yōu)解。而對于梯度優(yōu)化，從不同起點出發(fā)可能收斂于不同的局部最優(yōu)解。

圖3 二維Branin-Hoo函數(shù)示意圖Fig.3 Schematic of 2-D Branin-Hoo function

圖4 Branin-Hoo函數(shù)優(yōu)化收斂歷程Fig.4 Convergence history of Branin-Hoo function optimization

3.2 RAE2822翼型減阻優(yōu)化設(shè)計標(biāo)準(zhǔn)算例

RAE2822翼型減阻優(yōu)化設(shè)計算例為AIAA氣動設(shè)計優(yōu)化討論組(ADODG)定義的標(biāo)準(zhǔn)算例[24](Case 2)。翼型的設(shè)計狀態(tài)為：馬赫數(shù)Ma=0.734,雷諾數(shù)Re=6.5×106，其優(yōu)化數(shù)學(xué)模型為

(2)

式中：Cl、Cd和Cm分別為翼型的升力系數(shù)、阻力系數(shù)和力矩系數(shù)；A和A0分別為優(yōu)化翼型和初始翼型的面積。

采用CST(Class function/Shape function Transformation)[15]參數(shù)化方法，共18個設(shè)計變量。由于該算例升力系數(shù)約束為主動約束，為了避免處理等式約束，本文將其簡化為升力系數(shù)不減，并固定迎角為2.879 5°。

圖5給出了代理優(yōu)化過程中目標(biāo)函數(shù)的收斂歷程。共進(jìn)行了3輪優(yōu)化，通過第2輪和第3輪優(yōu)化，阻力系數(shù)分別進(jìn)一步下降了約0.5 counts和0.2 counts，說明優(yōu)化已經(jīng)充分收斂。表1給出了優(yōu)化翼型的氣動性能及面積等參數(shù)與基準(zhǔn)翼型的比較。在所有約束都得到精確滿足的前提下，阻力大大降低。與文獻(xiàn)[24]中梯度優(yōu)化結(jié)果的比較說明，代理優(yōu)化結(jié)果明顯優(yōu)于梯度優(yōu)化結(jié)果。圖6給出了優(yōu)化前后翼型的表面壓力系數(shù)Cp分布對比。可以看出，優(yōu)化后上表面的強激波被削弱成2道非常微弱的激波，優(yōu)化效果顯著。

圖5 RAE2822翼型多輪次減阻優(yōu)化目標(biāo)函數(shù)收斂歷程Fig.5 Convergence history of objective function for multi-round drag minimization of RAE2822 airfoil test case

表1 RAE2822翼型優(yōu)化前后的氣動性能及面積對比

圖6 RAE2822翼型優(yōu)化前后的表面壓力分布對比(Ma=0.734,Re=6.5×106,Cl=0.824)Fig.6 Comparison of pressure distributions for baseline RAE2822 and optimum airfoil (Ma=0.734,Re=6.5×106,Cl=0.824)

3.3 寬體客機翼身組合體機翼氣動優(yōu)化設(shè)計算例

寬體客機翼身組合體機翼氣動優(yōu)化設(shè)計算例為高亞聲速運輸機翼身組合體超臨界機翼優(yōu)化設(shè)計(該算例來自文獻(xiàn)[19])。機翼的設(shè)計狀態(tài)為：Ma=0.83,Re=4.34×107，其優(yōu)化數(shù)學(xué)模型為

分別采用代理優(yōu)化算法、梯度優(yōu)化算法和基于子空間的單純形搜索算法(Subplex)進(jìn)行優(yōu)化。其中代理優(yōu)化算法通過試驗設(shè)計產(chǎn)生80個樣本點來構(gòu)造初始代理模型，并采用改善期望(EI)優(yōu)化加點準(zhǔn)則進(jìn)行加點。梯度優(yōu)化算法則分別采用基準(zhǔn)外形、試驗設(shè)計中氣動性能最優(yōu)的外形作為初始外形進(jìn)行優(yōu)化。表2給出了3種優(yōu)化方法的結(jié)果和計算量對比，圖8給出了3種優(yōu)化方法的阻力系數(shù)收斂曲線?？梢钥闯?，采用代理優(yōu)化的結(jié)果最佳，而梯度優(yōu)化采用不同起點的優(yōu)化結(jié)果不同，且均收斂于局部最優(yōu)解。如果在代理優(yōu)化結(jié)果的基礎(chǔ)上進(jìn)一步開展梯度優(yōu)化，優(yōu)化外形的阻力系數(shù)幾乎不可能再下降，表明代理優(yōu)化已完全收斂。相比于Subplex優(yōu)化算法，代理優(yōu)化算法的效率提高了一倍以上。針對本算例40個設(shè)計變量，代理優(yōu)化在300個樣本點以內(nèi)基本收斂。對于梯度優(yōu)化，從基準(zhǔn)外形出發(fā)的優(yōu)化結(jié)果陷入局部最優(yōu)(橙色線)；如果從試驗設(shè)計樣本中的最優(yōu)外形出發(fā)，可以獲得更優(yōu)的結(jié)果，但阻力仍然比代理優(yōu)化結(jié)果更大。該算例針對一個接近實際飛機型號的復(fù)雜外形進(jìn)行優(yōu)化設(shè)計，首次展示了代理優(yōu)化算法在全局氣動優(yōu)化設(shè)計方面的巨大應(yīng)用潛力與研究價值。

圖7 翼身組合體表面網(wǎng)格和FFD控制框示意圖Fig.7 Schematic of surface grid of wing-body configuration and FFD control box

表2 某寬體客機翼身組合體外形采用不同優(yōu)化方法的優(yōu)化結(jié)果對比

圖8 某寬體客機翼身組合體外形減阻優(yōu)化收斂歷程Fig.8 Convergence history of objective function for drag minimization of wing-body configuration of a wide-body transport

4 大型寬體客機超臨界機翼氣動優(yōu)化設(shè)計

4.1 優(yōu)化設(shè)計流程及工具

4.1.1 SurroOpt

SurroOpt[13-22,25-28]是筆者所在課題組開發(fā)的基于代理模型的通用優(yōu)化程序，其可以求解任意單目標(biāo)、加權(quán)系數(shù)的多目標(biāo)、Pareto多目標(biāo)的無約束、多約束優(yōu)化問題。該程序包括拉丁超立方(LHS)、均勻設(shè)計(UD)、蒙特卡羅(MC)抽樣等試驗設(shè)計方法，二次響應(yīng)面(PRSM)、Kriging模型、梯度增強Kriging(GEK)模型[18]、分層Kriging(HK)模型[28]、徑向基函數(shù)(RBF)等多種代理模型方法，5種主要的優(yōu)化加點方法及其約束處理方法，Hooke-Jeeves模式搜索、擬牛頓梯度優(yōu)化、SQP法、單/多目標(biāo)遺傳算法等多種成熟優(yōu)化算法作為子優(yōu)化算法。該程序可以運行于Windows及Linux操作系統(tǒng)，可以進(jìn)行單核計算，也可進(jìn)行多核并行計算。該程序提供了一種通用的接口程序，用戶可以借助該接口程序搭建滿足自己需要的任意優(yōu)化設(shè)計問題。大量的標(biāo)準(zhǔn)測試算例和實際工程應(yīng)用的優(yōu)化算例表明，對于小規(guī)模設(shè)計變量的局部優(yōu)化問題，該程序的優(yōu)化效率(以目標(biāo)函數(shù)計算次數(shù)為評價標(biāo)準(zhǔn))與梯度優(yōu)化算法相當(dāng)；對于全局優(yōu)化問題，優(yōu)化效率比遺傳算法等其他全局優(yōu)化算法高1～2個數(shù)量級。此外，實踐表明，該程序具有良好的約束處理能力，具體約束處理方法見文獻(xiàn)[13，17]。

目前，SurroOpt可以應(yīng)用于基于高可信度CFD的氣動優(yōu)化設(shè)計，也可用于結(jié)構(gòu)優(yōu)化、多學(xué)科優(yōu)化及其他工程優(yōu)化設(shè)計問題。針對氣動優(yōu)化設(shè)計(見圖9)，該程序可首先通過試驗設(shè)計獲得設(shè)計空間內(nèi)一系列樣本點，使用參數(shù)化方法獲得樣本點所對應(yīng)的氣動外形，再采用Euler方程或雷諾平均Navier-Stokes(RANS)方程等CFD方法進(jìn)行繞流數(shù)值模擬得到每一個樣本點(即氣動外形)的氣動力特性。其次，根據(jù)這些樣本點及其氣動特性建立目標(biāo)函數(shù)和約束函數(shù)代理模型。根據(jù)一定的加點準(zhǔn)則，使用優(yōu)化算法在設(shè)計空間內(nèi)進(jìn)行尋優(yōu)，確定增加樣本點的位置。再采用Euler方程或RANS方程進(jìn)行數(shù)值模擬獲得新增樣本點的氣動數(shù)據(jù)。將新的樣本信息加入到樣本數(shù)據(jù)集中，重新建立代理模型，重復(fù)上述過程直至達(dá)到收斂標(biāo)準(zhǔn)。

圖9 通用代理優(yōu)化程序SurroOpt的優(yōu)化流程Fig.9 Framework of a generic SBO-type optimizer “SurroOpt”

4.1.2 幾何外形參數(shù)化

本文采用CST方法[29]和樣條函數(shù)分別對機翼剖面翼型和扭轉(zhuǎn)角分布進(jìn)行參數(shù)化。CST方法采用類函數(shù)和形函數(shù)相結(jié)合的方式描述幾何外形，其所需的控制參數(shù)較少，參數(shù)化精度高，并能夠自動保持翼型的光滑性。本文固定機翼平面形狀，對多個站位的剖面形狀和扭轉(zhuǎn)分布進(jìn)行優(yōu)化設(shè)計。

4.1.3 全速勢方程流動求解器

優(yōu)化過程中采用的是一種基于“全速勢方程”與黏性修正相結(jié)合的亞聲速和跨聲速繞流的快速計算程序“WBCAero”。該程序可以用于計算由機身+機翼+發(fā)動機短艙+尾翼組成的全機構(gòu)型的氣動特性，并考慮黏性和附面層影響。該程序能夠自動生成計算網(wǎng)格，針對全機構(gòu)型進(jìn)行一次完整計算僅需大約15 s(采用普通個人計算機)。本文在優(yōu)化過程中，未考慮機身黏性，所以阻力系數(shù)值偏小。優(yōu)化完成后，再加入機身黏性進(jìn)行評估。

4.2 優(yōu)化設(shè)計結(jié)果

為了驗證本文多輪次代理優(yōu)化算法的有效性，對某遠(yuǎn)程寬體客機超臨界機翼的高速巡航狀態(tài)進(jìn)行了多點設(shè)計(暫未考慮低速性能)。

4.2.1 寬體客機基準(zhǔn)機翼設(shè)計

基準(zhǔn)機翼是在CRM(Common Research Model)機翼[30]基礎(chǔ)上設(shè)計得到的，設(shè)計狀態(tài)為：Ma=0.85,Re=4.0×107,CL=0.5。機翼按展向位置等比例截取11個剖面的翼型，布置于基準(zhǔn)機翼相應(yīng)的展向位置，并采用CRM機翼的扭轉(zhuǎn)分布。該基準(zhǔn)機翼前緣后掠角約為36.5°，展弦比約為9。

為了保證更好的幾何兼容性和氣動兼容性，對展向11個站位處的剖面翼型進(jìn)行了手動修型設(shè)計，并使其滿足厚度約束要求。圖10為基準(zhǔn)機翼的上表面壓力云圖。圖11和圖12給出了其中9個站位處最終獲得的基準(zhǔn)機翼剖面翼型形狀和設(shè)計狀態(tài)的壓力分布，可以看出所設(shè)計的基準(zhǔn)機翼具有很好的幾何兼容性和氣動兼容性。

圖10 基準(zhǔn)機翼上表面壓力云圖Fig.10 Pressure contour of baseline wing on upper surface

圖11 基準(zhǔn)機翼沿展向9個站位處的剖面翼型形狀Fig.11 Sectional airfoil shapes of baseline wing at 9 span-wise locations

圖12 基準(zhǔn)機翼沿展向9個站位處的剖面壓力分布Fig.12 Sectional pressure distribution of baseline wing at 9 span-wise locations

4.2.2 寬體客機機翼第1輪氣動優(yōu)化設(shè)計

對本文所設(shè)計的基準(zhǔn)機翼進(jìn)行氣動特性分析發(fā)現(xiàn)，在Ma=0.83狀態(tài)機翼上表面中段出現(xiàn)較強的雙激波，導(dǎo)致阻力明顯增大，且設(shè)計機翼的誘導(dǎo)阻力也偏大。因此，針對翼身組合體構(gòu)型的機翼中部的5個剖面及扭轉(zhuǎn)分布開展了第1輪多點減阻優(yōu)化設(shè)計。剖面優(yōu)化設(shè)計問題可以描述為

(4)

該優(yōu)化問題共有5個設(shè)計目標(biāo)，分別對應(yīng)于5個計算狀態(tài)的阻力系數(shù)最小化。設(shè)計變量為5個剖面處的翼型形狀，每個剖面采用5階CST進(jìn)行參數(shù)化(即每個剖面具有12個設(shè)計變量)，共60個設(shè)計變量。優(yōu)化收斂曲線如圖13所示?？梢钥闯觯?個計算狀態(tài)下的阻力系數(shù)均有所減小。另外，對機翼扭轉(zhuǎn)分布進(jìn)行了優(yōu)化，優(yōu)化后的扭轉(zhuǎn)分布如圖14所示。優(yōu)化結(jié)果與基準(zhǔn)機翼的對比如表3所示。由于計算中沒有考慮機身黏性，計算獲得的阻力系數(shù)值偏小，但優(yōu)化前后阻力系數(shù)的差值反映了機翼優(yōu)化設(shè)計的減阻效果。圖15給出了第1輪優(yōu)化機翼的翼身組合體構(gòu)型在設(shè)計狀態(tài)與基準(zhǔn)外形的壓力云圖、升力分布、扭轉(zhuǎn)分布、厚度分布及典型站位壓力分布對比。圖中：b

為機翼半展長;α為迎角。可以看出，機翼表面的激波強度明顯減弱。從阻力系數(shù)分解來看，優(yōu)化外形的激波阻力和誘導(dǎo)阻力均顯著減小。

圖16給出了第一輪優(yōu)化機翼外翼段5個站位處的剖面翼型形狀。盡管優(yōu)化機翼的阻力系數(shù)明顯減小，但僅僅依靠數(shù)值優(yōu)化設(shè)計難以保證幾何兼容性。雖然優(yōu)化機翼可能并不實用，但第1輪優(yōu)化設(shè)計結(jié)果為設(shè)計者提供了一個可行的設(shè)計方向，即增加各設(shè)計剖面翼型的前后加載。

圖13 第1輪機翼氣動優(yōu)化設(shè)計阻力系數(shù)收斂歷程Fig.13 Convergence history of drag coefficient in the first optimization round

圖14 第1輪優(yōu)化機翼與基準(zhǔn)機翼的扭轉(zhuǎn)分布對比Fig.14 Comparison of twist distribution between optimal wing and baseline in the first round

表3 第1輪機翼氣動優(yōu)化設(shè)計結(jié)果Table 3 Optimization results of the first round wing aerodynamic design

圖15 寬體客機機翼第1輪優(yōu)化設(shè)計后的壓力云圖、升力分布、扭轉(zhuǎn)分布、厚度分布及典型站位壓力分布Fig.15 Pressure contour, lift, twist, thickness and pressure distribution of typical span-wise stations of optimal wing in the first round optimization design

圖16 第1輪優(yōu)化機翼的外翼段沿展向5個站位處剖面翼型形狀Fig.16 Sectional airfoil shapes of optimal wing at 5 span-wise positions in the first round

4.2.3 寬體客機機翼第2輪氣動優(yōu)化設(shè)計

觀察第1輪設(shè)計結(jié)果可以看出，由于對各剖面翼型厚度進(jìn)行了限制，導(dǎo)致在翼型下表面出現(xiàn)了不合理的凸起，從而影響了各剖面翼型的幾何兼容性。此外，發(fā)現(xiàn)當(dāng)考慮發(fā)動機短艙后，第1輪設(shè)計外形的氣動特性明顯惡化,因此，在第2輪優(yōu)化設(shè)計中，修改了厚度約束，并針對機身-機翼-發(fā)動機短艙的全機構(gòu)型開展外翼段7個剖面的多點優(yōu)化設(shè)計。

該優(yōu)化問題同樣具有5個優(yōu)化目標(biāo)，每個剖面采用5階CST參數(shù)化，共有84個設(shè)計變量。采用代理優(yōu)化算法進(jìn)行優(yōu)化，并結(jié)合人工修型使得優(yōu)化外形滿足幾何兼容性要求。該優(yōu)化問題可以描述為

minCD,1,CD,2,CD,3,CD,4,CD,5

1)Ma=0.83,CL=0.50,Re=4.0×107

2)Ma=0.85,CL=0.48,Re=4.0×107

3)Ma=0.85,CL=0.50,Re=4.0×107

4)Ma=0.85,CL=0.515,Re=4.0×107

5)Ma=0.87,CL=0.50,Re=4.0×107

Design sections: Kink, Mid1, Mid2, Mid3,

Mid4, Tip, Winglet

s.t. 1)t1>0.109

2)t2>0.095

3)t3>0.095

4)t4>0.095

5)t5>0.095

6)t6>0.095

7)t7>0.095

(5)

優(yōu)化機翼與基準(zhǔn)機翼的氣動特性對比如圖17所示?？梢钥闯?，第2輪優(yōu)化機翼的激波明顯減弱，并保證了較好的幾何兼容性和氣動兼容性，驗證了本文基于代理模型的多輪次氣動優(yōu)化設(shè)計方法的有效性。

圖17 寬體客機機翼第2輪優(yōu)化設(shè)計后的壓力云圖、升力分布、扭轉(zhuǎn)分布、厚度分布及典型站位壓力分布Fig.17 Pressure contour, lift, twist, thickness and pressure distribution of typical span-wise stations of optimal wing in the second round optimization design

5 寬體客機全機巡航構(gòu)型氣動特性分析

本節(jié)對基準(zhǔn)外形和優(yōu)化外形，采用RANS求解器進(jìn)行安裝優(yōu)化機翼的全機巡航構(gòu)型的氣動特性分析與對比。出于驗證本文方法的考慮，只評估了高速巡航狀態(tài)性能。

5.1 計算網(wǎng)格生成

采用ICEM六面體結(jié)構(gòu)網(wǎng)格，網(wǎng)格單元數(shù)約為5 700萬。遠(yuǎn)場距離取為參考弦長的100倍。網(wǎng)格劃分時，保證物面第1層網(wǎng)格的y+約為0.5，附面層內(nèi)網(wǎng)格法向增長率為1.25。同時，各翼面前后緣處的網(wǎng)格弦向尺度為當(dāng)?shù)叵议L的0.1%，翼根與翼梢處網(wǎng)格的展向尺度約為半展長的0.1%，機身上頭尾兩處網(wǎng)格的弦向尺度約為2%平均氣動弦長。全機半模的表面網(wǎng)格如圖18所示。

圖18 設(shè)計外形的全機巡航構(gòu)型表面網(wǎng)格示意圖Fig.18 Schematic of surface grid of cruise configuration of designed shape

5.2 RANS方程求解器及其驗證

為了獲得更準(zhǔn)確的機翼氣動性能，并驗證優(yōu)化結(jié)果的有效性，本文采用基于RANS方程的求解器“TNS”對優(yōu)化前后的全機巡航構(gòu)型進(jìn)行評估?？臻g離散采用中心格式，湍流模型采用SST湍流模型，全湍流計算。

在評估前，通過對DLR-F6翼身組合體對求解器精度進(jìn)行了驗證。采用的計算網(wǎng)格為AIAA阻力預(yù)測大會DPW提供的六面體結(jié)構(gòu)網(wǎng)格，網(wǎng)格量約為800萬。計算狀態(tài)為：Ma=0.75,Re=3.0×106,CL=0.5。圖19給出了采用TNS計算獲得的2個典型站位處剖面壓力分布與實驗結(jié)果的對比[31]。從評估結(jié)果可以看出，采用TNS求解器獲得的壓力分布與實驗結(jié)果吻合良好，能夠較準(zhǔn)確地捕捉激波位置，從而驗證了采用TNS求解器進(jìn)行評估的可靠性。

圖19 DLR-F6翼身組合體采用TNS求解器評估的剖面壓力分布與實驗結(jié)果對比[31]Fig.19 Comparison of sectional pressure distribution evaluated by TNS and experimental results of DLR-F6 wing-body configuration[31]

5.3 氣動特性分析評估

優(yōu)化前后全機巡航構(gòu)型的評估狀態(tài)為：Ma=0.85,Re=4.0×107,CL=0.5，計算結(jié)果如表4所示。作為參考，表4中也給出了WBCAero評估的全機巡航構(gòu)型且考慮機身黏性的計算結(jié)果?？梢钥闯?，采用TNS評估的優(yōu)化外形相對基準(zhǔn)外形阻力系數(shù)減小7.8 counts，充分驗證了本文優(yōu)化得到的氣動外形具有較優(yōu)的巡航氣動性能。除此之外，表4還給出了其他RANS求解器的評估結(jié)果，采用筆者所在課題組開發(fā)的PMNS3DR求解器評估優(yōu)化外形的阻力系數(shù)減小了9.2 counts。多種求解器評估的結(jié)果雖然略有差別，但都表明優(yōu)化外形相對于基準(zhǔn)外形的巡航氣動性能得到顯著改善。圖20和圖21分別給出了采用TNS求解器評估的基準(zhǔn)外形與優(yōu)化外形上表面壓力云圖和6個機翼展向位置處的壓力分布對比圖?？梢钥闯?，各站位處的激波強度均有不同程度減弱，驗證了本文基于代理模型的多輪次優(yōu)化設(shè)計方法能夠有效用于大型民機的機翼氣動優(yōu)化設(shè)計。需要說明的是，從驗證優(yōu)化設(shè)計方法的角度出發(fā)，本文只考慮了高速巡航狀態(tài)設(shè)計。而對于實際設(shè)計，還必須考慮低速設(shè)計點，以及考慮阻力發(fā)散馬赫數(shù)和抖振邊界的其他性能指標(biāo)。

以上算例研究結(jié)果，充分展示了本文多輪次代理優(yōu)化算法在解決復(fù)雜外形氣動優(yōu)化設(shè)計方面的有效性和高效性，為解決大型民機復(fù)雜氣動設(shè)計問題提供了一種有力途徑。但需要說明的是，本文首先采用全速勢方程加黏性修正的流動數(shù)值模擬方法進(jìn)行優(yōu)化設(shè)計，然后采用RANS方程求解器對優(yōu)化結(jié)果進(jìn)行評估。但由于求解器精度限制，如果要進(jìn)一步提高其氣動性能，則必須直接采用RANS方程求解器進(jìn)行優(yōu)化。

表4 采用TNS求解器對基準(zhǔn)外形與優(yōu)化外形評估的氣動性能對比Table 4 Comparison of aerodynamic performance of baseline and optimal configuration evaluated using TNS

圖20 RANS求解器計算得到的基準(zhǔn)外形與設(shè)計外形表面壓力云圖對比Fig.20 Comparison of surface pressure contour of baseline and optimal configuration using RANS solver

圖21 RANS求解器計算得到的設(shè)計外形與基準(zhǔn)外形剖面壓力分布對比Fig.21 Comparison of sectional pressure distribution of baseline and optimal configuration using RANS solver

6 結(jié) 論

本文介紹了大型民機超臨界機翼氣動優(yōu)化設(shè)計的基本思想與要點，研究提出了一種面向工程應(yīng)用的、基于代理模型的多輪次氣動優(yōu)化設(shè)計方法。與人工修型相結(jié)合，將本文方法應(yīng)用于寬體客機超臨界機翼氣動優(yōu)化設(shè)計，并對優(yōu)化前后全機巡航構(gòu)型采用RANS方程求解器進(jìn)行了典型狀態(tài)的氣動特性評估，結(jié)果表明氣動性能得到較顯著提升。

1) 對于寬體大型民機超臨界機翼設(shè)計，必須開展一體化優(yōu)化設(shè)計。在三維流動環(huán)境下進(jìn)行機翼剖面優(yōu)化設(shè)計，同時還需要考慮機身和發(fā)動機短艙等對機翼氣動性能的影響。

2) 代理優(yōu)化算法已發(fā)展至適用于100維以內(nèi)的優(yōu)化設(shè)計問題，能夠有效應(yīng)用于飛行器復(fù)雜外形的氣動優(yōu)化設(shè)計，具有較強的全局優(yōu)化能力、較高的優(yōu)化效率和很強的約束處理能力。

3) 采用基于代理模型的多輪次氣動優(yōu)化設(shè)計方法，通過不斷調(diào)整優(yōu)化問題的設(shè)計空間、優(yōu)化目標(biāo)、約束，并充分利用之前優(yōu)化設(shè)計中產(chǎn)生的分析數(shù)據(jù)加速優(yōu)化收斂，能快速獲得滿足設(shè)計要求的優(yōu)化設(shè)計結(jié)果。

4) 對于飛行器氣動設(shè)計問題，采用數(shù)值優(yōu)化與人工修型相結(jié)合是必要的。僅依靠數(shù)值優(yōu)化，難以保證優(yōu)化結(jié)果滿足幾何與氣動兼容性，需要通過人工經(jīng)驗對優(yōu)化外形進(jìn)行修型，以滿足設(shè)計要求；另一方面，數(shù)值優(yōu)化也能夠為人工經(jīng)驗提供指導(dǎo)，明確設(shè)計方向，有利于獲得更優(yōu)的氣動外形。

本文從驗證和展示優(yōu)化設(shè)計方法的角度，主要考慮了大型寬體民機超臨界機翼高速巡航狀態(tài)的多點氣動設(shè)計。在實際設(shè)計中，還需要考慮低速性能及阻力發(fā)散馬赫數(shù)和抖振邊界能綜合性能。所面臨的問題主要是所采用設(shè)計變量數(shù)還不能滿足寬體客機機翼精細(xì)化設(shè)計的需求，需要對優(yōu)化方法和技術(shù)進(jìn)行進(jìn)一步改進(jìn)。

致 謝

感謝德國宇航院Stefan G?ertz博士在代理優(yōu)化算法方面的建議和有益討論。感謝中國商飛上海飛機設(shè)計研究院張淼、程攀、馬涂亮、劉鐵軍等同志的指導(dǎo)和大力支持。