黃 亮,劉君強,貢英杰

(南京航空航天大學(xué)民航學(xué)院,江蘇 南京 211106)

0 引 言

有效地對航空發(fā)動機的剩余壽命進行預(yù)測,并在預(yù)測結(jié)果的基礎(chǔ)上提早地制定相應(yīng)的維修及備件計劃,能夠很大程度上減少維修費用。為降低航空發(fā)動機的維修保障費用,保證飛機的運行安全,亟需開展航空發(fā)動機剩余壽命預(yù)測的有關(guān)研究[1]。

壽命預(yù)測的研究方向主要可以分為兩類:基于數(shù)據(jù)驅(qū)動的方法和基于物理失效模型。發(fā)動機的失效是由于材料在高溫高壓的環(huán)境下逐漸退化的結(jié)果[2]。隨著各類產(chǎn)品變得越來越復(fù)雜化,基于物理失效模型的方法很難建立可靠的故障模型來對應(yīng)產(chǎn)品的失效過程[3]。基于數(shù)據(jù)驅(qū)動相關(guān)的方法主要依靠產(chǎn)品運行過程中監(jiān)測的退化數(shù)據(jù),建立相應(yīng)的退化模型,具有明顯的計算和建模優(yōu)勢[4-5]。

目前,壽命預(yù)測的研究重點是對監(jiān)測數(shù)據(jù)進行分析。文獻[6]等運用數(shù)據(jù)融合的方法,實現(xiàn)了基于可靠性和數(shù)據(jù)驅(qū)動的在線剩余壽命預(yù)測;文獻[7]則利用Kalman濾波算法融合傳感器監(jiān)測數(shù)據(jù),實時地對產(chǎn)品退化狀態(tài)進行估計。文獻[8]設(shè)計和分析改進反向傳播神經(jīng)網(wǎng)絡(luò)來預(yù)測發(fā)動機剩余壽命。由于Wiener過程可以用來表示退化趨勢的非單調(diào)性和獨立增量性,這引起了國內(nèi)外學(xué)者的極大的關(guān)注[9]。目前被大量地用于對高可靠性的復(fù)雜產(chǎn)品進行性能退化建模[10-13]。

上述研究豐富了剩余壽命預(yù)測理論的研究,但卻沒有考慮到性能退化的多階段性以及將維修前的退化監(jiān)測數(shù)據(jù)用于維修后的剩余壽命預(yù)測。不同樣本的信息或同一樣本的不同階段的信息是否可以融合,是基于退化模式的一致性,因此一致性問題在實際工程中是普遍存在的。當(dāng)退化達(dá)到給定閾值時發(fā)生故障[14-15],產(chǎn)品可以通過更換一些組件來繼續(xù)工作。針對現(xiàn)有發(fā)動機壽命預(yù)測存在的問題,提出了基于多階段Wiener過程的正態(tài)總體均值和變異系數(shù)的一致性檢測方法,根據(jù)一致性檢測結(jié)果,決定是否可以融合維修前的退化監(jiān)測數(shù)據(jù),從而提高參數(shù)估計與剩余壽命預(yù)測的準(zhǔn)確性。

1 剩余壽命預(yù)測模型

1.1 模型分析

若一個非單調(diào)的隨機過程{X(t),t≥0}同時滿足以下性質(zhì):X(t)為獨立隨機增量過程

ΔX=X(t+Δt)-X(t)～N(uΔt,σ2Δt),Δt>0

式中,X(t)是關(guān)于時間t的連續(xù)函數(shù)。那么則稱該隨機過程{X(t),t≥0}是Wiener過程。

基于Wiener隨機過程的性能退化模型表示為

X(t)=ut+σB(t)

(1)

式中,X(t)為t時的性能退化量;u表示W(wǎng)iener過程的漂移系數(shù);σ表示擴散系數(shù);B(t)是標(biāo)準(zhǔn)布朗運動。

由于實際航空發(fā)動機的衰退過程具有隨機性和多階段性,本文選用多階段的Wiener隨機過程構(gòu)建衰退模型。多階段的性能退化模型可以表示為

X(t)=[X(0)+X1(t)]·I(0,t1)(t)+

[X(t1)+X2(t-t1)]·I(t1,t2)(t)+…+

[X(tn-1)+Xn(t-tn-1)]·I(tn-1,∞)(t)

(2)

式中,I(t)是示性函數(shù);ti為衰退到退化分界值的時間;X(ti)為ti時的性能退化量;Xn(t-tn-1)為第n個退化階段的退化函數(shù)。

由式(1)和式(2)可知,航空發(fā)動機多階段退化函數(shù)X(t)表示為

X(t)=[u1t+σ1B(t)]·I(0,t1)(t)+

[W1+u2(t-t1)+σ2B(t-t1)]·I(t1,t2)(t)+…+

[Wn-1+un(t-tn-1)+σnB(t-tn-1)]·I(tn-1,∞)(t)

(3)

式中,Wk為第k衰退階段的分界值。

1.2 剩余壽命分布推導(dǎo)

剩余壽命Lt′是指從t′時的退化量X(t′)到其第一次超過失效閾值W所用的時間,Lt′表示為

Lt′=inf{X(t+t′)≥W|X(t′)≤W}=

inf{X(t+t′)-X(t′)≥W-X(t′)}

(4)

式中,Lt′表示t′時刻的剩余壽命;X(t′)和X(t+t′)分別表示在t′和t′+t時刻的性能退化量;W表示失效閾值。

隨機變量Lt′的分布是逆高斯分布,其概率密度函數(shù)為

(5)

1.3 剩余壽命預(yù)測的多階段性

令ξk表示為X(t)從初始時間到其首次超過第k個階段分界值Wk的時間,ξk表示為

ξk=inf{t:X(t)≥Wk,t>0}

(6)

根據(jù)Wiener過程的齊次馬爾可夫性質(zhì),則ξk可以改寫為

ξk=ξk-1+inf{t:X′(t)≥Wk-Wk-1,t>0}

(7)

式中,X′(t)表示第k個階段的性能退化模型。

令Δξk=ξk-ξk-1,ΔWk=Wk-Wk-1,則式(7)可以改寫為

Δξk=inf{t:X′(t)≥ΔWk,t>0}

(8)

根據(jù)式(5)和式(8)可得第k個退化階段的壽命分布的概率密度函數(shù)為

f(Δξk|uk,σk)=

(9)

由式(9)可以求得壽命T的期望表示為

(10)

當(dāng)t′時的X(t′)滿足X(ξk-1)≤X(t′)≤X(ξk)時,Lt′可以表示為

Lt′=(ξk-t′)+Δξk+1+…+Δξn

(11)

2 參數(shù)估計

假設(shè)共有m個同一種類型的發(fā)動機,第i(i≤m)個發(fā)動機在第k個階段的監(jiān)測時間序列為tik|1,tik|2,…,tik|n,相對應(yīng)的退化監(jiān)測值記為yik|1,yik|2,…,yik|n。未知參數(shù)uk,σk基于第i個發(fā)動機的第k個階段的監(jiān)測值的對數(shù)似然函數(shù)為

(12)

式中,Δtik|j=tik|j-tik|j-1；Δyik|j=yik|j-yik|j-1。

根據(jù)式(12)得m個產(chǎn)品的完全對數(shù)似然估計為

(13)

極大化l(uk,σk|X,m),對u,σ求偏導(dǎo)得

(14)

3 基于多階段Wiener過程的正態(tài)總體均值的一致性檢驗

3.1 正態(tài)總體的均值區(qū)間

如果樣本數(shù)據(jù)y1,y2,…,yn滿足正態(tài)分布N(u,σ2),其中u,σ2都是未知的,那么稱u在置信水平為1-a下的置信區(qū)間為正態(tài)總體的均值區(qū)間。

在等時間間隔下Δt=ti,j-ti,j-1,Wiener過程ΔX～N(uΔt,σ2Δt),ΔX1,ΔX2,…,ΔXn是來自第k個退化階段的維修前的監(jiān)測數(shù)據(jù),數(shù)據(jù)的均值和方差計算公式為

(15)

(16)

(17)

Uk是自由度為n-1的t分布,對于給定置信水平1-a下,可得

P{-ta/2(n-1)

(18)

將式(17)代入式(18)中,此時可以得到

(19)

這時可獲得ukΔt的以1-a為置信水平的置信區(qū)間

(20)

3.2 正態(tài)總體均值一致性檢驗

若對于維修后的監(jiān)測數(shù)據(jù)序列y(t)的每個階段的均值均落在維修前監(jiān)測數(shù)據(jù)的正態(tài)總體均值ukΔt的置信區(qū)間內(nèi),則認(rèn)為這兩組監(jiān)測數(shù)據(jù)的均值是一致的。該方法考慮了復(fù)雜系統(tǒng)性能退化的多階段性和退化數(shù)據(jù)均值的隨機性,求得每階段置信水平為1-a的置信區(qū)間,從而在一定范圍驗證維修前后的發(fā)動機監(jiān)測數(shù)據(jù)的均值是否具有一致性。

當(dāng)維修前后的監(jiān)測數(shù)據(jù)通過了正態(tài)總體均值的一致性檢測之后,還需進一步檢測兩組數(shù)據(jù)的波動是否具有一致性。

4 基于多階段Wiener過程的變異系數(shù)的一致性檢測

不同樣本或同一樣本不同階段的信息是否可以融合,是由退化模式的一致性決定。針對一致性驗證的問題,在本文中,提出了一種基于Wiener過程變異系數(shù)的一致性測試方法。

4.1 構(gòu)造含變異系數(shù)的隨機變量

變異系數(shù)定義為

(21)

(22)

證明

(23)

(24)

(25)

(26)

所以Hk的概率密度函數(shù)為

(27)

證畢

4.2 一致性測試

根據(jù)證明求得的Hk的概率密度函數(shù),對基于多階段Wiener過程的變異系數(shù)的一致性測試提出如下結(jié)論:

(4)根據(jù)樣本觀測值計算hk。如果hk處于第k個階段原假設(shè)的置信區(qū)間中,則接受原假設(shè)。否則拒絕原假設(shè),即表明維修前后的退化數(shù)據(jù)存在顯著差異。

5 算法步驟

先利用正態(tài)總體均值一致性檢驗方法,檢驗維修前后退化監(jiān)測數(shù)據(jù)的各階段正態(tài)總體均值是否一致,不一致則不融合維修前的數(shù)據(jù),若一致則繼續(xù)用基于變異系數(shù)的一致性檢測方法檢測維修前后各階段的退化數(shù)據(jù)的離散程度是否一致,不一致則不融合維修前的數(shù)據(jù),若一致則融合維修前的數(shù)據(jù)。具體步驟如圖1所示。

步驟1根據(jù)式(20)對發(fā)動機維修前后的數(shù)據(jù)進行正態(tài)總體均值的一致性檢測,若不通過檢驗則不融合維修前的退化檢測數(shù)據(jù)。

步驟2若通過正態(tài)總體均值檢測,則根據(jù)式(22)進行基于變異系數(shù)的的一致性檢測,檢測維修前后退化監(jiān)測數(shù)據(jù)的離散程度。若不通過，則不融合維修前的退化檢測數(shù)據(jù)。

步驟3若通過基于變異系數(shù)的的一致性檢測,則融合維修前后的退化監(jiān)測數(shù)據(jù)。

步驟4根據(jù)進行一致性檢測后的數(shù)據(jù),依據(jù)式(3)建立多階段的性能退化模型,進行維修后發(fā)動機的剩余壽命預(yù)測。

圖1 基于一致性檢驗的剩余壽命預(yù)測流程Fig.1 Prediction process of residual lifetime based on consistency test

該一致性檢驗算法考慮了復(fù)雜系統(tǒng)性能退化的多階段性和退化數(shù)據(jù)均值的隨機性,能夠有效檢驗維修前后的數(shù)據(jù)是否具有一致性,利用變異系數(shù)的倒數(shù)構(gòu)造隨機變量,因變異系數(shù)的值受到數(shù)據(jù)的方差和平均值的影響,它描述了測定結(jié)果的偏差程度。

6 實例驗證

航空發(fā)動機的排氣溫度裕度(exhaust gas temperature margin,EGTM)是用于監(jiān)測發(fā)動機運行狀態(tài)的關(guān)鍵指標(biāo)。由文獻[18]的研究結(jié)果可知,航空發(fā)動機在壽命末期的退化速度加快,退化過程具有多階段性。

表1為18臺某型號發(fā)動機的性能退化數(shù)據(jù),前12臺為維修前的EGTM退化監(jiān)測數(shù)據(jù),后6臺為維修后的EGTM退化監(jiān)測數(shù)據(jù)。根據(jù)發(fā)動機廠商設(shè)定,EGTM的失效閾值為0℃。

表1 航空發(fā)動機的EGTM數(shù)據(jù)Table 1 EGTM data of aeroengines

圖2為這18臺發(fā)動機的EGTM的退化曲線,根據(jù)文獻[19]中的啟發(fā)式分割算法,計算性能退化分界點為27.1 ℃。因此,用兩階段的Wiener過程建立退化建模。

圖2 EGTM退化過程Fig.2 Degradation process of EGTM

6.1 基于Wiener過程的正態(tài)總體均值的一致性檢驗

表估計結(jié)果Table 2 Estimated result of

表3 維修后EGTM數(shù)據(jù)的Table of EGTM data after maintenance

6.2 基于Wiener過程的變異系數(shù)的一致性檢測

根據(jù)1～12號發(fā)動機的數(shù)據(jù),利用式(14),進行性能退化量總體兩階段的u和σ2的估計,估計的結(jié)果如表4所示。

表4 維修前的u和σ2的估計值Table 4 Estimated value of u and σ2 before maintenance

表5 編號16發(fā)動機的EGTM數(shù)據(jù)的值Table 5 Values of and S2 of EGTM data for No.16 engine

表6 編號18發(fā)動機的EGTM數(shù)據(jù)的值Table 6 Values of and S2 of EGTM data for No.18 engine

編號16發(fā)動機的第一階段的h的值為

編號16發(fā)動機的第二階段的h的值為

編號18發(fā)動機的第一階段的h的值為

編號18發(fā)動機的第二階段的h的值為

根據(jù)式(27)的Hk的概率密度函數(shù),可求得編號為16發(fā)動機的一致性檢測的第一階段置信區(qū)間為[-2.603,2.00],-2.348 8在該置信區(qū)間內(nèi),一致性檢測的第二階段的置信區(qū)間為[-4.75,2.07],-1.182 5在該置信區(qū)間內(nèi),所以編號為16的發(fā)動機的維修前后的退化監(jiān)測數(shù)據(jù)的波動沒有顯著差異。

編號為18發(fā)動機的一致性檢測的第一階段置信區(qū)間為[-2.594,2.00],-0.673 5在該置信區(qū)間內(nèi),一致性檢測的第二階段的置信區(qū)間為[-4.75,2.07],-0.637 4在該置信區(qū)間內(nèi),所以編號為18的發(fā)動機的維修前后的退化監(jiān)測數(shù)據(jù)的波動沒有顯著差異。

6.3 剩余壽命預(yù)測

根據(jù)第6.1節(jié)和第6.2節(jié)的計算結(jié)果可知,編號為16和18的發(fā)動機的性能退化路徑與維修前的12臺發(fā)動機性能退化路徑?jīng)]有顯著差異,所以可以融合維修前的退化監(jiān)測數(shù)據(jù)。下面以16號發(fā)動機為例,進行融合維修前數(shù)據(jù)與不融合維修前數(shù)據(jù)的剩余壽命預(yù)測。

表7 16號發(fā)動機融合EGTM數(shù)據(jù)后的u和σ2估計值Table 7 Estimated value of u and σ2 for NO.18 engine after integrating EGTM data

融合維修前EGTM數(shù)據(jù),利用融合后的EGTM數(shù)據(jù)對衰退模型參數(shù)進行估計,并對發(fā)動機的剩余壽命進行預(yù)測,預(yù)測結(jié)果如圖3和圖4所示。

圖3 融合EGTM數(shù)據(jù)后的第一階段剩余壽命概率密度分布圖Fig.3 Probability density distribution of remaining life for the first stage after integrating EGTM data

圖4 融合EGTM數(shù)據(jù)后的第二階段剩余壽命概率密度分布圖Fig.4 Probability density distribution of remaining life for the second stage after integrating EGTM data

為直觀反映出融合了維修前的監(jiān)測數(shù)據(jù)與未融合維修前的退化監(jiān)測數(shù)據(jù)在航空發(fā)動機剩余壽命預(yù)測問題中的精確性,定義相對誤差(relative error,RE)指標(biāo)為

(28)

式中,RUL為預(yù)測的剩余壽命;T為實際的剩余壽命。

兩種模型的誤差計算結(jié)果如表8和圖5所示。

表8 預(yù)測誤差對比Table 8 Predictive error comparison

圖5 預(yù)測誤差對比Fig.5 Predictive error comparison

從表8可知,融合維修前退化數(shù)據(jù)的壽命預(yù)測的平均誤差小于未融合維修前退化數(shù)據(jù)的壽命預(yù)測。由圖5可知,隨著發(fā)動機運轉(zhuǎn)時間的增加,融合維修前監(jiān)測數(shù)據(jù)的發(fā)動機剩余壽命預(yù)測的結(jié)果更精確。

產(chǎn)生以上誤差對比結(jié)果的原因是當(dāng)通過退化數(shù)據(jù)一致性檢驗后,此時的剩余壽命預(yù)測能夠融合維修前的監(jiān)測數(shù)據(jù),使訓(xùn)練模型所用的數(shù)據(jù)量變多,排除數(shù)據(jù)的偶然性,求得的模型參數(shù)更精確,從而使剩余壽命預(yù)測的結(jié)果更加準(zhǔn)確。

7 結(jié) 論

(1) 提出了基于多階段隨機Wiener過程的正態(tài)總體均值的數(shù)據(jù)一致性檢驗方法,用于驗證維修前后退化監(jiān)測數(shù)據(jù)的正態(tài)總體均值是否具有一致性。

(2) 根據(jù)變異系數(shù)的倒數(shù)構(gòu)造H隨機變量,并導(dǎo)出抽樣分布,基于抽樣分布提出了一種通用的一致性檢驗方法。該方法用于驗證維修前后退化監(jiān)測數(shù)據(jù)的離散程度是否具有一致性。

(3) 當(dāng)維修前后的退化監(jiān)測數(shù)據(jù)通過本文提出的兩種一致性檢驗后,就可將維修前的退化監(jiān)測數(shù)據(jù)用于維修后發(fā)動機的剩余壽命預(yù)測,從而提高剩余壽命預(yù)測的精度。