侯成義
(西北工業大學國防科技研究院,710072西安,peggy.cao@163.com)
遺傳算法、粒子群算法等隨機優化算法同基于梯度的優化算法相比具有很強的全局性和魯棒性,然而這些隨機優化算法進行需要用CFD、FEM分析等耗時的目標值評估的工程優化時,存在著計算量大的問題,嚴重制約了在工程中的應用.用代理模型代替費時的精確模型計算進行優化的方法得到了發展[1-2],常用的代理模型有多項式響應面模型、人工神經網絡模型、Kriging模型等.首先由一定數量的樣本構造代理模型,用優化算法搜索代理模型的最優解并進行精確模型校驗,進而更新校驗點到代理模型直至滿足收斂要求.由于代理模型的計算量比精確模型小得多,該方法可大幅度提高優化效率,但是其優化精度嚴重依賴于代理模型預測的能力,容易陷入局部最優[3].M.Schonlau[4]提出的EGO(Efficient Global Optimization)算法在選取校驗點時綜合考慮了代理模型的預測值和預測精度,有效避免了陷入局部最優的風險,該方法在工程優化中得到了應用[5-6].
本文采用拉丁超立方實驗設計方法產生樣本點,用EI(Expected Improvement)最優策略綜合評價由Kriging代理模型獲取的預測適應值及預測精度,并結合響應最優策略,用粒子群算法搜索最優校驗點,結合Hicks-Henne型函數翼型參數化方法與雷諾平均N-S方程流場求解器進行了翼型優化設計.
選取合理分布的樣本點對構造代理模型具有重要的影響,一般要求有限的樣本應該盡可能全面的反映設計空間的特性.常用的實驗設計方法有全析因設計、正交設計、均勻設計等,本文采用了拉丁超立方設計方法產生樣本點.
粒子群優化算法(Particle Swarm Optimization)[7]是一種模擬鳥群覓食行為的群體智能算法,該算法實現簡單,參數設置方便,收斂速度快,能有效解決復雜優化問題.該算法以“粒子”作為個體,按照設定的規則進行群體的協作從而進行有效的搜索.首先隨機初始化粒子種群位置和初速度,每個粒子通過跟蹤自身所找的個體最優粒子pbest以及整個種群找到的全局最優粒子gbest,在粒子當前速度、pbest和gbest的位置的共同引導下,粒子在下一個時刻將飛行到新的位置.粒子的速度和位置更新公式分別為

式中:d=1,2,…,n(n為搜索空間維數);i=1,2,…,m(m為種群規模);t為當前進化代數;分別為對應粒子的歷史最優位置和全局最優位置;c1,c2分別為學習因子;R1,R2分別為(0,1)的隨機數;ω為慣性權重因子,慣性權重因子決定了粒子先前速度對當前速度的影響程度.Shi Y采用隨迭代進行慣性權重因子減小的方法,在優化初期提高全局搜索能力,而在后期加速收斂.
Kriging代理模型[8]起源于地理空間統計學,是一種估計方差最小的無偏估計模型,通過相關函數的作用,具有局部估計的特點,它可以較好的預估未知點處函數值的分布情況.Kriging模型假設目標函數值與設計變量之間的真實關系可以寫為

式中:f(x)為回歸模型,是對設計空間的全局近似,為確定性部分,可以一個常數β表示;Z(x)為均值為0,方差為σ2的統計過程,表示對全局近似的背離,2個插值點的協方差為

式中R為點x(i)和x(j)的相關函數,可用高斯函數表示為

未知點x0處的預測值^y(x0)和方差估計值^σ通過下列形式給出:

可以通過極大似然估計確定相關參數θk為

選擇代理模型預測值的最優點作為校驗點的方法稱為響應最優策略.但是代理模型的預測精度總是有限的,如果代理模型在真實的最優解處預測精度較差,就有可能在選擇校驗點時總是無法搜索到該區域,使得優化陷入局部最優.因此在選擇校驗點時有必要綜合考慮Kriging模型獲取的預測適應值及預測精度.本文采用了EI(Expected Improvement)方法,根據下式確定個體的EI值,EI值越大說明該個體越需要校正.

式中:fmin為樣本點中最小適應值;^y為預測適應值為預測標準差;Φ,Ψ分別為標準正態分布和標準正態分布密度函數.當個體預測值越小,預測標準差越大,則個體的EI值越大,個體越需要校正.
本文在選擇校驗點時,一部分校驗點用響應最優策略選擇,另一部分采用EI最優策略,2種方法交替進行.優化流程為:
1)由拉丁超立方設計方法在設計空間內生成一定數量的樣本.
2)由樣本點構建Kriging代理模型.
3)確定是采用響應最優策略還是EI最優策略,基于粒子群算法優化選擇最優響應或最優EI校驗點.
4)判斷是否滿足收斂條件,如果滿足則退出,不滿足則執行步驟5).
5)將校驗點加入樣本,執行步驟2).
翼型參數化方法采用的是Hicks-Henne型函數法[9],新翼型由基本翼型的上表面或下表面加擾動構成,表達為

式中n為控制上表面或下表面變量的個數,本文選取n為7,即共計14個設計變量.
本文以RAE2822翼型為初始翼型進行翼型氣動優化設計研究.采用雷諾平均N-S方程求解繞翼型流場,湍流模型為k-ω SST模型.設計狀態為巡航馬赫數0.73,雷諾數6.5E6,設計升力系數0.74,約束條件為最大相對厚度不減,低頭力矩特性不差于初始翼型,優化目標為降低阻力系數.
由拉丁超立方設計方法產生120個翼型作為代理模型初始樣本點,校驗次數為240,EI最優策略進行2次后響應最優策略進行1次,如此交替進行,即分別進行160次和80次,PSO的種群規模為30,迭代步數為60.為了比較本文優化方法,還采用響應最優法和粒子群算法進行翼型優化,其中響應最優法初始樣本數為120,校驗次數為240,而粒子群算法種群規模為30,迭代步數為40,設計結果如表1所示.單純的響應最優法優化陷入局部最優,其他2種方法得到了相似的優化結果,阻力降低了22.4%,粒子群算法用了1 200次流場計算,而本文方法僅用360次.圖1和圖2分別為本文方法優化前后翼型形狀和壓力分布,與初始翼型相比,設計翼型翼型完全消除了激波,阻力系數降低了0.003 4,而翼型最大相對厚度沒有降低.設計翼型整體彎度降低,而通過增加后加載來保持升力系數.同時設計翼型增加了前加載和頭部上表面吸力峰值,一方面增加了翼型升力,另一方面抵消了后加載所附加的低頭力矩,使設計翼型保持著良好的低頭力矩特性.

表1 不同優化方法計算次數及優化結果對比

圖1 設計翼型與初始翼型形狀
雖然構造Kriging模型以及預測氣動力需要耗費一定時間,但是與調用氣動力計算的耗時相比要小得多,上文算例優化耗時約比粒子群算法減少68%.

圖2 設計翼型與初始翼型壓力分布
1)本文構建了將Kriging代理模型和粒子群算法相結合的高效優化系統,采用粒子群算法搜索由EI最優策略和響應最優策略交替選擇的最優校驗點,前者綜合考慮代理模型的預測值及預測精度,可以有效避免陷入局部最優,后者可以促使快速收斂.
2)翼型優化設計算例表明本文方法克服了響應最優法容易陷入局部最優的缺點,保證優化精度,同時又有效地減少了計算量,大大提高了翼型氣動優化設計的工程實用性.
[1]DUVIGNEAU R,VISONNEAU M.Hybrid genetic algorithms and neural networks for fast CFD-based design[R].Atlanta:AIAA 2002-5465,2002.
[2]RAJAGOPAL S,GANGULI R.Multidisciplinary design optimization of an UAV wing using kriging based multi-objective genetic algorithm[R].AIAA 2009-2219,2009.
[3]JONES D L.A taxonomy of global optimization methods based on response surfaces[J].Journal of Global Optimization,2001,21(4):345-383.
[4]SCHONLAU M.Computer experiments and global optimization[D].Ontario Canada:Waterloo,1997.
[5]王紅濤,竺曉程,杜朝輝.基于Kriging代理模型的改進EGO算法研究[J].工程設計學報,2009,16(4): 266-270.
[6]蘇偉,高正紅,夏露.一種代理遺傳算法及其在氣動優化設計中的應用[J].西北工業大學學報,2008,26(3):303-307.
[7]KENNEDY J,EBERHART R C.Particle swarm optimization[C]//Proceedings of IEEE International Conference on Neural Networks.Washington,DC:IEEE Service Center,1995:1942-1948.
[8]CRESSIE N A C.Statistics for spatial data[M].New York:John Wiley&Sons,1993.
[9]HICKS R M,HENNE P A.Wing design by numerical optimization[J].Journal of Aircraft,1978,15(7): 407-413.