汪曉姍， 彭志科， 陳是扦

(上海交通大學(xué) 機(jī)械系統(tǒng)與振動國家重點實驗室，上海 200240)

顫振是金屬加工中出現(xiàn)的自激振動，嚴(yán)重制約了加工精度和加工效率。隨著機(jī)床加工柔性化的日益發(fā)展，要求加工能對不同的工件在不同的條件下進(jìn)行，因此不能從根本上杜絕顫振。近年來，于英華等[1-2]很多國內(nèi)外學(xué)者開始著力于研究顫振的在線監(jiān)測和預(yù)報控制。

銑削顫振的形成是一個孕育的過程，在顫振早期，監(jiān)測系統(tǒng)搶先把銑削過程即將發(fā)生顫振的預(yù)兆先通報給控制系統(tǒng)，是實現(xiàn)顫振預(yù)報的關(guān)鍵。各種各樣的傳感信號可以被用來采集顫振信號，比如動態(tài)銑削力信號[3-4]，加速度信號[5]，力矩信號[6]等。Kuljanic等[7]提出了利用多傳感器來監(jiān)測識別顫振，這種方法提高了顫振監(jiān)測的精確性和魯棒性；Tansel等[8]通過銑削位移信號的諧波分量來監(jiān)測顫振；Tarng等[9]在發(fā)生顫振時切削力信號在的頻譜上的窄帶特性來預(yù)測顫振。

工程中，直接根據(jù)傳感信號難以識別顫振，因此需要從采集來的信號中獲得對故障敏感的特征量。時域中主要采用的是概率統(tǒng)計分析，如互相關(guān)系數(shù)[10]、標(biāo)準(zhǔn)差[11]等；頻域中，遲玉倫等[12]基于通過識別銑削力的譜峰頻率來監(jiān)測顫振。然而，這些傳統(tǒng)的方法只適合分析平穩(wěn)信號，而且無法準(zhǔn)確預(yù)報顫振的發(fā)生。對于成分復(fù)雜的顫振信號，需要將其分解成單分量再進(jìn)行分析。Cao等[13]應(yīng)用EMD分解銑削力信號，提取顫振的特征，但是無法克服EMD分解的模態(tài)混疊。該方法易產(chǎn)生虛假的信號分量，干擾故障特征的診斷；Zhang等[14]提出用VMD分解銑削力信號，但是VMD本質(zhì)上通過自適應(yīng)的頻域濾波器組來提取信號，其分解銑削力信號得到的子分量沒有物理意義。Vold-Kalman濾波方法[15]基于數(shù)據(jù)方程誤差和最小化結(jié)構(gòu)，有效地進(jìn)行多分量信號分解。不僅可以避免EMD分解帶來的模式混淆，而且可以有效的從復(fù)雜的顫振信號中提取單分量。Vold-Kalman濾波提取單信號首先需要根據(jù)時頻脊線估計單分量信號的瞬時頻率。經(jīng)典的瞬時頻率估計方法中，過零法和相位差分法較為簡單，但是誤差大。HHT繼承了EMD的缺點。在參數(shù)化時頻分析的基礎(chǔ)上，Yang等[16]提出了基于優(yōu)化頻譜集中性指標(biāo)，對信號瞬時頻率參數(shù)估計的方法，該方法具有抗噪性強(qiáng)、精確度高等特點。

顫振產(chǎn)生時，信號的頻率成分和能量分布發(fā)生改變。穩(wěn)定銑削時，顫振的能量分布在噪聲、主軸的轉(zhuǎn)速頻率及其諧波處。發(fā)生顫振時，振動的能量逐漸轉(zhuǎn)移到顫振頻率處。因此本文提出基于頻域能量分布改變的顫振識別方法。針對復(fù)雜的銑削顫振信號，考慮Vold-Kalman濾波方法在單分量提取和優(yōu)化頻譜集中性指標(biāo)在瞬時頻率計算方面的獨特優(yōu)勢，本文提出了基于兩者的多分量信號分解方法。采用分解后的子信號在頻域上的能量分布作為顫振識別的特征。多次銑削實驗驗證了本文所提出的方法的有效可行。

1 多分量信號分解

1.1 基于優(yōu)化頻譜集中性指標(biāo)估計信號瞬時頻率

該過程可以描述為，選取一個目標(biāo)分量，采用合適的變換核參數(shù)構(gòu)造對應(yīng)的旋轉(zhuǎn)算子，將原信號旋轉(zhuǎn)變換到它的一個旋轉(zhuǎn)域。為了在上述旋轉(zhuǎn)域評估分量能量的集中性，提出頻譜集中性指標(biāo)。旋轉(zhuǎn)域中能量越集中，集中性指標(biāo)越大，可認(rèn)為變換核參數(shù)與分量模型的相位參數(shù)越相近。

單分量解析信號可用多項式相位信號形式表示[17]

(1)

(2)

sd(t;C)=s(t)Φ(t;C)

(3)

sd(t;C)=a(t)exp(j(2πc0t+φ0))

(4)

此時能量集中于c0處， 為衡量頻譜集中性，定義頻譜集中性指標(biāo)

(5)

式中：E(·)表示期望算子； F (·)表示傅里葉變換，以SCI最大為參數(shù)估計的指標(biāo)。

(6)

本文，估計參數(shù)采用粒子群優(yōu)化算法，通過找到解調(diào)信號的最大譜峰值對應(yīng)的頻率來估計初始頻率參數(shù)

(7)

通過式(6)、式(7)估計瞬時頻率的所有參數(shù)(c0,…,ck)。

1.2 Vold-Kalman濾波

Vold-Kalman濾波方法可以從復(fù)雜多分量信號中提取有效的頻率成分。與傳統(tǒng)的濾波方法相比，該方法避免了由時域至頻域變換帶來的相位偏差，可以將復(fù)雜多分量信號分解為單分量。Vold-Kalman濾波的中心頻率能根據(jù)瞬時頻率進(jìn)行自適應(yīng)調(diào)節(jié)，因此能夠有效分離在時頻域內(nèi)鄰近甚至交叉的信號分量[18]。

調(diào)制信號可以表示為

(8)

式中：Ak(t)為第k個分量的幅值包絡(luò)；Θk(t)為載波信號；k為階次。

(9)

振幅包絡(luò)線是一個低頻率調(diào)制的載波信號，可以表示為一個低階多項式。對于離散信號

(10)

Ak(n-1)-2Ak(n)+Ak(n+1)=εk(n)

(11)

對于實際信號，Ak(0)=0。 當(dāng)n=1時，Ak(0)-2Ak(1)+Ak(2)=εk(1)； 當(dāng)n=N時， 其中是信號的長度。 式(11)的矩陣的形式

(12)

矩陣形式為

MA=ε

(13)

實際測試信號

(14)

為噪聲或誤差項，式(14)可以寫成矩陣形式

X-BA=δ

(15)

優(yōu)化頻譜集中性指標(biāo)， 估計信號的瞬時頻率ωk(n)， 從而得到B， 根據(jù)(13)和(16)得到各信號分量的幅值包絡(luò)矩陣A， 由矩陣C和A可以重構(gòu)各信號分量

S=AC

(16)

信號的非線性調(diào)頻分量分解的步驟如下：

(1) 優(yōu)化頻譜集中性指標(biāo)，估計分量的瞬時頻率；

(2) 根據(jù)估計到的瞬時頻率應(yīng)用Vold-Kalman濾波器濾出分量ui；

(3) 從信號中減去復(fù)原的信號分量作為初始信號，重復(fù)以上步驟，得到信號的各個分量。

2 特征提取

為了顫振預(yù)報及時并且精確，正確提取特征是關(guān)鍵。所選擇的特征參數(shù)要充分反映銑削顫振過程的本質(zhì)，同時還必須考慮信號采集和數(shù)據(jù)處理的簡易可行，使監(jiān)測系統(tǒng)可以在幾秒時間內(nèi)完成計算工作，向控制系統(tǒng)發(fā)出預(yù)報信號。

頻帶熵、信息熵、近似熵等的引用使信號復(fù)雜性表征成為可能，但是它們都存在缺點：在單一尺度上對時間序列進(jìn)行分析、需要比較自身數(shù)據(jù)段且參數(shù)改變時結(jié)果一致性差。因此，本文采用多分量信號分解和熵結(jié)合的方法。在穩(wěn)定銑削時，除了主軸轉(zhuǎn)速頻率及其諧波處，其他頻率處具有能量小、分布相對平均、不確定的特點。當(dāng)發(fā)生顫振時，除了主軸轉(zhuǎn)速頻率及其諧波處，在相應(yīng)的頻帶內(nèi)就會出現(xiàn)相應(yīng)共振頻率，能量便會集中在此頻率帶內(nèi)，使能量分布的不確定性減少。因此可以通過計算分解得到的各子信號的能量分布以檢測顫振。信號經(jīng)過分解可以得到n個單分量u1,u2,…,un， 計算出分別對應(yīng)的能量

(17)

假設(shè)殘余分量可以忽略，由于本文的分解方法具有正交性，n個子信號的能量之應(yīng)該等于原始信號的能量。由于各個子信號包含不同的頻率成分和能量，從而形成在銑削力信號在頻域的能量分布。由此引入能量熵銑削過程中的振動信號包括：由于銑刀旋轉(zhuǎn)引起的周期信號部分、由于再生效應(yīng)引起的顫振和由于系統(tǒng)噪聲等引起的隨機(jī)擾動部分。從實驗信號中提取顫振分量對分析顫振有重要意義。本文提出的顫振監(jiān)測方法如圖1(b)所示。首先，動態(tài)銑削力信號經(jīng)過Vold-Kalman濾波器去除由于主軸旋轉(zhuǎn)產(chǎn)生的信號分量。然后優(yōu)化頻譜集中性指標(biāo)對信號進(jìn)行瞬時頻率估計。Vold-Kalman濾波將復(fù)雜的多分量信號分解為單分量，從時域中提取對分析顫振有意義的信號分量。最后計算分量信號的能量熵作為識別顫振的特征。

(18)

3 實驗設(shè)置

實驗裝置圖1(a)所示，所有的銑削實驗都是在數(shù)控銑削機(jī)床上完成的，動態(tài)銑削力信號測量系統(tǒng)由Kistler 9272測力計、Kistler 5070A電荷放大器、5697數(shù)據(jù)采集卡和PC組成，采樣頻率為20 000 Hz。工件材料為400 mm←100 mm←100 mm的鋁合金， 刀具為四齒硬質(zhì)合金球頭銑刀， 半徑4 mm， 懸伸長度為40 mm。 顫振的產(chǎn)生與銑削參數(shù)有直接的關(guān)系， 因此可以通過改變主軸轉(zhuǎn)速n、進(jìn)給速度f、軸向銑削深度ap獲取感興趣的實驗數(shù)據(jù)。

本文進(jìn)行了5次實驗，實驗條件如表1所示。 Test1、Test2、Test3、Test4構(gòu)成第一組實驗； Test5分別在穩(wěn)定、過渡、顫振時采樣，這3段采樣信號構(gòu)成第二組實驗。

圖1 實驗裝置及方法流程圖Fig.1 The experimental platform and the flow chart

狀態(tài)n/(r·min-1)f/(mm·min-1)ap/mm第一組Test1穩(wěn)定4 2005806.6Test2穩(wěn)定3 6002502.6Test3顫振4 0005806.6Test4顫振4 0006406.7第二組Test5過渡4 0005806.3

4 結(jié)果分析

對于第一組實驗，分別對Test1、Test2、Test3、Test4進(jìn)行多分量信號分解，計算子信號的能量在頻域上的分布。對比分析在不同實驗條件下，穩(wěn)定和顫振的子信號的能量熵。

圖2(a)為Test5采集到的進(jìn)給方向動態(tài)銑削力的時域波形，Test5經(jīng)歷了從穩(wěn)定銑削過渡到顫振的過程，實驗前0.6 s銑削穩(wěn)定，0.6 s～0.8 s信號振幅急劇增大， 0.8 s開始顫振。 分別自0.2 s、0.8 s、2 s起， 往后截取5 000個采樣點長的信號記作Sample1(穩(wěn)定)、Sample2(過渡)、Sample3(顫振)，作為第二組實驗。對比Sample1、Sample2、Sample3經(jīng)過信號分解得到的子信號的能量分布，分析在實驗條件不變的情況下，隨著銑削從穩(wěn)定過渡到顫振，子信號在頻域的能量熵的變化。

圖2 Test5及各采樣點時域波形Fig.2 Time domain waveform of test5 and sampling points

4.1 數(shù)據(jù)分析

第一組原始實驗數(shù)據(jù)的時域波形和頻譜如圖 3所示，Test1(圖3(a))、Test2(圖3(b))銑削穩(wěn)定，頻譜由主軸轉(zhuǎn)速頻率及其諧波、噪聲組成；改變銑削條件，動態(tài)銑削力振幅急劇增加，發(fā)生顫振(圖3(c))，頻譜上的出現(xiàn)能量遠(yuǎn)遠(yuǎn)高于其他頻率的顫振峰值；繼續(xù)增加銑削深度，顫振的更加劇烈，如圖3(d)。

圖3 Test1、Test2、Test3、Test4的波形和FFTFig.3 Vibration signal in time domain andfrequency domain of Tests

第二組原始實驗數(shù)據(jù)的頻譜如圖4(a)、圖4(b)、圖4(c)所示，f表示主軸轉(zhuǎn)速頻率，fc表示顫振頻率。銑削穩(wěn)定時(圖4(a))，頻譜由主軸轉(zhuǎn)速及其諧波成分和分布不確定的噪聲組成；過渡狀態(tài)(圖4(b))頻譜除了主軸轉(zhuǎn)速諧波成分和噪聲，還出現(xiàn)了顫振分量和它的各階諧波成份，此時工件表面還沒有明顯的顫振；隨著顫振程度增加(圖4(c))，顫振頻率及其諧波對應(yīng)的峰值也急劇增長，信號的能量除了集中在主軸轉(zhuǎn)速頻率處，都匯集在fc及其諧波上，且顫振峰值還再繼續(xù)增加。

在分解信號、分析子信號能量頻域分布之前，為了避免主軸轉(zhuǎn)速及其諧波分量對顫振識別產(chǎn)生的干擾，先估計這些分量的瞬時頻率，用Vold-Kalman濾波提取此分量信號并從原始信號中剔除這一部分。對預(yù)處理過的銑削信號進(jìn)行FFT，穩(wěn)定銑削時頻譜(圖4(d))只剩下在頻率域上分布均勻且峰值較低的隨機(jī)振動；過渡到顫振時(圖4(e)、圖4(f))， 頻譜在fc及2倍的fc處出現(xiàn)能量很高的峰值，除此之外還有能量比較低的噪聲。

通常動態(tài)銑削力信號具有非線性、非平穩(wěn)等特征，直接進(jìn)行傅里葉變換無法揭示頻率分量隨時間的變化規(guī)律；短時傅里葉變換、小波變換等能夠建立信號在時間和頻率上的分布，然而它們是先驗性的自適應(yīng)能力差。并且它們只能用來離線判斷銑削狀態(tài)，無法運用在實際中。預(yù)報銑削顫振要準(zhǔn)確、及時，一旦出現(xiàn)顫振發(fā)生的預(yù)兆，計算機(jī)可以識別并且迅速的啟動控制系統(tǒng)，調(diào)整轉(zhuǎn)速、銑削用量等，將顫振控制在孕育階段。故障振動信號的能量分布會隨著顫振發(fā)生相應(yīng)的改變，故采用分解后子信號在頻域分布的熵值的范圍可用來判斷是否顫振。

圖4 Sample1、Sample2、Sample3濾波前后的FFTFig.4 FFT of Sample1、Sample2、Sample3 beforeand after filtering

首先，Vold-Kalman濾波器波除去信號的轉(zhuǎn)速頻率及其諧波成分。優(yōu)化頻譜集中性指標(biāo)來估計發(fā)生顫振時預(yù)處理過的信號的瞬時頻率參數(shù)。用估計到的參數(shù)將非線性調(diào)頻信號解調(diào)為平穩(wěn)信號。最后Vold-Kalman濾波器提取解調(diào)單分量信號。本文對實驗信號精細(xì)分解，按照分解后單分量子信號能量從高到低，得到這些分量集中度較好的頻譜表示。由于主要的故障信息集中在前幾個子信號中，剩余信號能量遠(yuǎn)小于原始信號總能量，可以忽略，因此本文選用了前8個子信號，8個子信號的能量之和應(yīng)該恒等于原始信號的總能量。由于篇幅限制，本文僅展示預(yù)處理后的第二組信號的Sample1、Sample3的分解結(jié)果。Sample1的經(jīng)過分解后的各子信號按能量從高到低如圖5(a)所示，均為隨機(jī)分布的噪聲，能量分布相對平均和不確定；Sample3經(jīng)過分解后的各子信號能量從高到低如圖5(b)所示，其中u1，u2能量遠(yuǎn)高于其他各子信號能量，信號頻譜的幅值范圍均在0～4。u1的峰值正好在銑削系統(tǒng)的固有頻率附近，證明本實驗出現(xiàn)的是典型的再生顫振。u2為顫振頻率的2倍頻。剩余各階子信號是被齒頻調(diào)制的顫振分量或噪聲，其時域波形的振幅小，能量明顯低于u1和u2。

經(jīng)過預(yù)處理的動態(tài)銑削力信號經(jīng)過分解后得到多個分量，且各個分量的瞬時頻率在頻譜上不相交，子信號的能量隨分解階數(shù)的增加而降低。在信噪比低的情況下，基于頻譜集中性指標(biāo)估計分量的瞬時頻率較精確，抗噪性強(qiáng)。即使對于能量較小的單分量的瞬時頻率參數(shù)的估計也比較準(zhǔn)確。

圖5 本文方法分解結(jié)果Fig.5 Decomposition results of this method

4.2 顫振識別

隨著顫振的發(fā)生，振動的能量發(fā)生轉(zhuǎn)移，所以振動的概率分布也發(fā)生了變化，這些都可以用熵來進(jìn)行判斷。由于各階子信號含不同的頻率成分且具有不同的能量，分別求它們的能量E={E1,E2,…,En}。首先根據(jù)式(18)將各子信號的能量歸一化。分析第一、二組實驗經(jīng)過分解后的能量分布，分別如圖6(a)，圖6(b)所示。對于第一組實驗，將分解得到的子信號歸一化，u1～u8的能量從高到低，分布相對均勻(圖6(a))，熵值較大；而Test3能量集中在u1和u2，所以熵值相對較??；Test4顫振更加劇烈，隨著顫振的程度增加，u3～u8所占能量幾乎為零，所以計算得到的熵值最小。計算Test1、Test2、Test3、Test4的能量熵結(jié)果如圖7(a)所示。對于第二組實驗，短暫的穩(wěn)定銑削過程中，子信號能量分布較為均勻。隨著顫振的出現(xiàn)，能量急劇轉(zhuǎn)移到u1和u2，顫振劇烈時，u1和u2匯集了幾乎全部的能量。根據(jù)式(19)，計算Sample1、Sample2、Sample3的能量熵結(jié)果如圖7(b)所示，可以看出銑削從穩(wěn)定過渡到顫振，熵的值急劇減小。第二組的分析更加證明了本文所提出的顫振識別方法的可靠性。

圖6 第一組、第二組歸一化后各子信號能量分布Fig.6 Normalized energy ration of the sub-signals

穩(wěn)定銑削時，預(yù)處理后的銑削力信號經(jīng)分解，得到的子信號能量分布較平均和不確定，能量熵較高。當(dāng)出現(xiàn)顫振時，在相應(yīng)的頻帶內(nèi)就會出現(xiàn)相應(yīng)的顫振頻率。此時，能量集中在顫振頻率上，使能量分布的不確定性減小，從而熵減小。通過第一組、第二組實驗的分析，可見本文提出的方法可以有效識別顫振。然而由于實驗條件的限制，本文沒有給出發(fā)生顫振時能量熵的閾值。

圖8所示為Test5從穩(wěn)定銑削到顫振的過程，可見在0.6 s顫振發(fā)生時，能量熵明顯變小，證實了能量熵可以作為監(jiān)測顫振的指標(biāo)。本文的不足之處在于，由于實驗前后時間跨度大，無法提供工件表面加工質(zhì)量前后對比照片。

圖7 第一組、第二組能量熵Fig.7 Energy entropy of the tests

圖8 Test5能量熵隨銑削力的變化Fig.8 The vibration signal and the relatedenergy entropy in Test 5

5 結(jié) 論

本文提出了一種簡單有效的銑削顫振識別方法。動態(tài)銑削力信號具有不平穩(wěn)、頻率時變、噪聲大等特點，分解信號前先經(jīng)過Vold-Kalman濾波器去除信號中的主軸轉(zhuǎn)速及其諧波分量，再通過優(yōu)化頻譜集中性指標(biāo)來估計信號瞬時頻率參數(shù)，用估計到的參數(shù)將非線性調(diào)頻信號解調(diào)為平穩(wěn)信號，最后用Vold-Kalman濾波器提取單分量信號，并得到能量從高到低的各階子信號。該信號分解方法抗噪性好，對實驗信號分解精細(xì)，并得到這些分量集中度較好的頻譜表示。隨著銑削條件變化，經(jīng)過分解得到的各子信號，其能量分布也呈現(xiàn)變化。因此，本文提出表征能量分布特征的能量熵作為顫振識別的特征。在穩(wěn)定銑削時，實驗信號進(jìn)行預(yù)處理后，其分解得到的各階子信號主要為噪聲，因此能量熵小。隨著顫振越來越劇烈，振動的能量集中到含有顫振信息的子信號上，因此能量分布從分散變得集中起來，能量熵變大。經(jīng)過實驗驗證，該方法簡單有效，尤其在顫振發(fā)生的早期識別精度高，是實現(xiàn)顫振在線實時控制的關(guān)鍵。