航發葉片機器人磨削顫振檢測方法*

李鵬偉，戴士杰，張文華，母嘉恒，胡天榮

（河北工業大學，天津 300400）

隨著機器人智能控制與各個學科的交叉融合，在航空發動機葉片表面磨拋中高效率、高精度的機器人加工備受青睞。航空發動機葉片屬于典型的薄壁結構零件，在機器人高速磨削加工中很容易發生顫振。機器人加工系統的顫振是發生在刀具和工件之間的相互振動，是一種復雜的彈性動力學現象，主要由機器人加工系統整體結構從周期性切削力中吸收能量自激產生[1]。顫振不僅會造成磨削力波動，而且會嚴重降低加工表面質量[2]。為了防止顫振對磨削系統破壞，早期、快速、準確的顫振檢測至關重要。

顫振檢測的過程主要分為3部分：數據采集、特征提取、狀態識別[3]。其中，在顫振信號的數據采集過程中，可以利用各種加工信號，例如，聲音信號采集簡單可靠，但容易受到外界噪聲的影響；力信號對切削狀態的變化反應迅速，但傳感器安裝位置的選擇比較困難；加速度信號對振動的表現最為直觀，但采集加速度信號需要昂貴的儀器[4-7]。在顫振信號的特征提取過程中，可以利用信號處理的方法構造特征指標，該指標為后期定義顫振判據提供理論依據。Chen等[8]依靠數據處理算法，從提取的多個特征參數中自動選擇有效的顫振識別特征參數。Yao等[9]利用顫振信號的小波變換標準差和小波包變換能量比作為二維特征向量進行顫振檢測。在顫振的狀態識別過程中，實時性是顫振識別的一個重要要求，這意味著顫振應該事先被識別出來。Chen等[10]提出一種基于EEMD和特征選擇的多通道顫振智能檢測策略，但煩冗的計算量會導致系統實時性較差。Li等[11]通過經驗模態分解成功提取顫振識別的特征向量，該方法用時不足1 s。

目前，關于加工顫振檢測的研究主要圍繞銑削、車削等展開，而較少關注機器人磨削[12]，而且大部分檢測方法仍然不能在工件損壞之前識別顫振，達不到防止工件損壞的核心目標。因此，本文提出一種基于排列熵的機器人磨削顫振檢測方法：首先，采集機器人末端的顫振加速度信號；然后，利用滑動窗口算法從原始信號中提取連續的排列熵值；最后，通過排列熵經驗閾值判斷末端執行器是否發生顫振。在此過程中，采用基于序數模式的排列熵算法，利用序數模式及其分布規律，簡化了排列熵計算過程，提高了排列熵特征的計算效率。仿真和試驗結果表明，排列熵值的驟降時刻遠遠早于顫振完全產生時刻，說明該方法能夠在顫振孕育過程的初始階段檢測出顫振，為后續采取抑制措施提供更多的反應時間。

1 基于序數模式的排列熵

1.1 序數模式

序數模式描述了時間序列相鄰值之間的順序關系。可設時間序列 （xm，xm-t，…，xm-dt） 的d階序數模式為序數模式的計算表達式為

式中，m為時間序列中數據的數量；h=1，2，…，d；t為時間延遲，代表相鄰序數模式的中點間的距離；d為階數，代表一個序數模式中含有（d+ 1）個數據。

通過給每一種序數模式分配代號，可以更加方便地觀察序數模式的分布情況，其中代號按式 （2）進行一對一的分配。

例如，當階數d= 2時，會產生6種序數模式類型，分別對應數字代號“0、1、2、3、4、5”，如圖1所示。

圖1 d = 2時序數模式的類型Fig.1 Types of ordinal patterns at d = 2

1.2 基于序數模式的排列熵計算

排列熵是一種衡量時間序列復雜程度的特征，能夠準確地定位系統發生突變的時刻，可以對突變信號的微小變化時刻進行放大。

滑動窗口內的時間序列 （xm，xm-1，…，xm-M-dt+1） 在時間延遲t和階數d下的排列熵可以根據式 （3）進行計算。

式中，M為窗口中存在的序數模式的數量；qj= #[k∈{m，m-1，…，m-M+ 1}|（xk，xk-t，xk-dt）存在序數模式j]。

設由10個值組成的時間序列為 （x1，x2，…，x10），此時間序列中時間延遲t= 2、階數d= 2的序數模式如圖2所示。

圖2 滑動窗口中的序數模式Fig.2 Ordinal patterns in the sliding window

由圖2可知，兩個窗口存在4種序數模式的“重疊”部分。根據式 （3）可在窗口1和窗口2中相應地求得排列熵（9）和（10），即

1.3 高效計算連續序數模式

圖3 序數模式“0”以及其連續序數模式的可能類型Fig.3 Ordinal pattern “0” and its possible types of consecutive ordinal patterns

由圖3可知，對于序數模式“0”的下一個值有3個可能的位置，則下一個連續序數模式有3種可能的類型：“0、3、4”。

同理可得，可以對其他階數的給定序數模式的下一個連續序數模式的可能類型進行分析。當前序數模式已知，通過下一個數據與當前序數模式中數據的大小比較，可求得下一個序數模式的類型，如式 （4）所示。

式中，ωd表示根據當前序數模式來確定下一個序數模式的函數；I表示 （xm，xm-t，…，xm-(d–1)t）中有多少個值大于或等于xm+t。

根據式 （4）和 （5），可以求得所有階數序數模式的連續序數模式。d=2階的序數模式及其連續模式如圖4所示。例如，對于序數模式“0、1、2”，下一個值有3個可能的位置“I= 0、1、2”，也就有3個可能的連續序數模式“0，3，4”。對于序數模式“3、4、5”，下一個值有3個可能的位置“I= 0，1，2”，也就有3個可能的連續序數模式“1，2，5”。同理可得，每個代號為 （d+1）（k-1），…，（d+1）k-1的序數模式具有相同的連續序數模式，其中k=1，2，…，d!。根據上述方法，可以求得所有的序數模式以及其連續序數模式。

圖4 階數d =2的所有序數模式以及其連續序數模式Fig.4 All ordinal patterns of order d =2 and their consecutive ordinal patterns

排列熵是在具有最大重疊的大小為M的滑動窗口中計算的，即連續窗口的第1個值是前一個窗口的第2個值。連續的窗口 （xm-1，xm-2，…，xm-M-dt）和 （xm，xm-1，…，xm-M-dt+1）在值xm和xm-M-dt上有所不同，因此窗口中的序數分布在序數模式（m）和（m-M）的出現頻率上有所不同。所以在已知排列熵（m-1）的情況下，排列熵（m）的計算方法可表示為

1.4 基于序數模式的排列熵的計算效率

本文提出的基于序數模式的排列熵計算方法通過比較信號數據中相鄰數值的大小確定序數模式的類型（序數模式類型的定義已經給出）。利用圖4的序數模式分布規律可快速獲得連續的序數模式分布情況。最后通過連續序數模式的類型，以及式 （1）～（6）計算基于序數模式的排列熵值，并利用滑動窗口算法得到連續的基于序數模式的排列熵值。

分別利用傳統排列熵算法和所提出的基于序數模式排列熵算法對不同信號長度以及不同階數的隨機信號進行分析，然后對兩種算法排列熵提取的計算時間進行比較 （計算時間是多次運行的平均值），結果如表1所示。

表1 不同方法計算排列熵值的執行時間Table 1 Execution time for calculatating permutation entropy using different methods

由表1可知，當原始信號長度加倍時，兩種計算方法的執行時間都有所增加。然而，相比于傳統排列熵計算方法，基于序數模式排列熵計算方法的計算時間要低于一個數量級。因此，基于序數模式的排列熵計算方法效率更高，更有利于實現顫振的及時檢測。

基于序數模式的排列熵方法對上述有效計算進行總結，計算方案如圖5所示。其中，W為時間序列的長度。

圖5 基于序數模式的排列熵計算方案Fig.5 Permutation entropy calculation scheme based on ordinal pattern

2 仿真算例分析

顫振的發生一般有3個階段，即平穩階段、過渡階段和顫振階段[13]。顫振信號通常表現為單次諧波特征，因此，當顫振發生時，加速度信號可用單頻諧波函數描述。為了驗證所提出的磨削顫振檢測方法的有效性，現利用仿真軟件，直觀地表達出顫振從穩定到不穩定逐步發展的過程，且將存在平穩、過渡和顫振3個階段的顫振信號模擬出來，顫振加速度模擬信號如圖6所示。

圖6 經歷平穩、過渡和顫振階段的加速度模擬信號Fig.6 Simulated acceleration signal through smooth, transition and chatter phases

本研究選取階數d= 3，時間延遲t= 1，滑動窗口大小M= 512對顫振原始信號進行窗口滑動計算，每次對每個窗口范圍的信號進行一次排列熵值計算，最終得到連續的排列熵值數據，可以更直觀地觀察加速度信號變化與排列熵值變化的差異對比情況。滑動窗口算法計算出基于序數模式排列熵的動態變化情況如圖7所示。

圖7 顫振信號的排列熵特征Fig.7 Permutation entropy characteristics of chattering signal

由圖7可知，在平穩階段的排列熵特征曲線大致保持在水平穩定狀態。由于周期性振動的積累，顫振逐漸產生，排列熵特征曲線發生驟降，而且驟降發生在顫振過渡階段的前期。排列熵特征驟降可以作為顫振即將發生的征兆，而且排列熵驟降時刻要遠遠早于顫振完全產生時刻，因此可以利用排列熵驟降作為顫振檢測的特征。

3 試驗驗證

3.1 試驗平臺搭建

本節通過航空發動機葉片機器人磨削試驗來驗證本文所提出的磨削顫振檢測的可靠性，試驗平臺的搭建如圖8所示，加工工件是材質為鎳基合金的航空發動機壓氣機葉片，磨具為粒度500目的柔性砂輪。

圖8 航空發動機葉片機器人磨拋試驗平臺Fig.8 Robot grinding experimental platform for aero-engine blade

機器人磨削顫振信號采集裝置如圖9所示，相比于末端執行器，三軸加速度傳感器的體積及重量極小，在此可忽略加速度傳感器重量對末端顫振的影響。其他試驗設備及相關參數如表2所示。

表2 試驗設備及其主要參數Table 2 Experimental equipment and its main parameters

圖9 機器人磨削顫振信號采集裝置Fig.9 Acquisition device of chattering signal of robot grinding

3.2 顫振信號分析

通過改變磨削加工參數，采集機器人末端在不同加工狀態下的顫振信號，磨削加工參數如表3所示。

表3 試驗中的磨削加工參數Table 3 Grinding parameters in the experiment

試驗由安裝在磨削末端執行器上的三軸加速度傳感器測量顫振加速度信號，通過HY-DAS數據系統采集顫振加速度信號，采用上位機對顫振加速度信號進行實時排列熵提取。其中，利用本文所提出基于序數模式的排列熵算法進行磨削顫振檢測，即選取階數d= 3、時間延遲t= 1、滑動窗口大小M= 512提取排列熵特征。根據表3中的加工參數進行航空發動機葉片機器人磨削，在試驗結果中選取幾組代表性數據進行分析，得到顫振原始信號及其連續排列熵特征，如圖10所示。

圖10 顫振原始信號和排列熵Fig.10 Chattering raw signal and permutation entropy

圖10（a）中排列熵值在0.23 s時低于0.95，早于顫振完全發生階段0.41s；圖10（b）中排列熵值在0.2 s時低于0.95，早于顫振完全發生階段0.5 s；圖10（c）中排列熵值在0.88 s時低于0.95，早于顫振完全發生階段0.52 s。由此可知，排列熵值低于0.95時刻早于顫振完全發生時刻約0.48 s。

由圖10可知，當顫振未發生時，排列熵特征處于水平穩定的狀態。隨著顫振的產生，原始信號幅值隨之增加，在過渡階段的初始位置，排列熵值發生驟降 （虛線為排列熵驟降時刻），此時即可判斷出磨削過程中發生了顫振，而且排列熵值驟降時刻遠遠早于顫振完全產生時刻，即可以在較早階段檢測出顫振。因此以排列熵值的驟降作為顫振檢測的標準，能夠在顫振完全產生前0.48 s檢測出顫振，為磨削系統采取顫振抑制措施提供更多的反應時間。該試驗證明所提出的方法在機器人磨削顫振檢測方面是有效的。

3.3 特征對比

基于振動信號時域分析理論，選取方差、裕度因子作為顫振檢測的特征量來反映磨削顫振信號的變化趨勢，利用滑動窗口算法提取方差、裕度因子特征曲線，并與排列熵特征進行對比，結果如圖11所示。

圖11 顫振檢測特征：排列熵、方差、裕度因子Fig.11 Features of chatter detection: permutation entropy, variance,margin factor

由圖11可知，方差特征曲線會隨著信號幅值的增減而增減，其變化受原始信號幅值的影響較大，而且方差特征曲線的變化幅度小于原始信號變化幅度，所以其檢測準確性較差。裕度因子特征曲線變化不穩定，曲線波動范圍較大，所以其檢測精度較差。相比于方差和裕度因子，小波包熵特征驟降時刻要更早 （圖11中的虛線為各個特征的轉折時間點），而且其波動較小、變化穩定。

綜上所述，排列熵特征對顫振發生的反應速度更快，因此利用排列熵特征能夠更早地檢測出顫振，為磨削系統采取顫振抑制措施提供更多反應時間。

3.4 顫振狀態識別

顫振狀態識別最簡單的方法就是特征閾值法，通過對比某個特征量在平穩階段與顫振征兆出現階段的變化情況，設定一個經驗型的閾值，一旦超過或低于預先設定的閾值則判定為出現顫振，即可提醒系統采取顫振抑制措施。根據圖10中的排列熵特征的動態變化規律，可設置經驗閾值為0.95，即當排列熵低于0.95則判定此時加工系統中發生顫振。

綜上所述，機器人磨削末端顫振檢測方法的流程如圖12所示。首先設置磨削加工參數，利用機器人帶動末端執行器進行航空發動機葉片磨削，并通過安裝在末端執行器上的加速度傳感器對磨削顫振信號進行采集；然后利用序數模式進行排列熵特征連續提取；最后當排列熵特征值小于0.95時，即判定此時機器人末端正在產生顫振，可為系統發出警告，提醒系統采取顫振抑制措施。

圖12 顫振檢測方法流程圖Fig.12 Flow chart of chatter detection method

試驗結果表明，當顫振信號加速度幅值剛開始增加時排列熵值就會發生驟降，因此，在過渡階段的初始位置將顫振檢測出來 （排列熵值低于0.95時即可認定此時的磨削過程中發生了顫振）為進行抑制操作提供了更多反應時間。

4 結論

針對航空發動機葉片機器人磨削顫振檢測不及時導致顫振抑制滯后的問題，本文提出了一種基于排列熵的航空發動機葉片機器人磨削末端顫振檢測方法。該方法通過對原始信號進行排列熵特征提取，并利用該特征進行顫振檢測，當排列熵值低于0.95時，判定機器人末端正在產生顫振，立即向系統發出警告，提醒系統采取顫振抑制措施。

針對實時性要求，采用基于序數模式的排列熵算法對顫振原始信號實時計算排列熵，該方法相較于傳統排列熵算法計算時間減小了一個數量級，當信號長度為10000時，基于序數模式排列熵算法的計算時間約為0.25 s。

最后，試驗結果表明，通過排列熵特征提取算法的改進，基于排列熵的航空發動機葉片機器人磨削末端顫振檢測方法能夠在顫振爆發前約0.48 s檢測出顫振，可及時向系統發出警告，為后續采取措施提供了更多的反應時間。

目前本文中排列熵特征閾值的應用范圍有待進一步拓展，后續還需進一步研究狀態識別閾值的自適應選擇標準。