基于PCA-HicksHenne 方法的幾何不確定性穩健優化設計

張威, 王強, 路嘉晨, 閻超,*

(1. 北京航空航天大學 航空科學與工程學院, 北京 100191; 2. 中國航天空氣動力技術研究院, 北京 100074)

在飛行器的整個壽命周期中存在著各種各樣的不確定性,而考慮不確定因素的氣動優化研究在飛行器氣動外形設計中扮演著至關重要的角色。 考慮不確定性因素的優化設計主要包括穩健優化設計[1]及可靠性優化設計[2]。 穩健設計的優化目標是為了在性能提升的同時降低目標函數對不確定因素的敏感度,而可靠性設計則是為了降低故障概率。 本文的關注點主要在于穩健優化設計研究。 現有的穩健優化設計主要著重于針對環境擾動下的優化設計研究,而對于幾何擾動的關注則相對較少。

對于幾何不確定性的考量,主要有2 方面:①局部位置的凸起凹陷,主要是由于工作狀態下的磨損、侵蝕等導致的。 由于其大小、位置、形狀、分布形式等較為依賴對大量實驗數據的分析,而實驗數據難以獲取使得該類幾何不確定性量化也較難實現。 ②加工過程中存在不確定性。 現有對該類幾何不確定性的研究中,由于認知不足,或是根據經驗直接在如翼型前緣等關鍵位置添加滿足正態分布的幾何擾動[3],然后進行不確定性量化;或是采用區間法[4]直接進行幾何不確定性量化。 這些方法的使用,無論是對正態擾動添加位置的選擇與確定,還是區間法的直接應用,經驗性都較強,忽視了對于全局的考量,及潛在幾何模態的影響。

本文以RAE2822 翼型為例,采用主成分分析(principal component analysis, PCA)方法[5-6],對整個參數化過程進行分析,增強對全局范圍內主要幾何變形模態的認知,了解對應參數的分布規律。 應用敏感度分析方法,針對不同目標函數(以升力系數、阻力系數、翼型厚度等為例),找出主要變形控制模態,作為幾何擾動來源。 同時,進行了考慮幾何不確定性的穩健優化設計研究。

1 PCA-HicksHenne 方法

1.1 Hicks-Henne 參數化方法

Hicks-Henne 方法[7]是常用的參數化方法,本文以此對RAE2822 翼型進行參數化,具體的幾何外形表達式如下:

式中:y為變形后翼型的上(下)表面函數;yorigin為初始翼型的上(下)表面函數;fk(x)為擾動型函數,ck為對應參數。

原始的Hicks-Henne 方法中的擾動型函數定義如下:

研究表明[8-9],該方法存在后緣擾動不足的問題。 對此,本文將翼型后緣擾動函數做了式(3)的修改。 其導函數(式(4))在后緣處取值不為0,可以提供有效的后緣擾動。

1.2 PCA 方法

PCA 是一種常用的數據降維方法。 通過對高維特征數據進行降維,僅保留重要特征,從而實現數據處理效率的提升。

基于最近重構性或最大可分性[10],可將PCA的優化目標寫為

式中:X為翼型樣本幾何向量組成的矩陣;W為協方差矩陣XXT的特征向量矩陣;I為單位矩陣。

采用拉格朗日乘子法可得

式中:wi為矩陣W第i列的向量,也是XXT的特征向量;λi為對應的特征值。

因此,僅需對協方差矩陣XXT進行特征值分解,得到特征值如式(7)所示。 一般對所得到的特征值按照由大到小順序進行排序,特征值越大,表明對應的“成分”的貢獻越大。 隨后再選取前d′個特征值對應的特征向量W*= (w1,w2,…,wd′),即為主成分分析的解。 特征值由大到小如下:

1.3 PCA-HicksHenne 參數化方法

采用PCA 方法,對經由Hicks-Henne 方法得到的大量翼型進行分析,提取目標翼型變形的主要模態,并以此為基礎對翼型進行變形。 該方法可以提高對翼型變形內在幾何模態的認知,進一步降低參數化所需參數維度,提高效率,這就是本文所采用的PCA-HicksHenne 方法。

以RAE2822 翼型為例,設定Hicks-Henne 方法中的型函數個數為21(上下翼面合計變量維度為42)。 較多的型函數個數可以賦予參數化方法更好的設計空間表達能力。 在給定設計空間內隨機選取100 000 個樣本翼型數據,如圖1 所示,并對其進行PCA 分析。 給出前20 個主要變形模態所對應的能量分布,如圖2 和圖3 所示。 可以看出,前8 個變形模態表達了98% 以上的變形信息,這表明僅采用前8 個模態就能夠很好地表征原設計空間中的翼型變形。

圖1 樣本翼型及RAE2822 翼型示意圖Fig.1 Schematic of sample airfoils and RAE2822 airfoil

圖2 主要幾何變形模態對應模態能量分布Fig.2 Energy amplitude of main geometric transformation modes

圖3 主要幾何變形模態對應累積能量分布Fig.3 Accumulated energy of main geometric transformation modes

圖4 給出經PCA 分析得到的前8 個變形模態。 模態1 和模態2 分別為翼型的厚度變形模態及彎度變形模態,模態3 和模態4 分別為上下翼面的最大厚度位置的軸向位移模態,模態5 ～模態8 則分別為翼面的擠壓模態。

圖4 主要幾何變形模態示意圖Fig.4 Schematic of main geometric transformation modes

本文以這些模態為設計變形模態,通過線性疊加到基準翼型RAE2822 上,從而實現翼型的變形。 同時,觀察本節隨機選取的100 000 個翼型對應在這8 個模態構建的參數空間中的參數分布情況,如圖5 所示。 觀察可知,所有模態的參數分布與N(0,0. 012)的正態分布幾乎完全一致,因此,可以認為第3 節中所用到的不確定模態也滿足正態分布。 同時將各模態參數的設計范圍定義為[ -0.03,0.03],以滿足99.7% 以上的參數變化范圍。

圖5 主要幾何變形模態對應參數分布情況Fig.5 Parameters distribution of main geometric transformation modes

2 重要參數識別

本文的主要目的是為了發展基于PCA-Hicks-Henne 參數化方法的幾何不確定性穩健優化設計方法。 因此,為了選定特定參數模態作為幾何擾動模態,需要對各模態針對目標函數(如升力系數、阻力系數、翼型厚度等)的參數敏感度進行分析。

Sobol 指數方法是一種常用的全局參數敏感度分析方法[11-12]。 通過Sobol 指數方法,可以得到各參數對于目標函數的影響,主要包括參數的主指標、參數之間的交互指標及參數的總指標。參數的主指標是指完全由該參數產生的影響;交互指標是指參數之間相互作用,共同產生的影響,通常用來判斷參數之間是否存在關聯;總指標則是由參數主指標與交互指標共同組成,通常表示該參數產生的所有影響的總和。

本文以RAE2822 翼型的升力系數、阻力系數、翼型厚度為目標,研究第1 節中所得的各變形模態在設計空間內的參數敏感度,結果分別如圖6 ～圖8 所示,計算狀態如表1 所示。 可以看出,對于升力而言,模態2 是主要影響模態,這是由于跨聲速條件下翼型彎度是升力產生的主要來源;對于阻力而言,模態3 是主要影響模態,這是由于在模態3 下,上翼面的最大厚度位置在軸向上發生變化,導致上翼面激波產生位移,激波強度、激波阻力發生變化,從而阻力發生變化;對于翼型厚度而言,模態1 即厚度變形模態是主要影響模態。

圖6 升力系數的參數敏感度分析Fig.6 Parametric sensitivity analysis of lift coefficient

圖7 阻力系數的參數敏感度分析Fig.7 Parametric sensitivity analysis of drag coefficient

圖8 翼型厚度的參數敏感度分析Fig.8 Parametric sensitivity analysis of airfoil thickness

表1 計算狀態Table 1 Computational condition

因此,模態1 ～模態3 是翼型性能的主要影響模態。 這些模態的微小變化,會對翼型性能產生較大的影響。 因此,本文選取這3 個模態作為主要幾何擾動模態,進行穩健優化設計。

3 幾何擾動定義

在進行穩健優化設計之前,還需對幾何擾動模態的參數進行定義。 根據第1 節研究結果,可以假定3 個擾動模態仍滿足正態分布,其均值為0。 為更好展現本文方法的性能,分別設定擾動模態參數的標準差為設計參數的0.2 ～0.4 倍,在基準RAE2822 翼型上,探索這些不同標準差假設下幾何擾動的影響。

混沌多項式方法是一種廣泛使用的不確定度量化分析方法[13-14]。 本文選用二階混沌多項式方法來量化不同擾動量下的翼型升力系數、阻力系數、翼型厚度的不確定度,結果如表2 ～表4 所示。 可以看出,隨著幾何擾動量的增大,各個目標函數的不確定性都在增強。 同時也注意到,0.4 倍擾動量下,99.7% 置信區間的翼型厚度變化范圍達到10 倍風洞模型允許公差[15],遠遠偏離實際允許情況,不予考慮。

表2 不同擾動量下升力系數的不確定度對比Table 2 Comparisons of lift coefficient uncertainty under different perturbations

表3 不同擾動量下阻力系數的不確定度對比Table 3 Comparisons of drag coefficient uncertainty under different perturbations

表4 不同擾動量下翼型厚度的不確定度對比Table 4 Comparisons of airfoil thickness uncertainty under different perturbations

因此,在第4 節的穩健優化設計中,定義模態1 ～模態3 為幾何擾動模態,符合正態分布,均值為0,標準差為0. 3 倍的設計參數標準差(即0.003)。 擾動下的翼型變化情況如圖9 所示。

圖9 幾何擾動下的翼型變化情況Fig.9 Airfoil transformation under geometric perturbation

4 穩健優化設計

4.1 優化流程及計算說明

圖10 給出了幾何擾動下的穩健優化設計流程。 在無擾動的常規優化流程中,先對初始外形進行PCA-HicksHenne 參數化設計;其次對參數進行試驗設計,得到采樣點;再次進行配套網格變形及計算流體力學(computational fluid dynamic,CFD)解算,得到采樣點的響應;隨后在采樣點上構建代理模型,并進行優化算法尋優;最終得到優化外形。 而在考慮幾何擾動的優化流程中,還需要在參數化完成之后進行敏感度分析,尋找幾何擾動參數,并根據CFD 結果進行擾動下的不確定度量化,隨后再構建代理模型、進行算法尋優。

圖10 幾何擾動下穩健優化設計流程Fig.10 Flow chart of robust optimization design under geometric perturbation

本文選取RAE2822 翼型作為基礎翼型,在自研CFD 計算平臺MI-CFD[16]上,構建克里金模型[17],并利用多島遺傳算法[18]進行尋優。 計算狀態如表1 所示。 為展示穩健優化的性能,將其與確定優化的結果進行對比,具體見4.2 節和4.3 節。

4.2 案例1

本節案例是單目標優化案例。 確定優化目標及對應的穩健優化目標定義分別為

maxCl

式中:Cl和θ分別為翼型的升力系數和翼型厚度;μ和σ分別為對應下標變量在擾動下的均值和標準差。

優化結果如表5 所示。 本文中無擾動翼型的升阻力性能與擾動下平均性能接近。 可以看出,確定優化外形在升力均值方面提升很多,但同時也不可避免地伴隨著升力標準差的增長。 性能提升的同時,不確定度也在增大,這是不利的。 而對于穩健優化外形,因為在優化過程中已經考慮了對不確定性因素的約束,所以在其性能提升的同時,不確定度基本無變化。

表5 案例1 優化結果對比Table 5 Comparisons of optimal results in Case 1

4.3 案例2

本節案例是多目標優化案例。 確定優化目標及穩健優化目標定義分別為

式中:Cd為翼型的阻力系數。

圖11 和圖12 給出2 個優化過程的Pareto前沿面,其中穩健優化結果僅給出平均升力與平均阻力的對應關系。 在圖11 和圖12 中分別選出3 個設計點,并將結果在表6 和表7 中列出。

圖11 案例2 確定優化Pareto 前沿面Fig.11 Pareto frontier of deterministic optimization in Case 2

圖12 案例2 穩健優化Pareto 前沿面Fig.12 Pareto frontier of robust optimization in Case 2

表6 案例2 確定優化結果Table 6 Solutions of deterministic optimization in Case 2

表7 案例2 穩健優化結果Table 7 Solutions of robust optimization in Case 2

如表6 和表7 所示,所有得到的優化翼型均滿足翼型厚度約束。 其中,確定優化得到的3 個優化翼型在平均升阻力性能得到提升的同時,對阻力系數的不確定度恰好也有減弱的效果,而對升力系數的不確定度則恰恰相反。 3 個翼型的升力系數標準差均有不同幅度的增長,最高達到了近200%。 對于飛行器穩健性的設計要求,這是不利的。 而穩健優化得到的3 個優化翼型,則在目標性能相斥的情況下,兼顧了升阻力的均值及標準差等各方面的性能要求。 在提升某些方面性能的同時,也能保證其他性能的不惡化。

由表6 和表7 確定優化得到的翼型③(以下為確定優化翼型)和穩健優化得到的翼型⑤(以下稱穩健優化翼型)的平均性能相近,對其作對比分析。

圖13 給出了RAE2822 翼型、確定優化翼型及穩健優化翼型的形狀對比,并在圖14 ～圖16中給出了3 個翼型及其對應擾動翼型的壓力系數Cp的分布情況。

圖13 優化翼型對比Fig.13 Comparisons of optimal airfoils

圖14 RAE2822 翼型及擾動翼型的壓力分布Fig.14 Surface pressure distribution of RAE2822 airfoil and its corresponding airfoils under perturbation

圖15 確定優化翼型及擾動翼型的壓力分布Fig.15 Surface pressure distribution of deterministic optimal airfoil and its corresponding airfoils under perturbation

如圖13 ～圖16 所示,在設計狀態下,2 個優化外形在下翼面后緣附近都存在由于彎度增加的后緣加載,升力因此增加。 此外,穩健外形在上翼面前緣附近曲率更大,產生了一個較小的吸力峰,因此其升力相對更大。 而在阻力方面,2 個優化外形的上翼面相對RAE2822 翼型都更加平緩,確定優化外形分解為2 個小激波,穩健優化外形激波強度減弱,二者阻力均降低。

圖16 穩健優化翼型及擾動翼型的壓力分布Fig.16 Surface pressure distribution of robust optimal airfoil and its corresponding airfoils under perturbation

而在幾何擾動下,穩健優化外形的激波強度和位置的變化與RAE2822 翼型相似,因此兩者的升力系數變化幾乎一致,標準差相當。 另外穩健外形的激波強度相對較小,因此其阻力系數的標準差也相對較小。 對于確定優化外形,其激波位置在擾動下變化劇烈,由此導致升力變化較大。而其激波強度較弱,變化也不大,因此阻力變化不大。

5 結 論

本文構建了基于主成分分析的PCA-HicksHenne參數化方法,并在此基礎上進行了考慮幾何不確定性的穩健優化研究,結論如下:

1) 對原始Hicks-Henne 方法進行了修改,使其在翼型后緣處也能提供有效擾動,提高了參數化方法的空間表達能力。

2) 構建了基于主成分分析的PCA-HicksHenne參數化方法,揭示了翼型的主要幾何變形模態,驗證了對應參數的分布規律,同時實現了參數化方法的降維,提高了優化效率。

3) 利用參數敏感度分析方法,發現翼型厚度變形模態、彎度變形模態及上翼面最大厚度位置的軸向位移模態是翼型性能的主要影響模態,必須予以重視。

4) 發展了基于PCA-HicksHenne 方法幾何不確定性穩健優化設計方法,拓展了對幾何不確定性穩健優化的研究,優化得到的目標翼型在性能提升的同時,穩健性也得以保證。