Online chatter detection in robotic machining based on adaptive variational mode decomposition

Chatter is the main problem that limits the application of industrial robots in the field of machining process. It is critically important to establish an adaptive chatter detection solution for robot machining process and realize the online detection of chatter. However, different from machine tool chatter, the chatter in robotic machining process is more complex to be detected due to the variable stiffness characteristics and weaker stiffness of normal industrial robot, and the existing literature has less research on this problem. This paper presents a comprehensive solution for online chatter detection in robotic machining process. Firstly, in order to detect the chatter in robotic machining process and avoid mode mixing problem in variational mode decomposition (VMD) process, an adaptive variational mode decomposition (AVMD) method based on kurtosis and instantaneous frequency is proposed, which realizes the adaptive selection of the decomposition parameter. Secondly, optimal decomposition parameters are calculated by using genetic algorithm. By optimizing the discrete step length of decomposition parameter, it greatly reduces the optimization time. Last but not least, approximate entropy, energy entropy, and proposed entropy drift coefficient are extracted to distinguish chatter and stable machining state. Simulation and experimental results show that the proposed method can meet the real-time requirements of online detection and detect the occurrence of chatter effectively.


Introduction
During the last few decades, industrial robots have been widely used in many industries for various applications, including welding [1], material handling [2], painting [3], assembly [4], and machining [5]. Among them, robotic machining has attracted significant attention from both academics and industry because of its low cost, strong flexibility, and taking up small space [6,7]. However, the stiffness of an industrial robot is usually less than 1 N/μm, while the stiffness of machines tools is often greater than 50 N/μm [8], which makes chatter easier to happen in robotic machining than traditional machine tools.
Chatter in the machining process belongs to the selfexcited vibration, which seriously affects the surface quality of the workpiece and limits the machining efficiency [9,10]. It is critical that vibration needs to be controlled within an acceptable range [11]. Chatter vibration has become a major obstacle to achieve stable and efficient robotic machining and is generally caused by mode-coupling or regeneration of chip thickness [12]. There are many literatures addressing the chatter issues during the machining process including the chatter mechanism, chatter detection, and suppression strategies in the past decades, such as Hayati et al. in [13] proposed a boring bar with an internal frictional damper to reduce chatter vibrations of boring bars. Dong et al. in [14] presented an updated numerical integration method (UNIM) to analyze the stability of delay-differential equation (DDE) to avoid chatter vibrations in milling process. Yuan et al. in [15] developed a novel method to mitigate chatter vibrations in the thin-wall milling of the structures with half-opened side walls through designing a supplementary device. However, few literatures have been devoted to chatter detection in robotic machining process. Cen et al. [16] presented a discrete wavelet transform (DWT)-based online chatter detection algorithm and chatter suppression strategy applied to robot milling. Tao et al. [17] proposed a novel approach to identify the chatter in robotic drilling process based on multi-synchrosqueezing transform (MSST) and energy entropy. Sun et al. [18] used rotary ultrasonic milling (RUM) technology to suppress chatter in robotic milling process. Pan et al. [8] claimed that the mode coupling chatter was the dominant chatter in the robotic milling process, and the occurrence of the chatter during other machining applications was ignored. Through the mode analysis and milling experiments, Marcel Cordes et al. [19] found that only low-frequency robot's structural modes cause chatter in milling of titanium. Lin et al. [20] proposed a SC-based (spindle configuration) deformation model, which establishes a mapping between the spindle configuration and the deformation of the robot end effector. The loss of stability in robotic machining process often results in violent vibration of the robot itself, which will exacerbate the damage to the robot body and cutter. Early chatter detection allows operators to interfere in the robotic machining process and avoid damage to robot body. Significant research has been conducted on chatter detection in traditional machining process [21]. Various sensor signals have been used to monitor chatter, such as cutting force, vibration acceleration, motor current, sound, torque signal, or the combination of them. Wang et al. [22] proposed to segment vibration signal and current signal to create a multisensor fusion network for motor condition monitoring, which not only improves the effectiveness of feature extraction but is also adaptive to varying machining parameters. After years of research, scholars generally believe that acceleration signals and cutting forces are more capable of detecting chatter in machine tool processing for their comprehensive performance in price, installation, signal quality, and information usability.
In order to obtain more information in robotic machining dynamics, a tri-axial acceleration sensor is installed on the spindle, and the collected acceleration signal is used as monitoring vibration signal.
The statistical characteristics of vibration signal in the time domain and frequency domain varied with the change of vibration, and this change is usually disturbed by the background noise [23]. Therefore, advanced signal processing methods are critical for the subsequent analysis in the chatter detection process. Extensive research has been focused on researching signal processing methods, and significant progress has been achieved [24][25][26]. Traditional spectral methods based on fast Fourier transform (FFT) are excluded because of the nonlinear and non-stationary nature of chatter signals. Wavelet analysis and later wavelet transform (WT) and wavelet packet transform (WPT) are powerful analytical tools in both the time and frequency domains. But it is difficult to determine the suitable wavelet base and decomposition layers. Chen et al. [27] proposed a feature extraction and optimization method based on recurrence quantitative analysis and affine propagation clustering to monitor chatter. Shi et al. [28] introduced a new method using ordered-neuron long short-term memory (ON-LSTM) and population-based training (PBT) to detect chatter in high-speed milling processes. Dong et al. [29] used the complexity index of vibration signals to identify chatter in milling of thin-walled part and found the change laws of the complexity of vibration signals from stably milling to chatter. Especially, more and more machine learning methods have been used in signal feature extraction and achieved good results. Qiao et al. [30] proposed a novel endto-end deep learning network named adaptive weighted multiscale convolutional neural network (AWMSCNN) to adaptively extract robust and discriminative multiscale fusion features from raw vibration signals. Yang et al. [31,32] respectively implemented the FTNN (feature-based transfer neural network) and PK-MMD (polynomial kernels induced maximum mean discrepancy) to construct models for intelligent signal analysis and both methods presented higher classification accuracy. Empirical mode decomposition (EMD) was proposed by Huang et al. [33] as a self-adaptive analysis method. This method is based on the local characteristic time scales of a signal and may decompose a complicated signal function into a set of complete and almost orthogonal components named intrinsic mode functions (IMFs) and without requiring the selection of the wavelet base and the decomposition layers. A lot of studies in recent years have investigated efficient methods of detecting chatter with EMD in traditional machine tool processing. In order to detect the onset of chatter vibrations, the recurrence quantification analysis and the EMD are employed in the milling process of titanium alloy Ti6242 [34]. Peng [35] adopted EMD based on timefrequency analysis for the detection of tool breakage in milling process. Cao et al. [21] combined the benefits of WPT and EMD to extract features according to the Hilbert-Huang spectrum for chatter detection. However, EMD suffers from many issues, such as mode mixing, boundary effect, and sensitive to noise. To solve these problems, some alternative methods have been introduced. Cao et al. used a self-adaptive analysis method, termed ensemble empirical mode decomposition (EEMD), to analyze vibration signals and extract two nonlinear indices to indicate chatter [36]. Fu et al. [37] presented a comprehensive solution by combining EEMD and Hilbert spectral analysis method to extract features in milling process. Although EMD can effectively decompose the non-stationary signal, it lacks theory basis and cannot be evaluated from mathematical basis.
The chatter frequency generated during the robotic machining process may be changeable because of the constantly changing stiffness with different joint configuration [19], which makes the chatter detection more complicated in robotic machining process. Variational mode decomposition (VMD) was introduced by Dragomiretskiy and Zosso [38], which has solid theoretical basis and can decompose the signal into several sub-components concurrently. Meaningfully, VMD can overcome the mode mixing problem of EMD, and can exhibit good performance of signal processing. However, the problem of selecting the number of decomposed modes K and the penalty coefficient α still affects the decomposition results.
In order to determine the parameters of VMD, Yang et al. [39] searched for the optimal decomposition parameters based on the simulated annealing method (SA), and accurately detected the occurrence of chatter by calculating the approximate entropy and the sample entropy, but the SA method requires large amounts of time. Liu et al. [40] introduced an automatic selection method based on particle swarm optimization (PSO) and the maximum crest factor of the envelope spectrum (CE) to solve the problem of parameter selection. Liu et al. [41] proposed an automatic selection method of decomposition parameters based on kurtosis values. They selected the optimal decomposition parameters by comparing the kurtosis values of the reconstructed signal and achieved good decomposition results. However, choosing the optimal sub-component based on chatter frequency was not suitable for situation in robotic machining where the chatter frequency is uncertain. The complexity of chatter characteristics in robotic machining increases the difficulty of extracting chatter components from vibration signal. In order to prevent chatter from occurring and avoid severe chatter from damaging the robot and tool system, it is urgent to propose a chatter detection method suitable for robotic machining process.
In this study, a comprehensive solution for chatter detection in robotic machining process has been proposed. Firstly, in order to avoid possible mode mixing problem due to overdecomposition, the adaptive variational mode decomposition (AVMD) method based on kurtosis and instantaneous frequency is presented, which realizes the adaptive selection of the decomposition parameter. Secondly, optimal decomposition parameters are calculated by using genetic algorithm. By optimizing the discrete step length of decomposition parameter, it greatly reduces the optimization time. Last but not least, approximate entropy, energy entropy, and proposed novel chatter feature based on entropy are extracted to distinguish chatter and stable machining state.
The rest of paper is organized as follows. Section 2 briefly describes the methodology of chatter detection in robotic machining. Section 3 describes experimental setup and the chatter characteristics of the robotic milling. The analysis of simulation and experimental results are presented in Section 4, and conclusions are given in Section 5.
2 Methodology of chatter detection in robotic machining

Variational mode decomposition method
Variational mode decomposition (VMD) is a non-recursive method which has a rigorous mathematical derivation process, and overcomes the problem of mode mixing compared with empirical mode decomposition (EMD). The purpose of VMD is to decompose a real valued input signal into a discrete number of intrinsic mode functions (IMFs), u k , that have specific sparsity properties while reproducing the input. It is assumed each mode to be mostly compact around a center pulsation w k , which is to be determined along with the decomposition. VMD can be treated as a constrained variational problem, which is represented as follows: s:t:∑ k u k ¼ y ð2Þ where {u k } = {u 1 , …, u k } represents all modes and {w k }-= {w 1 , …, w k } represents their center frequencies, respectively. y is the original signal, δ is the Dirac distribution, and * denotes the convolution operator. Meanwhile, by introducing a quadratic penalty term α and Lagrange multipliers λ, the constrained problem is transformed into an unconstrained problem: By mean of alternate update r nþ1 k , w nþ1 k , λ n + 1 , Lagrange expression's saddle point can be looked for. Then the solution of this quadratic optimization problem is obtained by letting the first variation vanish for the positive frequencies: The center frequencies are solved as: where b r nþ1 k w ð Þ is identified as a wiener filtering of the current residual b y w ð Þ−∑b r i w ð Þ, and w nþ1 k is the center of gravity of the corresponding mode's power spectrum.

Adaptive variational mode decomposition method
The penalty coefficient α and the number of decomposed modes K are critically important in the decomposition process for variational mode decomposition (VMD) and must be chosen in advance. In this study, to evaluate the decomposition performance of K and α in signal processing, adaptive variational mode decomposition (AVMD) method based on kurtosis and instantaneous frequency is firstly established. On the one hand, in order to avoid the mode mixing problem caused by over decomposition in VMD, an adaptive search method for the maximum value of decomposed mode number K based on instantaneous frequency is proposed in this paper. On the other hand, the kurtosis reflects the index of the degree of deviation of the signal from the Gaussian distribution, and is particularly sensitive to the chatter component in the machining signal. Therefore, the kurtosis is used as an index to evaluate the enhancement performance of the chatter component. The steps of the proposed adaptive VMD method are as follows: Step (1). Automatic search of parameter K max . K max is the maximum value of K that guarantees no over-decomposition. The initial values of K and α are given as K = 2, α = 1/2fs (fs means the sampling frequency), and α remains unchanged in this step. The signal is decomposed into K sub-components using VMD, and the instantaneous frequencies of sub-components are calculated by the Hilbert method respectively. Then instantaneous frequencies of any two sub-components are compared with each other. If any two sub-components have close frequency (5%), the maximum value is K max = K and the process switches to Step (2). Else, K = K + 1 and do Step (1) again.
Step (2). Determine the parameter value range. It is critically important that α cannot be less than 1/2 of the sampling frequency for avoiding mode mixing problem [42]. Therefore, ranges of decomposition parameters are given by K ∈ [2, K max ] and α ∈ [1/2fs, 2fs].
Step (3). Calculate the kurtosis of sub-components. The signal is decomposed into K sub-components by VMD. Then solve the kurtosis of sub-components separately, and the kurtosis is calculated by where x i represents the discrete data in sub-components, x is the mean value, N is the data length, and σ i is the standard var.
Step (4). Obtain the decomposition parameter set. Under the ith set of decomposition parameters {K i , α i }, the local maximum kurtosis value can be calculated as Y i . Then the parameter set can be constructed

Optimization method of decomposition parameters
In this research, genetic algorithm (GA) was selected to search for the global optimal parameters based on the natural evolutionary processes that enable species to adapt to their environment [43]. In the process of GA optimization, the discrete step length of individual variables has a key impact on the efficiency of parameter optimization, especially for parameter α with a larger range of values. In order to shorten the calculation time and ensure the stability of the optimization parameters, the influence of the step length of α on the optimization result will be analyzed in the research. The optimization process of decomposition parameters is as follows: Step (1). Given the ranges of decomposition parameters and the parameters are binary coded.
Step (2). Select the step length of individual variables and generate initial population.
Step (4). Choose individual with high fitness value and achieve new population.
Step (5). If the termination condition is met, then the last generation is decoded in binary and the optimal results can be obtained. If no, proceed to Step (3).

Sampling strategy and chatter features
Researchers generally need to detect the occurrence of chatter by calculating changes in time-domain and frequency-domain characteristic parameters, and cannot directly judge from the time-domain signal trend, whether it is the original vibration signal or the sub-components obtained by variational mode decomposition (VMD) method. The occurrence of chatter is often accompanied by the increase of signal complexity and the change of energy. In view of good performance of quantifying the complexity of the signal, the energy entropy [41] and the approximate entropy [44], together with the proposed entropy drift coefficient, were used as features for chatter detection.

Sampling strategy
A rolling time window sample strategy proposed by Fu et al. [37] was used in this research. Frame length and frame shift are two key sampling parameters and the size of the frame length mainly affects the calculation time, while the size of the frame shift affects the sampling time. In order to meet the real-time requirements, it is very important to choose the appropriate sampling parameters. The rolling sampling strategy is described as Fig. 1.

Entropy drift coefficient
Using an absolute threshold of entropy would be problematic because the levels of the entropies change with the cutting conditions. In order to detect the changes of entropy quantitatively, a trend evaluation coefficient based on entropy has been proposed in this research, which is called entropy drift coefficient (EDC). The entropy drift coefficient can be calculated as follows: En(i) means the entropy value of the ith sampling segment (frame length), and m is the extension length of the segment. This coefficient can detect whether the entropy value continuously changes in a short time, and then reflect the change characteristics of the vibration signal.

Online adaptive chatter detection strategy based on AVMD
The chatter detection strategy proposed in this paper mainly includes three parts: (1) adaptive adjustment of ranges of decomposition parameters; (2) parameter optimization process using genetic algorithm (GA); (3) online chatter detection. The first two processes have to be conducted in advance because the optimization process is time-consuming, and the offline source signal comes from historical data with similar cutting conditions, which would be updated with the increase of cutting process. The online chatter detection strategy of robotic milling is described in Fig. 2. 3 Chatter characteristics of the robotic milling

Experiment setup
In order to test the effectiveness of the proposed method, a series of end milling tests were conducted. Experiments were performed on a MOTOMAN ES165D robot equipped with a ABG0311A1 spindle. Detailed parameters of machining robot are listed in Table 1. A tri-axial accelerometer (PCB356A15) was mounted on the spindle housing to measure the vibrations in three directions of X, Y, and Z during the milling process. Vibration signals were sampled with NI-9218 and then transmitted to a PC which was used for data storage and signal processing. The sampling frequency was set as fs = 5120Hz. Cutting experiments were conducted by straight cut. The chatter monitoring platform of robotic milling is depicted in Fig. 3.
In this paper, 6061 aluminum alloy is selected as the workpiece to be processed. During the machining process, no cutting fluid was used. The cutter was a solid carbide end mill with three teeth, and the detailed characteristics parameters are listed in Table 2.

Hammer tests
Chatter vibrations are usually caused by the most flexible, dominant structural modes of the machine-workpiece system. To analyze the source of chatter in robotic milling, frequency response functions (FRFs) were obtained by exciting the system with an impact hammer and measuring the vibration response with accelerometers. The results are illustrated in Fig. 4. H ij means the FRF in i direction when striking from the j direction. The highest compliance is observed at around 565 Hz mode in H yy and some lower compliance were found near 1338 Hz and 1805Hz. These high-frequency modes describe the highfrequency natural frequency of the robotic milling system. As shown in Fig. 4b, some low-frequency modes surely exist, while their amplitudes are very small compared with those of high-frequency mode. In order to prevent the failure in measuring low-frequency modes of accelerometers, Marcel Cordes et al. [19] used a laser Doppler vibrometer Polytec to test the FRFs. They found that the low-frequency modes did exist in robotic milling system (14Hz and 26Hz in their research) with higher amplitudes, and proved that the lowfrequency modes are mainly related to the robot's configuration, while the high-frequency modes (850Hz and 2325Hz in their research) belong to the spindle-tool assembly. In order to describe the chatter characteristics of the robotic milling system more completely, this paper will make an analysis based  Fig. 2 The flow diagram of online chatter detection method based on AVMD

Chatter tests in robotic milling
In traditional machining process, chatter vibration is closely related to spindle speed and axial cutting depth according to   Altintas [45]. In order to activate the chatter phenomenon in robotic milling, this research purposely chose spindle speeds n from 3000 to 12,000 rpm, the axial depth of cut from 1 to 6 mm, and the radial depth of cut from 2 to 10 mm. In addition, the same posture as the hammer tests was adopted in the milling tests and the frequency response functions are assumed to be constant because the milling distance is not large (the maximum is 50mm in the y direction). Since this article focuses on chatter detection in robotic milling, it is unreasonable to analyze all groups of signals and eight representative results were picked out. The measured signals and obtained spectrums using fast Fourier transform (FFT) are depicted in Fig. 5, where f represents the spindle rotation frequency. The parameters of these results are listed in Table 3.
In this research, the workpiece factor can be ignored when analyzing the source of chatter because the rigidity of workpiece clamping is far greater than that of robot and tool installation. By removing the rotation frequency of spindle and cutter teeth, abnormal frequency components (indicated by the red arrow and red words) are observed from Fig. 5. Based on the experimental results of hammer tests, abnormal frequencies (1367Hz and 1368Hz) in cases 2 and 4 are close to one of the natural frequencies of the robotic milling system, 1338Hz, obtained from the hammer test. These kinds of vibration are mainly caused by spindle-tool system and the chatter source could be the spindle shaft and slender tool, while abnormal frequencies (49Hz, 55Hz, 62Hz) in cases 6, 7, and 8 are closely related to the low-frequency modes, which are mainly caused by the structural vibration of the robot. In order to verify the effectiveness and feasibility of the proposed method, a simulation signal is constructed and analyzed. The simulation signal is made up of four components: (a) stable cutting signal, which represents the spindle rotation frequency and the angular frequency is 100 Hz; (b) lowfrequency chatter signal, whose angular frequency is 50 Hz; (c) high-frequency chatter signal, whose angular frequency is 300 Hz; (d) wns(t) describes the white noise signal. The simulation signal is constructed as follows:

Adaptive adjustment of ranges about decomposition parameter
The initial decomposition parameters are set as K = 2, α = 1/ 2fs. According to the calculation steps described in Section 2.2, the result of K max is 9. Therefore, the range of decomposition parameters is K ∈ [2, K max ] and α ∈ [1/2fs, 2fs].

Optimizing process using genetic algorithm
The genetic algorithm (GA) is applied to the optimization process of decomposition parameters. The optimization step of K is 1, and the optimization step of α is 0.02fs. The optimization calculation is carried out according to the steps described in Section 2.3. The results show that the optimal

Chatter detection based on proposed features
According to the decomposition parameters obtained above, the simulation signal x(t) is decomposed by variational mode decomposition (VMD). Set the sampling parameters as frame length = 64 (number of data) and frame shift = 32, then select the sub-component with maximum kurtosis value to calculate the approximate entropy, energy entropy, and entropy drift coefficient (EDC) value frame by frame. The results are shown in Fig. 6.
The simulation signal consists of two abnormal frequency components; the first one (300Hz) appears at 1s and the second (50Hz) appears at 2s. As shown in Fig. 6, the approximate entropy and energy entropy show numerical mutation at 1s and 2s, which is synchronized with the time-domain signal. Obvious changes in entropy value show that approximate entropy and energy entropy are sensitive enough to be used for detecting whether there is abnormal chatter component in the signal or not. Huge changes also appear in EDC parameters.
Empirical mode decomposition (EMD) method uses the signal's intrinsic characteristics and does not require the selection of the decomposition layers [33]. As a comparison, simulation signal is decomposed into several sub-components by using EMD and the calculation results of chatter features are illustrated in Fig. 7. The energy entropy obtained by EMD shows obvious mutation when the signal component changes, but the approximate entropy does not have a huge upward or downward trend. Similar results are found in the diagrams of EDC (Fig. 7(c)). The EDC value of energy entropy fluctuates greatly in the middle section (1-2s), which makes it difficult to set a stable threshold. Therefore, the calculation results of the simulated signal show that the VMD-based feature extraction method proposed in this paper is convenient for calculation and can effectively detect the chatter component in the signal.

Results of the AVMD method
To prove the feasibility of the proposed method, three different values of parameter α are selected to calculate K max by proposed adaptive variational mode decomposition (AVMD) method, and the results of eight groups of signals are shown in Fig. 8. K max has an obvious increasing trend with the increase of α. Considering that the range of α is [1/2fs, 2fs], it can be found that the selected K max , through utilizing α = 1/2fs, can reach a minimum value of K max . Thus, it is effective to avoid over-decomposition by using K max calculated from α = 1/2fs. The decomposition results of part signals (2-3s) in case 2 are shown in Fig. 9. The signal in case 2 is decomposed into eight sub-components from low frequency to high frequency, where each frequency band contains clear frequency components and no overlapping frequency component is found in adjacent sub-components.
As a comparison, part signals (2-3s) in case 2 are decomposed into twelve sub-components by using empirical mode decomposition (EMD) method and the results of the first four sub-components are illustrated in Fig. 10. Different from the filtered signal by AVMD, the sub-components decomposed by EMD present overlapping frequencies, which cause the chatter component to be scattered on different subcomponents, resulting in loss of signal energy, and greatly affecting the detection effect of abnormal frequency components. Therefore, proposed method can achieve subcomponents with individual frequency, which is essential to improve the accuracy of chatter detection.

Results of the decomposition parameter optimization
By setting different kinds of discrete step length of individual variable α, the maximum kurtosis value and the decomposition parameters are calculated according to the genetic algorithm (GA). The calculated maximum kurtosis values and decomposition parameters of group 1 signal are depicted in Fig. 11. The optimal parameters of the eight groups of steps are clearly marked in the figure (red points and red texts). The optimal results of K = 4 are calculated under eight kinds of step lengths, while the parameter α fluctuates around 4800 in cases of short step length (0.001fs, 0.002fs, 0.005fs, 0.01fs, 0.02fs, 0.05fs). The results show obvious convergence, that is, the change of step length of α will not affect the optimization result of parameter K, but only parameter α. When the step size reaches 0.1fs, the optimization result of α reaches 2560, which is far away from the optimization results in cases In order to better analyze the influence of discrete step length of parameter α on the parameter optimization results, the parameter K which is not affected by the step length is ignored. Eight kinds of step lengths are used to calculate the optimization parameters of each group of signals respectively and the results are shown in Fig. 12, in consideration of the maximum kurtosis value calculated by using 0.001fs which is better than others. The maximum kurtosis value of 0.001fs is taken as the benchmark. The kurtosis ratios are calculated as follows: It is not difficult to find that the kurtosis ratio is close to one when the step length is less than 0.01fs, and the calculated optimal penalty coefficient α has a slight fluctuation. But the kurtosis ratio calculated from the 2, 3, and 6 groups of signal data will drop significantly when the step length is greater than 0.02fs, and the penalty coefficient will change drastically. The kurtosis values calculated from the 1, 4, 5, 7, and 9 groups of signal data drop significantly when the step length is greater than 0.05fs, and the optimal penalty coefficient α begins to change significantly. Therefore, relatively stable and reliable optimized parameter of α can be obtained when the step length is less than 0.01fs. This discovery will provide a reference for scholars to optimize the decomposition parameters of VMD, and save the calculation time while ensuring the optimization performance.

Verification of real-time performance
It is critical to consider the calculation time of proposed adaptive variational mode decomposition method (AVMD) for an online detection system. A shorter calculation time helps to achieve a faster detection frequency, which can detect chatter as early as possible, so as to control the vibration and damage within an acceptable range [46]. In this section, through analyzing the calculation time of proposed method, a reasonable selection rules of sampling parameters would be provided according to the real-time requirements and calculation time in actual detection situation.
Liu et al. [47] has analyzed the time complexity of variational mode decomposition (VMD), but few scholars have analyzed the performance of VMD in the case of online detection. In the sampling strategy described in Section 2.4, frame length and frame shift are two key sampling parameters and the size of the frame length mainly affects the calculation time, while the size of the frame shift affects the sampling time. As shown in Fig. 13, the boundary line divides the graph into two regions: A and B. T-F curve describes the mean computing time (calculation time of proposed method using data from case 1) under different frame lengths. In order to ensure the continuity and stability of online detection, the selection of sampling parameters must follow two principles: (1) The mean calculation time meets the requirements of online detection time (assuming that the fluctuation of calculation time is very small). (2) The sampling time meets the requirements of online detection time, and in a detection cycle, sampling and calculation are carried out at the same time. The online detection cycle was set as T d ≤ 0.01s, then the first principle determines the upper bound of the frame length. The second principle determines the upper bound of frame shift length. In view that the frame length cannot be less than the sampling length, the sampling parameters must be selected in the semi-part B. Thus, the selection of sampling parameters should be limited in the red line area in Fig. 13. Considering the fluctuation of calculation time, the sampling length is set to 32, so that the sampling time is 0.0063s, less than the required 0.01s. For the frame length, Scholars believe that its value should be large enough to reduce energy loss and ensure signal quality [37], but too large length will lead to a drastic increase in computing time. In this study, the frame length was set as 64, and the mean calculation time was 5.2ms < 10ms.
In order to test the real-time performance of proposed method, some monitoring experiments were carried out on vibration signals under four groups of milling parameters. The monitoring process lasted for 1.5 s, and a total of 239 calculations were performed. The calculation time is shown in Fig. 14 and the red line represents the requirement of real time. For signals with small K, such as case 1, case 3, and case 4, the calculation time is almost less than 0.01s, which can ensure that the monitoring period would not be affected. For signals with large K, such as case 2, calculation time in some segments is more than 0.01s. The occurrence of this situation is related to the timeconsuming VMD process, and calculation result shows that VMD decomposition accounts for 95% of the calculation time. For the signal with large decomposition parameters, the real-time performance of chatter detection can be ensured by reducing the frame length.
As shown in Fig. 15, reducing the frame length can effectively reduce the detection time, and most of the calculations can be controlled within 0.01s. Part of the timeout point does not appear continuously; it will not have a destructive impact on the continuity of detection. Therefore, the proposed method can control the detection time in a certain range and meet the requirements of online detection.

Chatter features analysis
The sampling strategy mentioned in the Section 2.4 was used to divide the signals of each group (frame length = 64 (0.0125s), frame shift = 32 (0.00625s)) and the approximate entropy of each segment was extracted by using proposed adaptive variational mode decomposition (AVMD). As depicted in Figs. 16 and 17, the approximate entropy of case 1 decreases at 0.81s and the spectrum of the signal in this period shows an abnormal frequency of 2476 Hz, which is 0.14s ahead of the violent growth time (0.95s) of the signal amplitude. In addition, abnormal tool marks appeared on the workpiece of group 1, as shown in Figs. 18, 19, 20, and 21. For signal in case 4, the approximate entropy decreases in 0.42s, and the spectrum shows that there is an abnormal frequency of 1368hz with large amplitude. Similar situations also occurred in cases 2, 3, 6, and 8. The large decrease of approximate entropy is often accompanied by the emergence of abnormal frequency components and tool marks, which represents the generation of chatter phenomenon. It is worth mentioning that the lower chatter frequency of 56 Hz appears in the signals of case 7 with the decrease of approximate entropy. In the signals of case 5, there is no continuous decrease or increase of approximate entropy and no abnormal tool marks shown in Fig. 21, which indicates that it is a stable cutting state. The method proposed in this paper can effectively detect the chatter component in the signal, and the detection result is not affected by the chatter frequency. Thus, there is no need to consider whether the chatter is regenerative or coupled, or both, which would perform an important role in promoting chatter suppression in robot machining. In order to compare the advantages and disadvantages of approximate entropy and energy entropy in signal analysis, energy entropy is also extracted according to the same steps as approximate entropy. As shown in Fig. 18, the energy entropy of case 1 decreases at 0.93s, which is almost consistent with the time when the signal amplitude rises sharply, but 0.12s later than the time when the approximate entropy drops. For the signals of case 2, the energy entropy decreases at 0.52s, 0.15s later than the time when the approximate entropy decreases at 0.37s. The times of large changes of the approximate entropy and energy entropy are summarized in Table 4. Compared with energy entropy, approximate entropy is more sensitive and is able to detect the change of signal complexity earlier. Therefore, the approximate entropy was chosen to calculate the entropy drift coefficient (EDC), and the calculation results are shown in Fig. 19.
For the stable cutting vibration signals in case 5, the EDC value is relatively stable. EDC values fluctuate in a certain range, and there is no abnormal extremum point, while extremum points of EDC appeared in those signals with abnormal frequency components (cases 1-4 and 6-8). The appearing times of extremum points are also summarized in Table 4. Relative to the amplitude of the signal, the calculated entropy can detect abnormal components in the signal earlier. Although the amplitude of the signal will also change with the occurrence of chatter, it is obvious that the entropy value can provide an early warning of chatter so as to minimize the damage to the surface quality as much as possible.
As shown in Fig. 16e, the approximate entropy of signals in case 7 decreases in 0.47s, and presents a continuous and significant downward trend. Correspondingly, the entropy drift coefficient (EDC) has an extremum point in 0.50s, as shown in Fig. 19g. Researchers can only observe the downward trend of approximate entropy by naked eyes, but cannot detect it by the value of approximate entropy itself. Fortunately, EDC can reflect the drift trend of approximate entropy by numerical value. Table 4 shows that the occurrence times of EDC extremum points are close to the drift times of approximate Fig. 17 Approximate entropy obtained using AVMD and spectrogram of transitional area: a time series and approximate entropy for case 5; b spectrogram for case 5 (0-1s); c time series and approximate entropy for case 6; d spectrogram for case 6 (0.23-0.36s); e time series and approximate entropy for case 7; f spectrogram for case 7 (0.44-0.54s); g time series and approximate entropy for case 8; h spectrogram for case 8 (0.49-0.59s) entropy, which just verifies that EDC is consistent with approximate entropy in detecting abnormal signal. The emergence of these extremum points verifies that EDC can reflect the drift of entropy and provide conditions for setting detection threshold. Through combining the simulation and experimental results, the absolute value of detection threshold about EDC can be set to 0.2. EDC is calculated based on the approximate entropy value, by calculating the difference of continuous m + 1 points and obtaining a coefficient reflecting the trend of entropy drift (m is set according to the sampling frequency and frame length and it is set as 5 in this research). The larger the EDC, the more away from zero the entropy kept. Although the detection times of approximate entropy in Table 4 is earlier than EDC, it is difficult to determine a stable threshold value of approximate entropy. Therefore, the EDC proposed in this paper is more suitable for online detecting the chatter components in the signal.
For comparison, the signals were also decomposed by using empirical mode decomposition (EMD) and the results of energy entropy obtained from EMD are shown in Fig. 20.
From the energy entropy shown in Fig. 20a, c, d, and g, the trend of the entropy value is similar to the results obtained by AVMD, which can be used to distinguish between chatter state and stable state with difficulty. From the energy entropy shown in Fig. 20d, chatter occurs at about 0.51s, which is later than the onset time indicated by AVMD (Fig. 18d). From the graph shown in Fig. 20b, it is difficult to distinguish between chatter and stable state. Another strange thing is that there is an abrupt change in the energy entropy as shown in Fig. 20f, but it does not exist for a long time. That makes it hard to distinguish whether there is abnormal component or not. This kind of unclear sub-components may be caused by mode mixing in EMD decomposition, which will lead to multiple features in one sub-component, making the chatter information submerged in the mixing signal and difficult to be detected. A similar situation can also be seen in Fig. 20g. Compared with EMD, the filtered signal using AVMD can detect the abnormal components in the signal earlier and obtain subcomponents with more chatter information, which is proved to be more suitable for online detection of chatter in robot machining. Fig. 18 Energy entropy extracted by using AVMD: a energy entropy for case 1; b energy entropy for case 2; c energy entropy for case 3; d energy entropy for case 4; e energy entropy for case 5; f energy entropy for case 6; g energy entropy for case 7; h energy entropy for case 8 Aiming at the online chatter detection problem in robotic machining, a comprehensive solution based on variational mode decomposition (VMD) was put forward in this paper. Firstly, the adaptive variational mode decomposition (AVMD) based on instantaneous frequency and kurtosis can realize the adaptive selection of the decomposition parameter and avoid the phenomenon of mode mixing effectively in VMD process. Compared with the empirical mode decomposition (EMD), proposed method possesses better decomposition performance and can achieve sub-components with individual frequency, which is essential to improve the accuracy of chatter detection. Secondly, genetic algorithm was used to optimize the decomposition parameters of VMD and a selection strategy of discrete step length was given to reduce calculation time. This work will provide a reference for scholars to optimize the decomposition parameters of VMD, and save the calculation time while ensuring the optimization performance. Finally, approximate entropy, energy entropy, and proposed entropy drift coefficient are extracted to detect chatter information in the signal.
Compared with the approximate entropy and energy entropy, proposed entropy drift coefficient (EDC) can not only detect the chatter component in the signal, but also achieve reasonable threshold to distinguish between chatter and stable machining state. It is worth mentioning that this paper also provides the selection strategy of sampling parameters in real-time chatter detection application based on VMD method, which can provide reference for scholars. Simulation and experimental results showed that the methods proposed in this paper can meet the real-time requirements of online detection and detect the occurrence of chatter effectively.
In this research, two main factors affect the performance of online detection, one is the decomposition time of VMD, the other is the threshold setting of chatter feature, although the proposed scheme has realized the real-time monitoring in robot milling. It is still a challenging task for further reduction of computation time and acquisition of stable threshold. In addition, the distinction of chatter types is necessary for better understanding of chatter in robotic machining process. Therefore, the optimization of detection algorithm and distinction of chatter types are our next step work.  Fig. 21 The surface topography: a group 1; b group 2; c group 3; d group 4; e group 5; f group 6; g group 7; h group 8 Data availability The data used or analyzed during the current study are available from the corresponding author on reasonable request.

Declarations
Additional declarations for articles in life science journals that report the results of studies involving humans and/or animals Not applicable.
Ethics approval The authors consciously assure that for the manuscript has not been published and is not under consideration for publication elsewhere.
Consent to participate All the authors consent to participate in this research and contribute to the research.
Consent for publication All the authors consent to publish the research. There are no potential copyright/plagiarism issues involved in this research.

Competing interests
The authors declare no competing interests.