基于離散伴隨方法的三維機(jī)翼外形優(yōu)化設(shè)計(jì)

摘要：目前民用飛機(jī)燃油消耗問(wèn)題嚴(yán)重，因此開(kāi)展氣動(dòng)外形優(yōu)化對(duì)提高整機(jī)氣動(dòng)性能和燃油效率具有重要意義。針對(duì)全局優(yōu)化方法和局部?jī)?yōu)化方法的特點(diǎn)，本文對(duì)民機(jī)開(kāi)展基于離散伴隨方法多約束氣動(dòng)問(wèn)題進(jìn)行研究，極大降低了計(jì)算量并提升了整機(jī)氣動(dòng)性能。在保證升力大小不變、機(jī)翼厚度和俯仰力矩不低于約束值的條件下，以降低飛機(jī)機(jī)翼阻力、提高機(jī)翼的升阻比為設(shè)計(jì)目標(biāo)，對(duì)波音737-800機(jī)翼進(jìn)行了氣動(dòng)外形優(yōu)化設(shè)計(jì)。基于離散伴隨方法來(lái)獲得梯度信息，使用自由變形法進(jìn)行機(jī)翼表面參數(shù)化，數(shù)值模擬基于可壓縮歐拉方程進(jìn)行外形優(yōu)化。結(jié)果表明，自由變形法切實(shí)可行，在滿足上述約束條件的前提下，波音737-800機(jī)翼在跨聲速及小迎角條件下的阻力系數(shù)顯著下降，升阻比得到大幅提升，取得了較好的優(yōu)化效果。本文建立的優(yōu)化設(shè)計(jì)方法為民機(jī)氣動(dòng)優(yōu)化設(shè)計(jì)提供了參考并且可以廣泛應(yīng)用于工程實(shí)踐。

關(guān)鍵詞：離散伴隨方法；自由變形技術(shù)；氣動(dòng)外形優(yōu)化；升阻比；CFD技術(shù)

中圖分類號(hào)：V211.4文獻(xiàn)標(biāo)識(shí)碼：ADOI：10.19452/j.issn1007-5453.2024.11.003

當(dāng)前實(shí)現(xiàn)碳中和、碳達(dá)峰目標(biāo)已經(jīng)提升到國(guó)家戰(zhàn)略高度。在民用航空領(lǐng)域提高客機(jī)的經(jīng)濟(jì)性、減少碳排放主要通過(guò)降低飛機(jī)所受的空氣阻力、提高航空發(fā)動(dòng)機(jī)的性能、降低飛機(jī)的重量三個(gè)途徑來(lái)實(shí)現(xiàn)，除此以外，還有使用可持續(xù)航空燃料（SAF）與液氫這類低碳零碳甚至負(fù)碳的液體燃料來(lái)替代航空燃料等途徑。其中，降低飛機(jī)所受的空氣阻力主要通過(guò)氣動(dòng)外形設(shè)計(jì)實(shí)現(xiàn)。根據(jù)統(tǒng)計(jì)結(jié)果，一架中寬體民航客機(jī)提高1%的升阻比，每年就可以節(jié)約10萬(wàn)美元的燃料費(fèi)用[1]。隨著計(jì)算流體力學(xué)（CFD）技術(shù)的快速發(fā)展，基于CFD技術(shù)的氣動(dòng)外形優(yōu)化可分為全局優(yōu)化方法和局部?jī)?yōu)化方法兩類。全局優(yōu)化方法的計(jì)算量相對(duì)較大；局部?jī)?yōu)化方法則通過(guò)求解梯度信息，通過(guò)梯度下降法獲得局部最優(yōu)解。局部?jī)?yōu)化方法計(jì)算量相對(duì)較小，但在設(shè)計(jì)變量較多的情況下，求解梯度信息的計(jì)算量仍然很大。特別是在三維復(fù)雜外形的優(yōu)化設(shè)計(jì)過(guò)程中，所需的計(jì)算時(shí)間是難以承受的。為了解決這一問(wèn)題，Jameson[2]提出將伴隨方法應(yīng)用于氣動(dòng)外形優(yōu)化設(shè)計(jì)，該方法大大降低了求解梯度信息所需的計(jì)算量。伴隨方法在近幾十年得到了迅速發(fā)展，并廣泛應(yīng)用于各類飛行器的優(yōu)化設(shè)計(jì)過(guò)程中[3-6]。國(guó)內(nèi)研究人員也針對(duì)伴隨方法開(kāi)展了大量研究，楊旭東[7]在國(guó)內(nèi)首先實(shí)現(xiàn)了三維機(jī)翼的減阻優(yōu)化設(shè)計(jì)。吳文華等[8]發(fā)展了氣動(dòng)布局優(yōu)化設(shè)計(jì)軟件，并在升力不變、保證一定的機(jī)翼容積的前提下，將機(jī)翼阻力降低了4.72%；陳曉東[9]使用海克斯-亨尼型函數(shù)法（Hicks-Henne）和類別形狀函數(shù)變換（CST）參數(shù)化法進(jìn)行翼型的氣動(dòng)外形優(yōu)化，應(yīng)用遺傳算法探尋設(shè)計(jì)空間內(nèi)的最優(yōu)方案，搭建起翼型優(yōu)化設(shè)計(jì)的自動(dòng)化平臺(tái)。張科施等[10]研發(fā)了對(duì)運(yùn)輸機(jī)機(jī)翼進(jìn)行變可信度氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)的自主軟件平臺(tái)，實(shí)現(xiàn)了減重優(yōu)化設(shè)計(jì)；顧文婷等[11]采用CST參數(shù)化幾何外形，結(jié)合代理模型和遺傳算法，構(gòu)建動(dòng)力短艙多點(diǎn)優(yōu)化設(shè)計(jì)平臺(tái)。

盡管伴隨方法已經(jīng)提出近幾十年的時(shí)間，然而該方法在國(guó)內(nèi)的研究工作并不多，主要由西北工業(yè)大學(xué)及中國(guó)空氣動(dòng)力學(xué)發(fā)展研究中心對(duì)該方法的應(yīng)用進(jìn)行了一部分的研究工作。約束條件對(duì)于氣動(dòng)外形優(yōu)化設(shè)計(jì)具有重要的意義，無(wú)約束條件下的減阻可能會(huì)導(dǎo)致升力的下降或機(jī)翼截面厚度降低進(jìn)而影響其強(qiáng)度。本文的氣動(dòng)外形優(yōu)化設(shè)計(jì)利用斯坦福大學(xué)于2016年開(kāi)發(fā)的開(kāi)源SU2氣動(dòng)外形優(yōu)化設(shè)計(jì)代碼[12]，基于離散伴隨方法及無(wú)黏流動(dòng)控制方程（歐拉方程），在同時(shí)限定升力、機(jī)翼多個(gè)截面翼型厚度和俯仰力矩的條件下，對(duì)波音737-800機(jī)翼進(jìn)行了氣動(dòng)外形優(yōu)化設(shè)計(jì)，優(yōu)化設(shè)計(jì)結(jié)果顯示，在滿足約束的條件下，機(jī)翼在跨聲速狀態(tài)時(shí)的阻力系數(shù)下降明顯，升阻比得到大幅度提升。

1流動(dòng)控制方程及求解

本文的流場(chǎng)求解基于守恒形式的歐拉方程

式中，W為守恒變量，F(xiàn)和G為對(duì)流項(xiàng)通量，Ω為流場(chǎng)邊界。流場(chǎng)計(jì)算過(guò)程中對(duì)流項(xiàng)采用顯式中心差分有限體積法（JST）格式進(jìn)行離散［13］，時(shí)間項(xiàng)采用歐拉隱式格式進(jìn)行離散。

2伴隨方法

設(shè)置待優(yōu)化的氣動(dòng)外形邊界位置為X，流場(chǎng)變量為F，優(yōu)化設(shè)計(jì)的目標(biāo)函數(shù)為I，則其表達(dá)式為

I=（F，X）（2）

通過(guò)對(duì)式（2）進(jìn)行變分，進(jìn)而能夠得到優(yōu)化設(shè)計(jì)目標(biāo)函數(shù)對(duì)流場(chǎng)變量F及所設(shè)計(jì)外形邊界位置X的梯度

式中，δ為極小的變化量，定義偏導(dǎo)數(shù)過(guò)程。

流動(dòng)控制方程可表示為如下形式

R（W，S）=0（4）

式中，R為流動(dòng)控制方程的殘差，在控制方程完全收斂的情況下，殘差R的值近似為0，對(duì)式（4）進(jìn)行變分，可以得到

由于δR=0，因此可引入任意拉格朗日乘子與式（5）相乘，并與式（3）合并得

因拉格朗日乘子φ為任意值，因此可以假設(shè)其滿足以下條件

即可將式（6）轉(zhuǎn)化為以下形式

通過(guò)上述過(guò)程可以看出，通過(guò)構(gòu)造伴隨方程，即式（7），只需要求解流動(dòng)控制方程及其伴隨方程，即可得到目標(biāo)函數(shù)對(duì)設(shè)計(jì)變量的梯度信息，而伴隨方程與流動(dòng)控制方程的計(jì)算量較為接近，因此，梯度求解過(guò)程計(jì)算量與求解兩次流動(dòng)控制方程基本相當(dāng)，與傳統(tǒng)的有限差分法相比，計(jì)算效率能夠顯著提高。

3優(yōu)化設(shè)計(jì)流程

本文所采用的氣動(dòng)優(yōu)化設(shè)計(jì)流程主要由自由變形技術(shù)參數(shù)化表示及網(wǎng)格變形模塊[14]、CFD氣動(dòng)力計(jì)算模塊、優(yōu)化搜索模塊組成。幾何參數(shù)化表示及網(wǎng)格變形技術(shù)是氣動(dòng)外形優(yōu)化設(shè)計(jì)中的關(guān)鍵一環(huán)，自由變形法（FFD）技術(shù)起源于計(jì)算機(jī)圖形學(xué)，其基本原理是將網(wǎng)格與可變形的FFD控制體相關(guān)聯(lián)，通過(guò)改變控制體的形狀進(jìn)而產(chǎn)生網(wǎng)格的變形。因此本文的設(shè)計(jì)變量為FFD控制體上的控制點(diǎn)，通過(guò)求解目標(biāo)函數(shù)對(duì)各控制點(diǎn)的梯度，計(jì)算控制點(diǎn)在每一次優(yōu)化迭代過(guò)程中的變形量，最后由FFD算法得出機(jī)翼表面及內(nèi)部流場(chǎng)網(wǎng)格點(diǎn)的變形量。優(yōu)化設(shè)計(jì)目標(biāo)為機(jī)翼的阻力系數(shù)，并沿著機(jī)翼展向均勻分布的5個(gè)截面上施加了翼型厚度約束條件；流場(chǎng)控制方程為無(wú)黏可壓縮的歐拉方程，通過(guò)CFD計(jì)算獲取機(jī)翼外形在一定工況下的氣動(dòng)特性，通過(guò)離散伴隨方法獲得氣動(dòng)特性隨幾何外形參數(shù)變化梯度，結(jié)合優(yōu)化算法，最終得到滿足約束條件的優(yōu)化設(shè)計(jì)結(jié)果。本文的氣動(dòng)外形優(yōu)化流程如圖1所示。

4氣動(dòng)外形優(yōu)化設(shè)計(jì)

4.1CFD計(jì)算方法驗(yàn)證

CFD計(jì)算方法的可靠性直接關(guān)系到優(yōu)化設(shè)計(jì)結(jié)果的可信度。因此在優(yōu)化設(shè)計(jì)前，對(duì)RAE2822翼型在迎角為2.32°、Ma0.75的情況下，進(jìn)行了流場(chǎng)仿真計(jì)算。計(jì)算結(jié)果與試驗(yàn)結(jié)果的壓力系數(shù)Cp對(duì)比情況如圖2所示，流場(chǎng)馬赫數(shù)Ma分布情況如圖3所示。結(jié)果顯示，計(jì)算結(jié)果與試驗(yàn)結(jié)果非常接近，馬赫數(shù)分布情況也符合預(yù)期。由此可見(jiàn)，本文優(yōu)化設(shè)計(jì)中所采用的CFD計(jì)算方法具備較高的可靠度，能夠應(yīng)用于氣動(dòng)外形優(yōu)化設(shè)計(jì)中。

4.2波音737-800機(jī)翼氣動(dòng)外形優(yōu)化設(shè)計(jì)

本文對(duì)波音737-800機(jī)翼進(jìn)行了氣動(dòng)外形優(yōu)化設(shè)計(jì)，機(jī)翼模型繞流場(chǎng)網(wǎng)格采用結(jié)構(gòu)化網(wǎng)格，機(jī)翼翼面及對(duì)稱面上的網(wǎng)格分布情況如圖4所示，流場(chǎng)中網(wǎng)格數(shù)量為614804。工作條件設(shè)置為Ma0.8395。在優(yōu)化設(shè)計(jì)過(guò)程中約束升力系數(shù)CL=0.3，同時(shí)約束沿展向均勻分布的5個(gè)截面厚度不小于原厚度的90%，俯仰力矩不小于初始值，截面分布情況如圖5所示。機(jī)翼的外形采用FFD控制體進(jìn)行參數(shù)化表達(dá)，設(shè)置8個(gè)控制頂點(diǎn)后通過(guò)SU2_DEF模塊生成了由352（11×16×2）個(gè)控制點(diǎn)組成的控制體，其分布如圖6所示。

機(jī)翼的外形優(yōu)化設(shè)計(jì)迭代了10步后阻力系數(shù)CD趨向于收斂，優(yōu)化過(guò)程中升力系數(shù)CL及阻力系數(shù)CD隨迭代步數(shù)的變化情況如圖7所示。從圖7中可以看出，在優(yōu)化設(shè)計(jì)過(guò)程中升力系數(shù)基本保持不變，較好地滿足了約束條件，阻力系數(shù)則顯著降低。

優(yōu)化設(shè)計(jì)前后的升力及阻力系數(shù)的對(duì)比情況見(jiàn)表1，優(yōu)化前機(jī)翼的阻力系數(shù)為CD=0.0137，優(yōu)化后變?yōu)镃D=0.0078，較優(yōu)化前降低了43%，升阻比則提高了L/D=77.8%。圖8給出了5個(gè)截面優(yōu)化前后的幾何形狀對(duì)比情況。由圖8可以看出，機(jī)翼截面的形狀變化不大，滿足約束條件。

選取機(jī)翼上均勻分布的5個(gè)截面的上下翼面馬赫數(shù)分布、壓力系數(shù)分布作為對(duì)整體優(yōu)化效果的評(píng)判，如圖9、圖10所示。從圖中可見(jiàn)，在機(jī)翼同一位置處，優(yōu)化后的壓力系數(shù)大于優(yōu)化前的壓力系數(shù)，因而其產(chǎn)生的向前的推力大于優(yōu)化前的，可見(jiàn)壓差阻力減小，這與圖11、圖12所示云圖情況相符合。

圖11給出了優(yōu)化前后的壓力云圖對(duì)比情況，從圖11中可看出，優(yōu)化前上翼面的激波強(qiáng)度較大，低壓區(qū)分布不規(guī)則，壓差阻力較大；優(yōu)化后的機(jī)翼低壓區(qū)顯著減少，壓力分布較為均勻，壓差阻力較小。圖12給出了優(yōu)化前后的馬赫數(shù)云圖對(duì)比情況，從圖12中可以看出，優(yōu)化前機(jī)翼上表面高流速區(qū)域與低流速區(qū)域分隔明顯，因此具有較高的激波阻力，而優(yōu)化后的機(jī)翼在升力幾乎不變的情況下極大地降低了高流速區(qū)域的覆蓋面積，馬赫數(shù)沿弦向的變化趨于平緩。

5結(jié)論

通過(guò)研究，可以得出以下結(jié)論：

（1）本文基于離散伴隨方法，以降低飛機(jī)阻力為目標(biāo)，用機(jī)翼表面壓力分布和馬赫數(shù)分布作為評(píng)價(jià)指標(biāo)，開(kāi)展了多約束氣動(dòng)優(yōu)化設(shè)計(jì)研究，具有工程應(yīng)用價(jià)值。

（2）在跨聲速條件下，機(jī)翼表面激波強(qiáng)度大，通過(guò)本文建立的優(yōu)化設(shè)計(jì)平臺(tái)阻力降低43%、升阻比提高77.8%，證明了此方法的有效性。

（3）基于離散伴隨的梯度類優(yōu)化方法計(jì)算速度相對(duì)較快，計(jì)算量相對(duì)較小，可以應(yīng)用于其他設(shè)計(jì)方案。

參考文獻(xiàn)

［1］LynchFT.Commercialtransports-aerodynamicdesignforcruiseperformanceefficiency[J].TransonicAerodynamics，1981，30（5）：81-147.

［2］JamesonA.Aerodynamicdesignviacontroltheory[J].JournalofScientificComputing，1988，3（3）：233-260.

［3］MavriplisDJ.Formulationandmultigridsolutionofthediscreteadjointproblemonunstructuredmeshes[M].Berlin：Springer，2004.

［4］MavriplisDJ.Discreteadjoint-basedapproachforoptimizationproblemonthreedimensionalunstructuredmeshes[J].AIAAJournal，2007，45（4）：741-750.

［5］ChristakopoulosF，JonesD，MüllerJD.Pseudo-timesteppingandverificationforautomaticdifferentiationderivedCFDcodes[J].Computersamp;Fluids，2011，46（1）：174-179.

［6］LyuZ，KenwayGK，PaigeC，etal.AutomaticdifferentiationadjointoftheReynolds-averagedNavier-Stokesequationswithaturbulencemodel[C].21stAIAAComputationalFluidDynamicsConference，2013：2581.

［7］楊旭東.基于控制理論的氣動(dòng)優(yōu)化設(shè)計(jì)技術(shù)研究[D].西安：西北工業(yè)大學(xué)，2002.

YangXudong.Researchonpneumaticoptimizationdesigntechnologybasedoncontroltheory[D].Xi’an：NorthwesternPolytechnicalUniversity，2002.（inChinese）

［8］吳文華，范召林，陳德華，等.基于伴隨算子的大飛機(jī)氣動(dòng)布局精細(xì)優(yōu)化設(shè)計(jì)[J].空氣動(dòng)力學(xué)學(xué)報(bào)，2012，30（6）：719-724+760.

WuWenhua，F(xiàn)anZhaolin，ChenDehua，etal.Fineoptimizationdesignoflargeaircraftaerodynamiclayoutbasedonconcomitantoperator[J].ActaAerodynamicaSinica，2012，30（6）：719-724+760.（inChinese）

［9］陳曉東.基于遺傳算法的翼型優(yōu)化設(shè)計(jì)[D].哈爾濱：哈爾濱工程大學(xué)，2018.

ChenXiaodong.Optimizeddesignofairfoilbasedongeneticalgorithm[D].Harbin：HarbinEngineeringUniversity，2018.（inChinese）

［10］張科施，凌圣博，韓忠華.跨聲速運(yùn)輸機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化平臺(tái)AeroStruct的發(fā)展及應(yīng)用[J].航空科學(xué)技術(shù)，2022，33（4）：47-56.

ZhangKeshi，LingShengbo，HanZhonghua.DevelopmentandApplicationofAeroStruct，anaerodynamic/structuraloptimizationplatformfortransonictransportaircraftwings[J].AeronauticalScienceamp;Technology，2022，33（4）：47-56.（inChinese）

［11］顧文婷，馬坤，韓杰.動(dòng)力短艙多點(diǎn)氣動(dòng)優(yōu)化設(shè)計(jì)[J].航空科學(xué)技術(shù)，2024，35（2）：23-30.

GuWenting，MaKun，HanJie.Multipointaerodynamicdesignoptimizationofpowered-onnacelle[J].AeronauticalScienceamp;Technology，2024，35（2）：23-30.（inChinese）

［12］EconomonTD，PalaciosF，CopelandSR，etal.SU2：anopensourcesuiteformultiphysicssimulationanddesign[J].AIAAJoural，2016，54（3）：828-846.

［13］JamesonA，SchmidtW，TurkelE.NumericalsolutionoftheEulerequationsbyfinitevolumemethodsusingRungeKuttatimesteppingschemes[C].14thFluidandPlasmaDynamicsConference，1981：1259.

［14］SederbergTW，ParrySR.Free-formdeformationofsolidgeometricmodels[C]//Proceedingsofthe13thAnnualConferenceonComputerGraphicsandInteractiveTechniques，1986：151-160.

基金項(xiàng)目：內(nèi)蒙古自然科學(xué)基金（2024LHMS01012）