魏闖,楊龍,李春鵬,張鐵軍
中國航空工業空氣動力研究院,高速高雷諾數氣動力航空科技重點實驗室,沈陽 110034
傳統的依賴經驗、以試湊(Cut and Try)為主的人工修型氣動設計方法,設計效率低,很難獲得最佳的氣動外形。隨著計算流體力學(Computational Fluid Dynamics,CFD)技術、計算機輔助幾何造型設計技術(Computer-Aided Design,CAD)以及數值優化等技術的發展和日趨成熟,基于高可信度計算流體力的氣動分析與優化設計技術,目前已廣泛應用于航空、航天飛行器的設計與研制,在提高飛行器氣動與綜合性能方面正發揮著越來越重要的作用[1-2]。
完整的飛行器氣動外形優化流程包括氣動外形參數化、網格自動生成、氣動性能評估和高效優化算法等。根據優化算法的不同現有的優化設計方法主要分為3類[3]:①遺傳算法、粒子群算法等啟發式優化算法;②梯度優化算法;③代理優化算法。3種方法各有優劣,國內外研究者均開展了廣泛的研究。如斯坦福大學Jameson等在伴隨方法方面開展了卓有成效的先驅性研究工作[4],而后該校Alonso團隊開發了基于伴隨方法的SU2[5-6]開源氣動優化程序,波音公司研制了MDOPT[7]優化軟件用于氣動外形優化,美國NASA埃姆斯研究中心CFD求解器FUN3D[8-9]、德國宇航院非結構求解器TAU[10-11]以及法宇航的求解器e1sA[12-13]均集成發展了氣動外形優化設計功能。國內中國空氣動力研究與發展中心[14-16]、西北工業大學[17-19]、清華大學[20]、南京航空航天大學[21-23]、中國航空工業空氣動力研究院[24]等也開展了廣泛深入的研究,分別建立了各自的氣動優化設計工具。
本文簡要介紹了中國航空工業空氣動力研究院氣動優化軟件(Aerodynamics Research Institute Optimization Code,ARI_OPT))各個模塊的基本原理和方法,給出了基于數值優化算例的功能驗證,詳細闡述了其在翼型、多段翼型和飛翼布局機翼等典型單目標和多目標氣動外形優化問題的應用算例,旨在促進氣動外形優化設計領域的交流。
ARI_OPT是中國航空工業空氣動力研究院針對飛行器氣動外形優化的氣動優化軟件。它可以求解任意多約束的單目標、多目標優化問題,主要應用于基于高精度CFD的氣動優化設計,也可用于氣動/結構、氣動/隱身等多學科優化及其他工程優化設計問題,在軍、民用飛行器氣動外形優化設計上得到了較為廣泛的應用。優化軟件集成于中國航空工業空氣動力研究院大型計算集群,可進行大規模并行優化。該軟件包含氣動外形參數化、網格自動變形、高逼真度CFD數值模擬、代理模型和高效優化算法等模塊。
ARI_OPT包含多種氣動外形參數化方法,如Hicks-Henne[25]、類別形狀函數變換(Class-Shape Transformation,CST)方法[26-27]、自由曲面變形(Free Form Deform,FFD)[28]、直接操作自由曲面變形(Directly manipulated FFD,DFFD)[29]等參數化方法。同時發展了基于CATIA二次開發的氣動外形參數化方法,該方法能準確地生成飛行器三維外形,同時可以方便地輸出表面面積、內部容積、前后梁高度和位置等幾何信息。另外可直接輸出可供上下游設計階段直接利用的CAD模型,相比其他參數化方法如直接控制網格點的方法等,無需再根據優化結果重新建立外形CAD模型,避免了模型精度損失,提高了工程實用性。典型的機翼參數化流程如圖1所示。

圖1 基于CATIA民機機翼參數化建模Fig.1 Parametric modeling of wing based on CATIA automation
優化過程中利用徑向基函數(Radial Basis Function,RBF)網格變形技術實現網格自動生成。徑向基函數的基本形式為[30]
(1)

本文徑向基函數采用計算效率與網格變形質量都較好的Wendland’s C2函數[30-31],并通過貪心算法[32]縮減變形矩陣維數和基于MPI(Message Passing Interface)的并行求解技術提高變形效率。典型帶流板下偏的三段翼構型,8.1萬 網格條件下,30次網格變形平均耗時不超過1 s,網格變形結果如圖2所示;NASA的CRM翼身組合體構型,1 017萬網格條件下,30次網格變形耗時平均約9.2 s,如圖3所示,具有較高的網格變形質量。

圖2 多段翼型網格變形結果Fig.2 Grid deformation of multi-element airfoil

圖3 翼型組合體網格變形結果Fig.3 Grid deformation of wing-body
優化過程中氣動性能評估依托中國航空工業空氣動力研究院航空大規模CFD平臺ARI_CFD(如圖4所示),主要使用求解雷諾平均Naiver-Stoke(RANS)方程的結構網格模塊ARI_ENSMB和非結構動態重疊網格模塊ARI_Overset,2個模塊都采用有限體積法,包含多種常用空間差分格式和湍流模型,在長期應用中得到了大量算例的廣泛驗證[33-41]。對于ARI_ENSMB求解模塊,黏性項采用2階中心差分格式,空間離散無黏項采用Roe-FDS(Roe’s Flux Difference Splitting)分裂格式,隱式LU-SGS(Lower-Upper Symmetric-Gauss-Seidel)時間推進,物面給定絕熱壁、無滑移條件,遠場采用無反射壓力遠場邊界條件,湍流模型采用k-ωSST二方程模型。圖5給出了1 000萬計算網格下DLR-F6翼身組合體構型ARI_ENSMB計算馬赫數Ma=0.75、雷諾數Re=3.0×106時不同迎角α下升力系數CL、阻力系數CD與德宇航TAU、NASA FUN3D、商業軟件FLUENT等求解器計算以及風洞試驗值的對比結果[42],可見ARI_ENSMB與試驗值符合較好,具有較高的計算準度。

圖4 ARI_CFD航空大規模CFD平臺Fig.4 Aeronautical large-scale CFD platform of ARI_CFD

圖5 DLR-F6翼身組合體構型ARI_ENSMB計算與試驗值及其他求解器結果對比(Ma=0.75、Re=3.0×106)Fig.5 Comparison of results between test data and ARI_ENSMB and other solvers for DLR-F6 wingbody(Ma=0.75,Re=3.0×106)
代理模型技術采用的算法和程序模塊主要來自西北工業大學韓忠華教授團隊[3,19,43-47]。試驗設計方法包括拉丁超立方(LHS)、均勻設計(UD)、蒙特卡洛抽樣(MC)等,代理模型包含二次響應面(PRSM)、Kriging模型、梯度增強Kriging模型(GEK)、分層Kriging模型(HK)、徑向基函數(RBFs)等多種代理模型,改善期望(maximizing Expected Improvement,EI)、目標函數值最小(Minimizing Surrogate Prediction,MSP)、目標函數誤差最大 (maximizing Meansquared Error,MSE)、低置信邊界方法(minimizing Lower-Confidence Bounding,LCB)以及概率提升(Probability of Improvement,PI)等多種優化加點準則[19],既可單獨使用也可根據優化任務需求組合使用多種加點方法。
優化策略包含遺傳算法優化和基于代理模型的優化方法(簡稱代理優化(Surrogate-Based Optimization,SBO))2種,其中代理優化方法具有效率高、魯棒性好的優點,得到了較為廣泛的應用。典型代理優化方法流程如圖6所示,典型優化流程如下:① 對設計空間進行試驗設計,獲得初始樣本點并調用數值求解模塊獲得響應值,構建初始代理模型;② 基于代理模型,采用傳統優化算法求解相應的子優化問題,以很小的計算代價,對最優解進行預測,按照一定的優化加點準則獲得新樣本點;③ 調用數值求解模塊計算得到新樣本點響應值,并將結果添加到現有數據集中,不斷更新代理模型,直到所產生的樣本點序列收斂于局部或全局最優解。子優化算法包括Hooke-Jeeves模式搜索、擬牛頓梯度優化、序列二次規劃法(SQP)、單/多目標遺傳算法等多種成熟優化算法等。對于典型三維優化問題,在一般情況下,利用500 CPU核完成一輪優化耗時不超過24 h。

圖6 ARI _OPT軟件基于代理模型優化流程Fig.6 Flow chart of surrogate-based optimization of ARI _OPT software
六峰值駝背測試函數(Aix-hump Camel Function)具有6個局部極小點,其中有2個為全局極小點,優化問題描述為
(2)
分別進行了遺傳算法和代理優化方法測試,其中遺傳算法采用多島遺傳算法(Multi-Island GA,MIGA),每代種群總數為32,分4個島,優化50代,共調用1 600 次計算。代理優化方法初始采樣點20個,試驗設計方法為LHS,代理模型為Kriging模型,EI+MSP+LCB+PI混合加點方式,每次加16點,共加16輪,共調用計算276次。表1給出了30次測試的平均結果,可見代理優化方法在效率和優化結果上都優于多島遺傳算法。

表1 多島遺傳算法和代理優化方法測試結果Table 1 Comparison of test results of SBO with MIGA
ZDT3問題為多目標優化領域普遍使用的標準測試函數,真實解集Pareto凸、非連續。優化問題描述為
(3)
分別進行了遺傳算法和代理優化方法測試,其中遺傳算法采用快速非支配排序遺傳算法(NSGA-2),每代種群總數為160,經過70代優化達到收斂要求,共調用11 200次計算。代理優化方法初始采樣點20個,試驗設計方法為LHS,代理模型為Kriging模型,Pareto優化時的MSP加點方法,總計調用計算1 260次。圖7給了2種方法的優化結果,可見2種方法都得到的Pareto前沿與真實Pareto前沿基本重合,都獲得了較優的優化結果,但代理優化方法調用計算次數較少,效率較高。

圖7 NSGA-2和代理優化方法測試結果Fig.7 Test results of NSGA-2 and SBO
針對RAE2822翼型進行跨聲速減阻優化設計,設計狀態如下:
Ma=0.73,Re=6.5×106,CL=0.8
優化問題描述為

(4)
式中:Cm為俯仰力矩系數。
參數化方法采用CST翼型參數化方法,上下翼面總計18個設計變量,通過求解定常可壓RANS方程獲得優化過程中翼型氣動性能,湍流模型為k-ωSST兩方程模型。優化采用代理優化方法,初始采樣點36個,試驗設計方法為LHS,代理模型為Kriging模型,EI+MSP+LCB+PI混合加點方式,每次加點4個,共加43次,總計調用CFD計算208次。表2給出優化翼型和RAE2822氣動力結果對比,優化構型減阻約31.74%。

表2 RAE2822和優化翼型氣動性能對比Table 2 Comparison of aerodynamic performance of RAE2822 and optimal airfoil
圖8和圖9分別給出了優化翼型和RAE2822外形壓力系數Cp分布對比曲線,c為弦長。可見優化翼型顯著地減弱了激波強度,使阻力大幅降低。

圖8 RAE2822和優化翼型外形對比Fig.8 Comparison of shapes for RAE2822 and optimum airfoils

圖9 RAE2822翼型優化前后的表面壓力系數分布對比Fig.9 Comparison of surface pressure coefficient distributions for RAE2822 airfoil before and after optimum
針對某寬速域翼型,為提高其跨聲速和超聲速氣動性能進行多目標優化,優化目的是提高Ma=0.72和Ma=3.15時的最大升阻比K。
優化目標:
(5)
約束條件:
(6)
參數化方法采用改進的Hicks-Henne翼型參數化方法優化,上下翼面總計16個設計變量,通過求解定常可壓RANS方程獲得優化過程中翼型氣動性能,湍流模型為k-ωSST兩方程模型。優化算法采用NSGA-2多目標優化算法,每代種群有60個個體,總代數為50,總共完成3 000個翼型構型計算評估。
圖10給出優化得到的Pareto前沿,共計89個 優化構型,相比初始構型,Pareto前沿優化構型Ma=0.72、α=3.0°(KDP1)和Ma=3.15、α=6.0°(KDP2)下升阻比都有所增加,其中OPT1和OPT3分別為KDP2和KDP1最優構型。沿著Pareto前沿由構型OPT1到OPT3,優化構型Ma=0.72、α=3.0°升阻比越來越大,相比初始構型最大增加37.17%(OPT3);反之,Ma=3.15,α=6.0° 升阻比數越來越大,相比初始構型最大增加11.39%(OPT3),OPT2為KDP1和KDP2都增加較多優化構型,相比初始構型分別增加21.48% 和8.01%。
圖11給出了優化翼型外形對比結果對比,圖12 和圖13給出了優化構型和初始構型Ma=0.72 和Ma=3.15下升力和升阻比對比結果。可見Ma=0.72下優化構型升力性能略有增加,Ma= 3.15下升力性能基本不變,滿足約束,Ma= 0.72和Ma=3.15優化構型升阻比都獲得了較大提升,達到了優化目的。

圖10 寬速域翼型多目標優化Pareto前沿Fig.10 Pareto multi-optimal front of wide-speed range airfoil

圖11 優化前后翼型外形對比Fig.11 Comparison of shapes for initial and optimum airfoils

圖12 Ma=0.72時初始和優化翼型氣動性能對比Fig.12 Comparison of aerodynamic performance of initial and optimal airfoils at Ma=0.72

圖13 Ma=3.15時初始和優化翼型氣動性能對比Fig.13 Comparison of aerodynamic performance of initial and optimal airfoils at Ma=3.15
針對帶擾流板下偏先進三段翼型著陸構型,開展氣動外型和縫道參數優化,以提高其升力性能,計算條件為Ma=0.2、Re=45.6×106,優化問題描述如下:
(7)
通過橢圓方程控制生成多段翼型外型,包括6個外型設計參數,前緣縫翼和后緣襟翼縫道寬度、搭接量和偏角6個設計參數,以及擾流板下偏角,總計13個設計變量。優化采用NSGA-2每代種群有70個個體,總代數為50,總共完成3 500個 構型計算評估。圖14給出優化結束后,優化得到的Pareto前沿,相比初始構型,Pareto前沿優化構型迎角8°和27°下升力系數都有所增加,其中OPT1和OPT3分別為迎角27°和8°升力性能最優構型。沿著Pareto前沿由構型OPT1到OPT3,迎角8°時升力系數越來越大,反之,迎角27°升力系數越來越大,OPT2為綜合迎角8°和迎角27°升力性能的優化構型。
圖15給出了多段翼型典型優化構型外形對比,表3給出了典型優化構型設計點升力系數對比結果。相比基本構型,優化構型升力系數有明顯增加,其中8°迎角升力系數最大增加 9.24%(OPT3),27°迎角升力系數最大增加1.93%(OPT1);從圖16多段翼型典型優化構型升力系數隨迎角變化曲線可以看出,構型OPT1有最大的最大升力系數,線性段增加較小,而構型OPT3線性段升力系數增加較大,最大升力系數增加較小,OPT2構型8°迎角和27°升力系數都有一定的增加,分別增加7.74%和1.43%。

圖14 多段翼型多目標優化Pareto前沿Fig.14 Pareto multi-optimal front of multi-element airfoil

圖15 優化前后多段翼型外形對比Fig.15 Comparison of shapes for baseline and optimum multi-element airfoils

表3 基礎和優化多段翼型升力系數對比Table 3 Comparison of lift coefficient of baseline and optimal multi-element airfoils

圖16 基礎和優化多段翼型升力系數曲線Fig.16 Curves of lift coefficient of baseline and optimal multi-element airfoils
針對典型飛翼布局,開展跨聲速減阻多目標優化設計,優化目的是提高跨聲速設計點升阻比和阻力發散馬赫數,優化問題描述如下:
(8)
式中:y為機翼展向站位;B為機翼展長。
參數化方法采用翼型剖面CST+CATIA二次開發三維成型方法,共計5個優化剖面(圖17紅色曲線),每個剖面18個設計變量,總計90個設計變量。采用代理優化方法,初始采樣點270個,試驗設計方法為LHS,代理模型為Kriging模型,采用Pareto優化時的EI加點方法,總計調用CFD計算1 200次。
圖18給出了優化Pareto前沿,表4給出了典型優化構型與基礎構型阻力系數對比結果,相比基礎構型,Ma=0.72、CL=0.35時優化構型減阻13.17%,Ma=0.78、CL=0.30時優化構型減阻19.56%。從圖19和圖20不同馬赫數下優化構型和基礎構型上翼面壓力云圖對比結果可以看出,優化構型顯著地減小了上翼面激波強度,從而使阻力大幅降低。從圖21不同升力系數下阻力系數隨馬赫數變化結果可以看出,優化構型阻力發散馬赫數獲得了一定的增加,達到了優化目標。

圖17 基于CATIA飛翼布局參數化建模Fig.17 Parametric modeling of flying wing based on CATIA automation

圖18 飛翼布局多目標優化Pareto前沿Fig.18 Pareto multi-optimal front of flying wing

表4 基礎和優化構型氣動性能對比Table 4 Comparison of aerodynamic performance of baseline and optimal configuration
另外,此算例中樣本點總數達到了1 200個,進行代理模型訓練時耗時約2 h,超過了單個樣本CFD數值求解的時間,制約著優化效率的提高,如何提高大樣本點下代理模型訓練效率是ARI_OPT軟件后續的一個改進方向。

圖19 CL=0.35基礎和優化構型壓力系數云圖對比Fig.19 Comparison of pressure coefficient contour of baseline and optimal configuration atCL=0.35

圖20 CL=0.30基礎和優化構型壓力系數云圖對比Fig.20 Comparison of pressure coefficient contour of baseline and optimal configuration atCL=0.30

圖21 基礎和優化多段翼構型阻力隨馬赫數變化曲線Fig.21 Curvers of CD-Ma of baseline and optimal multi-element airfoils
本文簡要介紹了ARI_OPT氣動優化軟件各個模塊的基本原理和方法及單個模塊的驗證結果,給出了ARI_OPT針對數值函數算例的功能驗證結果和寬速域翼型、多段翼型和飛翼布局機翼等典型單目標和多目標氣動外形優化問題的應用算例,表明了其可靠性和適用性。得出以下結論:
1) 數值函數算例驗證結果表明,相比基于遺傳算法優化,代理優化方法計算較少的樣本點就能獲得較優優化結果,效率較高。
2) 經過數值函數驗證和寬速域翼型、多段翼型、飛翼布局機翼等典型氣動外形優化問題應用檢驗,ARI_OPT氣動優化設計工具是有效的且具有較高的優化效率。
3) 大樣本點(1 000以上)情況下進行代理模型訓練時耗時較長,制約著優化效率的提高,如何提高大樣本點下代理模型訓練效率是ARI_OPT軟件后續的一個改進方向。
4) 本文ARI_OPT軟件的多目標優化算例,僅僅是氣動性能的多目標優化,理論完全可以推廣至氣動/結構、氣動/隱身等多學科優化問題,如何進一步應用于多學科優化問題將是下一步工作的主要內容。