基于M估計的轉子不平衡量的識別

魏巍，張海鵬，談演

摘 要：轉子不平衡問題是旋轉機械故障中經常發生的問題，因此對不平衡量大小相位的識別至關重要。利用傳統的最小二乘法計算的轉子不平衡量受異常值影響較大。加權最小二乘法需要對振動數據的處理較為繁瑣并且權重大小難以確定。現提出一種基于M估計的轉子不平衡量的識別方法，該方法具有穩健性，受異常值的影響較小，且便于計算權重。通過實驗證明，該方法可以得到轉子的不平衡量。

關鍵詞：動平衡；影響系數法；最小二乘法；M估計

中圖分類號：V231.9 文獻標志碼：A 文章編號：2095-2945（2018）06-0087-03

Abstract： Rotor unbalance is a common problem in rotating machinery fault， so it is very important to identify the phase of unbalance. The rotor unbalance calculated by the traditional least square method is greatly affected by the outliers. The weighted least square method needs to deal with the vibration data more cumbersome and the weight is difficult to determine. In this paper， a method of rotor unbalance recognition based on M estimation is proposed. This method is robust， less affected by outliers， and easy to calculate weights. It is proved by experiments that the unbalance of rotor can be obtained by this method.

Keywords： dynamic balance influence coefficient method； least square estimation； M estimation

引言

轉子的不平衡量對轉轉子影響非常大，過大的不平衡量會直接損壞機器，所以對轉子不平衡量的準確識別具有很大的意義。傳統的影響系數法可以簡單方便識別，但存在病態方程。最小二乘法和影響系數法的結合能夠解決病態方程，但不具有穩健性；加權最小二乘法需要對轉子的振動數據進行深入分析。本文通過M估計對轉子的不平衡量進行識別，具有較好的穩健性，通過實驗證明，M估計可以有效識別轉子的不平衡量。

1 轉子回歸模型

影響系數法是轉子動平衡最常用的方法。但是影響系數法存在病態方程，所以現經常采用影響系數法和最小二乘法相結合。影響系數矩陣可采用有限元模型或歷史平衡記錄獲得。根據線性振動理論，在m個轉速，n個測振點，q個平衡平面下，轉子平衡方程為：

[L]=[L0]+[A]×[W]（1）

其中，[L]為各測點的殘余振動；[L0]為轉子在不平衡量下引起的各測點的初始振動；[A]為影響系數；[W]為校正質量。

矩陣形式表達為：

[L]、[L0]為p=m×n的列陣，[W]為q個平衡質量列陣。{A}是影響系數p×q的矩陣。當P=q時，方程有唯一解，而在實際平衡過程中經常出現P>q的情況，此時，需要對方程采用最小二乘法的原理，即最小化殘差平方和，用矩陣形式表示有：

[W]=-（[A]T[A]）-1[A]T[L]

式（2）又可以轉化為多元回歸問題：

li=l0i+ai1·w1+ai2·w2+…+aiq·wq（i=1，2，3，…，p）（3）

2 M估計

穩健估計基本可以分為三大類型：M估計、L估計、R估計。其中M估計（Maximum Likelihood Type Estimates）是經典極大似然估計的推廣，稱為廣義極大似然估計，本文主要討論M估計。

設觀測樣本L1L2…Ln，X為待估參數，觀測值的分布密度為f（li，■），極大似然估計準則為

■lnf（li，■）=max

用ρ（*）代替-lnf（*），則極大似然估計準則可以改為

■ρ（li，■）=min

對上式求導，得：

■φ（li，■）=0

其中：φ（li，■）=■

由以上可知，有一個ρ函數，就定義了一個M估計。

幾種常見的ρ函數

（1）huber函數

通常情況下k=1.345，既能保證穩健性，又具有較高的估計效率。

（2）tukey函數

通常k=4.685，估計方法既穩健又有較高的估計效率。

（3）Hampel函數

通常取a=1.5，b=3，c=8，Hampel函數更為復雜一些，但當|u|>c時權因子取到了0，因此它具有更好的穩健性。

其中u=■

■為尺度參數，取絕對離差中位數：■=median（|v-median（v）|）

各個函數的目標函數、導函數、權函數如圖1。

M估計的基本算法：

列誤差方程，令各個權因子初始值為1，即w1，w2…wn=1，W=I

利用最小二乘法求得參數■和殘差V的第一次估計

■（1）=（BTPB）-1BTPl

V（1）=B■（1）-l

（3）根據V（1）求出尺度參數，并標準化殘差u=■，將標準化的殘差帶入到不同的權函數中，得到權重矩陣。再求得參數■和殘差V的第二次估計

■（2）=（BTPB）-1BTPl

V（2）=B■（2）-l

（4）再由V（2）構造新的權數，再解算方程，類似迭代計算，直到前后差值符合限差要求。

M估計的目標曲線不同，但都是對殘差偏大的觀測值進行降權處理。M估計的崩潰點不是很高，但是其效率較高。

3 實驗驗證

通過對轉子進行動平衡實驗來驗證該方法。轉子實驗臺如圖2所示。

轉子影響系數矩陣通過多次實驗得到。在確定轉子無故障前提下，在平衡面A處0°位置添加已知不平衡量0.8g，測量轉子在此不平衡下兩個測點的振動信息如圖3所示。

現已知不平衡質量，對計算效率和穩健性等因素綜合考慮，采用M-huber。用最小二乘法和M-huber估計計算轉子不平衡量的大小和相位，結果如表1所示。采用不平衡量實際值w和估計值■的相對誤差比較兩種方法的效果，相對誤差為

M-huber估計在各個轉速下的權重如圖4所示

通過以上結果可知，利用M-Huber估計可以更好的識別轉子的不平衡量，相對誤差大大減少。根據M-huber估計得到的權重在轉子一階臨界轉速下減小，與轉子一階臨界轉速下轉子振動數據存在異常值一致。

4 結束語

在綜合考慮計算效率和穩健性等因素，采用M-huber估計。通過M估計得到的不平衡質量與傳統的最小二乘法影響系數法相比較，具有更高的精確度，相對誤差大大降低。應用M估計識別轉子不平衡量，降低了轉子一階臨界轉速附近的權重，與轉子一階臨界轉子下振動數據存在異常一致。

參考文獻：

[1]繆海林.基于模型的轉子系統不平衡量識別的魯棒方法研究[D].南京航空航天大學，2012.

[2]王海娜.線性回歸模型的若干穩健估計方法及應用實例[D].山東大學，2013.

[3]Pennacchi P，Chatterton S，Ricci R.Rotor balancing using high breakdown-point and bounded-influence estimators[J].Mechanical Systems and Signal Processing，2010，24（3）：860-872.

[4]楊慶華.用最小二乘法影響系數法對高速柔性轉子的動平衡研究[J].機電工程，1998.

[5]Pennacch P，Vania A，Bachschmid N.Increasing the robustness of fault identification in rotor dynamics by means of m-estimators[J].Mechanical Systems and Signal Processing.2007，21（8）：3003-3029.