NonlinearDynamics ofCuttingProcess consideringHigher-Order Deformation of Composite Cutting Tool

In this paper, the nonlinear dynamic analysis of the cutting process of composite cutting tool is performed. (e cutting tool is simplified to a nonplanar bending rotating shaft. (e higher-order bending deformation, structural damping, and gyroscopic effect of cutting tool are considered. It is assumed that cutting tool is subjected to a regenerative two-dimensional cutting force containing the first and second harmonic components. Based on the Hamilton principle, the motion equation of nonlinear chatter of the cutting system is derived. (e nonlinear ordinary differential equations in the generalized coordinates are obtained by Galerkin method. In order to analyze the nonlinear dynamic response of cutting process, the multiscale method is used to derive the analytical approximate solution of the forced response for the cutting system under periodic cutting forces. In the forced response analysis, four cases including primary resonance and superharmonic resonance, i.e., Ω � ω1, Ω � ω2, 2Ω � ω1, and 2Ω � ω2, are considered.(e influences of ratio of length to diameter, structural damping, cutting force, and ply angle on primary resonance and superharmonic resonance are investigated. (e results show that nonlinearity due to higher-order bending deformation significantly affects the dynamic behavior of the milling process and that the effective nonlinearity of the cutting system is of hard type. Multivalued resonance curves and jump phenomenon are presented.(e influences of various factors, such as ratio of length to diameter, ply angle, structural damping, cutting force, and rotating speed, are thoroughly discussed.


Introduction
As a high-efficiency, high-quality, low-cost machining method, high-speed cutting technology has been widely used in aerospace and mold manufacturing. However, chatter reduces the cutting stability during machining operations, causes a decrease in the machining quality and cutting efficiency of the workpiece, damages the tool, and shortens the service life of the machine. e passive chatter control methods are mainly based on various types of dynamic vibration absorbers [1,2] and impact dampers [3]. However, composite materials are known to have higher static stiffness, damping, and specific stiffness. It has been found that the dynamic stiffness and fundamental natural frequency of cutter bar may be improved simultaneously if composite is employed [4][5][6]. is is very beneficial for the stability of high rotational speed machining for deep holes.
Within the framework of linear chatter theories and based on different theories, several of chatter phenomena as well as their stability boundary in different cutting processes were discovered [7][8][9][10][11][12][13][14]. Tobias [8] introduced the time-delay instability in the cutting force and proposed the regeneration theory, which is considered to provide a complete explanation of the chatter phenomenon so far.
However, because linear theory cannot predict some important dynamic behaviors of cutting process, nonlinear modeling of cutting systems has received more attention. e nonlinearity in cutting systems is mainly caused by structural nonlinearity, delayed nonlinearity of cutting force, etc. [15][16][17][18]. Hanna and Tobias [15] first proposed a timedelay nonlinear model with quadratic and cubic structural stiffness and cutting forces, which has inspired great interest in analyzing the global dynamics of the cutting system. e effective mathematical methods for the nonlinear theory of cutting systems include center manifold theory, bifurcation theory, perturbation analysis, phase portraits, and Poincaré sections.
Pratt [19] used the multiscale method, harmonic balance method, and Floquet theory to study the models of Hanna and Tobias [18]. e results showed that subcritical Hopf bifurcation may occur due to the existence of cubic structural nonlinearity. Moradi et al. [20] studied the existence of different types of bifurcation in the cutting process considering the tool wear and process damping. ey developed a 2-dof linear model of the tool and a polynomial nonlinear model of the cutting force. In addition, the multiscale method was used to obtain the approximation analysis solution of primary resonance. Moradi et al. [21] used a model similar to that in literature [20] to analyze the forced vibration of the milling process. In the study, not only the primary resonance but also the superharmonic resonance and the internal resonance were investigated. Moradi et al. [22] investigated the internal resonance and regenerative chatter of the milling process considering both the cutting force and the structural nonlinearity.
However, in the above studies, the tool structure was simulated with 1-dof or 2-dof vibration system, in which the concentration parameters such as mass and stiffness were obtained through experiments or empirical method. ere is no obvious correlation between these simplified models and the continuous system equations of the tool.
In order to investigate the dynamic characteristics of cutting system and the stability mechanism of machining process, it is necessary to conduct a comprehensive analysis on various influencing factors. In this case, if the empirical method is used to establish the model of the tool structure, a large number of repetitive tests are needed to obtain the dynamic parameters of the tool structure of different sizes, geometries, and materials, which is very time consuming and of low effectiveness. erefore, as a more effective modeling method for cutting tool, the continuous parameter dynamics modeling of cutting tool based on analytical method, has received great attention [23][24][25][26][27][28][29][30][31][32][33].
However, in most of the existing continuous system dynamics models, the influence of the nonlinearity of the tool structure has not been considered. erefore, the existing dynamic model is only suitable for analyzing the linear dynamics of the natural frequency, chatter frequency, dynamic deformation, and chatter stability of the cutting system. In addition, the current continuous system dynamic studies are mainly focused on the tool structure made of isotropic metal materials. Nevertheless, there is no research about the nonlinear dynamics of the process with anisotropy composite tool structures. e structural nonlinearity is essentially arisen by the flexible nature of the cutting tools with small diameter or long extension part. e structural nonlinearity of the cutting tools may be described using higher-order large deformations [34] or geometric nonlinearity [35].
In this paper, the nonlinear dynamics of the cutting process with a composite cutting tool are investigated. e composite tool structure is simplified to a nonplanar cantilever rotating shaft. e structural nonlinearity is introduced by considering the higher-order deformation of cutting tool. It is assumed that the cutting tool is subjected to a regenerative cutting force containing harmonic components. Based on the Hamilton principle, the nonlinear chatter equation of the cutting system is derived. e Galerkin method is used to obtain the nonlinear ordinary differential equations in the generalized coordinates. e multiscale method is used to obtain the approximate solution of the forced vibration response of the cutting process subjected to periodic cutting forces. Nonlinear dynamics of the system are studied for four cases of primary and superharmonic resonances; i.e., Ω � ω 1 , Ω � ω 2 , 2Ω � ω 1 , and 2Ω � ω 2 are studied. e numerical calculation is conducted to investigate the effect of various parameters on the frequency response of the cutting system.

Kinetic Energy and Potential Energy.
e structural sketch of composite cutting tool with nonplanar bending is shown in Figure 1. e kinetic energy of the rotating cutting tool without considering torsional deformation is as follows: where v and w represent the displacement at any point on the neutral axis along the y and z directions, respectively. ψ z and ψ y represent bending rotation angles of the cross section around the y and z axes, respectively. ρA represents the mass per unit length, and ρI represents the mass moment of inertia of the cross section. e superposed dots represent derivatives with respect to the time t. Using the expression of the displacement fields for Euler-Bernoulli beam and the stress-strain displacement relations of the composite cutting tool, as shown in Appendix A, one can obtain the following expression for the potential energy of the rotating composite cutting tool with the higher-order deformation: where c is structural damping coefficient. e virtual work of the cutting force δW can be expressed as where is the delta function. e linear regenerative cutting forces F y and F z can be obtained by [20] where α 0 � 0.5ς 1 + 0.25πη 1 , β 0 � 0.5η 1 + 0.25πς 1 , in which τ � 2π/NΩ is the delay time of milling process; N is the number of teeth; (K tc , K te ) and (K rc , K re ) are the cutting force coefficients in the tangential and radial directions, respectively; a is the axial cutting depth; and c f is the feed per tooth per revolution.

Dynamic Model of Milling Process.
In order to obtain the equation of motion of the cutting system, the principle of Hamilton is used as follows: Substituting (1), (2), (4), and (5) into (8), the equations of motion in both y and z directions can be obtained: Here, (EI) equiv and (EA) equiv represent the equivalent bending and tensile stiffness of the cutting tool, respectively; (ρI) equiv and (m) equiv represent the equivalent diametrical mass moment of inertia and the equivalent mass per unit length, respectively; ρ k represents the density of the kth layer; and r k and r k+1 represent the inner diameter and outer diameter of the kth layer, respectively. e detailed expression of Q 11 can be found in Appendix B. In (9), the primes denote differentiation with respect to x. e following nondimensional quantities are defined: Using the above dimensionless variables, (9) is rewritten as where In (12), the primes denote differentiation with respect to x, and the superposed dots denote derivatives with respect to the time t. e solution of (12) can be written as For cantilever beams, ϕ 1 (x) has to meet the following boundary conditions: e mode function that satisfies the boundary condition in (12) can be expressed as follows: where β 1 L � 1.875. Substituting (14) into (9), simplifying the equation by the Galerkin method, and dropping the constant cutting force terms, the following ordinary differential equations can be obtained. where Here, the superscript (4) means the fourth-order partial derivative of x.

Perturbation Analysis of Milling Process Using Multiple
Scales Method. In order to solve (17) using the multiscale method, the structural damping and nonlinear and exciting force terms are scaled using small parameters ε as follows: V(t) and W(t) are expanded in the form where Taking the derivative of (20), one obtains the following equations: Substituting (20) and (21) into (17) and considering (19), one can obtain the following equations by comparing the coefficients of ε and ε 3 : O(ε): Shock and Vibration 5 O(ε 3 ): Assume the solution of (22) is as follows: represents the imaginary unit, F 1 (T 2 ) and F 2 (T 2 ) are the complex functions to be determined, and F 1 (T 2 ) and F 2 (T 2 ) represent the complex conjugate. ω 1 and ω 2 refer to the forward and backward whirling frequencies, respectively, which are expressed as follows: Substituting (24) into (23) and dropping the constant terms λ V and λ W , one obtains the following equations: where In this paper, four cases of primary and superharmonic resonances, i.e., case I: Ω � ω 1 , case II: Ω � ω 2 , case III: 2Ω � ω 1 , and case IV: 2Ω � ω 2 , are studied.

Case I: Primary Resonance
where σ represents the detuning parameter. Substituting (24) into (26), one can obtain the following equations: where where q 1 � (ς 2 − iη 2 )/2. e particular solutions of (29) are Substituting (31) into (29) and equating the coefficient of e iω 1 T 0 and e iω 2 T 0 in both sides of (29), one has Equations (32) and (33) After simplification, the solvability conditions are reduced to two independent governing equations in the following form: Assume that the solutions of (36) and (37) are where a j (T 2 ) and θ j (T 2 ) (j � 1,2) refer to the amplitudes and phase angles of the response, respectively. Substituting (38) into (36) and (37), separating the real and the imaginary parts, one can obtain the following equations: Shock and Vibration 7 where To determine the steady-state forced response, the time derivatives in (39) and (40) are set to zero. It can be immediately concluded from (40) that only solution for a 2 is zero. is shows that under the primary resonance Ω ≈ ω 1 , only the first mode, i.e., the forward mode, can be excited, while the second mode, i.e., the backward mode, does not participate in the primary resonance and remains stationary.
Substituting a 2 � 0 into (39) and solving for σ, one obtains To study the stability of the steady-state response of case I, the nature of steady-state response can be investigated numerically by linearizing (39) (with a 2 � 0) around (a 0 , ψ 0 ).

Case II: Primary
Resonance (Ω � ω 2 ). According to (26), under such conditions, one introduces detuning parameters as Using their solvability conditions (which are not shown) and (38) leads to By using steady-state condition (i.e., set a 1 ′ � 0, ψ 1 ′ � 0, a 2 ′ � 0, and ψ 2 ′ � 0), the following steady-state response of case II can be obtained: From (44), it can be seen that a 1 � 0 and a 2 ≠ 0, indicating that only the forced vibration response of backward mode is excited, while forward mode remains stationary. (2Ω � ω 1 ). For this case, the formulation of steady-state response is the same as case I, but in (41), ς 2 and η 2 are replaced with c fς 1 /2 and c fη 1 /2, respectively. (2Ω � ω 2 ). For this case, the formulation of steady-state response is the same as case II, but in (43), ς 2 and η 2 are replaced with c fς 1 /2 and c fη 1 /2, respectively.

Numerical Results and Discussions
In this study, the composite material of carbon fiber/epoxy resin is selected as the material of the cutting tool. e mechanical properties of the material are shown in Table 1. Coefficients of cutting forces for simulation are given in Table 2. e cutting tool has a hollow structure, the outer diameter of the cross section is D � 8 mm, the inner diameter is D � 4 mm, the thickness of the section is h � 2 mm, and the length L is determined by the given ratio of length to diameter.
e composite cutting tool has 16 layers with identical thickness, and the stack sequence is [ ± θ] 8 . In all cases, Ω � 200, except where other values are mentioned. Figure 2 shows the natural frequency versus rotating speed, which is generally known as a Campbell diagram. In vibration of gyroscopic systems, there are two natural frequencies associated with forward and backward whirling motions. In forward natural frequency, the natural frequency is measured when the rotating cutting tool whirls in direction of the rotation. However, in backward natural frequency, the natural frequency is measured when the cutting tool whirls in the opposite direction of the rotation. e forward natural frequency (black solid line) increases with the increase of the rotating speed, while the backward natural frequency (blue dashed line) decreases with the increase of the rotating speed. e intersections of the curves related to the natural frequencies with the straight line ω n � Ω determine the critical rotating speeds of the rotating cutting tool.

Stability Lobe Diagram.
By removing the nonlinear term and the harmonic cutting force in the right-side term of (17), the stability of the cutting system can be investigated. Figure 3 shows a stability lobe diagram where the Ω versus a lim curve separates the space into two regions. Any (Ω, a lim ) pair that appears above the collective boundary indicates unstable milling process where regenerative chatter or selfexcited vibration occurs, while any pair below the boundary is a stable milling process.

Case I.
e numerical solutions for the forced vibration responses of the cutting system with composite cutting tool are presented in Figure 4. As shown in Figure 4, by increasing detuning parameters at point A, the amplitude a 1 gradually increases until reaching point B. a 1 jumps downward from point B to point C. Afterward, a 1 gradually drops to point D while continuing increasing the detuning parameter. On the other hand, as the detuning parameter decreases from point D to point E, a 1 increases and jumps upward from point E to point F. As the detuning parameter further decreases, a 1 drops until arriving at point A. By considering the nonlinearity of higher-order bending deformation, the forced response curve of the cutting system deviates toward the right, suggesting hard spring vibration  Table 2: Coefficients of cutting forces [21]. behavior of Duffing type oscillator. As the detuning parameter varies between σ F and σ B , the cutting system has three steady-state vibration amplitudes, p 1 , p 2 , and p 3 , among which p 1 and p 3 are stable and p 2 is unstable. is means that the cutting system works in an unstable cutting state. erefore, cutting conditions and consequently initial conditions should be adjusted such that the stable branch with less vibration amplitudes is realized in practice. Figures 5-9 show the effects of ratio of length to diameter, structural damping, cutting force, ply angle, and rotating speed on frequency response curve, respectively, for the case Ω � ω 1 . As shown in Figure 5, when ratio of length to diameter increases, vibration amplitudes decrease and frequency response curve bends more strongly toward left.
is is physically expected, because the nonlinear stiffness coefficient λ � (EA) equiv L 2 /2(EI) equiv is proportional to ratio of length to diameter according to (12). As is observed and physically expected from Figure 6, the increase in the structural damping leads to the decrease in vibration amplitudes. In order to study the effect of cutting force coefficients, we define ς 2 � K f ς 2 , η 2 � K f η 2 , with the parameters ς 2 and η 2 given in Table 2. As is shown in Figure 7, by increasing cutting forces, vibration amplitudes increase. Figure 8 shows that the primary resonance amplitude increases with the increase of the ply angle because the elastic modulus E 11 along the longitudinal direction of the fiber is significantly larger than the transverse elastic modulus E 22 (as shown in Table 1). erefore, as the ply angle is greater, the equivalent bending stiffness (EI) equiv of the composite cutting tool is smaller; thus, the nondimensional cutting force coefficients are greater (as shown in (11)). As a result, the amplitude of the primary resonance response is larger. Figure 9 shows that vibration amplitudes decrease with the decrease of rotating speed. is is physically expected, because the equivalent damping of the cutting system (which will be introduced in the subsequent section) increases with the decrease of rotating speed due to the damping effect from regenerative chatter mechanism. Figures 10-14 show the amplitude versus damping coefficient with different ratios of length to diameter, detuning parameter values, cutting force coefficients, ply angles, and rotating speeds, respectively, for the case Ω � ω 1 . From these figures, it can be seen that for some values of L/d, σ, K f , or θ, there are multivalued curves. For example, when L/d � 10 and c ＜0.45, as is shown in Figure 10, the system has two stable and one unstable branches, but for L/d � 10 and c ＞0.45, there exists only   Figures 15-19 show the amplitude versus ply angle with different ratios of length to diameter, detuning parameter values, cutting force coefficients, damping coefficients, and rotating speeds, respectively (Ω � ω 1 ). Similar to the cases in Figure 10～14, again, for some values of L/d, σ, K f , or c, multivalued curves can be observed. Figures 20-24 show the amplitude versus cutting force coefficient with different ratios of length to diameter, detuning parameter values, ply angle, damping coefficients, and rotating speeds, respectively. Again, when cutting force coefficient is less than the specified value, multivalued curves can be obtained for some values of L/d, σ, θ, c, or Ω. Figures 25-27 show the effects of ply sequences (see Table 3) on the frequency response curve, the amplitude versus dumping coefficient, and the amplitude versus cutting force coefficient, respectively. It can be seen from Figure 25 that the configuration A leads to the largest vibration amplitude, while the configuration C leads to the smallest vibration amplitude. is phenomenon can also be found in Figure 8 where the vibration amplitude increases with the Stable Unstable   ply angle when using the same configuration. Again, it can be seen from Figures 26 and 27 that multivalued curves can be observed for some values of L/d, σ, K f , or c. Figures 28-33 show the effects of ratio of length to diameter, structural damping, cutting force, ply angle, ply sequence, and rotating speed on frequency response curve, respectively, for the case Ω � ω 2 . As is observed, the frequency responses and effect of various parameters for the cutting system are similar to those in case I; as ratio of length to diameter and structural damping increase, vibration amplitudes decrease. Furthermore, as cutting force coefficient, ply angle, and rotating speed increase, vibration amplitudes increase. However, for this case, vibration amplitudes are less than those of case I. It should be noted that the equivalent damping of the cutting system is composed of the structural damping and the damping from regenerative chatter mechanism. e equivalent damping coefficient of case I can be expressed as

Case II.
In addition, the equivalent damping coefficient of case II can be expressed as  Stable Unstable equivalent damping coefficients of case I are always larger than those of case II. is explains the reason why vibration amplitudes of case I are larger than those of case II. It can be also seen that the equivalent damping coefficients increase with ratio of length to diameter, ply angle, cutting force, and structural damping but decrease with rotating speed. erefore, the increase of rotating speed leads to large vibration amplitudes, as shown in Figures 9 and 32. When rotating speed approaches infinity, the equivalent damping coefficient approaches the structural damping c (see Figure 38). Moreover, for some values of L/d, σ, K f , c, θ, or Ω, multivalued solutions can be found from the amplitude versus damping coefficient with different ratios of length to diameter, detuning parameter values, cutting force coefficients, ply angles, and rotating speeds, respectively (not shown).   Figure 39 shows the frequency response curve for the four cases. As is observed, the curves of cases I and III are very similar, except that the amplitude in case I is larger than that in case III. Likewise, the curves of cases II and IV are very similar, except that the amplitude in case II is larger than that in case IV.

Case III.
For this case, the equivalent coefficient is identical with that in case I, but the amplitude of excitation force (c 2 fς 1 + c 2 fη 1 )/4 in this case is lower than (ς 2 2 + η 2 2 ) in case I (e.g., using the same parameters L/d � 10, θ � 90°, and K f � 3, the amplitudes of excitation force are 0.2695 and 0.2532 for cases I and III, respectively). erefore, generally,  Stable Unstable the behavior of nonlinear forced vibration in this case (see Figure 39) is similar to case I, despite having less vibration amplitudes compared with case I.

Case IV.
For this case, the equivalent coefficient is identical with that in case II, and the amplitude of excitation force in this case is lower than that in in case II. us, under the same conditions, vibration amplitudes are generally lower than those in case II. e results of case IV are not shown to save space.
In order to validate the calculated approximate solution using multiscale method, the numerical simulation results according to (17) are also displayed in Figures 40 and 41. ese two calculated results using different methods are in good consistency.

18
Shock and Vibration

Conclusions
In this paper, a nonlinear dynamic model of the milling process with composite rotational cutting tool is presented. e cutting tool is simplified to a nonplanar bending Euler-Bernoulli beam, and the structural nonlinearity is attributed to the higher-order bending deformation of the cutting tool. A linear cutting force model considering regenerative mechanism, which includes time-delay terms and the first and second harmonic terms, is used. In addition, structural damping and gyroscopic effect are also considered. Nonlinear ordinary differential equations in the generalized coordinates are derived using the Hamilton principle and the Galerkin method. e lobe diagram of the cutting system is obtained. e multiscale method is used to construct the analytical approximate solutions of the forced vibration frequency response. It is found that, for all cases of primary resonance and superharmonic resonance, excitation of the forward (backward) mode does not produce the vibration responses in the backward (forward) mode because of not taking into account internal resonance in the proposed model. For primary and superharmonic resonance conditions, the effective nonlinearity of the cutting system with higher-order bending deformation rotating composite cutting tool is of a hard type. Jump phenomenon and multivalued solutions can be observed. e influence of higher-order bending deformation causes the frequency response curve to bend more strongly toward right when ratio of length to diameter increases. e study also finds that the vibration amplitudes of the cutting tool increase with cutting force or ply angle. It can be seen that the damping of the cutting system has a significant influence on the vibration amplitude of the composite cutting tool. e damping capacities of the cutting system can be measured by the two different equivalent damping coefficients c equiv1 and c equiv2 which include the structural damping and the damping from regenerative chatter mechanism. Cases I and III have the same equivalent damping coefficient c equiv1 , and cases II and IV have the same equivalent damping coefficient c equiv2 . It has been found that c equiv2 is larger than c equiv1 , so the vibration amplitudes of case I are larger than those of case II. In addition, for all resonant cases, as rotating speed decreases or structural damping increases, the equivalent damping coefficient is strengthened, and consequently less vibration amplitudes are observed. erefore, from this point of view, low rotating speeds are preferable to keep vibration amplitudes small. [38].

T:
Kinetic energy L: Length of the cutting tool t: Time ρ: Density A: Cross-sectional area I: Second moment of area of the beam cross section v, w: Transverse deflections of the shaft element x, y, z: Variable coordinates ψ y , ψ z : Rotation angles of the cross section around the y and z axes Ω: Rotating speed of the shaft Q 11 : Off-axis stiffness coefficient of the kth layer for the composite cutting tool U: Potential energy W c : Rayleigh dissipative energy function of the cutting tool c: Structural damping coefficient δ: Variational operator

20
Shock and Vibration δW: Virtual work of the cutting force F y, F z : Linear regenerative cutting forces τ: Delay time of milling process N: e number of teeth K tc : Cutting force coefficients in the tangential directions K te : Edge force coefficients in the tangential directions K rc : Cutting force coefficients in the radial directions K re : Edge force coefficients in the radial directions a: Axial cutting depth c f : Feed per tooth per revolution (EI) equiv , (EA) equiv : Equivalent bending and tensile stiffness of the cutting tool (ρI) equiv , (m) equiv : Equivalent diametrical mass moment of inertia and equivalent mass per unit length ρ k : Density of the kth layer r k , r k+1 : Inner diameter and outer diameter of the kth layer ϕ 1 (x): Mode function ω 1 , ω 2 : Forward and backward whirling frequencies ε: Small parameter i:

B. The Expression of the Off-Axis Stiffness Coefficient
e off-axis stiffness coefficient of the kth layer for the composite cutting tool is defined as Q 11 � C 11 cos 4 θ (k) + C 22 sin 4 θ (k) + 2 C 12 + 2C 66 sin 2 θ (k) cos 2 θ (k) , where θ (k) is the layer angle of each layer of the material.
e expressions of C 11 , C 12 , C 22 , and C 66 are as follows:

Data Availability
e data used to support the findings of this study are included within the article.

Conflicts of Interest
e authors declare that they have no conflicts of interest.