董 慶 ，李本威

（1.海軍裝備部駐蘇州地區(qū)軍事代表室，江蘇 蘇州 215011；2.海軍航空大學(xué)航空基礎(chǔ)學(xué)院，山東 煙臺 264001）

0 引言

航空發(fā)動機是高度非線性的時變系統(tǒng)，發(fā)動機模型在其性能及仿真、故障診斷、控制規(guī)律設(shè)計及實時監(jiān)控等方面發(fā)揮著重要作用。針對渦軸發(fā)動機喘振試驗，建立準確的渦軸發(fā)動機模型可避免試驗中的某些風(fēng)險和經(jīng)濟損失。

目前，建立航空發(fā)動機模型的方法主要有解析法和系統(tǒng)辨識法。采用解析法建立發(fā)動機模型存在建模過程復(fù)雜且有一定的偏差和難以應(yīng)用到發(fā)動機實時監(jiān)控中等問題；而采用系統(tǒng)辨識的方法只需足夠的發(fā)動機輸入輸出數(shù)據(jù)就能簡單快速的建立較高精度的發(fā)動機模型。

傳統(tǒng)的辨識方法有最小二乘法、子空間狀態(tài)辨識法等。王磊等利用改進的子空間辨識法建立了某型渦扇發(fā)動機的小偏差狀態(tài)變量模型，但假設(shè)發(fā)動機是線性時不變系統(tǒng)，在一定程度上影響了模型性能。人工神經(jīng)網(wǎng)絡(luò)（Artificial Neural Networks，ANN）、支持向量機（Support Vector Machine，SVM）和極限學(xué)習(xí)機（Extreme Learning Machine，ELM）等具有很好非線性學(xué)習(xí)能力的機器學(xué)習(xí)算法已越來越多地應(yīng)用于航空發(fā)動機領(lǐng)域。丁凱峰等利用RBF 網(wǎng)絡(luò)建立了航空發(fā)動機辨識模型；陳超等利用人工神經(jīng)網(wǎng)絡(luò)對航空發(fā)動機起動模型進行了辨識與仿真。但以上人工神經(jīng)網(wǎng)絡(luò)都存在無反饋單元或者無延時單元等問題，而航空發(fā)動機是高度時變系統(tǒng)，傳統(tǒng)的神經(jīng)網(wǎng)絡(luò)建立的發(fā)動機模型精度有時不能滿足實際的精度要求。Asgari 等利用非線性自回歸網(wǎng)絡(luò)（Nonlinear Auto Regressive with Exogenous Inputs Model，NARX）建立了某型單軸燃氣渦輪起動階段模型并進行了仿真分析，獲得不錯的效果。但NARX 網(wǎng)絡(luò)參數(shù)由隨機函數(shù)產(chǎn)生，在一定程度上影響了其性能。

螢火蟲算法（Firefly Algorithm，F(xiàn)A）是一種新型智能優(yōu)化算法，具有參數(shù)少、尋優(yōu)速度快等特點，但由于步長因子的影響會導(dǎo)致前期收斂速度慢、出現(xiàn)局部最優(yōu)和出現(xiàn)峰值處不穩(wěn)定震蕩等問題，為此提出變步長螢火蟲算法（Change of Step of Firefly Algorithm，CSFA），利用該算法對NARX 網(wǎng)絡(luò)參數(shù)進行全局尋優(yōu)來改善網(wǎng)絡(luò)性能。

針對以上分析，結(jié)合已經(jīng)獲取的渦軸發(fā)動機臺架試車數(shù)據(jù)，本文提出以試車數(shù)據(jù)為驅(qū)動、利用CSFA優(yōu)化的NARX 網(wǎng)絡(luò)辨識渦軸發(fā)動機起動過程模型的方法。

1 CSFA-NARX網(wǎng)絡(luò)模型

1.1 NARX網(wǎng)絡(luò)

NARX 網(wǎng)絡(luò)為有外部輸入的非線性自回歸網(wǎng)絡(luò)。NARX 網(wǎng)絡(luò)不僅與當前時刻的輸入有關(guān)，也與過去的輸入輸出緊密相關(guān)，使用動態(tài)修正方法減少了模型修正的計算時間，樣本增加時矩陣階數(shù)不變，避免了矩陣求逆運算，提高了計算效率。

NARX 網(wǎng)絡(luò)通過引入延時單元及反饋來實現(xiàn)時間序列的動態(tài)預(yù)測，由帶有外部輸入的非線性差分方程表示為

式中：()為時刻的輸出值；()為時刻的輸入值；為自變量的非線性函數(shù)；為輸出延時階次；為輸入延時階次；(-)、(-)分別為輸出、輸入?yún)?shù)相應(yīng)延時的數(shù)值。

NARX 網(wǎng)絡(luò)主要由回歸分析器和非線性估計器2部分組成，包括輸入層、隱含層和輸出層，如圖1 所示?；貧w分析器（輸入層與隱含層）利用時刻及過去的輸入數(shù)據(jù)（），（-1），（-2），…，（-），和時刻過去的輸出數(shù)據(jù)（-1），（-2），…，（-），結(jié)合輸入權(quán)值和隱含層單元的偏置，由隱含層單元激活函數(shù)計算得到的值即為回歸分析器的輸出值，并傳遞給非線性估計器?；貧w分析器的輸出值作為非線性估計器的輸入值，非線性估計器利用其非線性函數(shù)模塊和線性函數(shù)模塊估算輸入輸出之間的映射關(guān)系。

圖1 NARX網(wǎng)絡(luò)結(jié)構(gòu)

某時刻第個隱含點的輸出h為

式中：為隱層節(jié)點使用的雙曲線正切sigmoid 函數(shù)；w()、w()分別為時刻第個隱層節(jié)點與輸入、輸出反饋之間的權(quán)值；b為第個隱層節(jié)點的偏置。

某時刻第個輸出層節(jié)點輸出()為

式中：λ()為時刻第個隱層節(jié)點到第個輸出層節(jié)點的權(quán)值；b為第個輸出層節(jié)點偏置；為隱含層節(jié)點數(shù)。

1.2 CSFA算法

FA 算法來源于大自然中螢火蟲的自然行為，群體中較亮的個體會吸引周圍的個體向這個方向移動來完成自身位置和亮度的更新，從而找出最優(yōu)位置來實現(xiàn)尋優(yōu)過程。FA 算法具有參數(shù)少、尋優(yōu)速度快、效率高等特點，從而得到廣泛應(yīng)用。但FA 算法迭代步長是固定的，導(dǎo)致若步長取較大值，則算法后期不利于求解局部最優(yōu)，易出現(xiàn)峰值處不穩(wěn)定震蕩問題；若步長取較小值，則算法前期不易快速找到全局最優(yōu)鄰域且收斂速度較慢。故提出改進的變步長螢火蟲算法。

螢火蟲的熒光亮度和相互吸引度是螢火蟲算法2 個重要參數(shù)，發(fā)光越亮的螢火蟲位置越優(yōu)，即目標函數(shù)值越優(yōu)，最亮螢火蟲的位置即是最優(yōu)解，發(fā)光越亮的螢火蟲對周圍的螢火蟲吸引度越強。

螢火蟲的相對熒光亮度和距離的平方成反比，相對熒光亮度為

式中：為螢火蟲最大亮度，即自身亮度（=0），目標函數(shù)值越優(yōu)，自身亮度越大；為光強吸收系數(shù)，因熒光會隨著距離的增加和傳播媒介的吸收逐漸減弱，設(shè)置光強吸收系數(shù)以體現(xiàn)此特性，可設(shè)置為常數(shù)。

2只螢火蟲之間的距離r為

相對吸引度為

式中：為最大吸引度，即=0處的吸引度。

最優(yōu)目標迭代，光亮較弱的個體會向光亮較強的個體移動，位置更新公式為

式中：x、x為螢火蟲和的空間位置；為步長因子；為[0 1]之間均勻分布的隨機因子；為可調(diào)整的正相關(guān)系數(shù)，即用kr代替固定步長。

CSFA算法尋優(yōu)基本流程：

（1）初始CSFA 算法特征參數(shù)。設(shè)置螢火蟲數(shù)目、、、、和最大迭代步數(shù)或者精度要求。

（2）隨機賦予螢火蟲初始位置，計算螢火蟲的目標函數(shù)值作為自身熒光亮度。

（3）計算每對螢火蟲之間的相對亮度和吸引度，個體移動的方向由相對亮度決定。

（4）更新每個個體的位置，隨機移動處于最佳位置的螢火蟲，并重新計算個體的自身熒光亮度。

（5）判斷是否滿足終止條件，若滿足進行下一步，否則返回步驟（3）繼續(xù)。

（6）輸出全局極點值和最優(yōu)個體值。

1.3 CSFA-NARX網(wǎng)絡(luò)模型

NARX 網(wǎng)絡(luò)初始參數(shù)初始化由隨機函數(shù)確定，大大影響模型的預(yù)測結(jié)果。因此，利用變步長螢火蟲優(yōu)化算法全局搜索來確定NARX 網(wǎng)絡(luò)的輸入權(quán)值、輸出權(quán)值、隱層節(jié)點的偏置和輸出層節(jié)點偏置b。將獲得的、、和b返回NARX網(wǎng)絡(luò)用于模型的辨識。

CSFA-NARX網(wǎng)絡(luò)模型流程如圖2所示。

圖2 CSFA-NARX網(wǎng)絡(luò)模型流程

目標函數(shù)選取均方誤差作為模型的性能評價標準，來衡量模型預(yù)測值與實測數(shù)據(jù)的吻合程度

式中：（）為模型的預(yù)測值；（）為實測值；為樣本數(shù)量。

2 渦軸發(fā)動機起動過程模型建立

某型渦軸發(fā)動機的起動由起動機、發(fā)動機、點火系統(tǒng)、供油系統(tǒng)以及旋翼負載系統(tǒng)等多個系統(tǒng)協(xié)同工作來完成。其中，不同大氣條件、起動機特性、點火和供油特性及其旋翼負載特性是發(fā)動機模型的輸入量；在起動過程中的狀態(tài)和性能參數(shù)是發(fā)動機模型的輸出量。因此渦軸發(fā)動起動過程影響因素有：大氣溫度、大氣壓力、起動機的輸出功率、負載力矩、燃油流量。在渦軸發(fā)動機的起動過程中，通常負載力矩在常規(guī)大氣溫度大氣壓力下變化較小可忽略不計。同時，選取渦軸發(fā)動機參數(shù)燃氣發(fā)生器轉(zhuǎn)速、輸出軸轉(zhuǎn)速和渦輪后溫度作為模型的輸出參數(shù)來衡量發(fā)動機的起動性能。

發(fā)動機的狀態(tài)不僅與當前時刻發(fā)動機的各參數(shù)密切相關(guān)，還與此時刻之前發(fā)動機的狀態(tài)密不可分。因此，需要設(shè)置輸入輸出延時階次提高模型的精度，故構(gòu)建的渦軸發(fā)動機起動過程數(shù)學(xué)模型表達式為

式中：、分別為模型的輸入、輸出量；為待辨識的渦軸發(fā)動機起動過程數(shù)學(xué)模型；、分別為輸入、輸出延時階次；根據(jù)發(fā)動機動態(tài)慣性，發(fā)動機模型輸入輸出參數(shù)延時階次通常設(shè)置為2。因此，將延時層的延時階次和均設(shè)置為2。

發(fā)動機轉(zhuǎn)子從靜止狀態(tài)逐漸加速到地面慢車狀態(tài)的過程稱為發(fā)動機的起動過程。根據(jù)起動系統(tǒng)和渦輪的功率變化情況，將起動過程劃分為3 個階段，建立發(fā)動機起動過程模型如圖3所示。

圖3 渦軸發(fā)動機起動過程不同階段

第1階段：按下起動按鈕1.25 s內(nèi)，由起動機帶動燃氣渦輪工作，發(fā)動機轉(zhuǎn)速在規(guī)定時間內(nèi)達到（0.08~0.12）×33400 r/min 時燃油系統(tǒng)開始工作點火。在此過程中發(fā)動機完全由起動機帶動工作，發(fā)動機的點火和燃油系統(tǒng)不工作，燃燒室無供油點火，渦輪后溫度約等于大氣環(huán)境溫度。

第2 階段：起動機和燃氣發(fā)生器共同工作帶動發(fā)動機加速工作，燃燒室開始供油點火燃燒。起動燃油調(diào)節(jié)器按起動燃油調(diào)節(jié)規(guī)律供油，直到燃氣發(fā)生器轉(zhuǎn)速達到（14000±1000）r/min時，起動機脫開。

第3 階段：起動機脫開后，燃氣發(fā)生器獨自帶動發(fā)動機加速工作至地面慢車狀態(tài)，起動過程完成。

整個起動過程需要限制在60 s 內(nèi)。渦軸發(fā)動機起動過程模型輸入輸出參數(shù)見表1。

表1 發(fā)動機起動過程模型輸入輸出參數(shù)

3 CSFA-NARX網(wǎng)絡(luò)模型辨識結(jié)果及分析

結(jié)合NARX 網(wǎng)絡(luò)模型和建立的渦軸發(fā)動機起動過程模型，根據(jù)表1 將渦軸發(fā)動機模型的輸入輸出參數(shù)分別作為NARX 網(wǎng)絡(luò)模型的輸入輸出參數(shù)，結(jié)合已經(jīng)獲取的試車數(shù)據(jù)，辨識得到渦軸發(fā)動機起動過程模型并進行驗證。

3.1 模型辨識數(shù)據(jù)的提取與預(yù)處理

本文采用的數(shù)據(jù)來自于某型渦軸發(fā)動機1000 h完整地面臺架試車試驗。該試驗的數(shù)據(jù)采集由測速發(fā)動機傳感器、不同類型動態(tài)壓力傳感器和不同類型溫度傳感器等配合數(shù)據(jù)采集系統(tǒng)完成，數(shù)據(jù)的采樣頻率為20 Hz，對發(fā)動機的進氣條件（如進氣溫度、進氣壓力和空氣流量等）、燃氣渦輪前后溫度、各轉(zhuǎn)速（燃氣發(fā)生器轉(zhuǎn)速和輸出軸轉(zhuǎn)速）、燃油消耗量和其他參數(shù)進行采集。

依據(jù)建立的渦軸發(fā)動機起動過程3 階段模型（表1），從提取的起動過程數(shù)據(jù)中選擇、、、、、和作為樣本參數(shù)。起動機的輸出功率由起動機的電流、電壓和效率計算得到。

本文所采用數(shù)據(jù)的4 組試車試驗起動過程各測試傳感器安裝位置相同且數(shù)據(jù)采集系統(tǒng)其它裝置均一致，進氣溫度設(shè)置為9.2、18.3、26.8、34.6 ℃，將4 組數(shù)據(jù)依次記為、、和。第1 階段每一溫度下選取25 個數(shù)據(jù)點，4 組共100 個數(shù)據(jù)點作為辨識樣本。第2、3 階段每一溫度下選取50 個數(shù)據(jù)點，4 組共200 個數(shù)據(jù)點分別作為辨識樣本，每組起動過程3 個階段共125個數(shù)據(jù)點。

對提取的辨識樣本進行剔除異常值、平滑、濾波和歸一化處理，在完成模型的訓(xùn)練和驗證后對歸一化的數(shù)據(jù)進行反歸一化還原。

3.2 CSFA及NARX網(wǎng)絡(luò)參數(shù)設(shè)置

關(guān)于CSFA 算法初始參數(shù)的設(shè)置，參考文獻[22]，選取螢火蟲數(shù)目=85，最大迭代次數(shù)取200，最大吸引度=1，步長因子=0.12，吸收系數(shù)=1.02，正相關(guān)系數(shù)=2。

NARX 網(wǎng)絡(luò)輸入層節(jié)點數(shù)、輸出層節(jié)點數(shù)分別與發(fā)動機模型的輸入、輸出參數(shù)數(shù)目相同。通常，將隱含層數(shù)設(shè)定為1，隱含層神經(jīng)元數(shù)由經(jīng)驗公式確定

式中：為10~14的常數(shù)。

因此，本NARX網(wǎng)絡(luò)隱含層節(jié)點數(shù)=15。延時層的延時階次與發(fā)動機起動過程數(shù)學(xué)模型（式（9））保持一致，輸入輸出延時階次均設(shè)置為2。

3.3 模型辨識結(jié)果及分析

CSFA 算法和FA 算法在對NARX 網(wǎng)絡(luò)參數(shù)全局搜索尋優(yōu)過程中，適應(yīng)度值的變化過程如圖4 所示。從圖中可見，CSFA 算法在迭代前期具有較高的全局搜索能力，后期收斂也較快。

圖4 適應(yīng)度值在尋優(yōu)過程中變化

用辨識模型預(yù)測值與實測值最大相對誤差以及最大相對誤差均值來評價建立的辨識模型性能。相對誤差為

式中：（）為模型的預(yù)測值；（）為實測值。

結(jié)合預(yù)處理后的發(fā)動機起動過程數(shù)據(jù)，利用CS?FA-NARX 網(wǎng)絡(luò)模型對渦軸發(fā)動機起動過程模型進行回歸辨識。先將發(fā)動機起動過程3 階段模型分別辨識得出輸出參數(shù)值，因?qū)嶋H起動過程是連續(xù)過程，所以將3 階段辨識模型得到的輸出參數(shù)數(shù)據(jù)點連接起來，得到完整起動過程的不同輸出參數(shù)的輸出數(shù)據(jù)。

綜合考慮所使用的計算機性能和模型樣本數(shù)量，本文采用留一交叉驗證法（Leave-One-Out Cross Val?idation，LOOCV）對模型的性能進行驗證。例如，以、和為訓(xùn)練樣本，以為驗證樣本得到辨識模型輸出參數(shù)、和的數(shù)值，與試車實測數(shù)據(jù)對比，結(jié)果和相對誤差分布如圖5~7 所示。4 組驗證樣本的驗證誤差見表2。

從圖5~7 和表2 中可見，輸出參數(shù)、和的辨識輸出結(jié)果都較好地逼近了實測數(shù)據(jù)，4 組驗證樣本中各出參數(shù)的最大相對誤差平均值分別為0.90%、1.51%和2.01%。表明本文提出的CSFA-NARX 網(wǎng)絡(luò)模型可用于渦軸發(fā)動機起動過程模型的辨識且模型精度滿足實際應(yīng)用的精度要求，取得了良好的效果。

圖5 M4為驗證樣本ng辨識結(jié)果

圖6 M4為驗證樣本nr辨識結(jié)果

圖7 M4為驗證樣本T4辨識結(jié)果

表2 4組驗證樣本的誤差

3.4 不同算法比較

在MATLAB R2018a 平臺上，利用相同的訓(xùn)練樣本和驗證樣本，通過對比CSFA-NARX、FA-NARX、NARX 和CSFA-BP 算法的辨識結(jié)果來驗證本文提出的CSFA-NARX算法的優(yōu)越性。

CSFA-NARX、FA-NARX 和NARX 3 種算法延時層階次均設(shè)置為2，4 種算法隱含層神經(jīng)元數(shù)根據(jù)式（10）設(shè)置為15，輸入輸出神經(jīng)元數(shù)根據(jù)表1 設(shè)置，NARX 網(wǎng)絡(luò)均采用（Levenberg-Marquardt，LM）訓(xùn)練算法，CSFA 和FA 算法螢火蟲數(shù)目85，最大迭代次數(shù)均取200，=1，=0.12，=1.02，CSFA算法中=2。

以、和為訓(xùn)練樣本和為驗證樣本，不同算法得到辨識模型輸出參數(shù)n值與實測數(shù)據(jù)對比結(jié)果和相對誤差分布如圖8 所示；采用不同算法4 組驗證樣本的驗證誤差對比如圖9所示。

圖8 M1為驗證樣本不同算法ng的辨識結(jié)果

圖9 不同算法辨識結(jié)果誤差對比

從圖8、9 中可見，改進的變步長螢火蟲算法優(yōu)化的NARX 網(wǎng)絡(luò)辨識模型精度高于螢火蟲算法優(yōu)化的NARX 網(wǎng)絡(luò)辨識模型的，更優(yōu)于直接采用NARX 網(wǎng)絡(luò)辨識模型的。而CSFA-NARX 網(wǎng)絡(luò)辨識模型精度與CSFA-BP 辨識模型的對比充分說明了具有反饋延時層的NARX網(wǎng)絡(luò)性能比傳統(tǒng)BP網(wǎng)絡(luò)更優(yōu)。

4 結(jié)論

（1）辨識模型輸出參數(shù)、和與試車實測數(shù)據(jù)的偏差都在合理范圍之內(nèi)，表明基于CSFA-NARX網(wǎng)絡(luò)模型對渦軸發(fā)動機起動過程模型辨識的方法具有可行性。

（2）采用留一交叉驗證方法對辨識模型的性能進行驗證，模型輸出參數(shù)、和的最大相對誤差均值分別為0.90%、1.51%和2.01%，表明基于CSFANARX 網(wǎng)絡(luò)模型對渦軸發(fā)動機起動過程辨識得到的模型精度滿足發(fā)動機實際要求。

（3）在相同的訓(xùn)練樣本和驗證樣本情況下，基于CSFA-NARX 網(wǎng)絡(luò)模型對渦軸發(fā)動機起動過程辨識模型的精度優(yōu)于FA-NARX、NARX 和CSFA-BP 算法辨識模型的，驗證了本文提出的CSFA-NARX 算法的優(yōu)越性。

（4）本文建立的渦軸發(fā)動機起動過程辨識模型可用于對其他不同環(huán)境條件下（如高溫、高壓等）渦軸發(fā)動機起動過程性能的遞推估算，實現(xiàn)發(fā)動機性能參數(shù)預(yù)測，也可為渦軸發(fā)動機后續(xù)的試車試驗提供指導(dǎo)和為其起動過程性能優(yōu)化奠定基礎(chǔ)。