蔡忠義, 郭建勝, 陳云翔, 董驍雄, 項(xiàng)華春

(空軍工程大學(xué)裝備管理與無人機(jī)工程學(xué)院, 陜西 西安 710051)

0 引 言

隨著狀態(tài)監(jiān)測技術(shù)的發(fā)展,以設(shè)備退化狀態(tài)監(jiān)控信息為依據(jù),及時(shí)準(zhǔn)確預(yù)測設(shè)備的剩余壽命是制定設(shè)備最優(yōu)維護(hù)方案、備件訂購策略和健康管理的重要依據(jù)[1]。剩余壽命預(yù)測主要是求解各監(jiān)測時(shí)刻上目標(biāo)產(chǎn)品的剩余壽命分布函數(shù)。基于Wiener過程的剩余壽命預(yù)測分為兩個(gè)環(huán)節(jié):

一是進(jìn)行退化建模,利用同類產(chǎn)品的退化信息,求解出反映總體退化特征的固定系數(shù)估計(jì)值和反映個(gè)體退化特征的隨機(jī)系數(shù)先驗(yàn)值。由于產(chǎn)品退化過程中受隨機(jī)不確性、測量不確定性等因素的影響,在退化建模時(shí)需要考慮以下退化特征:

(1)非線性特征。對(duì)具有非線性特征的退化數(shù)據(jù)處理,一般可采用兩種方法進(jìn)行處理。一是對(duì)非線性數(shù)據(jù)進(jìn)行線性化處理。Whtimore[2]首次提出時(shí)間尺度模型,將非線性退化數(shù)據(jù)轉(zhuǎn)化為線性退化數(shù)據(jù),給出產(chǎn)品概率密度函數(shù)(probability density function,PDF)表達(dá)式。由于計(jì)算過程簡單,使得該方法得以廣泛應(yīng)用[3]。二是直接建立非線性Wiener退化過程建模,給出退化產(chǎn)品壽命的PDF近似表達(dá)式。文獻(xiàn)[4-5]提出了非線性Wiener退化過程建模,給出了相應(yīng)的PDF近似表達(dá)式。相關(guān)的模型誤設(shè)研究[6],表明后者進(jìn)行建模評(píng)估的精度要高于前者,但后者建模和計(jì)算過程較為復(fù)雜。

(2)隨機(jī)效應(yīng)。隨機(jī)效應(yīng)是指同類產(chǎn)品中個(gè)體差異在退化過程中的具體表現(xiàn)。將描述隨機(jī)效應(yīng)的隨機(jī)變量稱為隨機(jī)系數(shù)。一般做法是將Wiener過程參數(shù)進(jìn)行隨機(jī)化處理。文獻(xiàn)[7]將漂移系數(shù)和擴(kuò)散系數(shù)都看作是隨機(jī)變量并作為先驗(yàn)分布,通過貝葉斯推斷,得到隨機(jī)系數(shù)的后驗(yàn)均值,但由于涉及復(fù)雜推導(dǎo),并未給出含有雙隨機(jī)系數(shù)的PDF表達(dá)式。文獻(xiàn)[8-9]認(rèn)為同類產(chǎn)品中個(gè)體差異主要受Wiener過程中漂移系數(shù)的影響,只將漂移系數(shù)看作為正態(tài)分布隨機(jī)變量。這種做法能較好地處理總體退化特征與個(gè)體退化差異之間的關(guān)系,所建的模型具有良好的模型擬合性。

(3)測量誤差。由于測量過程中受噪聲、擾動(dòng)、不穩(wěn)定測量設(shè)備的影響,測量數(shù)據(jù)中含有一定誤差。一般做法是將測量誤差看作為均值為零的正態(tài)隨機(jī)變量且與實(shí)際退化過程之間相互獨(dú)立[10]。文獻(xiàn)[11-12]建立了同時(shí)考慮隨機(jī)效應(yīng)和測量誤差的線性Wiener退化模型,給出了相應(yīng)的PDF近似表達(dá)式。文獻(xiàn)[13]綜合分析了帶非線性特征、隨機(jī)效應(yīng)和測量誤差的模型誤設(shè)問題,指出模型中考慮測量誤差可以有效提高評(píng)估精度。

二是進(jìn)行剩余壽命預(yù)測建模。剩余壽命建模實(shí)際上給出目標(biāo)產(chǎn)品剩余壽命的PDF。文獻(xiàn)[14]給出了含有測量誤差的線性Wiener過程的目標(biāo)產(chǎn)品剩余壽命的PDF近似表達(dá)式。對(duì)于帶測量誤差的剩余壽命預(yù)測的關(guān)鍵是對(duì)剩余壽命PDF中隨機(jī)系數(shù)或當(dāng)前實(shí)際退化狀態(tài)進(jìn)行實(shí)時(shí)更新,以體現(xiàn)目標(biāo)產(chǎn)品的個(gè)性特征。一般有以下兩種方法:

(1)隨機(jī)系數(shù)的貝葉斯更新方法。文獻(xiàn)[15]首次將貝葉斯推斷方法用于剩余壽命分布模型中隨機(jī)系數(shù)的更新,給出了先驗(yàn)信息與目標(biāo)產(chǎn)品退化信息之間的更新公式。文獻(xiàn)[13]應(yīng)用貝葉斯推斷方法,進(jìn)一步推導(dǎo)出了考慮測量誤差和隨機(jī)效應(yīng)的非線性Wiener過程的目標(biāo)產(chǎn)品剩余壽命中隨機(jī)系數(shù)的更新公式。這種方法充分利用了同類產(chǎn)品退化先驗(yàn)信息,使得目標(biāo)產(chǎn)品剩余壽命預(yù)測結(jié)果能反映其個(gè)體退化特征,但由于目標(biāo)產(chǎn)品實(shí)際退化狀態(tài)未知,可能存在一定誤差。

(2)隨機(jī)系數(shù)和實(shí)際狀態(tài)的聯(lián)合更新方法。常用的更新方法有隨機(jī)濾波模型[16]、隱馬爾可夫模型[17]等。其中,隨機(jī)濾波模型將隨機(jī)系數(shù)、當(dāng)前實(shí)際狀態(tài)都作為隱含狀態(tài),建立了觀測過程與隨機(jī)系數(shù)、當(dāng)前實(shí)際狀態(tài)之間的狀態(tài)空間模型,實(shí)現(xiàn)對(duì)隱含狀態(tài)的估計(jì)和更新。這種方法具有較高的預(yù)測精度,但計(jì)算過程較為復(fù)雜。

因此,本文針對(duì)SSADT下具有非線性、隨機(jī)效應(yīng)、測量誤差等特征的退化數(shù)據(jù),建立加速退化模型,給出總體觀測增量數(shù)據(jù)的輪廓對(duì)數(shù)似然函數(shù),得到模型中固定系數(shù)估計(jì)值和隨機(jī)系數(shù)先驗(yàn)值;建立與目標(biāo)產(chǎn)品個(gè)體退化特征相匹配的剩余壽命PDF,引入隨機(jī)系數(shù)更新貝葉斯方法,求出剩余壽命模型中隨機(jī)系數(shù)后驗(yàn)值。

1 步進(jìn)加速退化建模

1.1 模型假設(shè)

(1)受試產(chǎn)品的關(guān)鍵性能參數(shù)僅有一個(gè),其退化過程X(t)可用Wiener過程進(jìn)行描述

X(t)=λΛ(t;θ)+σBB(t)

(1)

式中,λ為漂移系數(shù);σB為擴(kuò)散系數(shù);B(t)為標(biāo)準(zhǔn)布朗運(yùn)動(dòng);Λ(t;θ)為時(shí)間t的連續(xù)非減函數(shù)(θ為未知參數(shù)向量),表征非線性特征。

(2)選擇溫度為加速應(yīng)力、Arrhenius模型為加速模型。

(3)受試產(chǎn)品在各加速應(yīng)力下發(fā)生的退化不可逆且退化機(jī)理不變,滿足累積退化模型要求;

(4)受試產(chǎn)品的性能退化數(shù)據(jù)在觀測中會(huì)產(chǎn)生測量誤差(在壽命預(yù)測性能約束的可行域內(nèi)),其觀測過程Y(t)表示為

Y(t)=X(t)+ε=λΛ(t;θ)+σBB(t)+ε

(2)

式中,ε表示測量誤差,一般認(rèn)為ε～N(0,σ2)且與λ,σB相互獨(dú)立。

1.2 步進(jìn)加速退化模型

由SSADT過程可知,當(dāng)前加速應(yīng)力下產(chǎn)品觀測初值是前一應(yīng)力下觀測末值。根據(jù)式(2),建立各加速應(yīng)力下觀測數(shù)據(jù)Y(t)的關(guān)系模型為

yi,k=

(3)

1.3 加速模型

一般認(rèn)為Wiener過程的漂移系數(shù)與應(yīng)力有關(guān)。采用Arrhenius模型來描述漂移系數(shù)與溫度應(yīng)力之間的關(guān)系，即

λi=aexp(-b/Si)

(4)

式中,a和b為待定參數(shù);Si為第i個(gè)應(yīng)力;λi為應(yīng)力Si下產(chǎn)品的漂移系數(shù)。

(5)

2 參數(shù)估計(jì)

根據(jù)多元Wiener過程性質(zhì),觀測增量數(shù)據(jù)向量Δyj服從n元正態(tài)分布[4]。

(1)當(dāng)給定λi時(shí),則有

(6)

(7)

(8)

(9)

根據(jù)式(9),基于m個(gè)產(chǎn)品觀測數(shù)據(jù)集Δy={Δy1,Δy2,…,Δym},建立未知參數(shù)集Θ的對(duì)數(shù)似然函數(shù)為

(10)

(11)

(12)

(13)

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

3.1 剩余壽命分布函數(shù)

將目標(biāo)產(chǎn)品在當(dāng)前觀測時(shí)刻ti,k處的剩余壽命Li,k定義為

Li,k=inf{li,k:X(ti,k+li,k)≥w|X(ti,k)

(14)

(15)

式中

β(li,k)=(li,k+ti,k)c-(ti,k)c-c(li,k+ti,k)c-1li,k

(16)

φ(li,k)=(li,k+ti,k)c-(ti,k)c

(17)

ε,ε～N(0,σ2)

(18)

?

(19)

當(dāng)前時(shí)刻ti,k處目標(biāo)產(chǎn)品剩余壽命的期望為

(20)

3.2 隨機(jī)系數(shù)在線更新

(21)

其中

(22)

(23)

(24)

(25)

4 仿真實(shí)例

本文依據(jù)文獻(xiàn)[18]給出的激光器實(shí)測退化數(shù)據(jù),仿真得到激光器在步進(jìn)加速應(yīng)力下的退化數(shù)據(jù),以驗(yàn)證所提出方法的正確性。文獻(xiàn)[8，19]也是基于該激光器數(shù)據(jù)仿真得到的加速退化數(shù)據(jù),研究加速退化建模問題。

已知該激光器的性能參數(shù)為工作電流(單位:mA),正常工作溫度為25 ℃。激光器的工作電流會(huì)隨著時(shí)間的延長而出現(xiàn)遞增的退化現(xiàn)象。圖1給出了15組激光器的實(shí)測退化數(shù)據(jù)。

圖1 15組激光器的實(shí)測退化數(shù)據(jù)Fig.1 15 sets of lasers measured degradation data

已知該激光器受溫度影響較為敏感,其性能退化速率與溫度之間的關(guān)系符合Arrhenius模型。引用文獻(xiàn)[8]中給出的激光器退化模型參數(shù)真值如表1所示。引入測量誤差,設(shè)定σ=0.08。采用蒙特卡羅仿真方法,設(shè)定步進(jìn)應(yīng)力為25 ℃、50 ℃、75 ℃,每個(gè)應(yīng)力下測量5次,測量間隔為150 h,得到8組仿真的步進(jìn)應(yīng)力下實(shí)際退化數(shù)據(jù)(不帶測量誤差,此時(shí)σ=0)和觀測數(shù)據(jù)見圖2、圖3(圖中藍(lán)線為目標(biāo)產(chǎn)品的退化軌跡)。

圖2 仿真的SSADT實(shí)際退化數(shù)據(jù)Fig.2 Simulated SSADT actual degradation data

圖3 仿真的SSADT觀測數(shù)據(jù)Fig.3 Simulated SSADT observed data

利用仿真得到的激光器SSADT觀測數(shù)據(jù)來驗(yàn)證本文所提出方法的正確性和優(yōu)勢。將文獻(xiàn)[19]中未考慮測量誤差的非線性加速退化建模方法記為M1;將帶測量誤差的線性加速退化建模方法(即c=1)記為M2;本文提出的帶測量誤差的非線性加速退化建模方法記為M3(由此可見,M1和M2是M3的特殊情況)。采用赤池信息準(zhǔn)則(AIC)和均方誤差(MSE)來判別各模型之間的擬合優(yōu)劣性。

4.1 退化模型的參數(shù)估計(jì)

表1 不同模型的參數(shù)估計(jì)結(jié)果

4.2 隨機(jī)系數(shù)在線更新

以SSADT中第8個(gè)樣本為目標(biāo)產(chǎn)品,其在第2 250 h處的實(shí)際退化數(shù)值為11.81 mA。為了驗(yàn)證本文在剩余壽命中引入隨機(jī)系數(shù)更新貝葉斯方法的有效性,假設(shè)該激光器的退化失效閾值為11.81 mA,在2 250 h處剛好失效。將M3的剩余壽命預(yù)測模型中加入隨機(jī)系數(shù)在線更新方法記為M4。通過對(duì)比M3、M4與真實(shí)值剩余壽命預(yù)測結(jié)果來驗(yàn)證M4的優(yōu)勢。

表2 觀測末期隨機(jī)系數(shù)后驗(yàn)值

由表2可知,M4計(jì)算出的隨機(jī)系數(shù)均值和方差的后驗(yàn)值介于其先驗(yàn)值與真實(shí)值之間;而M3計(jì)算出的隨機(jī)系數(shù)均值和方差的后驗(yàn)值離真值較遠(yuǎn)。這是因?yàn)镸4在M3的基礎(chǔ)上引入隨機(jī)系數(shù)更新方法,能更好地反映出目標(biāo)產(chǎn)品的個(gè)體退化特征,具有更好的估計(jì)優(yōu)勢。

4.3 目標(biāo)產(chǎn)品剩余壽命預(yù)測

分別將M3、M4的隨機(jī)系數(shù)后驗(yàn)值以及真值,代入式(19)和式(20),計(jì)算出目標(biāo)產(chǎn)品的剩余壽命PDF如圖4所示。

由圖4可知,M4的剩余壽命PDF比M3的剩余壽命PDF更窄且更接近真值的PDF;同時(shí)M4的目標(biāo)產(chǎn)品剩余壽命點(diǎn)估計(jì)值較M3更接近于真值的剩余壽命點(diǎn)估計(jì)值。這是因?yàn)镸4的剩余壽命PDF采用隨機(jī)系數(shù)更新方法,使得目標(biāo)產(chǎn)品在剩余壽命預(yù)測過程中不斷以自身觀測數(shù)據(jù)來更新隨機(jī)系數(shù),使得預(yù)測結(jié)果更逼近于其真實(shí)退化過程。說明本文所提出的帶測量誤差的退化模型和基于隨機(jī)系數(shù)貝葉斯更新的剩余壽命預(yù)測模型具有更好的預(yù)測精度。

圖4 不同方法下剩余壽命的PDFFig.4 Remaining lifetime PDF by different methods

5 結(jié) 論

(1)本文建立了帶測量誤差的非線性步進(jìn)加速退化模型,給出了基于Fminsearch函數(shù)的輪廓似然函數(shù)求解方法,得到了固定系數(shù)估計(jì)值和隨機(jī)系數(shù)先驗(yàn)值,為解決退化模型中多參數(shù)估計(jì)提供了可行途徑;

(2)引入隨機(jī)系數(shù)貝葉斯更新方法,推導(dǎo)出了目標(biāo)產(chǎn)品剩余壽命的PDF近似表達(dá)式,利用目標(biāo)產(chǎn)品退化信息來更新剩余壽命的隨機(jī)系數(shù)后驗(yàn)值,實(shí)現(xiàn)了目標(biāo)產(chǎn)品剩余壽命預(yù)測結(jié)果在線更新;

(3)結(jié)合仿真實(shí)例分析,與現(xiàn)有算法對(duì)比分析,驗(yàn)證了本文所建模型在估計(jì)誤差和預(yù)測精度上具有一定優(yōu)勢。