張 龍, 劉楊遠(yuǎn), 吳榮真, 王 良, 承志恒, 顏秋宏
(華東交通大學(xué) 載運(yùn)工具與裝備教育部重點(diǎn)實(shí)驗(yàn)室,南昌 330013)
軸承的性能從健康到失效是一個(gè)漸變的過程,性能退化評(píng)估不僅是對(duì)該過程實(shí)現(xiàn)定量評(píng)估的有效手段,而且是實(shí)現(xiàn)軸承故障預(yù)診斷的前提與基礎(chǔ),更是發(fā)揮狀態(tài)預(yù)測(cè)維修潛力的重要一環(huán)[1]。
軸承性能退化評(píng)估的一般流程為特征提取、特征選擇、建立劣化指標(biāo)、比較實(shí)測(cè)樣本與無故障基準(zhǔn)之間相似性[2]等幾個(gè)方面。其中比較度量相似性方法又可分為基于邊界的檢測(cè)方法、基于概率密度估計(jì)和基于模型重構(gòu)的檢測(cè)方法3個(gè)方面。典型的邊界檢測(cè)方法有支持向量機(jī)[3]和支持向量數(shù)據(jù)描述[4],該類方法通過無故障訓(xùn)練數(shù)據(jù)集來創(chuàng)建邊界,使用待評(píng)估數(shù)據(jù)與邊界的距離確定其相對(duì)于無故障狀態(tài)的類隸屬。但是基于邊界的異常檢測(cè)方法問題在于距離指標(biāo)的上限難以界定,從而難于從退化指標(biāo)中定量推斷性能退化程度?;诟怕拭芏裙烙?jì)的方法一般使用高斯混合模型[5]和隱馬爾科夫模型[6],該類方法首先對(duì)訓(xùn)練期間的無故障數(shù)據(jù)進(jìn)行概率密度建模,之后將待測(cè)數(shù)據(jù)輸入到訓(xùn)練好的模型中,樣本若落入低密度區(qū)域則不屬于無故障樣本。但是,高斯混合模型假設(shè)無故障數(shù)據(jù)的概率密度分布可由多個(gè)高斯分布組合得到,其獨(dú)立高斯函數(shù)數(shù)量難于精確估計(jì),不能保證與實(shí)際情況完全相符;而對(duì)于復(fù)雜的系統(tǒng)而言,隱馬爾科夫模型利用指數(shù)表達(dá)形式模擬健康狀態(tài)區(qū)間,該策略與現(xiàn)實(shí)情況相違背,往往導(dǎo)致精度較低。
基于模型重構(gòu)的方法通過映射目標(biāo)數(shù)據(jù)集,將待測(cè)數(shù)據(jù)輸入到經(jīng)過無故障數(shù)據(jù)訓(xùn)練的模型進(jìn)行映射,其重建誤差大小(例如歐幾里德距離)與退化性能程度密切相關(guān)。多元狀態(tài)估計(jì)(multivariate state estimation technique, MSET)就是一個(gè)典型的基于模型重構(gòu)的性能退化評(píng)估方法,本質(zhì)上是一種非線性的多元預(yù)測(cè)診斷技術(shù),也可稱為一種非參數(shù)建模方法。MSET通過對(duì)健康狀態(tài)基準(zhǔn)進(jìn)行相似性建模,然后輸入待測(cè)的觀測(cè)向量進(jìn)行模型重構(gòu),比較度量前后差異性從而判斷目標(biāo)的性能退化程度。例如Wang等[7]將最大條件互信息(conditional mutual information, CMI)特征選擇算法應(yīng)用于MSET的訓(xùn)練數(shù)據(jù),用K近鄰算法(K-nearest neighbor,KNN)動(dòng)態(tài)存儲(chǔ)歷史記憶矩陣,成功實(shí)現(xiàn)了風(fēng)電機(jī)組性能狀態(tài)檢測(cè)。MSET重構(gòu)模型相比于基于邊界和概率密度估計(jì)的性能退化方法優(yōu)勢(shì)明顯,其模型簡(jiǎn)單且物理意義明確,能及時(shí)準(zhǔn)確地反映設(shè)備性能退化趨勢(shì)。
在性能退化特征提取和選擇時(shí),常用的方法有時(shí)域、頻域特征提取、小波特征提取和時(shí)序模型分析等等。崔凱等[8]將加速度有效值、峭度值兩個(gè)時(shí)域指標(biāo)結(jié)合在一起聯(lián)合判斷風(fēng)機(jī)軸承運(yùn)行狀態(tài)。Zhou等[9]采用軸承信號(hào)小波包能量熵作為原始特征,然后將全壽命數(shù)據(jù)輸入到利用健康狀態(tài)特征建立好的RBF網(wǎng)絡(luò)模型中,從而得到性能退化指標(biāo)。Cong等[10]利用基于時(shí)序模型中的自回歸模型(AR模型)和Kolmogorov-Smirnov檢驗(yàn)相似概率值來反映滾動(dòng)軸承性能劣化的趨勢(shì)。但上述傳統(tǒng)的單個(gè)域指標(biāo)特征提取效果單一,對(duì)狀態(tài)變化感知不明顯,無法有效囊括軸承性能狀態(tài)且不利于早期故障的發(fā)現(xiàn)。
許多研究將多域特征指標(biāo)融合進(jìn)行性能退化評(píng)估,例如Liu等[11]從軸承振動(dòng)信號(hào)中提取多域統(tǒng)計(jì)特征,然后利用核聯(lián)合特征矩陣近似對(duì)角化算法將其融合為敏感特征,最后計(jì)算健康數(shù)據(jù)和當(dāng)前數(shù)據(jù)的類間和類內(nèi)離散度作為性能退化指標(biāo)。對(duì)于將多域特征組合的方法,雖然相較于單域指標(biāo)能更細(xì)致地描繪軸承退化過程,但是高維特征存在信息冗余且計(jì)算量大的問題,從而降低性能退化評(píng)估模型的性能。Dong等[12]將時(shí)域、頻域、時(shí)頻域特征利用主成分分析(principal component analysis,PCA)降維后作為性能退化指標(biāo)。叢華等[13]采用遺傳算法(genetic algorithm,GA)優(yōu)化提取特征,然后計(jì)算特征與支持向量數(shù)據(jù)描述超球體模型間的距離從而評(píng)估軸承的性能退化程度。然而,①PCA降維與其他優(yōu)化算法不同,并非直接從高維特征中挑選出具有代表性的低維特征,而是利用空間轉(zhuǎn)換將高維空間映射成低維空間,導(dǎo)致破壞了原有數(shù)據(jù)的結(jié)構(gòu),改變了原有數(shù)據(jù)的物理意義;②利用PCA或者GA單純對(duì)輸入特征進(jìn)行降維或優(yōu)化雖然一定程度上減少了冗余信息,但傳統(tǒng)的時(shí)頻域特征降維后的指標(biāo)不一定適配所需模型,匹配度有待進(jìn)一步考量。而目前做聯(lián)合降維性能退化評(píng)估的長(zhǎng)短期記憶網(wǎng)絡(luò)[14-15](long short-time memory, LSTM)等參數(shù)過于復(fù)雜,時(shí)間跨度大且網(wǎng)絡(luò)較深,導(dǎo)致計(jì)算量大且耗時(shí)費(fèi)力。
綜上,本文提出一種基于MSET重構(gòu)模型整體優(yōu)化的軸承性能退化評(píng)估方法。首先,基于軸承振動(dòng)信號(hào)提取m組時(shí)域和頻域特征、n階AR模型系數(shù)和k維三層小波包Renyi 熵歸一化后組成m+n+k維多域特征向量指標(biāo),并將得到健康狀態(tài)的高維特征向量構(gòu)建MSET模型的歷史觀測(cè)矩陣;然后,利用遺傳算法對(duì)軸承高維特征向量和多元狀態(tài)估計(jì)中的歷史記憶矩陣的觀測(cè)向量個(gè)數(shù)和維度進(jìn)行同步聯(lián)合優(yōu)化,從而得到低維特征向量和歷史記憶矩陣觀測(cè)向量個(gè)數(shù)與該低維特征向量最匹配的MSET 模型,實(shí)現(xiàn)了特征優(yōu)選和滾動(dòng)軸承性能評(píng)估模型的整體自適應(yīng)優(yōu)化;最后,構(gòu)建余弦相似度劣化指標(biāo)實(shí)現(xiàn)軸承性能退化評(píng)估。
小波包相較于小波而言,在信號(hào)高頻部分分解得更加精細(xì),能更好地刻畫信號(hào)的突變性和非平穩(wěn)性。其本質(zhì)就是集合數(shù)個(gè)高低通濾波器從而將頻譜劃分多個(gè)層次,得到相應(yīng)的高低頻成分,更好地表達(dá)信號(hào)。其分解所得的頻帶之間沒有疏漏和冗余,相互獨(dú)立。
小波包分解表達(dá)式[16]如下
(1)
小波包重構(gòu)表達(dá)式如下
(2)
小波包變換實(shí)際上就是小波基函數(shù)與信號(hào)進(jìn)行卷積的過程。卷積在時(shí)域中的表現(xiàn)形式為移動(dòng)的內(nèi)積,從而使基函數(shù)和信號(hào)之間相似部分放大,反之則被抑制,因此小波基函數(shù)的選擇將直接影響最后結(jié)果。信號(hào)進(jìn)行j層小波包分解后將分成2j個(gè)子帶信號(hào)。Nikolaou等[17]提出選擇3層或4層小波包分解時(shí),對(duì)所得到的信號(hào)分析結(jié)果影響不大。因此,為提高計(jì)算效率,本文采用db5小波基函數(shù),對(duì)信號(hào)進(jìn)行3層小波包分解。
當(dāng)設(shè)備發(fā)生故障后,信號(hào)的不確定信息及能量分布將較正常運(yùn)轉(zhuǎn)時(shí)發(fā)生改變,不同的故障類型以及不同程度的故障其能量分布更是不同。信息熵具有表征信號(hào)復(fù)雜程度的能力,信號(hào)隨機(jī)性越高、越復(fù)雜時(shí)信息熵越大,反之信號(hào)越規(guī)則信息熵越小。
信息熵計(jì)算公式如下
(3)
Renyi熵作為信息熵概念的擴(kuò)展,具有表征信號(hào)在時(shí)頻分布特征的能力[18]。Renyi 熵較信息熵而言對(duì)信號(hào)細(xì)微變化更加敏感,尤其是在端點(diǎn)處,能夠捕捉到隨機(jī)變量在概率分布中的差異并進(jìn)行放大。Renyi熵計(jì)算公式如下
(4)
式中,α≠1且α≥0。當(dāng)α=0時(shí)Hα(x)取最大值,Markel等[19]對(duì)Renyi熵的α進(jìn)行了研究,表明α=0.5為結(jié)果最優(yōu)值,故本文α選為0.5。
自回歸模型(autoregressive model, AR)是一種隨機(jī)信號(hào)參數(shù)化建模方法,常用于處理時(shí)間序列。AR模型利用有限的知識(shí)和自身的參數(shù)對(duì)前面多個(gè)觀測(cè)值進(jìn)行學(xué)習(xí),從而對(duì)下一時(shí)刻的輸出進(jìn)行預(yù)測(cè),其對(duì)系統(tǒng)狀態(tài)變化極其敏感,且能表征系統(tǒng)狀態(tài)特征。
取時(shí)間序列y(t),AR模型的分析階數(shù)為p,則關(guān)于時(shí)間序列的p階AR模型可以表示為
(5)
式中:L(t)為AR模型的殘差;αj為第j項(xiàng)的系數(shù)。
殘差L(t)為真實(shí)值和預(yù)測(cè)值兩者間的誤差,對(duì)于平穩(wěn)信號(hào)L(t)主要表征白噪聲。所以AR模型實(shí)際上就是調(diào)整模型參數(shù)使殘差更接近高斯白噪聲。本文選用最常用的最小二乘法來進(jìn)行AR譜估計(jì)以及模型參數(shù)的選擇,根據(jù)貝葉斯信息準(zhǔn)則(Bayesian information criterion,BIC)確定模型階數(shù)??纱笾路譃槿?
步驟1確定分析階數(shù)p的范圍繼而調(diào)整AR模型,本文p分別取1,2,…,300;
步驟2通過最小二乘法分別求得各階次下的自回歸參數(shù)αj(j=1,2,…,300),構(gòu)造式(5)所示的AR模型,進(jìn)而得到殘差L(t);
步驟3根據(jù)各階次殘差L(t)的BIC值最小原則,確定最優(yōu)階數(shù)。
多元狀態(tài)估計(jì)是一種非線性的非參數(shù)建模方法。MSET通過比較新輸入數(shù)據(jù)與無故障時(shí)的健康數(shù)據(jù),對(duì)各變量數(shù)據(jù)之間的內(nèi)在關(guān)聯(lián)做出判斷,根據(jù)歷史數(shù)據(jù)所獲得的權(quán)值向量等知識(shí)對(duì)實(shí)際運(yùn)行狀態(tài)做出估計(jì),通過殘差比較估計(jì)向量與新輸入的觀測(cè)向量的相似性,從而對(duì)設(shè)備進(jìn)行診斷。
假設(shè)在某一過程或系統(tǒng)中某個(gè)時(shí)刻tj有n個(gè)互相關(guān)聯(lián)的變量,將其記為觀測(cè)向量X(tj),如式(6)所示
X(tj)=[x1(tj)x2(tj) …xn(tj)]T
(6)
式中,X(tj)為tj時(shí)刻狀態(tài)變量xi的觀測(cè)值。
首先,通過一定的方法合理地挑選m列正常運(yùn)行狀態(tài)下的觀測(cè)變量即可組成歷史記憶矩陣Dn×m,如式(7)所示,從而完整地表示某一過程或系統(tǒng)正常運(yùn)行狀態(tài)下的整個(gè)動(dòng)態(tài)過程,換言之D即是對(duì)設(shè)備正常運(yùn)行數(shù)據(jù)的特征提取、選擇、學(xué)習(xí)和記憶的過程。
(7)
式中:n為在某一時(shí)刻互相關(guān)聯(lián)的狀態(tài)變量有n個(gè);m為取m個(gè)觀測(cè)向量X(tj),即取m個(gè)時(shí)刻。
然后,假設(shè)在某一過程或系統(tǒng)中某個(gè)時(shí)刻模型的輸入為觀測(cè)向量Xobs,將其輸入歷史記憶矩陣D之后的輸出為估計(jì)向量Xest。Xest是對(duì)輸入Xobs的預(yù)測(cè)值,即對(duì)當(dāng)前狀態(tài)的估計(jì)。對(duì)于任一輸入觀測(cè)向量Xobs,MSET均對(duì)應(yīng)生成一個(gè)m維的權(quán)值向量W
(8)
根據(jù)權(quán)值向量W對(duì)D中m個(gè)歷史觀測(cè)向量進(jìn)行線性組合即可得到估計(jì)向量Xest,如下所示
w1X(t1)+w2X(t2)+…+wmX(tm)
(9)
ε=Xest-Xobs
(10)
W=(DTD)-1DTXobs
(11)
在式(11)中由于Xobs的個(gè)數(shù)遠(yuǎn)大于監(jiān)測(cè)系統(tǒng)參數(shù)個(gè)數(shù),導(dǎo)致DTD不可逆,因此可采用非線性運(yùn)算符代替普通矩陣乘積運(yùn)算,即可得權(quán)重向量式(12)
W=(DT?D)-1·(DT?Xobs)
(12)
式中,?為非線性運(yùn)算符,避免因變量之間的相關(guān)性造成矩陣點(diǎn)乘之后的不可逆現(xiàn)象,從而擴(kuò)大關(guān)系式(12)的通用性。本文MSET重構(gòu)模型采用兩向量之間的歐式距離(Euclidean distance)作為非線性運(yùn)算符
(13)
式中,xi和yi均為變量。
將式(12)代入式(9)可得估計(jì)向量的表達(dá)式
Xest=D·(DT?D)-1·(DT?Xobs)
(14)
式(10)、式(13)和式(14)的物理意義為:當(dāng)新輸入的觀測(cè)向量Xobs仍然處于設(shè)備正常工作狀態(tài)區(qū)間范圍時(shí),其與歷史記憶矩陣D中某些觀測(cè)向量的組合會(huì)很相似,歐氏距離會(huì)很小,導(dǎo)致輸出的估計(jì)向量Xest預(yù)測(cè)值精度高,重構(gòu)精度高,殘差ε小;當(dāng)新輸入的觀測(cè)向量Xobs為設(shè)備出現(xiàn)故障時(shí),勢(shì)必導(dǎo)致其動(dòng)態(tài)特性發(fā)生變化,使Xobs偏離原始正常工作范圍,從而無法利用D中某些觀測(cè)向量的組合進(jìn)行重構(gòu),兩者相似度低,歐氏距離大,導(dǎo)致輸出的估計(jì)向量Xest預(yù)測(cè)值精度低,殘差ε大。由于Xobs和Xest兩向量間的異同性包含豐富的運(yùn)行狀態(tài)信息,因此常利用二者間的殘差進(jìn)行故障模式識(shí)別或性能退化評(píng)估。
由上可知MSET重構(gòu)模型能有效地獲取故障信息并對(duì)性能退化程度進(jìn)行評(píng)估。但應(yīng)進(jìn)一步考慮以下情況:①定量分析模型重構(gòu)前后的殘差從而對(duì)退化趨勢(shì)進(jìn)行量化描述;②確保高維多域特征經(jīng)降維后得到的低維特征與MSET模型歷史記憶矩陣中觀測(cè)向量個(gè)數(shù)是匹配的;③設(shè)置合理的預(yù)警閾值從而對(duì)設(shè)備故障進(jìn)行及時(shí)預(yù)警。這三點(diǎn)都是本文方法的關(guān)鍵所在,后續(xù)將進(jìn)一步進(jìn)行研究和分析。
軸承在實(shí)際運(yùn)行中,其故障的產(chǎn)生和發(fā)展是一個(gè)突發(fā)的、連續(xù)變化的過程,為了真實(shí)地反映故障連續(xù)變化的過程和驗(yàn)證本文所提方法有效性和可靠性,利用西安交大-昇陽科技聯(lián)合實(shí)驗(yàn)室(XJTU-SY)所提供的全壽命加速疲勞軸承試驗(yàn)數(shù)據(jù)[20]進(jìn)行分析。試驗(yàn)臺(tái)如圖1所示,包括驅(qū)動(dòng)部分、測(cè)試部分和液壓加載部分,具體器件如圖中所示。DT9837采集器通過兩個(gè)PCB 352C33單向加速度傳感器對(duì)型號(hào)為L(zhǎng)DK UER204的試驗(yàn)軸承進(jìn)行信號(hào)采集,轉(zhuǎn)軸轉(zhuǎn)頻為40 Hz,采樣頻率為25.6 kHz。試驗(yàn)設(shè)計(jì)了三種工況,每種工況有 5組測(cè)試軸承,工況類別見表1。
表1 加速疲勞試驗(yàn)工況Tab.1 Accelerated fatigue test condition
圖1 XJTU-SY試驗(yàn)臺(tái)Fig.1 The XJTU-SY test rig
本文以工況3第一個(gè)軸承Bearing 3_1為例進(jìn)行分析,該軸承從開始運(yùn)行到出現(xiàn)外圈故障直至失效累計(jì)時(shí)長(zhǎng)42 h 18 min,試驗(yàn)中每隔1 min采樣1.28 s數(shù)據(jù),共計(jì)2 538組數(shù)據(jù),每組數(shù)據(jù)32 768個(gè)點(diǎn),圖2為該軸承全壽命時(shí)域波形圖。利用軸承故障特征頻率計(jì)算公式所得的軸承各元件故障頻率如表2所示。
表2 故障特征頻率Tab.2 The fault characteristic frequency 單位:Hz
圖2 軸承3_1的時(shí)域振動(dòng)信號(hào)Fig.2 The temporal vibration signal of bearing 3_1
由前言可知,不同域的特征包含滾動(dòng)軸承振動(dòng)信號(hào)在不同狀態(tài)空間的信息,例如表3所示的11個(gè)時(shí)域特征中絕對(duì)平均值、標(biāo)準(zhǔn)差、均方根值和峰值4個(gè)指標(biāo)反映時(shí)域信號(hào)的能量大小,其余7個(gè)指標(biāo)反映信號(hào)在時(shí)域的分布特征。而為盡可能地反映軸承信號(hào)在頻域中所包涵的信息,取表4中13個(gè)常見的頻域統(tǒng)計(jì)特征,其中序號(hào)1的特征反映信號(hào)在頻域中振動(dòng)能量大小,序號(hào)2~9的特征反映頻譜的集中程度,序號(hào)10~13的特征反映主頻帶在頻譜中的位置變化。
表3 11個(gè)時(shí)域特征指標(biāo)Tab.3 11 time domain characteristic indexes
表4 13個(gè)頻域特征指標(biāo)Tab.4 13 frequency domain characteristic indexes
對(duì)各樣本提取AR 系數(shù)特征時(shí),AR 模型階數(shù)的確定至關(guān)重要,階數(shù)過低則無法準(zhǔn)確對(duì)模型進(jìn)行描述,反之則會(huì)導(dǎo)致過擬合,本文利用常見的BIC 準(zhǔn)則對(duì)階數(shù)進(jìn)行確定。Bearing 3_1軸承全壽命周期2 538 min,本文認(rèn)為試驗(yàn)剛開始時(shí)的前60 min為設(shè)備開機(jī)后的磨合階段,試驗(yàn)數(shù)據(jù)不準(zhǔn)確,于是取第61~第80共 20組測(cè)試軸承正常運(yùn)行下的壽命信號(hào)分別計(jì)算其AR 最優(yōu)階數(shù),由圖3可知,最優(yōu)AR 模型階數(shù)為26階,所以后續(xù)使用該軸承數(shù)據(jù)做故障診斷提取AR 系數(shù)特征時(shí),一律取26階AR 模型系數(shù)。
圖3 20個(gè)無故障信號(hào)的最優(yōu)階數(shù)Fig.3 Optimal orders of 20 bearing signals without faults
由1.1節(jié)可知,信號(hào)數(shù)據(jù)進(jìn)行j層小波包分解后將分成2j個(gè)子帶信號(hào),所以三層小波包分解將產(chǎn)生8個(gè)子帶信號(hào),依次對(duì)每個(gè)子帶信號(hào)計(jì)算其Renyi 熵,即可組成一個(gè)8維三層小波包Renyi 熵特征向量。
綜上,本文采用上述24組時(shí)域和頻域特征、26階AR模型系數(shù)和8維三層小波包Renyi 熵歸一化處理后共同組成58維多域特征向量指標(biāo)。
本文基于MSET重構(gòu)模型整體優(yōu)化的軸承性能退化評(píng)估模型如圖4所示,具體步驟如下:
圖4 性能退化評(píng)估模型Fig.4 Evaluation model of performance degradation
步驟1基于滾動(dòng)軸承振動(dòng)信號(hào)對(duì)原始時(shí)域數(shù)據(jù)提取24組時(shí)域和頻域特征、26階AR模型系數(shù)和8維三層小波包Renyi 熵歸一化后組成58維多域特征向量指標(biāo);
步驟2利用滾動(dòng)軸承正常運(yùn)行狀態(tài)下健康數(shù)據(jù)特征構(gòu)建多元狀態(tài)估計(jì)模型的歷史記憶矩陣;
步驟3利用遺傳算法對(duì)訓(xùn)練集高維特征向量和MSET中歷史記憶矩陣的觀測(cè)向量個(gè)數(shù)和維度進(jìn)行同步聯(lián)合優(yōu)化,實(shí)現(xiàn)特征降維和個(gè)數(shù)優(yōu)選;
步驟4對(duì)測(cè)試集數(shù)據(jù)提取低維特征,維度與GA優(yōu)化后一致,將得到的測(cè)試樣本低維特征向量作為觀測(cè)向量輸入MSET模型中;
步驟5輸出估計(jì)向量,計(jì)算觀測(cè)向量與估計(jì)向量之間相似性。
經(jīng)MSET模型重構(gòu)后的Xest為單重構(gòu)的低維向量,并不能直接反映劣化程度,本文選用重構(gòu)前后兩向量之間的余弦相似度(cosine similarity, CS)作為性能退化指標(biāo),從兩向量空間距離的角度評(píng)估滾動(dòng)軸承性能退化程度,CS表達(dá)式如下
(15)
GA相較于粒子群優(yōu)化[21]和模擬退火[22]等優(yōu)化算法而言,對(duì)初始條件的設(shè)置要求較低,只需設(shè)定相應(yīng)的適應(yīng)度函數(shù)即可,具有較強(qiáng)的魯棒性,能更快地搜索到全局最優(yōu)解。上述步驟3中利用GA對(duì)模型進(jìn)行同步優(yōu)化的具體操作如下:
步驟1采用二進(jìn)制編碼方式進(jìn)行編碼,對(duì)測(cè)試軸承前n個(gè)樣本和58維特征向量組成染色體長(zhǎng)度為n+58的基因組,其排列方式如圖5所示;
圖5 編碼排列方式Fig.5 Coding arrangement mode
步驟2初始化GA參數(shù),設(shè)置初始種群大小為50,遺傳代數(shù)為100,采用多點(diǎn)交叉方式,交叉率為0.7,變異率為0.05,代溝值為0.9;
步驟3設(shè)置訓(xùn)練樣本集合,主要由無故障階段和失效階段兩部分信號(hào)構(gòu)成;并采用所有正常樣本CS均值與所有故障樣本CS均值的比值設(shè)計(jì)適應(yīng)度函數(shù),從而計(jì)算個(gè)體適應(yīng)度,適應(yīng)度越高,個(gè)體越優(yōu);
(16)
式中:m為正常樣本數(shù)量;n為故障樣本數(shù)量;CSN為正常樣本的余弦相似度;CSF為故障樣本的余弦相似度。
步驟4根據(jù)步驟2中所確定的交叉和變異遺傳算子對(duì)當(dāng)前種群進(jìn)行處理,并判別其種群個(gè)體適應(yīng)度是否達(dá)到要求。若達(dá)到要求則進(jìn)行解碼,得到最優(yōu)低維向量的維數(shù)及其索引以及最優(yōu)歷史記憶矩陣觀測(cè)向量個(gè)數(shù)和維度,執(zhí)行步驟6,否則繼續(xù)下一步;
步驟5進(jìn)行選擇、交叉、變異過程,返回步驟3;
步驟6獲得最優(yōu)參數(shù)建立MSET單重構(gòu)模型,利用模型進(jìn)行性能退化評(píng)估。
余弦相似度CS是一條連續(xù)變化的曲線,將其作為劣化指標(biāo)能準(zhǔn)確地量化表征軸承目前運(yùn)行狀態(tài)與無故障運(yùn)行時(shí)的偏離程度。進(jìn)一步還需要對(duì)該曲線標(biāo)定一個(gè)恰當(dāng)?shù)念A(yù)警值從而及時(shí)對(duì)故障進(jìn)行示警。
由統(tǒng)計(jì)學(xué)知識(shí)可知,數(shù)學(xué)期望為μ,方差為σ2的正態(tài)分布(μ,σ2)數(shù)據(jù),樣本數(shù)值落在(μ-3σ,μ+3σ)區(qū)間內(nèi)的概率約為0.997 3,因此判定在樣本故障程度相近時(shí),CS值也符合正態(tài)分布要求。進(jìn)一步當(dāng)連續(xù)多個(gè)樣本的CS值超過由歷史樣本CS值所確定的(μ-3σ,μ+3σ)區(qū)間范圍時(shí),判定該時(shí)刻故障程度于前一刻已經(jīng)發(fā)生較大變化。CS值是連續(xù)單調(diào)遞減的,故只計(jì)算下限μ-3σ即可,如式(17)所示。本文利用3σ準(zhǔn)則對(duì)CS值進(jìn)行預(yù)警值標(biāo)定,未超過該值的視為正常運(yùn)行狀態(tài),同時(shí)利用當(dāng)前所有無故障樣本更新CS值的3σ區(qū)間,達(dá)到自適應(yīng)示警的目的。
式中:{D(i),i=1,2,…,n}為t=[1,n]時(shí)刻的CS值;T(t)為t時(shí)刻預(yù)警值。當(dāng)首次連續(xù)多個(gè)CS值小于te時(shí)刻的預(yù)警值T(te)時(shí),定義該時(shí)刻為早期故障時(shí)刻或故障程度加深時(shí)刻。為防止外界偶然因素的影響,本文早期故障時(shí)刻要求連續(xù)5個(gè)測(cè)試數(shù)據(jù)的CS值小于或等于te時(shí)刻的預(yù)警閾值。
對(duì)試驗(yàn)數(shù)據(jù)按照步驟1~步驟5進(jìn)行整體優(yōu)化性能退化評(píng)估操作,具體如下:首先從健康運(yùn)行軸承的采集信號(hào)中選擇300個(gè)樣本提取其58維特征作為觀測(cè)向量用來構(gòu)造歷史記憶矩陣D。然后額外選取100個(gè)正常運(yùn)行狀態(tài)的樣本和100個(gè)故障樣本提取其高維特征作為訓(xùn)練數(shù)據(jù),利用遺傳算法對(duì)訓(xùn)練集高維特征向量和MSET 模型中歷史記憶矩陣的觀測(cè)向量個(gè)數(shù)和維度進(jìn)行同步聯(lián)合優(yōu)化。最后將Bearing 3_1 軸承全壽命周期2 538個(gè)樣本作為測(cè)試數(shù)據(jù),提取其低維特征,維度與GA優(yōu)化后一致,將得到的測(cè)試樣本低維特征向量作為新輸入觀測(cè)向量輸入MSET模型中進(jìn)行后續(xù)分析。
對(duì)試驗(yàn)數(shù)據(jù)進(jìn)行樣本劃分后按上述過程進(jìn)行整體優(yōu)化性能退化評(píng)估,優(yōu)化結(jié)果如圖6所示。從圖6(a)中可以看出,經(jīng)優(yōu)化后歷史記憶矩陣所包含的觀測(cè)向量個(gè)數(shù)由原始的300個(gè)變?yōu)?33個(gè),且分布較為均勻,有效地剔除了部分冗余觀測(cè)向量,減少了后續(xù)性能退化評(píng)估的計(jì)算量;從圖6(b)中可以看出,在步驟3同步聯(lián)合優(yōu)化中,58維特征優(yōu)選成了9維特征,由小波Renyi 熵后三項(xiàng)和AR模型中6個(gè)系數(shù)組成,其中24組時(shí)域和頻域特征降維成了0維,表明時(shí)、頻域指標(biāo)對(duì)西安交大軸承數(shù)據(jù)中的Bearing 3_1試驗(yàn)數(shù)據(jù)分析結(jié)果較差。
圖6 優(yōu)化結(jié)果Fig.6 Optimization results
為進(jìn)一步分析該試驗(yàn)數(shù)據(jù)特征優(yōu)選后時(shí)域、頻域特征降維后變成了0維的原因,對(duì)3個(gè)時(shí)域指標(biāo)和3個(gè)頻域指標(biāo)進(jìn)行作圖,如圖7所示。從圖7中可以看出,脈沖因子和峭度值與軸承故障發(fā)展趨勢(shì)基本一致,但后期均出現(xiàn)先減小后增大的情況,與故障實(shí)際發(fā)展情況不相符,并且脈沖因子指標(biāo)前期波動(dòng)過大,不利于早期故障的發(fā)現(xiàn);加速度有效值波形后期振蕩過于嚴(yán)重,與故障實(shí)際發(fā)展情況不相符;其余3個(gè)頻域特征均不能準(zhǔn)確表征故障程度和性能退化趨勢(shì)。
圖7 軸承時(shí)域指標(biāo)Fig.7 Time domain index of bearing
對(duì)Bearing 3_1軸承全壽命2 538個(gè)樣本進(jìn)行分析,性能退化評(píng)估CS 值曲線經(jīng)過5點(diǎn)平滑處理后如圖8所示。由于故障程度具有漸變性,平滑處理不僅不會(huì)影響故障預(yù)警,更能削弱噪聲干擾,更好地表達(dá)信號(hào)。圖中虛線為自適應(yīng)預(yù)警值,實(shí)線為故障指標(biāo)CS 值,曲線整體呈現(xiàn)下降趨勢(shì),表明CS 值能準(zhǔn)確、一致地表征軸承故障嚴(yán)重度趨勢(shì),且相較于圖7的時(shí)域指標(biāo)對(duì)軸承早期故障更加敏感。
圖8 軸承全壽命CS值Fig.8 CS value of bearing life
圖8中性能退化評(píng)估曲線大致可劃分成無故障時(shí)期、初始故障時(shí)期、中度故障時(shí)期和惡化失效時(shí)期4個(gè)階段。由于該軸承的CS 值在534時(shí)刻起連續(xù)5個(gè)值小于自適應(yīng)預(yù)警值,故0~533時(shí)刻為無故障時(shí)期,534時(shí)刻為初始故障發(fā)生時(shí)刻。534~1 119時(shí)刻內(nèi)的CS 值呈明顯下降趨勢(shì),但在1 119時(shí)刻曲線出現(xiàn)階躍現(xiàn)象,CS 值急劇下降,這是因?yàn)樵摃r(shí)刻軸承故障程度持續(xù)加重,導(dǎo)致新輸入觀測(cè)向量與歷史記憶矩陣的偏離程度持續(xù)加深,最終在1 119時(shí)刻達(dá)到中度故障時(shí)期。1 119~1 378時(shí)刻內(nèi)的CS 值振蕩降低后垂直上升,這是因?yàn)楣收系陌l(fā)展是一個(gè)漸變的過程,從微小故障逐漸擴(kuò)大,然后磨平,隨著時(shí)間增加故障再度擴(kuò)大,所以導(dǎo)致CS 值出現(xiàn)了大范圍地反復(fù)波動(dòng)。1 378時(shí)刻,原始故障點(diǎn)已被磨平,因故障產(chǎn)生的突變沖擊減小,新輸入觀測(cè)向量與歷史記憶矩陣的偏離程度減小,CS 值突變?cè)龃蟆? 378~2 538時(shí)刻,軸承故障程度不斷加重,且CS 值波動(dòng)程度比之前各階段更大,波動(dòng)時(shí)間更短,表明故障點(diǎn)的“磨平- 惡化”情況更加明顯,且故障程度快速加深直至軸承失效。從圖中還可以看出,CS 值曲線整體毛刺較多,局部波動(dòng)性大,這是由于軸承在整體運(yùn)轉(zhuǎn)過程中微小的“故障磨平、再度劣化”現(xiàn)象一直存在。
為進(jìn)一步分析CS 值與故障程度之間的關(guān)聯(lián)性,選取無故障階段(0~533)、 輕度故障階段(534~1 118)、 中度故障階段(1 119~1 377)、重度故障階段(1 378~2 538) 四種狀態(tài)中各3個(gè)樣本,計(jì)算其低維特征后輸入聯(lián)合優(yōu)化模型,可得到輸出估計(jì)向量與觀測(cè)向量之間的余弦相似度差異,不同故障程度下CS 值的變化如圖9所示。從圖中可以看出,無故障樣本的CS 值最大,3個(gè)樣本的均值約為0.985,接近于1,這是因?yàn)闊o故障樣本提取的特征與歷史記憶矩陣中各向量相似程度高,輸出估計(jì)向量與觀測(cè)向量的偏斜度小,兩空間向量之間的夾角小,余弦相似度大。而隨著外圈故障尺寸的變大,估計(jì)向量與觀測(cè)向量之間的CS值越來越小,這是由于隨著故障嚴(yán)重度加深,軸承當(dāng)前狀態(tài)越來越偏離無故障運(yùn)行狀態(tài),輸入特征越來越偏離依靠健康樣本所建立的歷史記憶矩陣的觀測(cè)向量,兩者之間的偏斜度越來越大,在CS 值上的表現(xiàn)則呈現(xiàn)出遞減趨勢(shì),從而很好地區(qū)分軸承的故障程度。
圖9 不同程度外圈故障CS值Fig.9 CS value of outer ring faults in different degrees
為驗(yàn)證聯(lián)合優(yōu)化多域特征-MSET 重構(gòu)模型方法的有效性,對(duì)No.534和No.533樣本信號(hào)進(jìn)行以包絡(luò)譜譜峰因子為優(yōu)化指標(biāo)自適應(yīng)共振解調(diào)[23]處理,共振解調(diào)后對(duì)軸承信號(hào)進(jìn)行傅里葉變換,如圖10所示。
圖10 故障樣本包絡(luò)解調(diào)Fig.10 Envelope demodulation of fault samples
圖10(a) No.534樣本共振解調(diào)結(jié)果中出現(xiàn)明顯脈沖現(xiàn)象,且在頻譜的包絡(luò)譜中出現(xiàn)了與外圈故障頻率BPFO(123.3 Hz)相近的一倍頻124 Hz、二倍頻246 Hz和三倍頻370 Hz,此外轉(zhuǎn)軸轉(zhuǎn)頻40 Hz也清晰可見。但由于該時(shí)刻為早期故障點(diǎn),軸承外圈剛剛出現(xiàn)微小故障,導(dǎo)致二倍頻和三倍頻幾乎泯沒在噪聲當(dāng)中。而在圖10(b)No.533樣本包絡(luò)解調(diào)中沖激現(xiàn)象和123.3 Hz的外圈故障頻率成分并不明顯,且在534時(shí)刻的最優(yōu)濾波頻帶相較于533時(shí)刻的最優(yōu)濾波頻帶發(fā)生了明顯改變。因此,證明534時(shí)刻為軸承早期故障點(diǎn),驗(yàn)證了方法的有效性和可靠性。
為驗(yàn)證本文所提方法的優(yōu)越性,首先采用僅進(jìn)行特征維度優(yōu)化的MSET軸承性能退化模型進(jìn)行對(duì)比分析。利用GA對(duì)58維特征向量進(jìn)行降維,如圖11所示,從圖中可以看出僅對(duì)特征進(jìn)行優(yōu)化時(shí),時(shí)域、頻域的特征向量仍有所保留,與圖6的分析不相符,導(dǎo)致輸入將含有許多冗余信息。進(jìn)一步分析,將降維后的34維特征向量輸入未經(jīng)優(yōu)化的MSET模型,其性能退化曲線和自適應(yīng)預(yù)警值如圖12(a)所示,其早期故障點(diǎn)為No.306樣本,遠(yuǎn)遠(yuǎn)早于本文方法所判斷的早期故障點(diǎn),原因是因?yàn)槲唇?jīng)優(yōu)化的MSET模型歷史記憶中的觀測(cè)向量個(gè)數(shù)過多,而輸入的特征向量維度也過多,導(dǎo)致含有許多冗余信息從而提前誤判故障點(diǎn)。為驗(yàn)證該時(shí)刻軸承并未發(fā)生故障,對(duì)No.306樣本信號(hào)進(jìn)行包絡(luò)解調(diào),如圖12(b)所示,其頻譜中并無外圈故障頻率成分,證明此刻軸承仍處于正常運(yùn)行狀態(tài),并未發(fā)生故障,早期故障點(diǎn)判斷失誤。
圖11 多域特征向量?jī)?yōu)化結(jié)果Fig.11 Multi-domain eigenvector optimization result
圖12 對(duì)比試驗(yàn)一分析Fig.12 Analysis of the comparative experiment one
其次,采用僅優(yōu)化MSET中歷史記憶矩陣的軸承性能退化模型進(jìn)行對(duì)比分析。利用GA對(duì)觀測(cè)向量的個(gè)數(shù)進(jìn)行優(yōu)選后,如圖13所示。從圖13中可以看出,歷史記憶矩陣觀測(cè)向量個(gè)數(shù)由原來的300個(gè)變成了156個(gè),分布較為均勻,有效剔除了冗余觀測(cè)向量。進(jìn)一步分析,將提取的高維特征向量直接輸入優(yōu)化后的MSET模型,其性能退化曲線和自適應(yīng)預(yù)警值如圖14所示。從圖14中可以看出,在No.320樣本處CS值曲線與自適應(yīng)預(yù)警值曲線已有切點(diǎn),預(yù)示軸承已經(jīng)開始出現(xiàn)故障,但從圖14(b)的包絡(luò)解調(diào)結(jié)果可以看到該樣本僅包含轉(zhuǎn)頻40 Hz,因此該預(yù)測(cè)結(jié)果不準(zhǔn)確。而由于未對(duì)輸入高維特征進(jìn)行優(yōu)化,導(dǎo)致性能退化曲線并不能很好地區(qū)分軸承的故障程度,自適應(yīng)預(yù)警值也不能準(zhǔn)確有效地反映出早期故障發(fā)生時(shí)刻,致使長(zhǎng)時(shí)間后在No.2296樣本處兩曲線再次出現(xiàn)交點(diǎn),而該時(shí)刻雖已發(fā)生故障,但遠(yuǎn)遠(yuǎn)晚于本文方法所判斷的早期故障點(diǎn),時(shí)效性不強(qiáng),不具有使用意義。
圖13 觀測(cè)向量?jī)?yōu)化結(jié)果Fig.13 Optimization result of observation vector
圖14 對(duì)比試驗(yàn)二分析Fig.14 Analysis of the comparative experiment two
由以上兩個(gè)對(duì)比試驗(yàn)綜合分析可知,相比于本文所提出的聯(lián)合優(yōu)化多域特征-MSET模型,單優(yōu)化輸入特征維度和單優(yōu)化歷史記憶矩陣的模型效果均一般,與軸承的實(shí)際發(fā)展趨勢(shì)不一致,不能及時(shí)地、準(zhǔn)確地反映軸承早期故障,也無法很好地區(qū)分故障程度,由此說明本文所提的整體優(yōu)化MSET單重構(gòu)模型方法的準(zhǔn)確有效性和優(yōu)越性。
進(jìn)一步為驗(yàn)證所提模型的泛化性和優(yōu)越性,采用多個(gè)軸承全壽命數(shù)據(jù)進(jìn)行性能退化評(píng)估,并與其他文獻(xiàn)中結(jié)果進(jìn)行對(duì)比。從表5中可以看到,所用方法在不同軸承數(shù)據(jù)集下均能更準(zhǔn)確地判別軸承性能退化趨勢(shì),模型具有一定的泛化性。
表5 不同軸承數(shù)據(jù)間的試驗(yàn)結(jié)果Tab.5 The results between different bearing data
針對(duì)傳統(tǒng)單域特征指標(biāo)無法充分表征軸承狀態(tài)信息,而目前性能退化評(píng)估中基于多域高維特征向量的重構(gòu)評(píng)估模型存在信息冗余且易受到不一致優(yōu)化目標(biāo)的影響,不能在低維空間中保留必要的信息,而導(dǎo)致模型次優(yōu)性能。提出一種基于MSET重構(gòu)模型整體優(yōu)化的軸承性能退化評(píng)估方法,通過同步聯(lián)合優(yōu)化多域特征和MSET 模型歷史記憶矩陣從而實(shí)現(xiàn)特征優(yōu)選和性能評(píng)估模型的整體自適應(yīng)優(yōu)化。
(1)XJTU-SY 軸承全壽命數(shù)據(jù)分析表明經(jīng)聯(lián)合優(yōu)化多域特征-MSET 重構(gòu)模型后的所得CS 值能很好地表達(dá)軸承故障嚴(yán)重度變化,且隨著程度加深具有單調(diào)遞減趨勢(shì);并通過對(duì)傳統(tǒng)時(shí)域指標(biāo)進(jìn)行比較,表明本文所提方法對(duì)軸承早期故障更加敏感,且與故障發(fā)展趨勢(shì)的一致性更好。
(2)通過對(duì)兩組單優(yōu)化對(duì)比試驗(yàn)和不同方法進(jìn)行比較分析,并利用復(fù)小波共振解調(diào)法對(duì)早期故障點(diǎn)進(jìn)行驗(yàn)證,表明本文基于MSET重構(gòu)模型整體優(yōu)化的軸承性能退化評(píng)估方法更具有效性、可靠性和優(yōu)越性。