基于代理模型的轉子系統不平衡定量辨識

顧煜炯， 陳東超， 徐 婧， 何成兵

（華北電力大學 國家火力發電工程技術研究中心，北京102206）

旋轉機械是電廠、化工廠等生產單位的關鍵設 備，轉子系統不平衡是旋轉機械較為常見的一種故障，危及設備的安全運行.這種不平衡可能來自轉子制造階段，也可能在機組投運后由于運行或檢修原因產生.準確識別和診斷出不平衡量的位置、大小和相角，以便快速、及時地采取有效措施消除故障，對于旋轉設備的安全穩定運行具有重要意義.

對于振動信號的時域、頻域以及時頻域分析已廣泛應用于不平衡、碰磨和裂紋等故障的診斷［1］.這類方法一般被稱為基于信號的故障診斷方法，但也只能診斷出有無不平衡故障，不能給出不平衡的定量信息.另一種方法是基于模型的診斷方法，把測量的振動信號與結構模型聯合起來，對故障程度進行定量估計［2－6］.在基于模型的不平衡診斷和辨識方法中，主要有2種思路：（1）在每次假設一組不平衡參數后，將該組不平衡參數對應的不平衡力作用在轉子系統的有限元模型上，通過有限元計算得到測點處的振動響應.通過多次冗繁的有限元模型計算，尋找計算出的振動響應與實測振動響應最接近的不平衡參數，該參數便為所求.實際旋轉機械的轉子結構大多較為復雜，這種方法計算量相當繁雜、耗時；（2）轉子系統采用有限元方法進行建模，故障轉化為施加在模型上的等效外力，將實際系統中由于故障產生的等效力和理論故障力的絕對偏差值作為目標函數，采用最小二乘法等求解出使偏差最小時的故障參數，即得到故障的定量辨識結果.由于實際轉子系統中只有少數自由度布置振動測點，一般需要根據系統的振型，采用模態擴展的方法求得所有自由度上的等效力，但是模態擴展過程中會出現模態截斷現象，使故障參數辨識誤差增大；另外，在模態擴展過程中所使用的模態矩陣有可能是個比較嚴重的病態矩陣，無法準確辨識不平衡參數.

近年來，代理模型開始廣泛應用于工業設計等領域［7］，為了克服以往不平衡定量辨識方法的缺點，筆者提出了一種基于代理模型的轉子系統不平衡定量辨識方法.支持向量回歸（SVR）的泛化能力強，并且對小樣本具有很好的學習及預測能力，因此采用粒子群算法（PSO）優化的SVR（PSO－SVR）構造代理模型，建立不平衡參數與各測點振動響應的對應關系，降低反演優化迭代過程中冗繁的有限元計算；構建用于求解不平衡參數的目標函數，并借助建立的代理模型，采用PSO 算法尋找符合目標函數的全局最優解，從而達到不平衡定量辨識的目的.

1 不平衡參數辨識流程

所提出的不平衡辨識策略將拉丁超立方抽樣（LHS）算 法［8］、有 限 元 模 型 計 算［9］、SVR［10－11］和

式中：wi為權重參數，取值范圍為（0，1），和值為1；qn、ψn 分別為位置n 處不平衡的偏心量和相角；uic（qn，ψn）為在第i個測點處，由代理模型預測出的振動響應；uim為實測的振動響應.

如圖1所示，基于代理模型的不平衡參數辨識步驟可歸結為：

（1）利用LHS算法，生成不同偏心量和相角的不平衡參數初始樣本點，并通過有限元模型計算得到不同位置上施加以上不平衡力后各測點處的振動響應.PSO［12］算法等多種方法集成在一起，每種方法在整個辨識流程中所承擔的任務不同.集成思路及各方法的作用主要有：

（1）該不平衡辨識方法的主要目的是降低反演優化迭代過程中冗繁的有限元計算，然而這并不代表不需要任何的有限元計算.這是因為實際轉子系統運行過程中的不平衡故障樣本很少，需要采用有限元模型模擬實際轉子的運行過程，在有限元模型上施加不同不平衡工況下的故障力，并通過一定量的有限元模型計算得到構建代理模型所需的訓練樣本.

（2）LHS算法的主要任務是用來得到不平衡的大小、相角等參數樣本，是有限元計算產生訓練樣本的前一步工作.LHS算法可避免重復抽樣，同人為的選取不平衡參數樣本相比，其選取的樣本能更好地反映總體變異規律，更有利于獲得計算準確、覆蓋全面的代理模型.

（3）通過有限元模型獲得用于代理模型訓練的樣本后，選取SVR 進行代理模型的訓練.這是因為與響應面法、神經網絡等代理模型相比，SVR 泛化能力強，且對小樣本具有很好的學習及預測能力，模型所需訓練樣本的數量小，盡可能地降低了有限元模型在生成訓練樣本過程中的計算量.

（4）PSO 算法是一種容易實現、精度高、收斂快的尋優算法.PSO 算法在所提出的不平衡定量辨識策略中主要有2個作用：一方面，優化選取SVR 模型中的懲罰因子C、徑向基核函數的寬度參數g 和不敏感系數ε，避免憑經驗選取造成較大的回歸誤差；另一方面，不平衡識別問題可以歸結為優化問題，即在給定的搜索域內，借助建立的代理模型，通過PSO 算法尋找符合目標函數的不平衡參數，從而達到不平衡定量辨識的目的.假設已經測得轉子系統振動響應，目標函數可以表示為：

圖1 基于PSO－SVR代理模型的不平衡識別步驟Fig.1 Flowchart of unbalance identification based on PSO－SVR surrogate model

（2）以不平衡的偏心量和相角為輸入、振動響應為輸出，依次建立不平衡位于不同位置時，各測點振動響應的PSO－SVR 代理模型.

（3）根據各測點振動實測值和構建好的各PSO－SVR代理模型構建式（1）所示目標函數.

（4）采用PSO 算法求解步驟（3）中的多目標優化問題.依次在可能出現不平衡的位置進行不平衡辨識，當代理模型預測值與實測值誤差小于設定的閾值時，停止計算，此時的位置為不平衡所在位置，相應的最優解取為不平衡的定量辨識結果.

2 模型建立

2.1 拉丁超立方抽樣

Mckay等人在1979年提出了LHS算法［8］，其主要思想是基于逆函數轉換法，假設抽樣次數為N，輸入K 個隨機變量X1，X2，…，XK，變量Xk（1≤k≤K）的累積分布函數（CDF）為Yk＝Fk（Xk），此分布函數的取值范圍［0，1］被等分成N 個子區間，每個子區間內抽取一個隨機數Yk，然后由CDF的逆函數求得Xk，每個子區間都不重復抽取.

2.2 轉子系統建模

通常認為轉子系統是由一系列的剛性圓盤、帶有質量和彈性的軸段及軸承組成.依次建立各部分的運動方程后，可得轉子系統的組合運動方程［9］：

式中：M 為質量矩陣，包括軸系單元和圓盤的平動和轉動質量矩陣；C 為阻尼矩陣，包括陀螺力矩和軸承的阻尼；K 為剛度矩陣，包括軸段單元剛度和軸承剛度；Q 為激勵矩陣，包括圓盤的重量、不平衡激振力、軸承力、由于轉子故障引起的外力和其他外部力.

建立有限元模型時，每個節點考慮4個自由度：水平和垂直2個獨立方向的平移運動，圍繞水平軸和垂直軸2個方向的旋轉運動.關于轉子系統的詳細模型可參考文獻［5］，這里不再贅述.

2.3 不平衡故障模型

設轉子的速度為ω，作用在單個節點n、偏心量為qn、相角為的不平衡轉化為在水平和垂直方向的等效力，其表達式分別為：

2.4 PSO－SVR代理模型

2.4.1 支持向量回歸

SVR 模型就是采用支持向量回歸技術得到輸入和輸出間的映射關系，以代替計算量大且復雜的有限元模型等.

支持向量回歸的基本思想如下［10－11］：觀測樣本集為｛（xi，yi），i＝1，2，…，n，xi∈RN，yi∈R｝.設回歸函數為

式中：b為偏差量.

根據結構風險最小化原則，SVR 將學習過程轉化為優化問題：

利用拉格朗日函數和對偶原理，可以得到式（6）的對偶問題：

式中：α、α＊為拉格朗日乘子；Q 中的元素Qij＝ΦT（xi）Φ（xj）；I為單位矩陣.

求解此二次型規劃可求得α 的值，同時求得

式中：K（x，xj）＝ΦT（x）Φ（xj），為一個滿足Mercer條件的核函數.該函數可在不知非線性變換的具體形式下實現算法的非線性化，這是支持向量機的一

根據上述推導可得回歸函數f（x）的表達式為個顯著特點.

選擇徑向基核函數進行代理模型的建立，其表達式如下：

憑經驗選取SVR 模型中的參數C、g 和ε 會引起較大的回歸誤差，因此采用PSO 算法優化選取這幾個參數.

2.4.2 粒子群算法

PSO 算法是一種群智能計算方法［12］.該算法是將優化問題的解當成是搜索空間中無質量無體積的粒子，通過迭代找到最優解.在每一次迭代過程中，粒子通過跟蹤個體最優位置和全局最優位置來更新自己.粒子在搜索空間的速度和位置由下式確定：

式中：t為當前進化代數；r1和r2為均勻分布在［0，1］內的隨機數；c1和c2為加速常數；β為慣性權重；vid為第i個粒子速度向量的第d 個分量；pid為第i個粒子最優位置向量中的第d 個分量；pgd為群體最優位置向量中第d 個分量；xid為第i 個粒子位置向量的第d 個分量；d＝1，2，…，n，i＝1，2，…，m，m為種群規模.

2.4.3 PSO－SVR代理模型的構建

建立PSO－SVR 代理模型的具體步驟如下：

（1）讀取樣本數據，將樣本數據分割為訓練集和測試集.

（2）PSO 算法的初始化.初始化PSO 算法的相關參數，初始化慣性權重、加速常數、最大進化代數Tmax以及適應度閾值Acc.將需要優化的SVR 模型中的相關參數初始化為一群粒子，初始化粒子速度和粒子位置.

（3）根據適應度函數計算每個粒子的適應值，將粒子中最小的適應值所對應的初始位置設為初始群體最好位置.每個粒子的適應度函數定義為

式中：e為平均絕對百分比誤差；yi為實際值；f（xi）為SVR 模型輸出值；n為所選擇的樣本數目.

（4）更新粒子的速度和位置，產生新一代種群.

（5）根據設定的適應度函數，評價每個粒子的適應值.粒子的適應值越小，粒子的位置越好.

（6）檢查終止條件，如果滿足，則將全局最優粒子映射為SVR模型中的參數，并以此為優化結果得到一個建立好的SVR模型；否則，轉步驟（4）.尋優條件為進化代數達到最大或者適應度小于設定閾值.

3 仿真分析

3.1 分析對象

待分析的轉子系統結構簡圖見圖2.進行有限元模化后，共有13個節點，2個質量圓盤分別位于節點4和節點8，支撐軸承位于節點2和節點12，2個支撐軸承的支撐特性相同，圓盤2處有一定量的初始不平衡.轉子系統的特性參數如表1所示.采用有限元模型計算出系統的固有頻率和振型，其中水平方向前兩階的振型如圖3所示.

圖2 某轉子系統結構簡圖Fig.2 Schematic diagram of the rotor system

表1 轉子系統相關參數Tab.1 Main parameters of the rotor system

3.2 構建代理模型

設不平衡偏心量的范圍為0.06～0.6g·m，相角范圍為0°～360°，擬采用40 個樣本作為PSOSVR 代理模型的訓練樣本.采用LHS算法選取出的40個樣本如圖4所示.

將不平衡力分別加在圓盤1和圓盤2處，采用有限元模型求解以上40組不平衡工況下的振動響應，假設節點1和節點13處的水平和垂直方向的振動響應可以通過測量得到，采用PSO－SVR 構建這4個測點振動響應的代理模型.其中，當不平衡力施加在圓盤2處時，同一時刻下由PSO－SVR 代理模型得到的4個測點振動響應與不平衡工況的對應關系見圖5.

3.3 不平衡參數辨識

圖3 水平方向前兩階振型Fig.3 The first two vibration modes in horizontal direction of the rotor system

圖4 不平衡參數樣本Fig.4 Parameter samples of mass unbalance

在圓盤2處分別設置偏心量為0.2g·m和0.4g·m、相角為50°和200°時的4種不平衡工況，進行不平衡參數的辨識，得出不平衡位置在圓盤2處時，不平衡偏心量和相角的辨識結果如表2所示.由表2可知，采用PSO－SVR 代理模型準確地辨識出了轉子不平衡的偏心量和相角，辨識結果與理論值吻合得非常好，辨識誤差很小.

圖5 振動響應與不平衡參數的對應關系Fig.5 Correspondence between vibration response and unbalance parameters

表2 不同偏心量和相角下的不平衡參數辨識結果Tab.2 Identified unbalance parameters for different unbalance eccentricities and phases

汽輪發電機組等復雜旋轉機械的轉子系統建模過程中，由于其結構復雜，實際運行情況多變，建模過程中不可避免地會產生一定誤差.在圓盤2上施加偏心量為0.4g·m、相角為50°產生的不平衡力，在轉子系統模型具有5%支撐剛度誤差、10%支撐剛度誤差、5%圓盤質量誤差和10%圓盤質量誤差的情況下，進行不平衡參數的定量辨識，結果見表3.由表3可知，隨著模型誤差的增大，不平衡定量辨識結果的誤差也會增加，但辨識值與理論值吻合情況仍然較好，證明所提出的不平衡定量辨識方法具有較好的魯棒性.這是因為該不平衡辨識策略除有限元建模外的各個環節所引入的辨識誤差均較小，各方法的綜合集成比較成功，在較小的有限元模型誤差及其他各個環節引入誤差的綜合作用下，最終的辨識誤差并不是很大，辨識誤差并未因為較小的有限元模型誤差而發生特別劇烈的變化.但值得注意的是，無論是基于模型的還是本文的不平衡辨識策略，盡量提高有限元模型的精度仍是不平衡參數辨識結果準確性的重要保障.

表3 模型具有誤差時的不平衡參數辨識結果Tab.3 Identified unbalance parameters with model errors

4 結 論

（1）提出了一種基于PSO－SVR 代理模型的轉子系統不平衡定量辨識方法.算例中的不平衡偏心量和相角的最大辨識誤差在5%左右，驗證了該方法的準確性.

（2）當轉子系統模型具有10%支撐剛度誤差、10%圓盤質量誤差的情況下，不平衡偏心量和相角辨識誤差的最大增幅未到3%，表明該方法的魯棒性較好.

應該指出，本文基于較為簡化的轉子系統對所提出的方法進行了檢驗，尚缺少工程驗證.如何將所提方法應用于復雜旋轉機械實際運行中的不平衡定量辨識，是需要繼續開展的工作；另外，本文的方法只能進行單點不平衡的辨識，實際過程中也有可能出現兩點或者多點不平衡，多點不平衡的定量識別是一個難題，有待進一步深入研究.

［1］ 胡三高，安宏文，馬志勇，等.基于小波奇異值分析的汽輪機碰磨特征提取［J］.動力工程學報，2013，33（3）：184－188.HU Sangao，AN Hongwen，MA Zhiyong，et al.Feature extraction of rubbing fault for steam turbines based on wavelet singularity analysis［J］.Journal of Chinese Society of Power Engineering，2013，33（3）：184－188.

［2］ MARKERT R，PLATZ R，SEIDLER M.Model based fault identification in rotor systems by least squares fitting［J］.International Journal of Rotating Machinery，2001，7（5）：311－321.

［3］ BACHSCHMID N，PENNACCHI P，VANIA A.Identification of multiple faults in rotor systems［J］.Journal of Sound and Vibration，2002，254（2）：327－366.

［4］ LI Changyou，XU Mingqiang，GUO Song，et al.Model－based degree estimation of unbalance and misalignment in flexible coupling－rotor system［J］.Chinese Journal of Mechanical Engineering，2009，22（4）：550－556.

［5］ JALAN K A，MOHANTY A R.Model based fault diagnosis of a rotor－bearing system for misalignment and unbalance under steady－state condition［J］.Journal of Sound and Vibration，2009，327（3／4／5）：604－622.

［6］ SUDHAKAR G N D S，SEKHAR A S.Identification of unbalance in a rotor bearing system［J］.Journal of Sound and Vibration，2011，330（10）：2299－2313.

［7］ HUANG Zhangjun，WANG Chengen，CHEN Jian，et al.Optimal design of aeroengine turbine disc based on Kriging surrogate models［J］.Computers ＆Structures，2011，89（1／2）：27－37.

［8］ MCKAY M D，BECKMAN R J，CONOVER W J.A comparison of three methods for selecting values of input variables in the analysis of output from a computer code［J］.Technometrics，1979，21（2）：239－245.

［9］ 王勖成.有限單元法［M］.北京：清華大學出版社，2013.

［10］ VAPNIK V N.Statistical learning theory［M］.New York，USA：John Wiley，1998：2－5.

［11］ 牛培峰，麻洪波，李國強，等.基于支持向量機和果蠅優化算法的循環流化床鍋爐NOx排放特性研究［J］.動力工程學報，2013，33（4）：267－271.NIU Peifeng，MA Hongbo，LI Guoqiang，et al.Study on NOxemission from CFB boilers based on support vector machine and fruit fly optimization algorithm［J］.Journal of Chinese Society of Power Engineering，2013，33（4）：267－271.

［12］ 王曦，李興源，趙睿.基于相對增益和改進粒子群算法的PSS與直流調制協調策略［J］.中國電機工程學報，2014，34（34）：6177－6184.WANG Xi，LI Xingyuan，ZHAO Rui.Coordination strategy of PSS and DCM based on relative gain and improved PSO［J］.Proceedings of the CSEE，2014，34（34）：6177－6184.