運(yùn)行工況下轉(zhuǎn)子系統(tǒng)模態(tài)參數(shù)識(shí)別技術(shù)研究

臧廷朋， 王鳳仁， 溫廣瑞,2， 張志芬(. 西安交通大學(xué) 機(jī)械工程學(xué)院 智能儀器與監(jiān)測診斷研究所， 西安 70049; 2. 新疆大學(xué) 機(jī)械工程學(xué)院， 烏魯木齊 830047)

轉(zhuǎn)子系統(tǒng)的模態(tài)參數(shù)綜合反映了旋轉(zhuǎn)機(jī)械的動(dòng)力特性，準(zhǔn)確的識(shí)別模態(tài)參數(shù)對旋轉(zhuǎn)機(jī)械的狀態(tài)監(jiān)測和故障診斷有著重要意義[1]。隨著旋轉(zhuǎn)機(jī)械不斷朝著精密化、大型化、高速化發(fā)展，現(xiàn)代旋轉(zhuǎn)機(jī)械工作轉(zhuǎn)速可能接近第二階臨界轉(zhuǎn)速，因此必須精確的識(shí)別系統(tǒng)模態(tài)參數(shù)，從而保障機(jī)組高效、平穩(wěn)持續(xù)運(yùn)行[2]。

工程實(shí)際中常用傳遞矩陣法和有限元法識(shí)別系統(tǒng)的模態(tài)參數(shù)[3-5]，考慮到機(jī)組結(jié)構(gòu)的復(fù)雜性和機(jī)組相互之間關(guān)聯(lián)影響，建立準(zhǔn)確的系統(tǒng)模型比較困難。同時(shí)上述兩種方法沒有考慮實(shí)際運(yùn)行工況下系統(tǒng)剛度、阻尼變化對模態(tài)參數(shù)的影響，模態(tài)參數(shù)識(shí)別結(jié)果與真實(shí)值往往有較大差異[6]。另一方面，為了獲得系統(tǒng)的瞬態(tài)響應(yīng)需要施加人工激勵(lì)，對大型高速工業(yè)用轉(zhuǎn)子系統(tǒng)施加人工激勵(lì)費(fèi)時(shí)費(fèi)力，降低了模態(tài)參數(shù)識(shí)別的精度和效率。相對于廣泛應(yīng)用的頻域方法，模態(tài)參數(shù)時(shí)域識(shí)別方法直接使用不需要進(jìn)行時(shí)頻變換的時(shí)域信號(hào)進(jìn)行識(shí)別，避免了時(shí)頻變換導(dǎo)致的能量泄漏，同時(shí)不影響機(jī)組的連續(xù)正常運(yùn)行，因而模態(tài)參數(shù)的時(shí)域識(shí)別方法受到越來越多的關(guān)注[7-9]。鄭鋼鐵等[10]提出了一種基于時(shí)域振動(dòng)參數(shù)識(shí)別理論的模態(tài)參數(shù)識(shí)別方法，該方法通過獲取振動(dòng)響應(yīng)并結(jié)合系統(tǒng)先驗(yàn)知識(shí)，實(shí)現(xiàn)了模態(tài)參數(shù)的現(xiàn)場識(shí)別；宋志強(qiáng)[11]提出了基于遺傳算法的模態(tài)參數(shù)時(shí)域識(shí)別方法，在使用多信號(hào)分類確定系統(tǒng)階次的基礎(chǔ)上，利用遺傳算法尋優(yōu)獲得系統(tǒng)的自振頻率和阻尼比。目前基于運(yùn)行工況的系統(tǒng)模態(tài)參數(shù)識(shí)別尚缺少有效實(shí)用的方法。

本文針對運(yùn)行工況下系統(tǒng)模態(tài)參數(shù)識(shí)別問題，提出一種基于隨機(jī)子空間[12-15]和經(jīng)驗(yàn)?zāi)Ｊ椒纸庀嘟Y(jié)合的系統(tǒng)模態(tài)參數(shù)識(shí)別方法。該方法首先對轉(zhuǎn)子系統(tǒng)振動(dòng)響應(yīng)進(jìn)行窄帶濾波，獲取轉(zhuǎn)子振動(dòng)響應(yīng)中各階自由衰減分量，然后借助于EMD降低模態(tài)混疊和其他信號(hào)成分的影響，最后利用SSI識(shí)別各階模態(tài)的固有頻率。該方法直接處理系統(tǒng)時(shí)域響應(yīng)數(shù)據(jù)，完成保留了數(shù)據(jù)信息，同時(shí)避免了各階模態(tài)之間的相互影響，模態(tài)參數(shù)識(shí)別結(jié)果更具有可靠性和實(shí)際應(yīng)用價(jià)值。

1 經(jīng)驗(yàn)?zāi)Ｊ椒纸?/h2>

經(jīng)驗(yàn)?zāi)Ｊ椒纸鈁16-18]是一種基于數(shù)據(jù)本身特征把復(fù)雜的信號(hào)分解為有限個(gè)固有模式分量(IMF)的方法，可以自適應(yīng)處理非線性、非平穩(wěn)信號(hào)。其本質(zhì)是利用信號(hào)特征時(shí)間尺度來獲得信號(hào)的固有模式分量。EMD分解建立在三個(gè)假設(shè)的基礎(chǔ)上：信號(hào)至少有兩個(gè)極值點(diǎn)(一個(gè)極大值和一個(gè)極小值)；特征時(shí)間尺度由極點(diǎn)之間的時(shí)間間隔定義；如果缺少極值點(diǎn)但有變形點(diǎn)，可通過微分獲得極值點(diǎn)。

EMD分解得到一個(gè)固有模式函數(shù)主要思想為，通過提取信號(hào)x(t)的極大值點(diǎn)和極小值點(diǎn)，利用三次樣條得到上包絡(luò)和下包絡(luò)，獲取包絡(luò)均值為m(t)，用原始信號(hào)減去包絡(luò)均值得到h(t)；檢驗(yàn)h(t)是否為一個(gè)IMF分量，如果不是，令原信號(hào)是h(t)，重復(fù)上述過程直到找到第一個(gè)IMF分量c1(t)；將第一個(gè)IMF分量從x(t)分離出來作為新的數(shù)據(jù)，重復(fù)上面過程，可以得到一系列固有模式分量和一個(gè)殘余項(xiàng)。原始信號(hào)可以表示為

(1)

2 隨機(jī)子空間理論

2.1 隨機(jī)子空間模型

隨機(jī)子空間方法是一種適用于平穩(wěn)隨機(jī)激勵(lì)條件下模態(tài)參數(shù)識(shí)別的時(shí)域方法。其模型是由結(jié)構(gòu)動(dòng)力學(xué)MCK方程推導(dǎo)出來的，并假設(shè)噪聲信號(hào)是零均值的隨機(jī)平穩(wěn)噪聲。模型表達(dá)式為

(2)

式中:A為狀態(tài)矩陣;C為輸出矩陣;wk為零均值過程噪聲;vk為零均值測量噪聲;xk為狀態(tài)向量;yk為輸出的測量值。

2.2 隨機(jī)子空間原理

SSI方法的思想是把將來輸入的行空間投影到過去輸入的行空間上，然后利用投影的結(jié)果預(yù)測未來。根據(jù)傳感器輸出，將輸出信號(hào)構(gòu)建出Hankel矩陣

(3)

式(3)是由N個(gè)采樣數(shù)據(jù)構(gòu)成的，N=2i+j-1。其中p表示“過去”的數(shù)據(jù)，f表示“未來”的數(shù)據(jù)，要求j>>i。

對Hankele矩陣QR分解

(4)

由空間投影性質(zhì)可以進(jìn)行正交投影

(5)

式中：?表示偽逆。經(jīng)過QR分解可以使得數(shù)據(jù)量極大的縮減，有利于程序運(yùn)行。

對投影矩陣SVD分解

(6)

式中：U1∈Rli×n；S1∈Rn×n；V1∈Rj×n，系統(tǒng)的階次由S1的秩來確定。

對投影矩陣進(jìn)行卡爾曼狀態(tài)濾波

(7)

由式(6)和式(7)，可令

(8)

(9)

正交投影矩陣還可以定義為

(10)

(11)

類似上面的推理，可以得到

(12)

將卡爾曼狀態(tài)濾波與輸出代入系統(tǒng)狀態(tài)方程

(13)

式中：Yi,i表示只有一個(gè)行塊的輸出矩陣。

由于干擾與狀態(tài)不相關(guān)，可以采用最小二乘法來計(jì)算矩陣的漸近無偏估計(jì)

(14)

在確定系統(tǒng)的狀態(tài)矩陣和輸出矩陣之后，對系統(tǒng)的狀態(tài)矩陣A進(jìn)行特征值分解，求其特征值λi。離散時(shí)間系統(tǒng)與連續(xù)時(shí)間系統(tǒng)的特征值關(guān)系為

(15)

則信號(hào)的頻率為

(16)

3 基于SSI和EMD的系統(tǒng)模態(tài)參數(shù)識(shí)別

運(yùn)行工況下的轉(zhuǎn)子系統(tǒng)，其振動(dòng)激勵(lì)包括不平衡等自身激勵(lì)以及外部隨機(jī)環(huán)境激勵(lì)，前者會(huì)產(chǎn)生明顯的各階振蕩諧波，后者持續(xù)存在并激發(fā)產(chǎn)生振動(dòng)響應(yīng)中的自由衰減分量?；赟SI和EMD的系統(tǒng)模態(tài)參數(shù)識(shí)別方法，側(cè)重于通過隨機(jī)環(huán)境激勵(lì)產(chǎn)生的微弱自由衰減分量中實(shí)時(shí)的識(shí)別系統(tǒng)模態(tài)參數(shù)。首先根據(jù)先驗(yàn)知識(shí)和轉(zhuǎn)子系統(tǒng)振動(dòng)響應(yīng)設(shè)置濾波器參數(shù)，采用窄帶濾波方法提取轉(zhuǎn)子振動(dòng)響應(yīng)中的各階自由衰減分量，考慮到響應(yīng)信號(hào)中各階模態(tài)存在混疊的情況，利用EMD自適應(yīng)的進(jìn)一步提取更加精確的各階自由衰減分量信號(hào)，然后利用SSI求解系統(tǒng)狀態(tài)空間方程的狀態(tài)矩陣和輸出矩陣得到系統(tǒng)的模態(tài)參數(shù)。具體識(shí)別過程，如圖1所示。

圖1 模態(tài)參數(shù)識(shí)別流程圖

4 仿真分析

針對運(yùn)轉(zhuǎn)工況下轉(zhuǎn)子系統(tǒng)振動(dòng)響應(yīng)的特點(diǎn)，將各階振蕩諧波表示為幅值與頻率固定的余弦信號(hào)。將各階次自由衰減分量表示為頻率固定、幅值按指數(shù)形式變化的余弦信號(hào)?？紤]到轉(zhuǎn)子運(yùn)行時(shí)需要遠(yuǎn)離各階臨界轉(zhuǎn)速，因而各階自由衰減分量的幅值應(yīng)遠(yuǎn)小于諧波的振動(dòng)幅值。構(gòu)建如下的仿真信號(hào)

x(t)=30cos(2×π×50t+π/6)+

10cos(2×π×5t+π/8)+

0.1exp(-0.15t)cos(2×π×25t+π/4)+

0.05exp(-0.14t)cos(2×π×75t)

(17)

式中：第一項(xiàng)模擬一倍頻分量，一倍頻分量頻率為50 Hz；第二項(xiàng)模擬分倍頻分量，分倍頻分量頻率為5 Hz；第三項(xiàng)模擬系統(tǒng)一階模態(tài)自由衰減分量，頻率為25 Hz；第四項(xiàng)模擬系統(tǒng)二階模態(tài)自由衰減分量，頻率為75 Hz。添加方差為1的白噪聲，采樣頻率為1 kHz，采樣時(shí)間為2 s。仿真信號(hào)的時(shí)域波形如圖2所示。對仿真信號(hào)進(jìn)行窄帶濾波，并用EMD提取前兩階模態(tài)自由衰減分量，分別如圖3和圖4所示。

圖2 仿真信號(hào)時(shí)域波形

圖3 一階模態(tài)自由衰減分量

圖4 二階模態(tài)自由衰減分量

使用本文方法對提取的各階模態(tài)自由衰減分量進(jìn)行識(shí)別，結(jié)果如表1所示。通過對比，可以看出本方法可以很好的從模擬信號(hào)中提取出一階模態(tài)與二階模態(tài)固有頻率。

表1 模態(tài)參數(shù)識(shí)別結(jié)果

5 實(shí)驗(yàn)研究

在Bently RK4轉(zhuǎn)子實(shí)驗(yàn)臺(tái)上對本文提出的方法進(jìn)行驗(yàn)證，實(shí)驗(yàn)臺(tái)布置如圖5所示。設(shè)置采樣頻率為2 048 Hz，采樣時(shí)長為1 s，用DT9837B數(shù)采卡對轉(zhuǎn)子運(yùn)行工況下穩(wěn)定轉(zhuǎn)速過程進(jìn)行連續(xù)數(shù)據(jù)采集。有4個(gè)通道的數(shù)據(jù)，可以任意選取一組數(shù)據(jù)作為監(jiān)測數(shù)據(jù)，本文選用1_X(第一通道)通道的數(shù)據(jù)進(jìn)行驗(yàn)證分析。

轉(zhuǎn)子實(shí)驗(yàn)臺(tái)的工作轉(zhuǎn)速為4 000 r/min，介于一階與二階臨界轉(zhuǎn)速之間。在求解系統(tǒng)各階模態(tài)參數(shù)時(shí)，首先通過濾波器濾波，提取各階自由衰減分量。提取一階衰減分量濾波器中心頻率為30 Hz，帶寬是60 Hz；提取二階衰減分量濾波器中心頻率為110 Hz，帶寬是50 Hz。由于模態(tài)混疊及其他信號(hào)成分對模態(tài)參數(shù)識(shí)別精度有較大影響，對濾波后的信號(hào)進(jìn)行EMD分解，選取第一個(gè)IMF分量進(jìn)行隨機(jī)子空間模態(tài)識(shí)別。濾波和EMD提取的一、二階自由衰減分量分別如圖6～圖9所示。

圖5 實(shí)驗(yàn)臺(tái)結(jié)構(gòu)示意圖

圖6 濾波后的一階自由衰減分量

圖7 EMD提取的一階自由衰減分量

圖8 濾波后的二階自由衰減分量

圖9 EMD提取的二階自由衰減分量

對EMD提取后的各階模態(tài)衰減分量信號(hào)用SSI進(jìn)行模態(tài)參數(shù)識(shí)別，識(shí)別結(jié)果如表2所示。為了驗(yàn)證識(shí)別效果，通過一次起車過程獲取轉(zhuǎn)子系統(tǒng)的前兩階臨界轉(zhuǎn)速，起車過程bode圖，如圖10所示。其中1點(diǎn)處轉(zhuǎn)速為1 966 r/min，對應(yīng)的頻率是32.77 Hz，2點(diǎn)處轉(zhuǎn)速為6 496 r/min，對應(yīng)頻率是108.27 Hz。通過對比，可以看出本方法求解的一階、二階固有頻率與系統(tǒng)實(shí)際的特性相吻合；同時(shí)，考慮到起車過程數(shù)據(jù)采集具有滯后性，識(shí)別的各階模態(tài)頻率小于起車過程臨界轉(zhuǎn)速對應(yīng)的頻率。

表2 系統(tǒng)模態(tài)參數(shù)識(shí)別結(jié)果

為了分析起車升速過程中模態(tài)參數(shù)的變化規(guī)律，選取3 900 r/min至4 900 r/min之間11個(gè)不同轉(zhuǎn)速，對各個(gè)轉(zhuǎn)速時(shí)轉(zhuǎn)子系統(tǒng)進(jìn)行模態(tài)參數(shù)識(shí)別，識(shí)別結(jié)果如表3和圖11所示。圖中，兩條曲線分別是一階模態(tài)頻率隨轉(zhuǎn)速變化的情況和二階模態(tài)頻率隨轉(zhuǎn)速變化的情況。從表3和圖11中可以看出，系統(tǒng)模態(tài)頻率隨轉(zhuǎn)速變化會(huì)有輕微的波動(dòng)，整體變化是隨著轉(zhuǎn)速的上升，一二階模態(tài)頻率有降低的趨勢。

圖10 Bode圖

轉(zhuǎn)速/(r·min-1)各階模態(tài)頻率/Hz一階二階390032.5512108.4917400031.1822107.7785410030.8316106.8193420031.3412107.5438430030.1422107.1453440032.5127106.9571450031.0994106.2167460031.7180107.1076470031.8012106.4699480030.8582105.5993490030.4812105.6045

圖11 模態(tài)參數(shù)隨轉(zhuǎn)速變化趨勢

6 結(jié) 論

針對運(yùn)行工況下系統(tǒng)模態(tài)參數(shù)識(shí)別問題，本文提出一種基于隨機(jī)子空間和經(jīng)驗(yàn)?zāi)Ｊ椒纸庀嘟Y(jié)合的系統(tǒng)模態(tài)參數(shù)識(shí)別的方法。仿真數(shù)據(jù)和實(shí)驗(yàn)研究表明，該方法可以準(zhǔn)確識(shí)別系統(tǒng)一、二階模態(tài)參數(shù)。得出以下結(jié)論：

(1) EMD可以有效減少各階模態(tài)之間的混疊以及其他成分對模態(tài)參數(shù)識(shí)別的干擾，提高模態(tài)參數(shù)識(shí)別的精度。

(2) 該方法實(shí)施起來無需復(fù)雜設(shè)備也無需繁雜的計(jì)算，具有良好的抗噪性、通用性。

(3) 該方法直接處理系統(tǒng)運(yùn)行工況下的時(shí)域響應(yīng)數(shù)據(jù)，完成保留了數(shù)據(jù)信息，同時(shí)不影響機(jī)組的正常運(yùn)行，模態(tài)參數(shù)識(shí)別結(jié)果更具有可靠性和實(shí)際應(yīng)用價(jià)值。

[1] 段吉安,虞烈,謝友柏,等. 軸承-轉(zhuǎn)子系統(tǒng)模態(tài)參數(shù)時(shí)域識(shí)別的研究[J]. 機(jī)械科學(xué)與技術(shù),1997, 16(2): 11-15.

DUAN Ji’an, YU Lie, XIE Youbai, et al.The study on identification of model parameters of rotor-bearing systems[J]. Mechanical Science and Technology for Aerospace Engineering, 1997, 16(2): 11-15.

[2] 林勝,林春庭. 高速數(shù)控機(jī)床現(xiàn)狀和發(fā)展趨勢[J]. 精密制造與自動(dòng)化,2004(1):5-8.

LIN Sheng, LIN Chunting. Status and developing trend of high speed CNC machine tools[J]. Precise Manufaction & Automation,2004(1):5-8.

[3] 曹樹謙,陳予恕,丁千,等. 高速轉(zhuǎn)子動(dòng)平衡的傳遞函數(shù)法[J]. 機(jī)械強(qiáng)度,2002, 24(4):500-504.

CAO Shuqian, CHEN Yushu, DING Qian, et al. Transfer function technique of dynamic balancing for high-speed rotors[J]. Journal of Mechanical Strength, 2002, 24(4):500-504.

[4] 韓慶華, 張鵬, 蘆燕. 基于傳遞函數(shù)法的大跨建筑非結(jié)構(gòu)件動(dòng)力性能研究[J]. 土木工程學(xué)報(bào),2014, 47(增刊2):79-84.

HAN Qinghua, ZHANG Peng, LU Yan. Dynamic performance of non-structural components in large span buildings based on transfer function[J].China Civil Engineering Journal,2014, 47(Sup2):79-84.

[5] YANG J N, LEI Y. Identification of natural frequencies and dampings of in situ tall buildings using ambient wind vibration data[J]. Journal of Structural Engineering, 2004, 130(5):570-577.

[6] 孫偉,汪博,聞邦椿. 高速主軸系統(tǒng)靜止及運(yùn)轉(zhuǎn)狀態(tài)下動(dòng)力學(xué)特性對比分析[J]. 機(jī)械工程學(xué)報(bào),2012, 48(11):146-152.

SUN Wei, WANG Bo, WEN Bangchun. Comparative analysis of dynamics characteristics for static and operation state of high-speed spindle system[J]. Journal of Mechanical Engineering, 2012, 48(11):146-152.

[7] 任偉新.環(huán)境振動(dòng)系統(tǒng)識(shí)別方法的比較分析[J].福州大學(xué)學(xué)報(bào)(自然科學(xué)版).2001, 29(6): 80-86.

REN Weixin. Comparison of system identification methods using ambient vibration measurements[J]. Journal of Fuzhou University(Natural Science), 2011, 29(6): 80-86.

[8] BRINCKER R,ZHANG L M, ANDERSEN P. Output-only modal analysis frequency domain decomposition[C]∥Proc of 25th Int Seminar on Modal Analysis (ISMA25) Vol 2. Leuven, Belgium, 2000.

[9] 王彤,張令彌.運(yùn)行模態(tài)分析的頻域空間域分解法及其應(yīng)用[J].航空學(xué)報(bào).2006, 27(1): 62-66.

WANG Tong, ZHANG Lingmi. Frequency and spatial domain decomposition for operational modal analysis and its application[J]. Acta Aeronautica et Astronautica Sinica, 2006, 27(1): 62-66.

[10] 鄭鋼鐵,黃文虎. 轉(zhuǎn)子系統(tǒng)模態(tài)參數(shù)現(xiàn)場識(shí)別技術(shù)研究[J]. 力學(xué)學(xué)報(bào),1991, 23(5):634-640.

ZHENG Gangtie, HUANG Wenhu. An investigation into the field identification technique of modal parameters of rotating system[J]. Acta Mechanica Sinica,1991, 23(5):634-640.

[11] 宋志強(qiáng). 水電機(jī)組軸系統(tǒng)模態(tài)參數(shù)的時(shí)域識(shí)別研究[J]. 西北農(nóng)林科技大學(xué)學(xué)報(bào)(自然科學(xué)版),2011(11):229-234.

SONG Zhiqiang. Study on modal parameters time domain identification of generator set shaft system[J].Journal of Northwest A & F University (Natural Science Edition), 2011(11):229-234.

[12] 陳永高, 鐘振宇. 基于CEEMD分解和Data-SSI算法的斜拉橋模態(tài)參數(shù)識(shí)別[J]. 振動(dòng)與沖擊, 2016, 35(8):166-172.

CHEN Yonggao, ZHONG Zhenyu. Modal parameter identification of a cable-stayed bridge based on CEEMD and DATA-SSI algorithm[J]. Journal of Vibration and Shock, 2016, 35(8):166-172.

[13] 胡異丁, 李丹, 任偉新,等. 基于延時(shí)隨機(jī)子空間方法的非白噪聲環(huán)境激勵(lì)結(jié)構(gòu)模態(tài)參數(shù)識(shí)別[J]. 振動(dòng)與沖擊, 2015,34(8):71-76.

HU Yiding, LI Dan, REN Weixin, et al. Modal parameter identification of structures under non-white noise ambient excitations using delay-index-based stochastic subspace method[J]. Journal of Vibration and Shock, 2015, 34(8):71-76.

[14] REN W X, ZHAO T, HARIK I E. Experimental and analytical modal analysis of a steel arch bridge[J]. Journal of Structural Engineering, 2004, 130(7):1022-1031.

[15] 趙書強(qiáng), 張一,馬燕峰. 基于數(shù)據(jù)驅(qū)動(dòng)隨機(jī)子空間方法在低頻振蕩辨識(shí)中的應(yīng)用[J]. 電力系統(tǒng)保護(hù)與控制,2013(8):81-86.

ZHAO Shuqiang, ZHANG Yi, MA Yanfeng. Data-driven based stochastic subspace identification for low frequency oscillation analysis[J]. Power System Protection and Control, 2013(8):81-86.

[16] RILLING G, FLANDRIN P, GONCALVES P. On empirical mode decomposition and its algorithms[C]∥IEEE-EURASIP Workshop on Nonlinear Signal and Image Processing NSIP-03. Grado, Italy, 2003.

[17] GAI G H. The processing of rotor startup signals based on empirical mode decomposition[J]. Mechanical Systems and Signal Processing, 2006, 20(1): 222-235.

[18] 李天云,袁明哲,李軍強(qiáng),等. 基于EMD和SSI的電力系統(tǒng)低頻振蕩模態(tài)參數(shù)識(shí)別方法[J]. 電力系統(tǒng)保護(hù)與控制,2011, 39(8):6-10.

LI Tianyun, YUAN Mingzhe, LI Junqiang, et al. Method of modal parameter identification of power system low frequency oscillation based on EMD and SSI[J]. Power System Protection and Control, 2011, 39(8): 6-10.