基于參數優化VMD和多尺度排列熵的齒輪故障診斷

謝棕，張麗萍

(福州大學 機械工程及自動化學院，福建 福州 350108)

0 引言

齒輪及齒輪箱作為機械設備常用的調節轉速和傳遞轉矩的旋轉機械設備，不僅能夠傳遞較大的功率和載荷，而且具有較好的可靠性[1]。但是在高精度的切削加工中，當齒輪在變轉速、變載荷等復雜工況下工作，極易受到損傷、產生磨損、斷齒等情況，使得加工精度大打折扣。因此，齒輪的狀態監測和故障診斷變得尤為重要。

齒輪的故障信號多為非線性、非平穩信號。經驗模態分解(EMD)和集合經驗模態分解(EEMD)在分析非線性、非平穩信號時具有良好的優越性，但是存在嚴重的端點效應和模態混疊問題[2]。基于此，DRAGOMIRETSKIY K[3]在2014年提出一種自適應、非遞歸模式的信號分解方法——變分模態分解(variational mode decomposition,VMD)，它很好地解決了EMD算法存在的問題，并且具有更堅實的數學理論基礎。但VMD算法在分解過程中會同時受分解參數K和懲罰因子α這兩個參數的影響。因此，如何快速、準確地確定這兩個參數是研究的重點。劉長良等[4]提出一種觀察中心頻率的方法確定分解參數K值，當開始出現中心頻率相近時，認為出現了過分解，以此來確定K值，但卻忽略了懲罰因子α對分解結果的影響。瞿紅春等[5]利用手動尋優對VMD參數進行優化，先假定α不變，不斷尋優K值使得包絡熵最小，再將此K值設為定值，不斷尋優α值使得包絡熵最小，得出分解參數K和懲罰因子α。但這樣得出的K和α并不是最優值，且工作量巨大。基于以上分析，本文提出以IMF分量局部包絡熵最小為優化目標，利用改進的蝙蝠算法優化VMD分解參數K和懲罰因子α，進而得到最優的分解參數，通過計算K個IMF分量的多尺度排列熵(multiscale permutation entropy,MPE)構造故障特征向量，輸入到極限學習機(extreme learning machine,ELM)進行訓練和識別，進而實現齒輪的故障診斷。

1 參數優化VMD

1.1 VMD原理

VMD分解是一種非遞歸模式的信號分解方法，通過迭代搜尋約束變分模型最優解，將輸入信號f分解成k個離散的模態uk(t)，使得每個模態的估計帶寬之和最小，約束條件為各模態之和等于輸入信號f。

構造相應的約束變分模型如下：

其中：{uk}={u1,…,uk}，表示分解得到的k個模態函數；{wk}={w1,…,wk}，表示各模態分量的中心頻率[6]。

為求解式(1)的最優解，引入拉格朗日乘法算子λ(t)以及二次懲罰因子α，將約束性變分問題變為非約束性變分問題。λ(t)讓約束條件保持嚴格性，α能充分減少高斯噪聲對信號的干擾：

(2)

利用交替方向乘子算法不斷迭代更新各個模態及其中心頻率[7]，尋求式(2)的靶點。具體實現算法如下(n為迭代次數)：

2)迭代k=1∶K并更新各模態信號uk和中心頻率ωk，i∈[1,K]。

更新uk：

(3)

更新ωk：

(4)

3)根據下式更新λ，Γ為更新因子

(5)

4)若滿足下方收斂條件則停止迭代，否則重復1)、2)，ε為>0的正數。

(6)

結束迭代，得到K個IMF分量。

1.2 改進蝙蝠算法

蝙蝠算法(bat algorith,BA)是由英國學者YANG X S[8]在2010年提出的一種群智能優化算法。BA算法把求解問題目標函數的適應度值度量成蝙蝠個體所處位置的優劣。最優解類比蝙蝠的獵物，脈沖響度和頻率的變化在一定程度上表明與最優解的靠近程度。

在BA算法中，搜索頻率、位置和速度的更新公式如下：

fi=fmin+(fmax-fmin)β

(7)

(8)

(9)

其中：fi、fmax、fmin分別表示第i只蝙蝠當前時刻發出的聲波頻率以及聲波頻率的最大值與最小值；β∈[0,1]的隨機數；x*表示當前最優位置。

在局部搜索中，如果一個蝙蝠選擇了一個最優解，那么將在此最優解附近隨機產生一個新解：

Xnew=Xold+εAt

(10)

其中：Xold為從當前最優解中選擇的一個解；At為t時刻所有蝙蝠響度的平均值；ε是[-1,1]的一個隨機數。

蝙蝠脈沖的響度Ai和速率ri通過下式計算：

(11)

(12)

其中：α和γ是常數，并且0<α<1，γ>0。響度Ai和速率ri的持續迭代更新表明與最優解的靠近程度。

蝙蝠算法能夠在算法運行的前期通過將全局優化轉換到局部優化來實現算法的快速收斂，這同時也會導致該算法過早地處于停滯階段。為改善這樣的問題，本文在速度更新公式(8)中引入了線性遞減慣性權重ω(k)的概念：

(13)

ω(k)=wmax(wmax-wmin)(Tmax-k)/Tmax

(14)

其中：ωmax為初始慣性權重；ωmin為迭代至最大次數時的慣性權重；k為當前迭代次數；Tmax為最大迭代次數。查閱相關文獻知，ωmax=0.9，ωmin=0.4時算法性能最好[9]。迭代初期較大的初始慣性權重使算法保持了較強的全局搜索能力，而迭代后期較小的慣性權重有利于算法進行更精確的局部搜索。

1.3 VMD參數優化

VMD分解前，需要預先設定分量個數K和懲罰因子α的值。K和α的取值不同，會有不同的分解結果。K值過小，會造成模態混疊；K值過大，會產生虛假分量。α值過小，各IMF分量的帶寬越大；α值過大，各IMF分量的帶寬越小。因此，如何確定合適的K和α的值是研究的關鍵所在。

利用蝙蝠算法對VMD進行參數優化時，需要選定一個適應度函數。本文采用VMD分解后的局部極小包絡熵作為蝙蝠算法的適應度函數，搜尋最優的VMD參數組合。信號pj的包絡熵Ep可以表示成：

(15)

(16)

(17)

式中：pj是a(j)的概率分布序列；a(j)為信號x(j)經過Hilbert解調后得到的包絡信號。

改進BA算法優化VMD參數步驟如下：

步驟1 BA算法基本參數初始化。主要有種群大小sizepop、最大迭代次數maxgen、參數α和K的搜索范圍、最大脈沖響度A0、最大脈沖速率r0、脈沖頻率范圍、音量的衰減系數α以及搜索頻率的增強系數γ。

步驟2 初始化種群，以局部極小包絡熵作為適應度函數求出一組最優的α和K。

步驟3 根據式(7)、式(13)、式(9)不斷迭代更新種群的搜索脈沖頻率、速度和位置。

步驟4 生成均勻分布隨機數rand，如果rand>ri，則對當前最優解進行隨機擾動，產生一個新的解。

步驟5 如果計算出來的包絡熵小于所得出的極小包絡熵，對最優的α和K以及最小的包絡熵進行更新，然后按式(11)、式(12)對響度和頻率進行更新。

步驟6 重復步驟3-步驟5，直至達到最大迭代次數以及確定最小的包絡熵值，得出最佳參數α和K。

2 多尺度排列熵

排列熵是由BANDT C等[10]提出的一種衡量一維時間序列復雜度的平均熵參數，常用來提取機械故障的特征。但齒輪及齒輪箱運轉過程中的故障特征信息分布在多尺度中，僅用單一尺度的排列熵進行分析，會遺漏其余尺度上的故障特征信息。

針對排列熵的不足，COSTA M等[11]提出了多尺度排列熵的概念。多尺度排列熵是在排列熵的基礎上將時間序列進行多尺度粗粒化，然后計算不同尺度下粗粒化序列的排列熵，具體計算過程如下：

(18)

式中：s為尺度因子；[N/s]表示取整。

(19)

式中：l表示第l個重構的分量；τ表示延遲時間；m表示嵌入維數。

3)將式(19)得到的重構時間序列按升序排列可得

S(g)=(j1,j2,…,jm)

(20)

(21)

4)當Pg=1/m!時，每種符號序列的概率都相等，此時時間序列的復雜程度最高，排列熵最大，為lnm!。為方便表示，通常會將Hp(m)進行歸一化處理，即

Hp=Hp(m)/ln(m!)

(22)

Hp的取值范圍為[0,1]，Hp的大小反映了時間序列的復雜和隨機程度。Hp越大，說明時間序列越隨機，反之，說明時間序列越規則。

3 算法流程

本文提出基于參數優化VMD和多尺度排列熵的齒輪故障診斷方法，具體步驟如下：

1)利用改進的蝙蝠算法優化VMD分解參數K和懲罰因子α，進而得到最優的分解參數，得到K個模態分量IMF1、IMF2、…、IMFK。

2)計算被選模態分量的多尺度排列熵，選取合適尺度的排列熵值作為故障特征向量。

3)將故障特征向量輸入到已經訓練好的極限學習機中進行訓練和識別。

故障診斷流程如圖1所示。

圖1 齒輪故障診斷流程

4 實驗分析

4.1 實驗說明

為了驗證本文所提方法的有效性，在實驗平臺下(圖2)，采集試驗臺減速箱齒輪不同狀態下的實時振動信號。齒輪減速箱型號為JZQ200，采樣頻率為8 200 Hz。選擇齒輪正常、不平衡、松動、錯位、斷齒、裂痕這6種不同狀態的數據進行實驗驗證。

圖2 齒輪實驗平臺

4.2 故障特征提取

以減速箱齒輪出現裂痕為分析對象，本文在以下的實驗中選取的蝙蝠算法參數均為：sizepop=36，參數α和K的搜索范圍low=[1 500 3]，high=[2 500 8]，速度v的范圍為[-100 100]，A0=0.5，r0=0.25，脈沖頻率范圍[0 1]，α=0.8，γ=0.5。蝙蝠算法不同參數組合的優化如表1所示。

表1 蝙蝠算法不同參數組合的優化

圖3為局部極小包絡熵值隨尋優迭代次數變化曲線，優化結果如表2所示。

圖3 尋優迭代圖

表2 優化結果

作為對比，同時采用手動尋優的方法來求解最優的(K,α)的值。VMD參數設置如下：初始化alpha=2 000，tau=0，DC=0，init=1，tol=le-7。

表3 不同K值下的包絡熵值

由于篇幅原因，上述表格每個K值僅陳列5個極小包絡熵值，觀察到最小的包絡熵值僅為5.554 2。對比蝙蝠算法的迭代尋優結果，通過手動尋優求得的最小包絡熵值并不理想。

綜合分析，改進蝙蝠算法迭代尋優的最優參數K=6，α=247 4。選取裂痕齒輪一組數據中的4 500個樣本點進行VMD分解，得到6個IMF分量，分解結果如圖4所示。

圖4 VMD分解結果

MPE的計算需要設置嵌入維數m，延遲時間τ和尺度因子s。其中維數m取值范圍通常是3～7，m取值太小，算法的檢測性能降低，若m取值太大，將無法反映時間序列的細微變化。延遲時間τ對時間序列的計算影響較小，尺度因子s最大值一般取>10即可，本文選取m=3，τ=1，s=20。

于是就能求出前面VMD分解后6個IMF分量的排列熵，如圖5所示。

圖5 裂痕齒輪各IMF分量的MPE值

4.3 齒輪的故障診斷

對測得的正常、不平衡、松動、錯位、斷齒、裂痕這6種不同狀態的數據，每種數據取50組，每組4 500個點。對這50組數據進行VMD分解后并計算模態分量的MPE值，6個MPE分量，每個分量提取20尺度的因子，所以每個樣本就得到維度為120的特征。將故障特征向量輸入到已經訓練好的ELM中，取前40組數據作為已知故障樣本，后10組數據作為待識別故障樣本，診斷情況如圖6所示。

圖6 ELM故障識別

圖6中，橫坐標分別代表訓練樣本和測試樣本，縱坐標代表6種故障類別的標簽。測試樣本中的空心點代表預測類別。實心點代表實際類別。當測試的數據符合同類故障的標準時，空心點和實心點會重合，代表識別正確。本文測試集總的準確率達到了96.7%，對于絕大多數的故障情況，都能達到很好的識別效果。

5 結語

針對齒輪故障特征在單一尺度難以全面提取的問題，提出一種基于參數優化的VMD和多尺度排列熵的齒輪故障診斷方法。相比較手動尋優，該方法可以搜尋更優的VMD參數組合，更有效地提取出不同故障狀態下的信號特征參數。分析結果表明，通過參數優化VMD和多尺度排列熵的齒輪故障診斷方法，極大地提高了齒輪故障診斷的準確性。