秦海勤，趙杰，任立坤，李邊疆

海軍航空大學青島校區(qū) 力學工程系，青島 266000

航空發(fā)動機運行過程的可靠性對于飛行安全起決定性作用［1］。發(fā)動機長期工作在高溫、高壓、高轉速等惡劣環(huán)境下，氣路部件不可避免地會出現(xiàn)葉尖磨損、葉片腐蝕等現(xiàn)象［2-3］，進而導致氣路部件性能出現(xiàn)退化，使得超溫、失速、喘振等故障時有發(fā)生［4］。因此，為提高工作可靠性，有必要開展航空發(fā)動機的性能退化評估研究。

航空發(fā)動機的性能退化評估主要通過估計發(fā)動機健康因子（各部件的效率和流通能力）而實現(xiàn)［3］。目前發(fā)動機健康因子估計方法主要分為基于物理模型的方法和基于數(shù)據(jù)驅動的方法［1］。其中基于數(shù)據(jù)驅動的方法通過訓練大量“測量參數(shù)-健康因子”的映射樣本，使模型具備由測量參數(shù)計算發(fā)動機性能衰退趨勢的能力，該方法在退化數(shù)據(jù)集［5］上展現(xiàn)了出色的預測能力［6-8］。但實際上發(fā)動機的性能健康因子無法測量，因此基于數(shù)據(jù)驅動的方法難以應用于實際發(fā)動機性能退化評估中。而基于物理模型的方法則不需要發(fā)動機的退化先驗知識，而是通過發(fā)動機性能計算模型和發(fā)動機傳感器數(shù)據(jù)的偏差情況對發(fā)動機健康因子進行估計［9-11］，本質(zhì)是一個參數(shù)尋優(yōu)問題。

由于航空發(fā)動機是一個復雜的非線性系統(tǒng)，傳統(tǒng)的遺傳算法在進行發(fā)動機健康因子求解時存在精度不高、過早收斂等問題。差分進化（Differential Evolution，DE）算法具有結構簡單、控制參數(shù)少、計算精度高等優(yōu)點，在最近20 多年內(nèi)取得了迅猛的發(fā)展［12］。2013 年，Tanabe 和Fukunaga［13］提出了一種基于成功歷史的參數(shù)自適應差分進化算法SHADE（Success-History based Adaptive DE）算法。之后，又引入了線性種群縮減策略，得到了L-SHADE（Linear-SHADE）算法［14］，該算法在CEC2014 實參單目標優(yōu)化競賽中非混合算法獲得了最好的結果。L-SHADE 算法在眾多復雜問題上都具有很強的競爭力，因此被應用于許多實際工程問題中［15-16］。然而該算法的線性種群縮小策略在初期種群快速縮減中導致種群的多樣性受到影響，易陷入局部最優(yōu)，在面對復雜問題時效果不佳。

基于差分進化算法在參數(shù)尋優(yōu)問題求解中的有效性，同時針對現(xiàn)有L-SHADE 算法的不足，本文提出了一種改進的L-SHADE 算法，通過引入一種非線性種群縮減策略來解決種群快速縮減的問題，通過調(diào)整迭代不同階段的變異能力增強算法的魯棒性。利用30 個經(jīng)典基準函數(shù)驗證改進算法的有效性。在此基礎上將改進算法應用于某型發(fā)動機性能退化評估中。

1 算法理論

1.1 發(fā)動機退化估計

研究對象為某型小涵道比渦扇發(fā)動機，主要由風扇、壓氣機、燃燒室、高壓渦輪和低壓渦輪等部件組成，各個部件間存在著復雜的耦合關系。研究表明［17］，發(fā)動機的性能退化主要表現(xiàn)在旋轉部件氣路特性（流通能力和效率特性）的變化上。因此，通常定義發(fā)動機各旋轉部件流量和效率的健康因子分別為

式中：W、η分別為部件的換算流量和效率；下標d、c 分別表示退化發(fā)動機和未退化發(fā)動機的數(shù)據(jù)。

考慮部件性能退化后，發(fā)動機的非線性性能模型可表示為

式中：u是發(fā)動機輸入向量，z是發(fā)動機測量參數(shù)向量，θ為發(fā)動機各旋轉部件性能健康因子（流通能力SW 和效率特性SE）組成的向量，t為時間變量。

由于部件的健康因子θ無法直接測量或直接由傳感器測量參數(shù)計算，因此航空發(fā)動機氣路部件性能退化評估本質(zhì)是一個不可測參數(shù)的估計問題。通常根據(jù)測量參數(shù)逆向求解式（3）對部件健康因子進行估計：

式（4）的求解需要滿足測量參數(shù)的數(shù)量p不小于發(fā)動機健康因子數(shù)量n。受安裝空間和重量等因素影響，實際發(fā)動機傳感器數(shù)量有限，很難滿足測量參數(shù)數(shù)量不小于發(fā)動機健康因子數(shù)量這一要求，從而使式（4）求解成為欠定問題，理論上無法求解。

針對上述欠定問題，Stamatis 等［18］提出MOPA（Multiple Operating Points Analysis）方法，假設單個飛行循環(huán)過程中各部件健康因子不變，通過選取多個穩(wěn)態(tài)工作點拓展健康因子的求解方程組，以解決方程組求解中存在的測量參數(shù)不足的問題。該方法要求所選取的穩(wěn)態(tài)工作點間具有一定的獨立性，緩解健康因子求解中存在的多重共線性問題。但是Diakunchak［19］指出，發(fā)動機部件的流量和效率的退化程度會隨發(fā)動機工作狀態(tài)的改變而變化，這與MOPA 方法各穩(wěn)態(tài)工作點健康因子相等的假設相矛盾，為減小MOPA 方法假設帶來的誤差，多工作點選取時要求各工作點差異盡可能小。

通過對所研究發(fā)動機多架次的實際飛參數(shù)據(jù)分析，發(fā)現(xiàn)該型發(fā)動機穩(wěn)態(tài)工作點多集中于一個較小的高壓轉子轉速區(qū)間范圍（86%～92%）。由于該轉速范圍較小，因此近似認為在該區(qū)域內(nèi)部件性能健康因子相等。同時在多工作點選取中盡可能選擇單個飛行架次中該轉速范圍內(nèi)工作差異大的點，以緩解方程組求解中的多重共線性問題，從而獲得更好的方程適定性。

因此在單個區(qū)域范圍內(nèi)采用MOPA 方法，對式（4）進行拓展后，得

式中：q為MOPA 選擇的工作點數(shù)量，滿足pq≥n。

將式（5）的方程組求解問題轉化為優(yōu)化問題，而轉化為優(yōu)化問題的關鍵是確定適應度函數(shù)。綜合考慮測量參數(shù)的誤差范圍、測量參數(shù)的關系以及建模誤差和測量參數(shù)的比值，合理確定測量參數(shù)在適應度函數(shù)中的權值ωj，則得到的適應度函數(shù)為

式中：Ym為實際測量參數(shù)值；Ycal為計算參數(shù)值；a為工作點編號，a=1，2，…，q；b為參數(shù)編號，b=1，2，…，p。通過最小化適應度函數(shù)OF 即可估計得到健康因子。

1.2 改進L-SHADE 算法

1.2.1 標準L-SHADE 算法

L-SHADE 算法是基于線性種群規(guī)?？s減和成功歷史的參數(shù)自適應差分進化算法，其對差分進化算法的3 個控制參數(shù)（比例因子F、交叉率CR、種群規(guī)模N）進行了自適應調(diào)整。作為DE算法變體中最有潛力的算法之一，L-SHADE 算法的改進成為了近年來的研究熱點［20-21］。

標準L-SHADE 算法的主要步驟分為初始化、變異、交叉、選擇，在進化過程中通過參數(shù)自適應和線性種群縮減策略對控制參數(shù)進行調(diào)整。

初始化差分進化算法首先在種群向量的限制范圍內(nèi)隨機生成候選解的初始種群，其中第i個種群向量的第j維分量xi，j為

式中：rand［0，1］為0 和1 之間均勻分布的隨機數(shù)；下標max 和min 為種群限制范圍的上下限；i=1，2，…，N；j=1，2，…，d，其中d為種群向量的維度。

變異在差分進化算法中，第G代第i個種群父代向量xi，G會執(zhí)行變異操作，產(chǎn)生變異向量vi，G。L-SHADE 算法采用Zhang 和Sanderson［22］提出的current-to-pbest/1變異策略：

式中：Fi為縮放因子；xpbest，G是從第G代前pN（p∈[0，1]，乘積四舍五入到整數(shù)）個精英個體中隨機選擇的一個個體；xr1，G∈PG是從第G代N個種群向量集合PG中隨機選擇的一個個體；xr2，G∈PG∪A，其中A為一個外部存儲單元，存儲前G-1 代中選擇失敗的個體，當種規(guī)模超過N時，則隨機選擇一個A中的舊個體替換為新個體。通過引入存儲單元A，將劣解納入變異的過程，增加了種群的多樣性。

變異操作可能會使個體數(shù)值向量的某個值超出限制范圍，通過式（9）對數(shù)值進行修正：

交叉實驗向量個體vi，j，G由變異向量個體vi，j，G和父向量個體xi，j，G交叉產(chǎn)生：

式中：jrand為［1：d］中的隨機整數(shù)。通過引入jrand，使每個父向量中至少一維進行交叉操作，從而提高算法的收斂速度。

選擇分別計算實驗向量和父向量的適應度函數(shù)值，遵循適者生存的原則確定下一代個體，將表現(xiàn)更好的向量保留到下一代，數(shù)學表達式為

參數(shù)自適應控制參數(shù)決定了差分進化算法的精度，而人工設置控制參數(shù)無法使所設參數(shù)滿足所有優(yōu)化問題的需求，因此需要對控制參數(shù)進行自適應調(diào)整。

在L-SHADE 算法中，采用基于歷史信息存檔的自適應參數(shù)更新策略。在經(jīng)過選擇操作后，每一代選擇成功種群的交叉率SCR、縮放因子SF的Lehmer 均值MCR、MF存儲在大小為H的歷史存儲器中。最初MCR、MF初始化為0.5，當進化到新一代時，索引k=k+1，當k＞H時，則令k=1，該策略的歷史存儲器如表1 所示。

表1 歷史存儲器Table 1 History memory

在每一代中，每個種群個體xi采用各自的控制參數(shù)CRi、Fi，具體由式（12）和式（13）生成：

式中：randn 為正態(tài)分布，randc 為柯西分布，ri為［1，H］中生成的隨機整數(shù)。Lehmer 均值MCR、MF更新方法如下：

式中：Δfk=|f(uk，G)-f(xk，G)|，為成功選擇個體和父代個體適應度函數(shù)之差；Sk分別代指SCR和SF，為優(yōu)勝縮放因子和交叉率的集合。通過引入Lehmer 均值，使貢獻大的參數(shù)具有更大的影響力，較算術均值而言，Lehmer 均值的數(shù)值較大，可以避免控制參數(shù)過小而造成早熟的現(xiàn)象。

當進化過程中，某一代成功選擇個體的數(shù)量為0 時，則使CRi的值強制鎖定為0，使得每次的交叉過程只更改一個參數(shù)，以減緩算法的收斂速度，提高局部搜所能力。

種群規(guī)模縮減由于較大的種群規(guī)模會使優(yōu)化算法具備更好的種群多樣性，能夠擁有更好勘探能力，而較小的種群規(guī)模能夠增加對有希望區(qū)域的利用，使算法具有更好的開發(fā)能力，因此種群縮減技術受到了廣泛的關注［23］。

在L-SHADE 算法中，通過引入線性種群縮減技術（Linear Population Size Reduction，LPSR）提高算法的效率和收斂性，其中種群縮減遵循

式中：NG+1為下一代的種群規(guī)模；round［］為取整函數(shù)；Nmin=4 表示進化過程最終的種群大小；Ninit為初始種群大小；Nfes為當前的評估次數(shù)；Nfesmax為最大評估次數(shù)。

1.2.2 改進L-SHADE 算法

標準L-SHADE 算法在解決無約束優(yōu)化問題時表現(xiàn)較好，然而在處理復雜的非線性問題時易出現(xiàn)早熟的問題，過分依賴隨機生成的初始向量。在解決航空發(fā)動機性能退化計算問題時，如果發(fā)動機多工作點工作狀態(tài)偏差較小，進化會出現(xiàn)過早停滯的現(xiàn)象，導致收斂精度不高。由于LSHADE 算法引入了線性種群縮減策略，該方法在迭代初期種群就快速縮減，且在迭代過程中，比例因子過早達到了設定的最大值，這些都增加了過早收斂的可能。針對上述問題，提出了一種改進的L-SHADE 算法。

標準L-SHADE 算法的線性種群規(guī)?？s減策略種群變化過程單一，難以滿足不同優(yōu)化問題的需求。文獻［24］提出了一種非線性種群縮減策略（NLPSR），其策略如下：

式中：nfesr=Nfes/Nfesmax，為當前評估次數(shù)與最大評估次數(shù)的比值。

這種種群縮減策略較LPSR 而言減少了總體的評估次數(shù)，提高了算法的收斂速度，在針對簡單優(yōu)化問題時是非常有效的。然而與LPSR 技術相同，從第1 代種群就開始減少種群，并沒有給個體足夠的進化時間，意味著需要浪費部分計算資源來初始化個體并迅速刪除它們，造成了種群資源的浪費。且這種方法在進化初期種群規(guī)模減小太快，影響了種群的多樣性，增加了局部收斂的可能性。

由于航空發(fā)動機熱力學計算過程極為復雜，其健康因子求解過程并不是簡單的單峰值問題，因此需要減緩進化前期種群縮減速度，保持種群的多樣性，增強算法全局搜索的能力。同時為提高算法效率，需要在進化后期大幅降低種群規(guī)模。針對上述需求，提出了一種新的種群縮減方法，將當前評估次數(shù)與最大評估次數(shù)的比值nfesr作為指數(shù)，1-nfesr為底數(shù)，使種群規(guī)模在進化初期幾乎保持不變，進化后期大幅降低，其表達式為

圖1 為不同種群縮減策略種群規(guī)模隨評估次數(shù)變化的對比，其中NLPSR 為文獻［24］中的種群縮減策略，改進NLPSR 為提出的種群縮減策略。由圖1 可見，NLPSR 和LPSR 在迭代前期種群快速減少，而提出的種群縮減策略在迭代前期種群數(shù)量基本保持不變，能夠給予初代個體充分的進化時間，增強算法的多樣性。

圖1 種群規(guī)模隨評估次數(shù)的變化曲線Fig.1 Variation curves of population size with number of assessments

在對種群規(guī)?？s減策略進行改進后，發(fā)動機性能退化評估求解過程中面臨的另一個問題是選擇合適的變異策略。

理論上要求所選擇的變異策略在迭代初期使用較小的比例因子，以降低貪婪算子在變異過程中的權值，增強算法的全局搜索能力；在迭代后期使用較大的比例因子，以提高算法的局部搜索能力，為此引入Brest 等［25］在jSO 算法中提出的變異策略current-to-pBest-w/1：

式中：Fw為加權比例因子。在發(fā)動機性能退化迭代過程中，發(fā)現(xiàn)當評估次數(shù)大于0.4Nfesmax時，F(xiàn)w的值會過早達到設定的最大值，導致算法過早收斂，出現(xiàn)早熟的情況。因此提出對Fw的表達式進行平滑處理：

通過對Fw的線性變化，減少了其達到最大值的速度，增強了種群的多樣性。這種加權的變異策略與非線性種群縮減技術的目的相同，能夠提高算法初期的勘探能力和末期的開發(fā)能力。

1.3 改進L-SHADE 算法的流程框架

改進L-SHADE 算法的偽代碼如算法1所示。

?

改進L-SHADE 算法求解發(fā)動機健康因子的基本流程如圖2 所示。

圖2 改進L-SHADE 算法流程圖Fig.2 Flow chart of improved L-SHADE algorithm

2 改進算法的驗證與應用

2.1 改進算法的有效性驗證

為證明所提算法的有效性，采用30 個經(jīng)典的基準函數(shù)［26］進行實驗，其中30 個函數(shù)可分為單峰可分離函數(shù)（F1～F5）、單峰不可分離函數(shù)（F6～F13）、多峰可分離函數(shù)（F14～F18）、多峰不可分離函數(shù)（F19～F25）和其他組合函數(shù)（F26～F30），具體函數(shù)信息可通過文獻［26］的網(wǎng)站查尋。將提出的改進L-SHADE 算法與標準LSHADE［14］、L-SHADE-cnEpSin （Ensemble Sinusoidal Differential Covariance Matrix Adaptation with Euclidean Neighborhood with L-SHADE）［27］和L-SHADE-SPACMA（L-SHADE with Semi-Parameter Adaptation Hybrid with CMA-ES）［28］進行對比，其中L-SHADE-cnEpSin 與L-SHADESPACMA 均為CEC2017優(yōu)化競賽的優(yōu)勝算法［29］。

具體的驗證細節(jié)如下：

1） 各算法均采用已有文獻中的默認參數(shù)設置。

2） 各算法的種群大小為100，最大迭代次數(shù)為100，每個算法均進行50次獨立的重復實驗。

3） 各算法的優(yōu)化結果以優(yōu)化得到解對應的函數(shù)值f(x)與函數(shù)理論最優(yōu)解對應的函數(shù)值f(x*)差值的絕對值來比較，因此0為優(yōu)化的最好結果。

4） 各算法的性能指標通過50 次實驗優(yōu)化結果|f(x)-f(x*)|的平均值、標準差、最小值、最大值和平均運行時間來體現(xiàn)。

5） 30 個經(jīng)典基準函數(shù)的維度（自變量個數(shù)）均設置為30。

6） 在Windows7 系統(tǒng)的MATLAB2021a 開展實驗，基于Intel?CoreTMi7-9700CPU @ 3.00 GHz 的個人電腦。

表2 為4 種算法在30 個經(jīng)典基準函數(shù)的驗證結果，加粗的數(shù)據(jù)為4 種算法中結果最好的值，空白處表示L-SHADE-cnEpSin 和L-SHADESPACMA 算法沒有相應的機制來控制其在指定范圍內(nèi)的搜索，進一步導致優(yōu)化不收斂的情況。30 個函數(shù)中，改進L-SHADE 算法取得最優(yōu)解的次數(shù)為14 次，L-SHADE 算法為3 次，L-SHADE-cnEpSin 算法為4 次，L-SHADESPACMA 算法為9 次。

表2 基準函數(shù)實驗結果Table 2 Experimental results of benchmark functions

其中在單峰函數(shù)（F1～F13）中，改進L-SHADE 算法明顯優(yōu)于其他3 種算法。在多峰函數(shù)和組合函數(shù)（F14～F30）中，改進L-SHADE 算法較標準L-SHADE 算法也具有明顯的優(yōu)勢，盡管L-SHADE-cnEpSin 算法和L-SHADE-SPACMA 算法在部分函數(shù)中取得了最優(yōu)解，但L-SHADE-cnEpSin 算法在F18、F20 函數(shù)中不收斂、在F16、F23 函數(shù)中收斂效果不好；L-SHADE-SPACMA 算法在F24 函數(shù)中不收斂。而改進L-SHADE 算法具有更好的魯棒性，能夠適用絕大多數(shù)函數(shù)。另外，改進L-SHADE 算法的運行時間較標準L-SHADE算法增加不大，但比其他2 種算法運行時間明顯減少，降低了算法的復雜度。綜上，提出的改進L-SHADE 算法較其他3 種算法具有明顯優(yōu)勢。

2.2 改進L-SHADE 算法在發(fā)動機性能退化中的應用

研究對象共布置有5 個氣路傳感器，分別為燃油流量Wf，高壓轉子轉速N2、低壓轉子轉速N1、高壓壓氣機后總壓Pt3、低壓渦輪后總溫Tt5；其中Wf為發(fā)動機模型的輸入，模型的輸出量p為4 個。研究對象的健康因子為該型發(fā)動機旋轉部件（風扇、低壓壓氣機、高壓壓氣機、高壓渦輪、低壓渦輪）的效率和流通能力，分別表示為SE12、SW12、SE2、SW2、SE26、SW26、SE41、SW41、SE46、SW46，數(shù)量n=10。

為使發(fā)動機健康因子估計方程組為非欠定方程，即pq≥n，多工作點數(shù)量應滿足q≥3。理論上q的值越大，方程組越冗余，健康因子估計效果越好。但在實際情況下，隨著q的增加計算成本顯著增加，且受發(fā)動機傳感器測量噪聲的影響，數(shù)據(jù)越多，發(fā)動機健康因子的求解越不易收斂，且更容易出現(xiàn)“拖尾”現(xiàn)象，即一個部件的性能退化可能會反映到另一個部件上，使得健康因子向噪聲干擾的方向變化。因此q選擇為3。

為驗證改進L-SHADE 算法估計航空發(fā)動機健康因子性能的有效性，以確保輸入在發(fā)動機的真實飛行包線中，在某型小涵道比渦扇發(fā)動機的實際飛參數(shù)據(jù)中提取高壓轉子特定轉速范圍（86%～92%）的穩(wěn)態(tài)數(shù)據(jù)點，在發(fā)動機性能計算程序中預先植入假定的發(fā)動機各部件性能退化值，在程序中再輸入提取的飛參數(shù)據(jù)中的發(fā)動機的運行條件（飛行高度、馬赫數(shù)、大氣溫度）和控制變量（燃油流量Wf），經(jīng)計算，輸出3 型發(fā)動機傳感器參數(shù)（高壓轉子轉速N2、低壓轉子轉速N1、高壓壓氣機后總壓Pt3、低壓渦輪后總溫Tt5）的數(shù)值。選用3 組穩(wěn)態(tài)工作點數(shù)據(jù)代入式（6）中，通過最小化適應度函數(shù)求解發(fā)動機健康因子。

隨機選擇2 組預先植入的性能退化量的數(shù)據(jù)分別為算例1 和算例2，使用L-SHADE、L-SHADE-cnEpSin、L-SHADE-SPACMA 和改進L-SHADE 算法進行健康因子估計，其中各算法的種群大小為100，最大迭代次數(shù)為500，健康因子搜索范圍為［0.95，1.05］，迭代過程的對比如圖3 所示，計算結果如表3 所示。

圖3 不同算法迭代過程比較Fig.3 Comparison of iterative process of different algorithms

表3 航空發(fā)動機性能退化植入和各算法估計結果Table 3 Engine performance degradation implantation and estimation results of each algorithm

由表3 的對比結果可見，改進L-SHADE 算法的計算結果更好，能夠滿足發(fā)動機健康因子估計的工程精度需求。

由圖3 的迭代過程比較可見，L-SHADE 算法 、L-SHADE-cnEpSin 算 法 和 L-SHADESPACMA 算法受其種群縮減策略的影響，沒有給初代個體充分的進化時間，影響了種群的多樣性，在迭代后期易陷入局部最優(yōu)；而改進LSHADE 算法增強了算法迭代前期的種群多樣性和算法后期的開發(fā)能力，計算精度較其他3 種算法精度更高。

為進一步驗證算法的魯棒性，對上述2 個算例進行了20 次的重復實驗，其計算結果如表4 所示，結果表明所提的改進L-SHADE 算法在航空發(fā)動機性能退化評估中具有較強的魯棒性，計算精度較標準L-SHADE 算法平均提高了65.5%。

表4 各算法平均計算誤差Table 4 Average calculation errors of each algorithm

3 結 論

針對航空發(fā)動機的特點，通過引入改進種群縮減策略和加權變異策略，提出了一種改進的L-SHADE 算法，在對該算法有效性驗證的基礎上，開展了航空發(fā)動機性能退化評估的應用研究，得到的主要結論如下。

1） 改進L-SHADE 算法在航空發(fā)動機健康因子搜索過程中能夠增強迭代前期的勘探能力和后期的開發(fā)能力。

2） 通過30 個經(jīng)典基準函數(shù)對該算法進行了驗證，與其他DE 算法相比，該算法具有精度高、魯棒性強的優(yōu)點。

3） 航空發(fā)動機性能退化評估的應用表明，改進L-SHADE 算法具有較強的工程適應性，計算精度較標準L-SHADE 算法平均提高了65.5%。