Intermittency and multiscale dynamics in milling of fiber reinforced composites

We have analyzed the variations in cutting force during milling of a fiber-reinforced composite material. In particular, we have investigated the multiscale dynamics of the cutting force measured at different spindle speeds using multifractals and wavelets. The multifractal analysis revealed the changes in complexity with varying spindle speeds. The wavelet analysis identified the coexistence of important periodicities related to the natural frequency of the system and its multiple harmonics. Their nonlinear superposition leads to the specific intermittent behavior. The workpiece used in the experiment was prepared from an epoxy-polymer matrix composite reinforced by carbon fibers.


Introduction
Milling is one of the most common machining operations performed in the manufacturing industry. It is often used as a process for material removal, edge finishing and other functions. The cutting force in a milling operation evolves on multiple timescales and exhibits complex dynamics. Over the past few years, new technological development in milling provided a reliable high-speed cutting procedure. However, in spite of recent progress in the understanding of the nonlinear mechanisms leading to vibrations [1][2][3][4][5], the strategies for controlling conditions for stable machining are not clearly understood. Consequently, it would of great interest to gain a deeper understanding of the complex dynamics of the milling process. Some progress has been made in this direction using the adaptive control concept, based on identification of relatively short time series [6,7].
Fiber-reinforced composites, in view of their high specific strength and stiffness, are now widely used in various industrial applications. However, due to material discontinuity, non-homogeneity and anisotropy, machining of composites is more challenging than machining of simple metals and their alloys. Because of various possible damage mechanisms such as fiber pullout, fiber fragmentation and delamination, matrix burning, and matrix cracking, poor surface quality can occur [8,9]. Workpiece-tool vibrations appearing during machining increase considerably the temperature of contact, and this effect cannot be minimized by cooling fluid as the material can easily absorb it. Recently, Rusinek [10] and Litak et al. [11] investigated the cutting dynamics in milling of a fiber-reinforced composite material using nonlinear time series analysis techniques. In this study, we analyze the dynamics of the cutting force variations in milling of fiberreinforced composites using wavelets and multifractals.
The present paper is composed of five sections. Following this section which provides an introduction to the main topics of the paper, we describe the experimental set up and the measurement procedure in Sect. 2. Section 3 is devoted to the multifractal approach, whereas Sect. 4 presents the wavelet analysis. Finally in Sect. 5, a few concluding remarks are given.

Experimental set up and force measurement
The cutting force measurements were performed using the experimental set up shown schematically in Fig. 1a. The set up (Fig. 1b) consists of a CNC milling machine, a piezoelectric dynamometer, a charge amplifier, a sample-and-hold unit, an analog-to-digital (A/D) converter, and a computer. The force signals are transmitted from the dynamometer to the A/D converter and stored in the computer. Machining of epoxy-polymer matrix composite reinforced by carbon fibers (EPMC) was performed at various spindle speeds ranging from 2000 rpm to 8000 rpm, with the feed rate fixed at 520 mm/min with a depth of cut equal to 0.5 mm. The mill itself is made of diamond-coated steel with a diameter of 12 mm. The largest force component F x (oriented in the direction of machining) was measured. The time series of this cutting force at five different speeds ranging from 2000 to 8000 rpm are presented in Fig. 2. The total time of recording is 5 seconds and the sampling frequency is 4000 Hz. At a first glance, one can observe that the amplitude of oscillations is initially growing with increasing speed to reach the highest value at N = 6500 rpm; beyond

Multifractal analysis
A characteristic feature of the complex dynamics of the milling process is that they occur on multiple time scales. A convenient way to describe the dynamics of such multiscale processes is to use a multifractal formalism. The multifractal approach is based on the spectrum of Hölder exponents which can be used as a measure of complexity [12][13][14]. A fractal (or a monofractal) process is self-similar in the sense that its dynamics can be described in terms of a single powerlaw scaling exponent such as the Hurst exponent [15], and may be considered a homogeneous process. Accordingly, its complexity can be described by means of a single fractal dimension. In contrast, a multifractal process is heterogeneous and evolves on different time scales with different scaling exponents. It is therefore necessary to use several scaling exponents or fractal dimensions to describe the multiscale features of a multifractal process. This can be done by calculating the singularity spectrum in terms of the so-called Hölder exponent. For an infinitely long monofractal process, the singularity spectrum reduces to a single point. On the other hand, if the singularity spectrum does not reduce to a single point, it is indicative of multifractal behavior. A multifractal process may be considered to be locally self-similar and the Hölder exponent may be treated as a local Hurst exponent.
The broadness of the singularity spectrum is a measure of complexity of the multifractal process. Multifractal processes are known to occur in a wide variety of applications (see, for example, [16][17][18]).
Consider a real-valued function F (t). The Hölder exponent of this function at a point t = t 0 is defined as follows [13]: where P n (t) is a polynomial of degree less than α, and C is a constant. Note that for a function to be differentiable at the point t 0 , the value of α must be 1 at that point. On the other hand, if α < 1 at some point t 0 , the function will not be differentiable there and α will characterize the strength of the singularity at that point. Therefore the Hölder exponent of a function is a measure of the smoothness (or lack thereof) of a function at a point. The singularity spectrum, f (α), is defined as the Hausdorff dimension of a set of singularities of strength α.
We use the following two attributes of a multifractal spectrum: (i) the value of the Hölder exponent, α = α 0 , corresponding to the spectral peak, and (ii) the broadness, α, which is the distance between the (extrapolated) points of intersection of the spectral curve with the α-axis. The parameter, α 0 , represents the most dominant fractal exponent, and it reflects the degree of persistence or correlation in a time series. In particular, the value α 0 = 0 corresponds to Gaussian white noise, α 0 = 0.5 to Brownian walk, both of which indicate an uncorrelated process. On the other hand, the values of α 0 < 0.5, and α 0 > 0.5 indicate anti-persistent and persistent walks, implying positive and negative correlations, respectively, between the events in the time series [16][17][18]. The broadness of a singularity spectrum describes the range of possible fractal exponents and thus gives a measure of multifractality or complexity of the time series. A large value of broadness describes a richer multifractal structure whereas a small value approaches a monofractal limit.
Using the software provided in [12], we performed a multifractal analysis of the milling force (F x ) variations represented by the time series shown in Fig. 2.
The singularity spectra, f (α), for all five speeds are depicted in Fig. 3. From this figure, the values of the most probable Hölder exponent, α 0 related to the correlation of the time series, and the broadness of the singularity spectrum, α related to the degree of complexity, were obtained. These values are listed in Ta-  (Table 1) are not changing monotonically with increasing N . In particular, the average value is changing considerably and reaching the maximum of F x = 12.34 N for N = 5000 rpm, whereas the standard deviation is changing around 22 N. Notably, the case N = 5000 rpm is accompanied by the local increases of both multifractal parameters: correlation α 0 and complexity α. The increasing force could be related to the specific properties of the machined material as the distribution of carbon fibers may influence the resistance of the material.

Wavelet analysis
Wavelets have been used for time series analysis in a wide variety of applications. Wavelet analysis provides a spectral-temporal approach to identify the dominant modes of variability in a time series and to delineate how these modes vary over time. It is particularly useful for analyzing transient and intermit-  . 4 Tendencies in the change of multifractal parameters α 0 (a) and α, (b) versus N tent phenomena. A wavelet-based approach has advantages over the more traditional methods such as the Fourier transform or the windowed Fourier transform. The Fourier transform is a purely frequency domain technique which seeks to determine the periodicities in a signal through spectral peaks, but it cannot delineate the time spans over which the periodicities may persist. A windowed Fourier transform also known as a short-time Fourier transform (STFT) circumvents this limitation by applying the Fourier transform on a short segment of the signal using a window of fixed size and then sliding the window in time. The temporal variations of the periodicities, if any, can thus be determined. However, because a fixed-size window is used in STFT, the frequency resolution as well as the time resolution is fixed. As a consequence, for a given signal either the frequency resolution may be poor or the time localization may be less precise, depending on the size of the chosen window. In contrast, using variablesize windows, wavelet analysis provides an elegant way of adjusting the time and frequency resolutions in an adaptive fashion. A wavelet transform uses a window that narrows when focusing on small-scale or high-frequency features of the signal, and widens on large-scale or low-frequency features, analogous to a zoom lens [19]. Recently we have used wavelet anal-ysis in our studies of pressure fluctuations in internal combustion engines [16], and other applications [20][21][22] including the turning process [17]. We present below a brief description of the wavelet analysis methodology and then apply it to the cutting force time series shown in Fig. 2. A wavelet is a small wave with zero mean and finite energy. Consider a time series {x i } with i = 1, 2, 3, . . . , N. The continuous wavelet transform (CWT) of this time series with respect to a wavelet ψ(t) is given by the convolution of the time series with a scaled and translated version of ψ(t). The function ψ(t) is referred to as an analyzing wavelet or a mother wavelet. The convolution is expressed by [23] Here δt is the sampling interval, and an asterisk on ψ denotes its complex conjugate. The symbols s and n are called scale and time index, respectively. The scale parameter controls the dilation (s > 1) and contraction (s < 1) of the mother wavelet. The time index, n, indicates the location of the wavelet in time; in other words, as n varies, the signal is analyzed in the vicinity of this point. The amount of signal energy contained at a specific scale s and location n is given by the squared modulus of the CWT, and is referred to as the wavelet power spectrum |W n (s)| 2 . The wavelet power spectrum (WPS) is a measure of the variance at different scales or frequencies. The WPS which depends on both scale and time is represented by a surface. By taking contours of this surface and plotting them on a plane, a time-scale representation of the spectrum may be derived. A time-scale representation is found to be useful for extracting important features of signals arising in many applications. In our analysis we used a complex Morlet wavelet as the mother wavelet. A complex Morlet wavelet consists of a plane wave modulated by a Gaussian function and is described by [23] ψ(η) = π −1/4 e iω 0 η e −η 2 /2 .
Here ω 0 = 2πf 0 is the order of the wavelet, with f 0 being the center frequency. The value of ω 0 controls the number of oscillations that is present in the mother wavelet and thus influences the frequency and time resolutions of the corresponding wavelet transform. A larger value of ω 0 provides a higher frequency resolution whereas a smaller value improves the time resolution. In our computations, we have used a Morlet wavelet of order 6 as the mother wavelet. This choice provides a good balance between time and frequency localizations. In addition, for this choice, the scale is approximately equal to the Fourier period and thus the terms scale and period can be used interchangeably for interpreting the results. Using a continuous wavelet transform (CWT), we calculated the wavelet power spectra (WPS) of the various force F x time series (depicted in Fig. 2). The results are illustrated in Fig. 3. In this figure, the horizontal axis represents the number related to the sampling points, and the vertical axis denotes the period in sampling intervals. The colors red and blue represent the highest and lowest power levels, respectively, with the other colors denoting intermediate power levels. Figures 5a-c illustrate the results for particular cases of spindle speed N = 3500, 5000, and 6500 rpm. The main harmonics represented by dark red color are of about 69, 48, and 37 (measured in terms of the sampling interval) in Figs. 5ac, respectively. They are related to the corresponding spindle speeds. In every case we observe also the multiple of that frequency but interestingly the periodicity composed of about 10 sampling intervals is a fixed and present in all power spectra. In Fig. 5a (N = 3500 rpm), one can clearly see an intermittent (as in [24,25]) cutting behavior After increasing the speed, in Fig. 5b (N = 5000 rpm) the process becomes weaker and more homogeneous. In the last case (Fig. 5c, N = 6500 rpm), we observe that the intensity increases are grouped in a periodic pattern Additionally, Fig. 5c shows intermittent increases of intensity for very short periods of 2 sampling intervals. The above analysis clearly indicates that at a certain speed (N ≈ 5000 rpm), the dynamics of a milling process is changing. A competition between two time scales of system response leads to more complex behavior which can directly influence the final surface of the designed product.

Concluding remarks
We have examined the dynamics of cutting force variations in milling of a fiber-reinforced composite material by analyzing the experimental time series of the largest force component (F x ) using wavelets, and multifractals. The wavelet analysis revealed that the short or high-frequency periodicities of 2 and 10 sampling periods are intermittent. Interestingly, the largest detected intermittencies coincided with the maximum of fluctuations measured by the standard deviation. The above conclusions have been confirmed also by the multifractal measures which showed that there is a noticeable deviation of the changing trend in the correlation α 0 and complexity measure α at the particular spindle speed (ω = 3500 rpm). This peculiar behavior was formerly investigated by the multiscaled entropy showing its nonhomogeneous behavior against the scaling and similarity factors which confirmed the main changes in short period (intermittencies) as well as long-period (modulation) changes [11].
Note, in the standard approach to milling one can follow the uniform material which could fit the mathe-matical model and the corresponding regions of stability [3,[26][27][28][29][30]. On the other hand, the results obtained here, for a composite material, can be used further to study the stability of the milling process and to develop new methods of control. To get further insight, it would be of interest to analyze the cutting force variations during milling of different types of composite materials. Results of such analysis will be reported in a future publication.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.