基于變可信度代理模型的氣動優化

黃禮鏗，高正紅，張德虎

（西北工業大學 翼型葉柵空氣動力學國家重點試驗室，陜西 西安 710072）

0 引 言

隨著飛行器性能要求的不斷提高，氣動外形設計面臨更大的挑戰，出現了越來越多的高可信度氣動分析方法。高可信度方法通常需要很大的計算量，雖然在過去的幾十年間計算能力發生了巨大的變化，但高可信度模型仍然無法直接應用于優化設計。為了提高優化過程的效率，一種有效的方法是構造基于代理模型的優化方法［1－2］。

代理模型方法最基本的思路是通過構建目標函數的近似模型即代理模型，替代優化設計中復雜的氣動特性分析方法，從而減少優化問題求解的計算量。目前可行的代理模型技術包括數據擬合模型、降階模型、以及啟發式模型，數據擬合模型，通常通過對高可信度模型的樣本數據進行擬合或者退化得到高可信度模型輸入與輸出之間的關系，構造簡化的近似模型，常用的方法有BP神經網絡（Backpropogation Neural Net，BPNN）［3］、徑向基函數（Radial Basis Function，RBF）［4］、Kriging 模型［5］、高斯過程（Gaussian Process，GP）［6］、多項式響應面（Polynomial Response Surface Method，PRSM）［7］等等；降階模型，降階模型的構建依賴于控制方程信息，通常由偏微分方程或者常微分方程描述的系統可以構建降階模型，降階模型不適合于控制方程未知或基于經驗關系的情形，如使用正規正交分解（Proper Orthogonal Decomposition，POD）［8］和模型分析一類的方法得到的降階模型；啟發式模型，也被稱為多可信度、變可信度、變復雜度模型，變可信度模型利用兩個不同可信度的計算工具對相同的物理過程進行數值模擬，一般高可信度的方法精度較高，但需要更大的計算代價，低可信度的方法精度較低，但是需要的計算代價很低，利用大量低可信度數據來提高少量高可信度數據構造的代理模型精度，減少得到合理精度代理模型所需要的計算代價［9］。

對于復雜的多設計變量氣動外形優化問題，構造精確的代理模型需要大量的樣本計算量。如何在確保預測精度的前提下，降低代理模型樣本計算量是代理模型研究的重要方向。變可信度代理模型因引入兩種不同可信度目標函數計算模型的樣本，通過利用計算量較少的低可信度模型樣本，減少計算量較大的高可信度樣本的計算，得到了廣泛的關注［10］。本文采用Co－Kriging方法［11－12］開展變可信度代理模型研究，傳統的變可信度代理模型簡單的使用不同可信度模型的差值或比值實現，但是Co－Kriging方法使用兩組獨立的樣本數據，能夠處理更為復雜的關系，因此用其構建代理模型具有更大的靈活性，并且可以提供預測結果的誤差估計，給樣本補充提供了依據。文中通過典型函數問題研究了如何通過選樣方案設計，構造高效的變可信度代理模型。在此基礎上，針對RAE2822翼型氣動設計問題，基于不同代理模型的翼型優化設計對比證實了該方法的工程實用性。

1 Co－Kriging變可信度代理模型

代理模型可以使用有限數量的數值模擬結果來代表真實模型。為了使用有限數量的樣本數據得到高精度的高可信度模型近似，為了提高代理模型的計算效率和預測精度，人們對代理模型方法的結構和學習方法進行了大量研究，但是代理模型的預測精度與樣本數量直接相關，高精度氣動特性分析樣本的計算代價是制約代理模型精度提高的主要因素。變可信度代理模型可以利用大量低可信度模型信息來提高少量高可信度模型數據構造的代理模型的精度，在大大減少計算量的情況下得到較高的代理模型預測精度。傳統的變可信度代理模型方法利用高低可信度模型之間的差值，利用代理模型方法對高低精度模型的差值進行建模，用差值代理模型對低可信度模型的誤差進行修正，提高其預測精度。由于對同一樣本需分別采用高、低可信度模型進行分析，并且在使用中需要計算低可信度模型，其計算量仍然較大，而且其預測精度很大程度上依賴于低可信度模型的分析精度。Co－Kriging方法利用兩組獨立的高、低可信度樣本數據，充分利用樣本數據挖掘低可信度模型信息，降低了樣本數據和代理模型使用中的計算代價。Co－Kriging方法是Kriging方法的擴展，以使用于雙層甚至是多層模型的代理模型構造，其構造的代理模型的近似估計可以表示如下：

其中Zc（）代表對低可信度樣本數據構造的Kriging代理模型，Zd（）表示對高、低可信度樣本數據差值構造的Kriging代理模型。考慮兩組獨立的高、低可信度數據，高可信度數據有ne個樣本Xe，低可信度模型有nc個樣本Xc，樣本點之間的協方差可以表示為：

其中樣本之間的相關性矩陣R是由一個空間相關函數（Spatial Correlation Function，SCF）決定的：

若Rc（Xc，Xe）表示樣本數據Xc與Xe之間的空間相關性矩陣，則Co－Kriging方法的協方差矩陣即可表示為

因此Co－Kriging方法有兩種不同的空間相關函數，有更多的相關參數（θc、θd、pc、pd以及放大因子ρ）需要估計。由于低可信度模型數據與高可信度模型數據相互獨立，相關參數θc、pc的估計求解與標準Kriging方法的求解方法相同。為了估計θd、pd以及ρ，首先定義：

其中yc（Xe）是在點Xe處低可信度模型的值yc。這樣就可以采用與標準Kriging類似的優化過程來求得對θd、pd以及ρ的估計。Co－Kriging模型的預測估計式可表達如下：

其中＝（fTC－1f）－1fTC－1y，f表示值全為1的ns維列向量

同Kriging方法一樣，Co－Kriging方法也可以給出估計值的均方誤差估計，其表達式為：

2 選樣方案研究

代理模型的近似精度與樣本數量和樣本可信度直接相關，Co－Kriging變可信度代理模型可以使用兩組獨立的高、低可信度樣本數據，合理分配兩種可信度模型的選樣數量對Co－Kriging方法構造的代理模型的近似精度以及計算代價至關重要。為此使用兩個函數算例開展選樣方案研究，針對不同數量的高、低可信度樣本數據構造了三種代理模型進行對比，a）使用Kriging方法和高可信度樣本數據直接構造高可信度模型的近似（Kriging）；b）使用Kriging方法和高、低可信度樣本數據的差值構造差值的近似，再加上低可信度模型的值作為高可信度模型的近似（Delta－Kriging）；c）使用Co－Kriging方法和變可信度模型樣本數據作為高可信度模型的近似（Co－Kriging）。為了進行代理模型的近似誤差對比，利用拉丁超立方方法選取獨立的檢驗數據，誤差選為均方相對誤差（Mean Relative Square Error，MRSE）

及最大相對誤差（Max Relative Error，MRE）作為檢驗代理模型預測誤差的標準

2.1 單變量函數分析

假設高可信度模型為函數fe（x）＝（7x－1）2sin（14x－4），x∈｛0，1）。低可信度模型可以表達為fc（x）＝Afe＋B（x－0.5）－C其中A＝0.4，B＝10，C＝－5。由于問題非常簡單，樣本點的選取是認為指定的，通過測試得知nsh＝4，nsl＝11時，Co－Kriging方法就可以得出與真實函數吻合的近似結果（如圖1）。當高可信度模型的選樣數量nsh＝4即3nv＋1時，不同低可信度選樣數量nsl的近似結果如圖1，可見Co－Kriging對低可信度樣本數量非常敏感，較少的低可信度樣本數量的近似結果反映不了真實函數的基本趨勢，而低可信度模型的樣本數量nsl＝11即達到10nv＋1時，又可以得到與真實函數完全吻合的近似結果。當低可信度模型的選樣數量足夠11nv＋1時，不同高可信度選樣數量nsh的近似結果如圖2，Co－Kriging變可信度代理模型只需要較少的高可信度樣本（3nv＋1）就可以得到與真實函數吻合的近似結果，Kriging模型和Delta－krigng模型則需要更多的高可信度樣本（10nv＋1）。在高可信度模型樣本不足時，Delta－Kriging模型的近似結果要比Kriging模型的更精確。

圖1 低可信度選樣數量對代理模型預測結果的影響Fig.1 The influence of number of low－fidelity samples

圖2 高可信度選樣數量對代理模型預測結果的影響Fig.2 The influence of number of high－fidelity samples

2.2 雙變量函數

高可信度模型為函數fe（X）＝20＋＋x－10cos（2πx1）－10cos（2πx2），x1，x2∈［0，2］。該函數是一個具有多峰值的雙變量函數，構造的低可信度模型則表達為fc（X）＝Afe＋B（x1＋x2－0.5）－C。其中A＝0.6，B＝5，C＝1。由于問題的復雜度增加，樣本的選取采用了拉丁超立方方法（Latin Hypercube Method），給定較多高可信度選樣數量時nsh＝5nv＋1，不同低可信度選樣數量nsl的近似結果如圖3，可見當問題復雜度提高時Co－Kriging方法對低可信度樣本數量要求很高，需要低可信度樣本數量nsl達到35nv＋1，才能得到與真實函數吻合的近似結果。

圖3 低可信度選樣數量對代理模型預測結果的影響Fig.3 The influence of number of low－fidelity samples

在低可信度選樣數量足夠時nsl＝35nv＋1，不同高可信度選樣數量nsh的近似結果如圖4，使用Co－Kriging方法只需要3nv＋1的高可信度樣本就可以得到滿意的近似結果。可見問題的復雜程度決定Co－Kriging方法對低可信度模型的選樣數量，而對高可信度模型需要的樣本數量影響不大，由于低可信度模型的計算時間和資源消耗都很廉價，Co－Kriging方法可以大大減少構造高精度代理模型的代價。

3 氣動優化算例

圖4 高可信度選樣數量對代理模型預測結果的影響Fig.4 The influence of number of high－fidelity samples

由于測試函數體現不出復雜工程優化問題的非線性特性，因此進行一個RAE2822翼型減阻優化算例來進一步證實Co－Kriging方法的特性，計算狀態為M＝0.729以及攻角α＝2.31°。為了進行翼型的氣動特性計算，將全速度勢方程解作為低可信度模型，而把RANS方程的數值解作為高可信度模型。翼型采用Hicks－Henne方法進行參數化，總的設計變量個數為12個。根據函數算例的結果，使用拉丁超立方方法（Latin Hypercube Method）選取了150個低可信度樣本。

圖5為代理模型對升力系數預測結果的對比，基于差值的Delta－Kriging變可信度模型在樣本數量較少時其近似精度要比Kriging模型高，但隨著樣本數量的增加，其預測精度與Kriging模型的越來越一致，Co－Kriging變可信度模型在大量低可信度數據的輔助下使用很少的高可信度樣本就得到了很高的近似精度，這說明基于高、低可信度差值的Kriging模型可以提取的高低可信度模型間的關系信息很有限，而Co－Kriging方法可以充分挖掘低可信度模型的信息，大大提高了代理模型的效率，在基于高精度復雜數值模擬的優化設計時帶來極大的好處。

利用遺傳算法進行了基于代理模型和RANS方程數值模擬的翼型減阻優化，優化中對翼型的最大厚度以及升力系數進行約束。遺傳算法的種群規模為100，優化進行了30代，優化模型為：

圖5 代理模型預測誤差對比Fig.5 Comparison of lift coefficient approximation

CD、CL、CM分別為翼型的阻力系數、升力系數和俯仰力矩系數。目標函數的設置是期望升力系數盡量靠近0.743，同時最小化阻力，并且對控制翼型的最大厚度max和俯仰力矩系數進行約束。基于代理模型的優化中對每一代的最優解進行驗證并更新代理模型，優化后的翼型形狀和壓力分布如圖6，可見基于Co－Kriging變可信度模型的優化結果與基于RANS方程數值模擬（即高可信度模型）的優化結果基本一致，優化后翼型的激波得到了極大地削弱。而基于Kriging模型的優化結果雖然也削弱了激波，但是比另外兩種方法的要強一些。這說明利用足夠的低可信度樣本數據（大于10nv＋1）和少量高可信度樣本數據構造的Co－Kriging變可信度代理模型可以給出更好的近似結果和正確的優化方向，同時大大減少了優化過程的計算代價，可以滿足工程設計優化的需要并提高優化設計過程的效率。

圖6 優化結果翼型和壓力分布Fig.6 The optimized airfoil and pressure coefficient distribution

4 結 論

本文利用Co－Kriging方法構造了一種高效的變可信度代理模型，可以利用大量低可信度模型數據來增強模型的預測能力，提高代理模型的效率，在高可信度數據難以獲取時可以帶來極大的好處。代理模型的近似精度與樣本數量和樣本可信度直接相關，使用單變量函數和雙變量函數算例對代理模型進行的選樣方案研究表明，Co－Kriging方法比Kriging方法更加高效，只要超過3nv＋1就可以得到合理近似精度，其所需低可信度模型數據則要隨著問題的復雜程度不斷增加。用遺傳算法進行的基于代理模型的翼型減阻優化結果證實Co－Kriging變可信度代理模型更加高效，其優化結果和使用高可信度模型直接進行優化的結果基本一致，能夠滿足工程設計優化的需要。

［1］IKJIN LEE，CHOI K K，ZHAO L.Sampling－based RBDO using the dynamic Kriging（D－Kriging）method and stochastic sensitivity analysis［A］／／13th AIAA／ISSMO Multidisciplinary Analysis Optimization Conference［C］.Fort Worth，Texas，2010.AIAA 2010－9036.

［2］ALLAIRE D，WILLCOX K.Surrogate modeling for uncertainty assessment with application to aviation environmental system models［J］.AIAAJournal，2010，48（8）：1791－1803.

［3］閻平凡，張長水.人工神經網絡與模擬進化計算［M］.清華大學出版社，2005.

［4］BUHMAN N，MARTIN D.Radial basis functions：theory and implementations［M］.Cambridge university press，2003.

［5］SACKS J，WELCH W J，MICHELL T J，WYNN H P.Design and analysis of computer experiments［J］.StatisticalScience，1989，4（4）：409－435.

［6］GIBBS M N.Bayesian Gaussian processes for regression and classification［D］.［PhD Thesis］.University of Cambridge，1997.

［7］MYERS R，MONTGOMERY D.Response surface methodology［M］.New York：John Wiley ＆ Sons，Inc.，1995.

［8］ZIMMERMANNL R，GORTZ S.Non－linear reduced order models for steady aerodynamics［J］.ProcediaComputerScience，2010，1：165－174.

［9］LEIFSSON L，KOZIEL S.Multi－fidelity design optimization of transonic airfoils using physics－based surrogate modeling and shape－preserving response prediction［J］.JournalofComputationalScience，2010（1）：98－106.

［10］DANIEL M J，GEOFFREY T P，WILLIAM N D，et al.Robust multi－fidelity aerodynamic design optimization using surrogate models［R］.AIAA 2008－6052.

［11］FORRESTER A I J，So′BESTER A，KEANE A J，et al.Multi－fidelity optimization via surrogate modeling［J］.ProceedingsoftheRoyalSocietyA.，2007，463（11）：3251－3269.

［12］MYERS D E.Matrix formulation of Co－Kriging［J］.MathematicalGeology，1982，14（3）：249－257.

［13］SHAWN E G，HAROLD K，DON E B II.Comparison of three surrogate modeling techniques：datascaper，Kriging，and second order regression［R］.AIAA 2006－7048.

［14］MARTIN J D，SIMPSON T W.Use of Kriging models to approximate deterministic computer models［J］.AIAAJournal，2005，43（4）：853－863.

［15］張德虎，高正紅，李焦贊，等.基于雙層代理模型的無人機氣動隱身綜合設計［J］.空氣動力學學報，2013，31（3）：394－400.