本發(fā)明涉及設(shè)備剩余壽命預(yù)測(cè)技術(shù)領(lǐng)域，具體涉及一種梯形噪聲分布的指數(shù)模型機(jī)械設(shè)備剩余壽命預(yù)測(cè)方法。

背景技術(shù)：

隨著生產(chǎn)制造技術(shù)的快速發(fā)展和人類探索自然領(lǐng)域的不斷擴(kuò)展，許多設(shè)備變的越來(lái)越復(fù)雜。這些設(shè)備由于機(jī)械的復(fù)雜性和各種運(yùn)行因素(如磨損、外部沖擊、負(fù)載、運(yùn)行環(huán)境)的影響，其性能及健康狀態(tài)將會(huì)發(fā)生不可避免的退化，進(jìn)而造成最終的失效。對(duì)于實(shí)際工程設(shè)備，一旦發(fā)生由失效引起的事故，所造成的人員財(cái)產(chǎn)損失甚至環(huán)境破壞往往是不可估量的。因而如何有效評(píng)估機(jī)械設(shè)備的運(yùn)行狀態(tài)，防止因機(jī)械設(shè)備失效造成事故是當(dāng)前迫切需要解決的問(wèn)題。從機(jī)械設(shè)備開(kāi)始退化到完全失效，要經(jīng)歷一個(gè)逐步退化的過(guò)程。因而如果能在退化初期，根據(jù)監(jiān)測(cè)信息，及時(shí)發(fā)現(xiàn)異常、預(yù)測(cè)其剩余使用壽命，并據(jù)此確定對(duì)設(shè)備實(shí)施維護(hù)的最佳時(shí)機(jī)，對(duì)于切實(shí)保障復(fù)雜設(shè)備的運(yùn)行安全性、可靠性與經(jīng)濟(jì)性具有重要的意義。

指數(shù)模型剩余壽命預(yù)測(cè)方法由美國(guó)普渡大學(xué)gebraeel等人提出，試圖采用指數(shù)模型對(duì)衰退趨勢(shì)進(jìn)行描述，并根據(jù)觀測(cè)數(shù)據(jù)對(duì)模型參數(shù)進(jìn)行評(píng)估，以預(yù)測(cè)設(shè)備未來(lái)時(shí)刻健康衰退趨勢(shì)和剩余壽命。為了提高參數(shù)評(píng)估的準(zhǔn)確性，清華大學(xué)司小勝等人對(duì)指數(shù)模型進(jìn)行改進(jìn)，采用期望最大化和貝葉斯更新相結(jié)合的方法對(duì)模型參數(shù)進(jìn)行評(píng)估，得到了較好的參數(shù)評(píng)估效果。在指數(shù)預(yù)測(cè)模型中，模型準(zhǔn)確性是影響模型預(yù)測(cè)精度的關(guān)鍵因素，以上研究工作都假定模型噪聲項(xiàng)服從0初值的標(biāo)準(zhǔn)布朗運(yùn)動(dòng)分布，即為三角形噪聲分布，而根據(jù)大量的數(shù)據(jù)觀測(cè)，模型的噪聲項(xiàng)在壽命預(yù)測(cè)初始時(shí)刻都有著較大的波動(dòng)，并且隨著時(shí)間的增加，波動(dòng)幅度越來(lái)越大。因此傳統(tǒng)指數(shù)模型的噪聲項(xiàng)假定與實(shí)際情況不相符，導(dǎo)致模型預(yù)測(cè)精度和可靠性降低。

技術(shù)實(shí)現(xiàn)要素：

為了克服上述現(xiàn)有技術(shù)的缺點(diǎn)，本發(fā)明的目的在于提供了一種梯形噪聲分布的指數(shù)模型機(jī)械設(shè)備剩余壽命預(yù)測(cè)方法，提高指數(shù)模型的預(yù)測(cè)精度和可靠性。

為了達(dá)到上述目的，本發(fā)明采取的技術(shù)方案為：

一種梯形噪聲分布的指數(shù)模型機(jī)械設(shè)備剩余壽命預(yù)測(cè)方法，包括以下步驟：

第1步，建立梯形噪聲分布的指數(shù)模型:

y＝ae^bt+ε(1)

其中，y為指數(shù)模型健康狀態(tài)指標(biāo)，a,b為指數(shù)模型趨勢(shì)項(xiàng)參數(shù)，ε為指數(shù)模型噪聲項(xiàng)，服從均值為0、方差為的正態(tài)分布，t為從擬合起始時(shí)刻t0開(kāi)始的工作時(shí)間，為指數(shù)模型噪聲項(xiàng)初值，σ²為指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)；

第2步，實(shí)時(shí)監(jiān)測(cè)并采集機(jī)械設(shè)備中軸承、齒輪或轉(zhuǎn)子的振動(dòng)信號(hào)，計(jì)算健康狀態(tài)指標(biāo)序列，確定擬合起始時(shí)刻t0；

第3步，將健康狀態(tài)指標(biāo)序列根據(jù)擬合起始時(shí)刻t0分為兩個(gè)子序列：平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列y′＝[y1′,y′2,y3′,…,y′k]和指數(shù)模型健康狀態(tài)指標(biāo)觀測(cè)值序列y＝[y1,y2,y3,…,yk]；平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列y′＝[y1′,y′2,y3′,…,y′k]為擬合起始時(shí)刻t0前的健康狀態(tài)指標(biāo)子序列，yi′為ti′時(shí)刻的健康狀態(tài)指標(biāo)，ti′為工作時(shí)間，其中ti′＜t0；指數(shù)模型健康狀態(tài)指標(biāo)觀測(cè)值序列y＝[y1,y2,y3,…,yk]為擬合起始時(shí)刻t0后的健康狀態(tài)指標(biāo)子序列，yi為ti+t0時(shí)刻的健康狀態(tài)指標(biāo)，ti為從擬合起始時(shí)刻t0開(kāi)始的工作時(shí)間，其中ti＜tk，tk為從擬合起始時(shí)刻t0開(kāi)始的最后工作時(shí)刻，k為觀測(cè)時(shí)刻數(shù)；

第4步，估計(jì)梯形噪聲分布的指數(shù)模型中的參數(shù)，進(jìn)而得出梯形噪聲分布的指數(shù)模型；

第5步，由梯形噪聲分布的指數(shù)模型，得出在最后工作時(shí)刻tk后，各時(shí)刻預(yù)測(cè)健康狀態(tài)指標(biāo)yp服從均值為ae^bt，方差為σ²(t-tk)的正態(tài)分布，即

yp～n(ae^bt,σ²(t-tk))(12)

其中t＞tk；

第6步，令t＝tk+δt,tk+2δt,…，根據(jù)yp的概率密度函數(shù)進(jìn)行ns次隨機(jī)采樣，得到ns個(gè)指數(shù)模型健康狀態(tài)指標(biāo)預(yù)測(cè)序列yp＝[yp1,yp2,yp3,…]，即ns條退化軌跡，每條退化軌跡的失效時(shí)刻eol為指數(shù)模型健康狀態(tài)指標(biāo)預(yù)測(cè)序列內(nèi)第一次預(yù)測(cè)指標(biāo)yp超過(guò)失效閾值yh的時(shí)刻，即

eol＝inf(yp＞yh)(13)

剩余壽命rul為失效時(shí)刻eol與從擬合起始時(shí)刻t0開(kāi)始的最后工作時(shí)刻tk的差，即

rul＝eol-tk(14)

對(duì)ns條退化軌跡的剩余壽命rul進(jìn)行統(tǒng)計(jì)，繪制剩余壽命rul的頻率分布直方圖并且擬合出對(duì)應(yīng)的概率密度函數(shù)，確定95％剩余壽命置信區(qū)間，給出剩余壽命估計(jì)值。

所述的第4步的具體步驟為：

4.1)估計(jì)指數(shù)模型噪聲項(xiàng)初值

首先對(duì)平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列進(jìn)行中值濾波，得到平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列趨勢(shì)；平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列減去中值濾波得到的平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列趨勢(shì)得到平穩(wěn)運(yùn)行階段的健康狀態(tài)指標(biāo)序列噪聲并計(jì)算其方差作為噪聲項(xiàng)初值的估計(jì)；

4.2)用極大似然法估計(jì)指數(shù)模型趨勢(shì)項(xiàng)參數(shù)a,b：

因?yàn)橹笖?shù)模型噪聲項(xiàng)

所以指數(shù)模型健康狀態(tài)指標(biāo)

某時(shí)刻指數(shù)模型健康狀態(tài)指標(biāo)的概率密度函數(shù):

其中，θ＝(a,b)；

似然函數(shù)l(θ)：

對(duì)數(shù)似然函數(shù)l(θ)：

所以最大化似然函數(shù)即為最小化殘差和：

通過(guò)數(shù)值優(yōu)化算法求得指數(shù)模型趨勢(shì)項(xiàng)參數(shù)a,b的估計(jì)值；

4.3)估計(jì)指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)σ²：

將式(2)指數(shù)模型噪聲項(xiàng)變換為

ε～n(0,σ²(t+t0))(8)

其中

因此

根據(jù)(9)和(11)兩式，求解得出指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)σ²的估計(jì)值，進(jìn)而得出梯形噪聲分布的指數(shù)模型。

本發(fā)明的有益效果：本發(fā)明解決了傳統(tǒng)指數(shù)預(yù)測(cè)模型中噪聲項(xiàng)假定與實(shí)際情況不相符的問(wèn)題，將原有的三角形噪聲分布改為梯形噪聲分布，即增加了指數(shù)模型噪聲項(xiàng)的初始值，并且以平穩(wěn)運(yùn)行階段的健康狀態(tài)指標(biāo)序列噪聲作為指數(shù)模型噪聲項(xiàng)的初始值。相比傳統(tǒng)指數(shù)模型，該模型能更好地反映機(jī)械設(shè)備的真實(shí)衰退趨勢(shì)，提高了剩余壽命的預(yù)測(cè)精度和可靠性。

附圖說(shuō)明

圖1為本發(fā)明方法的流程圖。

圖2為測(cè)試軸承全壽命期內(nèi)的振動(dòng)信號(hào)。

圖3為測(cè)試軸承有效值指標(biāo)和回歸起始時(shí)刻選擇結(jié)果圖。

圖4為在時(shí)刻tk＝5000s時(shí)，梯形噪聲分布的指數(shù)模型與傳統(tǒng)指數(shù)模型噪聲項(xiàng)擬合對(duì)比，圖(a)為梯形噪聲分布的指數(shù)模型噪聲項(xiàng)擬合，圖(b)為傳統(tǒng)指數(shù)模型噪聲項(xiàng)擬合。

圖5為梯形噪聲分布的指數(shù)模型在時(shí)刻tk＝5000s時(shí)的預(yù)測(cè)結(jié)果圖，圖(a)為有效值預(yù)測(cè)趨勢(shì)圖，圖(b)為剩余壽命預(yù)測(cè)分布圖。

圖6為傳統(tǒng)指數(shù)模型在時(shí)刻tk＝5000s時(shí)的預(yù)測(cè)結(jié)果圖,圖(a)為有效值預(yù)測(cè)趨勢(shì)圖，圖(b)為剩余壽命預(yù)測(cè)分布圖。

圖7為梯形噪聲分布的指數(shù)模型從預(yù)測(cè)起始時(shí)刻開(kāi)始的剩余壽命估計(jì)和預(yù)測(cè)95％置信區(qū)間結(jié)果圖。

圖8為傳統(tǒng)指數(shù)模型從預(yù)測(cè)起始時(shí)刻開(kāi)始的剩余壽命估計(jì)和預(yù)測(cè)95％置信區(qū)間結(jié)果圖。

具體實(shí)施方式

下面結(jié)合附圖對(duì)本發(fā)明做進(jìn)一步詳細(xì)描述。

如圖1所示，一種梯形噪聲分布的指數(shù)模型機(jī)械設(shè)備剩余壽命預(yù)測(cè)方法，包括以下步驟：

第1步，建立梯形噪聲分布的指數(shù)模型:

y＝ae^bt+ε(1)

其中，y為指數(shù)模型健康狀態(tài)指標(biāo)，a,b為指數(shù)模型趨勢(shì)項(xiàng)參數(shù)，ε為指數(shù)模型噪聲項(xiàng)，服從均值為0、方差為的正態(tài)分布，t為從擬合起始時(shí)刻t0開(kāi)始的工作時(shí)間，為指數(shù)模型噪聲項(xiàng)初值，σ²為指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)，因?yàn)橹笖?shù)模型噪聲項(xiàng)的方差以為初值，隨時(shí)間t線性增加，呈梯形分布，因此稱該指數(shù)模型為梯形噪聲分布的指數(shù)模型；

第2步，實(shí)時(shí)監(jiān)測(cè)并采集機(jī)械設(shè)備中軸承、齒輪或轉(zhuǎn)子的振動(dòng)信號(hào)，計(jì)算健康狀態(tài)指標(biāo)序列，確定擬合起始時(shí)刻t0；

第3步，將健康狀態(tài)指標(biāo)序列根據(jù)擬合起始時(shí)刻t0分為兩個(gè)子序列：平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列y′＝[y1′,y′2,y3′,…,y′k]和指數(shù)模型健康狀態(tài)指標(biāo)觀測(cè)值序列y＝[y1,y2,y3,…,yk]；平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列y′＝[y1′,y′2,y3′,…,y′k]為擬合起始時(shí)刻t0前的健康狀態(tài)指標(biāo)子序列，yi′為ti′時(shí)刻的健康狀態(tài)指標(biāo)，ti′為工作時(shí)間，其中ti′＜t0；指數(shù)模型健康狀態(tài)指標(biāo)觀測(cè)值序列y＝[y1,y2,y3,…,yk]為擬合起始時(shí)刻t0后的健康狀態(tài)指標(biāo)子序列，yi為ti+t0時(shí)刻的健康狀態(tài)指標(biāo)，ti為從擬合起始時(shí)刻t0開(kāi)始的工作時(shí)間，其中ti＜tk，tk為從擬合起始時(shí)刻t0開(kāi)始的最后工作時(shí)刻，k為觀測(cè)時(shí)刻數(shù)；

第4步，估計(jì)梯形噪聲分布的指數(shù)模型中的參數(shù)，進(jìn)而得出梯形噪聲分布的指數(shù)模型，具體步驟為：

4.1)估計(jì)指數(shù)模型噪聲項(xiàng)初值

為了消除機(jī)械設(shè)備平穩(wěn)運(yùn)行階段的工況、負(fù)載等因素變化對(duì)噪聲項(xiàng)初值估計(jì)產(chǎn)生的干擾，首先對(duì)平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列進(jìn)行中值濾波，得到平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列趨勢(shì)；平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列減去中值濾波得到的平穩(wěn)運(yùn)行階段健康狀態(tài)指標(biāo)序列趨勢(shì)得到平穩(wěn)運(yùn)行階段的健康狀態(tài)指標(biāo)序列噪聲并計(jì)算其方差作為噪聲項(xiàng)初值的估計(jì)；

4.2)用極大似然法估計(jì)指數(shù)模型趨勢(shì)項(xiàng)參數(shù)a,b：

因?yàn)橹笖?shù)模型噪聲項(xiàng)

所以指數(shù)模型健康狀態(tài)指標(biāo)

某時(shí)刻指數(shù)模型健康狀態(tài)指標(biāo)的概率密度函數(shù):

其中，θ＝(a,b)；

似然函數(shù)l(θ)：

對(duì)數(shù)似然函數(shù)l(θ)：

所以最大化似然函數(shù)即為最小化殘差和：

通過(guò)數(shù)值優(yōu)化算法求得指數(shù)模型趨勢(shì)項(xiàng)參數(shù)a,b的估計(jì)值；

4.3估計(jì)指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)σ²：

將式(2)指數(shù)模型噪聲項(xiàng)變換為

ε～n(0,σ²(t+t0))(8)

其中

因此

根據(jù)(9)和(11)兩式，求解得出指數(shù)模型噪聲項(xiàng)方差擴(kuò)散系數(shù)σ²的估計(jì)值，進(jìn)而得出梯形噪聲分布的指數(shù)模型；

第5步，由梯形噪聲分布的指數(shù)模型，得出在最后工作時(shí)刻tk后，各時(shí)刻預(yù)測(cè)健康狀態(tài)指標(biāo)yp服從均值為ae^bt，方差為σ²(t-tk)的正態(tài)分布，即

yp～n(ae^bt,σ²(t-tk))(12)

其中t＞tk；

第6步，令t＝tk+δt,tk+2δt,…，根據(jù)yp的概率密度函數(shù)進(jìn)行ns次隨機(jī)采樣，得到ns個(gè)指數(shù)模型健康狀態(tài)指標(biāo)預(yù)測(cè)序列yp＝[yp1,yp2,yp3,…]，即ns條退化軌跡，每條退化軌跡的失效時(shí)刻eol為指數(shù)模型健康狀態(tài)指標(biāo)預(yù)測(cè)序列內(nèi)第一次預(yù)測(cè)指標(biāo)yp超過(guò)失效閾值yh的時(shí)刻，即

eol＝inf(yp＞yh)(13)

剩余壽命rul為失效時(shí)刻eol與從擬合起始時(shí)刻t0開(kāi)始的最后工作時(shí)刻tk的差，即

rul＝eol-tk(14)

對(duì)ns條退化軌跡的剩余壽命rul進(jìn)行統(tǒng)計(jì)，繪制剩余壽命rul的頻率分布直方圖并且擬合出對(duì)應(yīng)的概率密度函數(shù)，確定95％剩余壽命置信區(qū)間，給出剩余壽命估計(jì)值。

下面結(jié)合實(shí)施例對(duì)本發(fā)明做詳細(xì)描述。

實(shí)施例：以滾動(dòng)軸承加速壽命實(shí)驗(yàn)數(shù)據(jù)為例對(duì)本發(fā)明進(jìn)行驗(yàn)證。

滾動(dòng)軸承加速壽命實(shí)驗(yàn)在pronostia實(shí)驗(yàn)臺(tái)上操作完成，通過(guò)對(duì)軸承氣壓加載，使軸承在高負(fù)荷條件下工作，可以在數(shù)小時(shí)內(nèi)實(shí)現(xiàn)軸承從正常狀態(tài)退化到完全失效。實(shí)驗(yàn)過(guò)程中，軸承轉(zhuǎn)速為1800rpm，負(fù)載為4kn。采用加速度傳感器對(duì)軸承振動(dòng)信號(hào)進(jìn)行采樣，采樣頻率為25.6khz，數(shù)據(jù)長(zhǎng)度為2560，每次采樣持續(xù)時(shí)間為0.1s，采樣間隔為10s。當(dāng)振動(dòng)幅值超過(guò)20g時(shí)，軸承完全失效。實(shí)驗(yàn)軸承的全壽命振動(dòng)信號(hào)如圖2所示。

從振動(dòng)信號(hào)中提取有效值作為健康狀態(tài)指標(biāo)，并確定擬合起始時(shí)刻，如圖3所示，

選擇時(shí)刻tk＝5000s時(shí)對(duì)模型進(jìn)行評(píng)估，噪聲分布如圖4中實(shí)線所示，分別用梯形噪聲分布指數(shù)模型和傳統(tǒng)指數(shù)模型對(duì)參數(shù)進(jìn)行估計(jì)，得到的擬合噪聲項(xiàng)95％置信區(qū)間分別如圖4(a)和圖4(b)所示，從圖中可以直觀的看出，梯形噪聲分布指數(shù)模型相比于傳統(tǒng)指數(shù)模型，更符合數(shù)據(jù)分布的規(guī)律。

梯形噪聲分布指數(shù)模型的預(yù)測(cè)結(jié)果如圖5所示。圖5(a)為有效值預(yù)測(cè)期望及區(qū)間，圖5(b)為預(yù)測(cè)剩余壽命分布圖，剩余壽命估計(jì)為7840s，95％置信區(qū)間為[7160,8310]s，而真實(shí)剩余壽命為8000s。圖6為傳統(tǒng)的指數(shù)模型的預(yù)測(cè)結(jié)果，剩余壽命估計(jì)為6570s，95％置信區(qū)間為[5530,7320]s，而真實(shí)剩余壽命為8000s。對(duì)比結(jié)果可以看出，梯形噪聲分布指數(shù)模型預(yù)測(cè)更準(zhǔn)確，精度更高。

圖7顯示了梯形噪聲分布指數(shù)模型從tk＝1000s時(shí)開(kāi)始的剩余壽命估計(jì)和95％置信區(qū)間的預(yù)測(cè)結(jié)果，并且與傳統(tǒng)的指數(shù)模型預(yù)測(cè)結(jié)果(圖8)相對(duì)比。從圖中可以看出，梯形噪聲分布指數(shù)模型與傳統(tǒng)指數(shù)模型在預(yù)測(cè)前期都預(yù)測(cè)的波動(dòng)都不理想，隨著時(shí)間的推移，兩種方法都逐漸收斂到真實(shí)值。但是梯形噪聲分布指數(shù)模型的收斂速度更快，預(yù)測(cè)精度更高；且梯形噪聲分布指數(shù)模型得到的95％置信區(qū)間更小，剩余壽命估計(jì)更可靠

通過(guò)采用軸承加速壽命實(shí)驗(yàn)數(shù)據(jù)驗(yàn)證了本發(fā)明在機(jī)械設(shè)備壽命預(yù)測(cè)中的優(yōu)勢(shì)。本發(fā)明對(duì)傳統(tǒng)指數(shù)模型進(jìn)行改進(jìn)，對(duì)指數(shù)模型噪聲項(xiàng)的分布方差增加了初始值，使其更加符合機(jī)械設(shè)備退化數(shù)據(jù)的規(guī)律，從而提高了指數(shù)模型對(duì)機(jī)械設(shè)備壽命預(yù)測(cè)的精度和可靠性。

本發(fā)明所提出的一種梯形噪聲分布的指數(shù)模型機(jī)械設(shè)備剩余壽命預(yù)測(cè)方法，并不只局限于機(jī)械設(shè)備的剩余壽命預(yù)測(cè)，還可以應(yīng)用于其他電子元件的剩余壽命預(yù)測(cè)問(wèn)題。大量研究工作證明，本方法適用于各類具有指數(shù)衰退形式的機(jī)電產(chǎn)品的剩余壽命預(yù)測(cè)。實(shí)施者只需對(duì)本方法相應(yīng)步驟進(jìn)行適當(dāng)調(diào)整，以適應(yīng)不同產(chǎn)品的應(yīng)用需求，應(yīng)當(dāng)指出，在不脫離本發(fā)明構(gòu)思的前提下，所做的調(diào)整和變形，也應(yīng)視為本發(fā)明的保護(hù)范圍。