張銘光 駱學理 賈 登 張 易 馬 波*

(1.北京化工大學機電工程學院, 北京 100029;2.北京化工大學高端機械裝備健康監(jiān)控及自愈化北京市重點實驗室,北京 100029;3.中國石油集團工程技術研究院有限公司, 北京 102206)

引 言

往復式壓縮機、柴油機等復雜往復機械作為流程工業(yè)中的關鍵設備,其運行可靠性直接影響企業(yè)的安全生產,建立可靠的故障預警系統是保證設備安全運行的必要舉措。 由于往復機械結構復雜,其運行信號具有非平穩(wěn)性的特點,致使實現故障的早期預警較為困難[1]。

如何綜合利用多特征參數實現對往復機械運行狀態(tài)的有效監(jiān)測,已成為設備故障預警領域研究的熱點。 李一青等[2]提出基于自組織神經網絡的多特征融合方法,構建準確反映軋機自激振動趨勢的特征指標,通過判斷特征指標是否超過閾值來實現故障預警。 馬波等[3]基于主題模型技術建立相空間的主題分布,根據往復壓縮機不同工況相空間主題分布的差異性,實現設備的異常檢測預警及故障診斷。 吳定海等[4]提出基于雙樹復小波變換的振動特征提取方法,利用特征向量訓練得到的后向傳播神經網絡對柴油機進行故障預警。 馬波等[5]提出基于變分自編碼器的機械故障智能預警方法,自學習機械振動信號的高維特征統計分布模型,通過分析模型的變化實現智能預警;他們還提出基于降噪自編碼器的機械故障智能預警方法[6],將高維監(jiān)測數據編碼成低維特征,通過度量待測樣本編碼特征與基準的距離實現故障預警;之后他們又構建了高維特征相空間并采用無線t 混合模型對其進行擬合[7],依據模型的變化實現智能預警。 郭鵬飛等[8]將特征變量和標簽數據導入隨機森林算法,隨后建立了一種監(jiān)測變槳軸承磨損的預警模型。

基于神經網絡的故障預警方法需要大量的故障樣本用于訓練神經網絡,然而工業(yè)現場難以獲取故障數據,導致該類型方法較難實際應用。 基于統計學習的故障預警方法利用設備正常數據訓練模型,但目前該類型方法每次生成混合模型都需要重新估計先驗參數,未建立設備不同時刻運行狀態(tài)之間的有效聯系,同時增加了生成混合模型耗時,對故障預警的準確性和時效性均造成影響。

主題模型是文本主題歸納中常用的方法,動態(tài)主題模型(dynamic topic model, DTM)則是在傳統靜態(tài)主題模型基礎上引入時間維度,使得模型能夠顯示主題在特定時間間隔之間的變化。 趙美玲等[9]利用基于動態(tài)主題模型的輿情本體概念抽取方法,成功發(fā)掘出依時間變動的互聯網輿情主題。 蔣卓人等[10]采用一種結合有監(jiān)督學習的動態(tài)主題模型(該模型可準確反映文檔的主題結構),精確捕捉到了主題-詞匯概率分布的動態(tài)演化。 鑒于其優(yōu)異的動態(tài)建模能力,DTM 在文本主題挖掘和圖像分析等領域取得了非常優(yōu)異的應用效果,然而在機械故障預警領域尚未有所應用。 本文首次將DTM 引入機械故障預警領域,對基于統計學習的故障預警方法進行改進,提出一種DTM 結合學生t 分布的往復機械預警模型建模方法(DtMM)。 該方法固定混合模型內子成分的順序結構,依據DTM 建模原理建立混合模型與基準混合模型之間的參數演化關系,新的混合模型由基準混合模型結合輸入數據演化而來,通過計算混合模型與基準混合模型的差異來表征設備實時工況狀態(tài)與正常狀態(tài)的差異,以此實現對設備的狀態(tài)監(jiān)測。 最后分別采用往復壓縮機工程案例數據和故障試驗數據對提出方法的有效性進行驗證。

1 動態(tài)主題模型與變分推斷

1.1 動態(tài)主題模型

DTM 在傳統隱含狄利克雷分布(latent Dirichlet allocation,LDA) 方法的基礎上進行了延伸,在底層主題多項式的自然參數空間上使用狀態(tài)空間模型,使得主題可以在時間序列下演變。 對于樣本-主題概率分布,采用均值為α的邏輯正態(tài)分布來表示概率分布的不確定性,使用簡單的動態(tài)模型建立模型之間的順序演化關系[11]。 如式(1)所示,將主題和主題概率分布相關聯,得到一組主題模型。 一組有序樣本在時間切片t上的生成過程表述如下。

1)生成狄利克雷分布βk,t

式中,βk,t表示時間切片t所對應主題k的數據分布,N代表高斯分布,σ2I表示高斯噪聲過程。

2)描述主題隨時間的變化

式中,α表示每個樣本可能的主題分布,σα表示邏輯正態(tài)分布的標準差,I表示參數傳遞過程。

3)對每個輸入樣本

(a)生成樣本-主題分布參數η

(b)對于樣本中的每個數據點

i.生成樣本-主題分布Z

ii.生成各主題下數據Wt,d,n

式中,Mult表示映射過程,π(x)將多項自然參數映射為平均參數,可表示為

式中,w表示生成數據點的數量。

1.2 變分推斷

對于混合模型參數推斷,傳統的靜態(tài)主題模型采用Gibbs 采樣進行參數估計[12],但由于DTM 中高斯模型和多項式模型存在非共軛性,導致Gibbs 采樣收斂性差,因此采用收斂較快的變分推斷對混合模型參數進行估計。 通過不斷更新變分變量,來不斷減小變分分布與真實后驗分布之間的距離,當距離值小于預先設置的閾值后,即可將變分分布作為真實后驗分布的近似替代[13-14]。 為了求取統計變量Φ={Z,μ,Λ,u,V,α}的真實后驗分布p(Φ|X),引入一個變分分布q(Φ),其分解形式表示為

數據樣本X與全部統計變量Φ= {Z,μ,Λ,u,V,α}的聯合分布表示為

為了實現推斷目標,引入一個對數邊緣似然的等式,如式(9)所示。

式中,F表示lnp(X)的變分下界;KL(q‖p)為q(Φ)與p(Φ|X)之間的Kullback-Leibler 散度。 由于KL(q‖p)越小,q(Φ)與p(Φ|X)越相近,根據式(9),可通過最大化變分下界F來達到優(yōu)化q(Φ)的目的。

2 往復機械故障早期預警方法

2.1 方法流程

往復機械作為一種復雜機械,其振動響應信號往往由多個激勵源信號疊加而成,因此呈現出較強的非平穩(wěn)性。 若設備發(fā)生故障,往往會引起單個或多個激勵源信號分布發(fā)生變化,最終引起振動響應信號分布發(fā)生變化[15-16],其振動響應信號如式(10)所示。

式中,ri表示振源到傳感器的傳遞距離;Vk表示w(ri,t)間相互獨立的沖擊強度;w(ri,t)表示t時刻的振動沖擊波形;N(ri,t)即泊松過程,表示脈沖沖擊計數變量。 當考慮時變傳遞路徑h(ri,t)的影響時,非平穩(wěn)振動響應信號如式(11)所示。

式中,f(ri,τ)表示激勵源響應函數,i=1,2,…,m,j=1,2,…,n;N(r,t)表示高斯噪聲。

在設備故障早期,損傷部件與其他部件相互作用會使響應信號的分布產生變化。 DtMM 方法可自學習其統計分布,在基準模型基礎上不斷演化生成新的混合模型,通過計算設備實時運行工況下混合模型和基準混合模型的差異來判斷設備狀態(tài)是否出現異常,以此實現設備運行狀態(tài)的監(jiān)測。 響應信號的分布的變化可能由單個或多個部件故障引起[17],該變化可由某些特征反映,若設備運行狀態(tài)發(fā)生改變,其特征值必然發(fā)生變化。 因此,選取往復機械監(jiān)測參數的敏感特征構造特征矩陣,將特征矩陣作為DtMM 方法模型訓練和測試的輸入,可準確識別設備狀態(tài)的異常變化。

所提方法的預警流程如圖1 所示,主要分為特征矩陣構建、基準混合模型和實時混合模型訓練、報警閾值自學習與故障預警應用3 個部分。

圖1 DtMM 故障預警流程圖Fig.1 Fault early warning process of DtMM

2.2 DtMM 模型

DtMM 模型由多個t-分布線性疊加而成,t-分布具有一個控制拖尾長度的自由參數v,可以很好地擬合往復機械拖尾性振動信號u(ri,t),具備較強的魯棒性[7]。 混合模型概率密度函數表示為

式中,Θ= {πj,μj,Λj,vj}表示混合模型參數集,其中j為子模型序號,μ為均值矩陣,Λ為協方差矩陣,v為自由度;π=(π1,π2,…,πk)T為子模型權重向量;K為子模型的數量。 DtMM 概率圖模型如圖2 所示。

圖2 DtMM 模型生成示意圖Fig.2 DtMM model generation schematic diagram

對于數據集Di,t=(d1,n,d2,n,…,dm,n)T,di,n是由n個特征構成的特征向量的表示,DtMM 基于Di,t生成混合模型的步驟分為兩部分。

1)在初始時刻(t=0),初始化混合模型參數集Θ,依據振動數據自學習響應信號中激勵源的數量以及各激勵源信號的統計分布。 混合模型中子模型參數的先驗服從狄利克雷過程,在xn的生成過程中,引入一個隱統計變量zn,用于確定xn所屬的子模型, 生成過程如下。

式中,π(V) ={πj(V)},zn服從多項分布。2)設備故障的發(fā)生往往伴隨著緩慢的劣化過程,當前時間片的激勵源分布與前一時間片的激勵源分布存在著緊密的聯系。 為了準確描述這種關系,在新的混合模型生成過程中,將激勵源信號所服從的統計分布作為節(jié)點映射時序混合模型,通過邏輯正態(tài)分布鏈接前一時間片混合模型的樣本-激勵源分布,根據當前時間新增的數據樣本更新激勵源-數據分布,生成過程如下。

式中,N(Θt-1,σ2I)為正態(tài)分布。

2.3 構造特征矩陣

獲取往復壓縮機K組正常工況振動數據作為訓練集,M組不同類型故障的振動數據作為測試集,經特征提取后構建特征向量矩陣,記為Di,t。

式中,Di,t表示第i個樣本,t表示該樣本對應時刻,di,j表示第i組數據中第j個特征的特征值,i∈1,2,…,m,j∈1,2,…,n,m表示數據集的長度,n表示故障特征的數量。

2.4 訓練模型

訓練樣本Di,t0用于訓練基準混合模型,訓練樣本Di,ti,i∈[1,k]用于訓練正常狀態(tài)混合模型。 具體實施步驟如下。

1)輸入特征矩陣Di,t0,初始化基準混合模型,利用狄利克雷過程生成混合模型先驗參數,即基準混合模型中的子模型數量及權重參數。

2)基于步驟1)計算結果,結合變分推斷計算基準混合模型后驗參數,即子模型分布參數,得到基準混合模型mstandard并將其保存。

3) DtMM 基于mstandard分布參數演化生成正常工況混合模型先驗參數,結合訓練樣本D1,t1,采用變分推斷計算混合模型后驗參數,得到第一個正常狀態(tài)混合模型m1,并計算m1與mstandard的模型差異值。

4)依次輸入訓練樣本Di,t,i∈2,…,k,重復步驟3),得到正常狀態(tài)混合模型與基準混合模型差異集,用于自學習報警閾值。

2.5 自學習預警閾值

統計模型距離評價方法包括KL 散度、余弦距離、JS 距離和杰卡德相似系數等,其中KL 散度通過計算兩個模型的相對應成分間的差異來綜合表征兩個模型間的距離[18]。 DtMM 生成的模型由多個t-分布組成,KL 散度逼近方法通過計算t-分布的差異來表征不同模型之間的距離,故結果更加精確,因此本文采用KL 散度逼近方法計算混合模型mi與基準混合模型mstandard的距離,計算過程如式(23)所示。

其中,

式中,KL(mi(j)‖mstandard(j))表示混合模型中第j個子模型與基準混合模型中第j個子模型的KL 散度,即兩個t-分布的KL 散度。

預警閾值依據設備正常狀態(tài)數據自學習獲得,具有較強的自適應性和高的可靠性。 訓練k個正常狀態(tài)混合模型mi,統計mi與mstandard的距離,并計算其均值μkd和方差σkl。 設備處于正常運行狀態(tài)時,正常狀態(tài)混合模型mi與基準混合模型mstandard差異較小,可認為兩者距離近似服從正態(tài)分布。 因此,將報警閾值T設定為μkd±3σkl。 獲得報警閾值T后,當實時工況混合模型monline與基準混合模型間mstandard的KL 散度超過T時,認為設備發(fā)生故障,觸發(fā)報警。

3 試驗驗證與結果分析

為全面、準確地評估提出方法的有效性,分別利用往復壓縮機故障模擬試驗數據(數據集A)和工程案例數據(數據集B)對方法進行驗證。 為驗證提出方法的準確性,結合數據分析和測試結果,統計報警結果中存在的誤報、 漏報次數并計算預警準確率。 誤報、漏報均認定為錯誤報警,預警準確率以單個測試集為最小計算單位,每個測試集包含244 個測試樣本。 預警準確率A計算公式表示為

式中,DL表示單測試集包含的樣本容量,本文中DL均為244,FP表示誤報次數,UP表示漏報次數。 分別統計數據集A 和數據集B 中測試集的預警準確率,計算其平均預警準確率,以此對提出方法的準確性進行定量評估。 為驗證提出方法的時效性,從提前預警時間和生成混合模型耗時兩方面進行評估,預警時刻與故障時刻時間差越大,生成混合模型耗時越短,表明預警方法的預警時效性越好。 為驗證提出方法的可靠性,從預警指標變化趨勢方面進行評估,通過預警指標的平均標準差來表征其平穩(wěn)性。為反映提出方法的先進性,將變分自編碼器(variational auto-encoder,VAE)方法[5]、無限學生t 混合模型(infinite student's t-mixture model,iSMM)方法[7]和基于數據驅動隨機森林(random forest,RF)的風電機組變槳軸承磨損預警建模方法[8]作為對照試驗進行方法對比。

3.1 數據描述

(1)數據集A

數據集A 來自往復壓縮機試驗臺,包含8 種故障類型數據,即吸氣閥閥片斷裂、吸氣閥彈簧失效、排氣閥閥片斷裂、排氣閥彈簧失效、十字頭銷磨損、活塞環(huán)斷裂、活塞銷松動和基座松動,共計21 個測試集;設備正常狀態(tài)數據作為訓練集,數量為2,具體介紹如表1 所示。 為避免設備連續(xù)狀態(tài)數據對方法驗證造成誤差,訓練集與測試集的正常狀態(tài)數據分別來源于故障模擬試驗前后。 往復壓縮機試驗臺的監(jiān)測傳感器測點布置如圖3 所示,分別在缸頭、十字頭上方安裝加速度傳感器,曲軸箱殼體兩側安裝速度傳感器,傳感器采樣頻率為10 240 Hz,采樣點數為6 144,往復壓縮機轉速為500 r/min,采樣長度為往復壓縮機兩個運轉周期。

表1 數據集A 情況Table 1 Information for data set A

圖3 測點布局及傳感器類型Fig.3 Layout of the measuring point and sensor types

(2)數據集B

數據集B 來自某石化企業(yè)生產現場的往復壓縮機,故障案例數據包括活塞組件磨損、氣閥泄露和液擊3 種類型故障,每種機組類型分別選取1 組正常狀態(tài)數據作為訓練集,具體介紹如表2 所示。 數據采樣頻率為10 kHz,采樣點數為5 000,不同型號往復壓縮機轉速不同,轉速在300 ～350 r/min 范圍內,采樣間隔為往復壓縮機的兩個運轉周期。

相關測試表明,樣本容量的大小對于預警方法的報警準確率具有重要影響[7]。 隨著樣本容量的增大,預警方法的準確率隨之升高,但模型訓練耗時也隨之增加;當樣本容量超過100 時,預警方法的報警準確率受樣本容量影響較小。 因此,綜合考慮預警準確率和模型的訓練效率,樣本容量設定為100。

3.2 特征描述

為全面、準確地反映往復壓縮機的運行狀態(tài),通過分析往復壓縮機旋轉和往復兩種運動的形式特點,依據往復壓縮機不同類型故障所關聯的振動信號敏感特征,分別從時頻域和角度域選取了20 種特征構成模型訓練和測試的特征集,如表3 所示。

表3 特征集說明Table 3 Details of the feature set

3.3 結果與分析

分別采用試驗數據和工程案例數據對提出方法和對比方法進行測試,測試結果如表4 所示。

表4 方法測試結果統計Table 4 Method test result statistics

(1)數據集A 驗證結果

由表4 可見,對于試驗模擬的8 種故障,所提方法的平均預警準確率達到99%以上,相較于VAE方法高出1.28%,相較于iSMM 方法高出7.12%,相較于RF 方法高出24.28%,說明提出方法的預警準確率顯著優(yōu)于3 種對比方法。 在生成混合模型耗時方面,提出方法的平均生成混合模型耗時最短,僅為0.314 s, iSMM 方法耗時為提出方法的5 倍左右,VAE 方法耗時為提出方法的3 倍左右,RF 方法耗時為提出方法的2 倍左右,說明所提方法的時效性顯著提高。 當設備處于正常狀態(tài)時,提出方法的KL 散度即模型差異的平均標準差為1.02,較3 種對比方法顯著降低。

(2)數據集B 驗證結果

由表4 可見,對于27 組案例數據所包含的3 種故障,提出方法對27 組不同故障案例的平均識別準確率為100%,與VAE 方法相當,相比之下,RF 方法的預警準確率最低。 提出方法平均生成混合模型耗時僅為0.198 s,iSMM 方法耗時約為提出方法的7倍,VAE 方法耗時約為提出方法的3 倍,RF 方法約為提出方法的2 倍。 當設備處于正常狀態(tài)時,提出方法的KL 散度的平均標準差為1.35,較VAE 方法下降75%,較iSMM 方法下降83.3%,較RF 方法下降88.9%。 數據集B 中27 組案例數據來自不同的往復壓縮機機組設備,包含多種機組結構形式,提出方法均能有效實現預警,說明提出方法對于往復壓縮機具有良好的適應能力。

驗證結果表明:對于多種往復壓縮機故障,提出方法的預警準確率優(yōu)于VAE 方法、iSMM 方法和RF方法,生成混合模型的耗時顯著減少,預警時效性顯著提高。 在設備處于正常狀態(tài)時,提出方法的KL散度平均標準差最小,預警指標隨時間變化趨勢的平穩(wěn)性最強,說明該方法對正常運行狀態(tài)產生的非平穩(wěn)信號具備較強的魯棒性。 分析發(fā)現,iSMM 方法未能對4 缸-M 型機組活塞組件磨損故障的2 組故障樣本實現報警,這是由于當活塞組件磨損程度較輕時,因故障產生的數據離群點數量較少,致使iSMM 方法無法識別設備活塞組件磨損故障。 VAE方法和iSMM 方法生成的混合模型之間無有效關聯信息,在每次生成混合模型時都需要對混合模型先驗參數進行重新計算,該過程不僅導致生成混合模型耗時增加,同時在計算模型差異時,需要匹配混合模型之間子模型的對應關系,該過程存在不可避免的計算誤差,造成預警指標變化趨勢的平穩(wěn)性較差。RF 方法采用隨機森林算法,擅長處理高維數據和特征遺失數據,對于低維數據(特征較少的數據)的分類表現較差。 因此,對于往復機械活塞環(huán)磨損、拉缸等故障特征較少的機械故障,RF 方法的識別準確率較低。

以工程案例中2 缸D 型機組氣閥泄漏故障為例,通過將預警指標可視化,以圖像的形式直觀地將提出方法和3 種對比方法的驗證結果進行展示說明。 設備正常狀態(tài)和氣閥泄漏故障振動響應信號波形圖如圖4 所示,由圖4(a)可看出,往復壓縮機處于正常狀態(tài)時,特定相位的沖擊信號符合往復壓縮機在運轉過程中氣閥落座、彈簧釋放等部件的動作狀態(tài)。 由圖4(b)可看出,發(fā)生氣閥泄漏故障時,氣閥開閉相位附近出現多種頻率信號成分,信號幅值整體減小,表明氣閥泄漏故障造成振動激勵。

圖4 氣閥不同狀態(tài)的角域振動信號對比Fig.4 Comparison of valve vibration signals in different states

提出方法和3 種對比方法對于氣閥泄漏故障的測試結果如圖5 所示。 在氣閥泄漏故障案例測試集中,氣閥泄漏故障發(fā)生時刻為第148 h。 如圖5(a)所示,VAE 方法在第58 h 開始報警,相較故障發(fā)生時刻提前90 h;從預警指標(KL 散度)變化趨勢看,當設備處于正常運行狀態(tài)時,VAE 方法報警指標趨于平穩(wěn),但波動較為明顯,存在穿越報警線現象;在報警起始時刻至故障發(fā)生階段,報警指標呈緩慢上升趨勢,但仍在報警線上下波動,部分時刻預警指標低于報警線。 如圖5(b)所示,iSMM 方法在第67 h開始報警,相較故障發(fā)生時刻提前81 h;當設備處于正常運行狀態(tài)時,iSMM 方法報警指標波動范圍較大,可明顯觀測到穿越報警線現象。 如圖5(c)所示,RF 方法在第57 h 時開始報警,相較故障發(fā)生時刻提前91 h,從圖中可看到在報警時刻后,預警指標表現十分不穩(wěn)定,且在故障發(fā)生時刻之后預警指標回落到報警線下,存在較為明顯的漏報問題。 如圖5(d)所示,DtMM 方法在第43 h 開始報警,相較故障發(fā)生時刻提前105 h,預警時間比iSMM 方法提前24 h,比VAE 方法提前15 h,比RF 方法提前14 h;設備處于正常運行狀態(tài)時,DtMM 方法報警指標變化趨于平穩(wěn),且未發(fā)現較明顯波動;在設備故障早期階段和故障發(fā)生階段報警指標均高于報警線。 通過以上對比可以看出,提出方法可準確識別往復壓縮機此類復雜機械的異常狀態(tài),設備故障早期階段預警時間明顯增加,同時預警指標變化趨勢的平穩(wěn)性顯著優(yōu)于上述3 種方法。

圖5 4 種方法KL 散度測試結果Fig.5 KL dispersion test results of the four methods

驗證結果表明:在往復壓縮機運行狀態(tài)監(jiān)測過程中,當設備處于正常運行階段時,提出方法能夠有效克服非平穩(wěn)信號對預警指標造成的干擾,預警指標變化趨勢的平穩(wěn)性顯著增強,有效降低了誤報率;在設備故障早期階段,提出方法預警時間明顯增加,能夠在故障早期階段準確識別往復壓縮機的故障信號,預警時效性顯著提高。 提出方法的預警準確率和時效性明顯優(yōu)于3 種對比方法。 VAE 方法深入振動響應信號組成成分,對激勵源信號進行分析,但未考慮前一時間段的激勵源信號分布對后續(xù)信號變化的影響,每次生成混合模型都需要對混合模型先驗參數重新估計,生成的混合模型之間無關聯信息,因此計算實時工況混合模型與基準混合模型的差異時,需事先匹配子模型的對應關系,這一過程不僅增加了計算耗時,同時造成模型魯棒性較差,在設備正常狀態(tài)下易發(fā)生誤報。 iSMM 方法準確擬合了機械振動響應信號中的激勵源信號統計分布,但其生成混合模型機制同VAE 方法一致,造成在設備正常狀態(tài)時,報警指標易受非平穩(wěn)信號影響而產生大幅波動,導致預警模型可靠性較差。 RF 方法對低維數據的分類效果較差,對于往復機械部分敏感關聯特征較少的故障,該方法無法有效識別,導致預警模型的往復壓縮機監(jiān)測能力均存在一定的局限性。

4 結論

(1)本文提出一種基于動態(tài)主題模型的往復機械故障早期預警方法。 該方法在基準混合模型基礎上生成新的混合模型,而非完全基于數據,利用狀態(tài)空間轉移鏈接混合模型間的參數,利用輸入數據更新模型部分參數。 這一過程確立了混合模型中的子成分位置,依據輸入數據確定子成分分布,在計算混合模型與基準混合模型差異時省去了子成分匹配計算過程。 因此,提出方法預警模型的魯棒性更強,對于非平穩(wěn)性信號的適應性更好。

(2)分別利用往復壓縮機工程案例數據和試驗臺故障數據對提出方法的有效性進行了驗證。 結果表明,對于往復壓縮機故障早期預警,所提方法生成混合模型耗時減少、預警時間提前、預警誤報率低,取得了很好的效果。