- Email: [email protected]
Study of the gyroscopic effect of the spindle on the stability characteristics of the milling system
Journal of Materials Processing Technology 138 (2003) 379–384
Study of the gyroscopic effect of the spindle on the stability characteristics of the milling system G.L. Xiong*, J.M. Yi, C. Zeng, H.K. Guo, L.X. Li School of Mechanical Engineering, East China Jiaotong University, Nanchang 330013, PR China
Abstract A new dynamic milling model of a rotating spindle is developed and the gyroscopic effect of the spindle on the stability characteristics of the milling system is investigated for the first time. The results show that although the gyroscopic effect of the rotating spindle does not change the instability regions in milling, it increases the real parts of the eigenvalues of the system or reduces the critical axial depth of cut. In other words, it makes the stability prediction less conservative. # 2003 Elsevier Science B.V. All rights reserved. Keywords: Milling chatter; Rotating spindle; Gyroscopic effect
1. Introduction Machine tool chatter is a self-excited vibration produced by the interaction between the machine tool cutter and the workpiece. In the early research, almost all the researchers voluntarily or involuntarily employed Nyquist stability criterion and the concept of feedback control based on Tobias, Tlusty and Merritt’s work [1–3] to conduct a stability analysis in the frequency domain. According to the authors’ knowledge the gyroscopic effects of the rotating spindle on the chatter frequencies, the mode shapes and the stability characteristics of the milling system have not been considered and incorporated in the milling model in the previous literature. The stability analysis of the milling cutter-rotating spindle system becomes more complicated than the case of stationary milling cutter. The gyroscopic effects of the spindle will make the natural frequency responses become rotating speed-dependent and make one bending vibration mode split into two modes (i.e. the forward-wave and the backward-wave modes). In this case, the frequency response function (FRF) or the system matrices of the milling system become speed-dependent. Therefore, the chatter frequencies and the real parts of the eigenvalues (i.e. stability features) must be solved simultaneously at a specific rotating speed. Although the various methods [4,5] of solving the characteristic equation of the milling system have been well established for the stationary milling models, they cannot be used in the milling model with gyroscopic effects involved. * Corresponding author. E-mail address: [email protected] (G.L. Xiong).
In this paper, a new dynamic milling model including a rotating spindle is developed for the first time. The gyroscopic effect of the rotating spindle is investigated to explore how it relates to the stability characteristics of the milling system.
2. Dynamic milling forces The milling force model for the cutter with a zero helix angle has been well investigated for many decades [6–9]. The milling force model used in this study follows the one in Ref. [5], which is simple but capable of describing the primary characteristics of dynamic milling forces. Fig. 1 shows a cross-section of an end mill cutter, which interacts with a workpiece. The tangential and radial cutting forces applied on the jth tooth in the coordinates R–T rotating with the spindle are given by [5]: Ftj ðtÞ ¼ Kt ðyj Þahðyj Þ;
Frj ðtÞ ¼ Kr ðyj ÞFtj ðtÞ
(1)
where a is the axial depth of cut; Kt(yj) and Kr(yj) the cutting coefficients determined by experiment, which are assumed to be constant for an isotropic material workpiece; and h(yj) the chip thickness that can be expressed in the following form: hðyj Þ ¼ gðyj Þ½Vf sin yj þ ðRj0 ðyj Þ Rj ðyj ÞÞ
(2)
where Vf is the feed per tooth; Rj0 and Rj the radial displacement of the cutter at the previous and present cuts at the position yj, respectively; g(yj) a switch function which determine whether the jth tooth is in or out of cut, namely
0924-0136/03/$ – see front matter # 2003 Elsevier Science B.V. All rights reserved. doi:10.1016/S0924-0136(03)00102-X
380
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
avu ¼ a
Nt X 1 j¼0
2
gðyj Þ½Kt ðyj Þð1 cos 2yj Þ Kt ðyj ÞKr ðyj Þ sin 2yj Þ (7)
auv ¼ a
Nt X
1 gðyj Þ½Kt ðyj Þð1 þ cos2yj Þ Kt ðyj ÞKr ðyj Þsin2yj 2 j¼0 (8)
auu ¼ a
Nt X 1 gðyj Þ½Kt ðyj Þsin2yj Kt ðyj ÞKr ðyj Þð1 cos2yj Þ 2 j¼0
(9) TD
Fig. 1. Dynamic milling force model.
gðyj Þ ¼ 1 when yex yj yst (where yst and yex are the start and exit angular immersions of the cutter, respectively). Otherwise, gðyj Þ ¼ 0. In Eq. (2), the term gðyj ÞVf sin yj represents the static component of the chip thickness. It can be shown easily that the summation of the cutting forces due to these static components is periodic with period 1/ft, where ft is the tooth passing frequency. Therefore, the response to forced vibration may be amplified when ft is near to one of the natural frequencies of the system, which also plays an important role in the initialization of self-excited vibration. However, the periodic cutting force produced by these static components is independent of the response of the system, which does not give any contribution to the regenerative instability mechanism. Thus, it can be removed from Eq. (2). The total dynamic cutting force Fc(t) applied on the cutter teeth can be resolved into two components in the Y and Z directions in the stationary coordinates (XYZ): Fv ðtÞ ¼
N j 1 X j¼0
Fu ðtÞ ¼
N j 1 X j¼0
Fvj ¼
Nj1 X
ðFt j sin yj Fr j cos yj Þ
(3)
j¼0 Nj1 X Fu j ¼ ðFt j cos yj Fr j sin yj Þ
(4)
j¼0
where Nt is the total number of the teeth. Substituting Eqs. (1) and (2) into Eqs. (3) and (4) with Rj ðyj Þ ¼ u sin yj v cos yj and regrouping, gives v v0 Fv ðtÞ avv avu ¼ Fu ðtÞ auv auu u u0 v ¼ ð1 eTD Þ½aðtÞ (5) u
is a time delay operator (i.e. eTD vðtÞ ¼ vðt TÞ ¼ v0 ). e v and u represent the displacements of the spindle in the Y and Z directions, respectively, which are observed from the stationary coordinates. [a(t)] is referred to as the directional dynamic milling coefficient matrix which is periodic with the tooth passing frequency ft.
3. Dynamic finite element model of rotating spindle A finite rotating shaft element using Timoshenko beam theory with gyroscopic terms is shown in Fig. 2. There are two coordinate reference systems used to describe the motion of the element: a fixed coordinates (XYZ) and a coordinate (xyz) rotating with the shaft. The element is considered to be initially straight and is modeled as an eight DOF element. At a given node the element has four degrees of freedom: two displacements v and u, and two slopes about the Z and Y axes which are c and j, respectively. Shear, rotational inertia and gyroscopic effects are included in this element. The elemental kinetic energy consists of both translational and rotational terms, and the potential energy of the element includes elastic bending, shear and axial contributions. The total elemental potential and kinetic energies are developed, then the elemental translational and rotational mass matrices, elemental stiffness matrix and gyroscopic matrix are obtained by applying the extended Hamilton’s principle, which are expressed in the following forms [10,11] (the details of [M0], [M1], [M2], [N0], [N1], [N2], [G0], [G1], [G2], [K0] and [K1] matrices are given in Appendix A).
where avv ¼ a
Nt X 1 j¼0
2
gðyj Þ½Kt ðyj Þ sin 2yj Kt ðyj ÞKr ðyj Þð1 þ cos 2yj Þ (6)
Fig. 2. A typical rotating shaft element.
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
381
(1) Translational mass matrix. ½M ¼ ½M0 þ F½M1 þ F2 ½M2
(10)
where F is the transverse shear effect (F ¼ 12EI=kAGL2 ), EI the bending stiffness per unit of curvature, G the shear modulus, A the cross-sectional area of the shaft, k a constant determined by the shape of the cross-section of the shaft. (2) Rotational mass matrix. ½N ¼ ½N0 þ F½N1 þ F2 ½N2
(11)
Therefore, the total elemental mass matrix is ½M þ ½N. (3) Gyroscopic matrix. ½G ¼ ½G0 þ F½G1 þ F2 ½G2
(12)
(4) Stiffness matrix. ½K ¼ ½K0 þ F½K1
(13)
(5) Disc mass element. In the case of small displacements, the expression for the kinetic energy of a symmetric disc is given by [12]: TD ¼ 12 ðmD ðu2 þ v2 Þ þ IDz ðc2 þ j2 Þ þ IDx ðO2 þ 2OcjÞÞ
(14)
where mD is the mass of the disc, IDx and IDz the moments of inertia of the disc around the X and Z (or Y) axes, respectively. The application of LaGrange’s equations to Eq. (14) gives the mass and gyroscopic matrices in the space-fixed coordinates [12]: 2 38 9 mD 0 0 0 > v > > > > = 0 7 d @T @T 6 6 0 mD 0 7 ¼6 7 > dt @ d_ @d 4 0 0 IDz 0 5> > > >c> ; : j 0 0 0 IDz 2 38 9 0 0 0 0 > v_ > > > > = < > 60 0 0 0 7 6 7 u_ þ O6 7 _ 40 0 0 IDx 5> c> > > > ; : > j_ 0 0 IDx 0 (15)
Fig. 3. A dynamic finite model with a rotating spindle in milling.
displacements at the ith node of the spindle in the Y and Z directions, respectively; vb and ub the corresponding displacements of the bearing kyy, kyz, kzy; and kzz (or cyy, cyz, czy and czz) represent the stiffness (or viscous damping) coefficients of the bearing. The dynamic model shown in Fig. 3 is based on the acceptance of the assumption that the milling machine structure holding the spindle and workpiece clamping system are sufficiently rigid (or they will move together if they move) so that only the relative motion between the tool and the workpiece needs to be considered. Although a continuous spindle is modeled here, it is assumed that the whole axial depth of cut zone deflects the same amount because the stability characteristics of the cutting system are not changed by the very small slope of the spindle tip. Assembling the elemental rotating shaft matrices into the global matrices of the spindle and combining the spindle, the bearing and the milling cutting force models generates the whole dynamic milling model which can be expressed in the following form: s g þ ½GfU_ s g þ ½KfUs g þ ð1 eTD Þ½AðtÞfUs g ¼ f0g ½MfU (17) ^ s g=fU s gg represents the whole diswhere fUs g ¼ ffU ^ sg placement vector of the spindle including the vector fU associated with the degrees of freedom of the spindle tip sg which interacts with the workpiece and the vector fU associated with the remaining degrees of freedom of the spindle. ½M ¼ ½Ms
(18)
where [Ms] is the global mass matrix of the spindle.
T
where d ¼ ðv; u; c; jÞ . (6) Bearing stiffness and damping model. The elastic and viscous damping forces provided by the bearings supporting the spindle can be expressed in the form: kyy kyz FBYi vi vb ¼ kzy kzz FBZi ui u b cyy cyz v_ i v_ b (16) czy czz u_ i u_ b where FBYi and FBZi are the elastic and damping forces of the bearing support at the ith node of the spindle in the Y and Z directions, respectively; vi and ui the
½K ¼ ½Ks þ ½Kb
(19)
where [Ks] is the global stiffness matrix of the spindle, [Kb] the global bearing stiffness matrix associated with the stiffness matrix in Eq. (16). ½G ¼ ½Gs þ ½Cs þ ½Cb
(20)
where [Gs] is the global gyroscopic matrix of the spindle, [Cs] the global structural damping matrix of the spindle, which will be assumed in the form of ½Cs ¼ a½Ms or ½Cs ¼ a½Ms þ b½Ks (a and b are damping ratios), [Cb] the global bearing-damping matrix associated with the damping matrix in Eq. (16).
382
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
Table 1 The parameters of the finite element model of the spindlea Element no.
Length (m) Outer diameter (m) Inner diameter (m) E (Pa) a
1
2
3
4
5
6
0.0450 0.0190 0.0000 4.5Eþ11
0.0400 0.0745 0.0000 2.1Eþ11
0.0500 0.0400 0.0000 2.1Eþ11
0.0200 0.0600 0.0000 2.1Eþ11
0.2000 0.0600 0.0400 2.1Eþ11
0.2000 0.0600 0.0400 2.1Eþ11
All the stiffness coefficients of the bearings located at nodes 5 and 7 are kyy ¼ kzz ¼ 5:0E þ 8 (N/m) and kyz ¼ kzy ¼ 0 (N/m), respectively.
[A(t)] is a time-varying matrix associated with the cutting forces, which can be written in the form: ½aðtÞ 0 ½ AðtÞ ¼ (21) 0 0 where [a(t)] is given by Eqs. (6)–(9). It can be verified easily that [A(t)] is periodic with the tooth passing period T. eTD is a time delay operator.
4. Application of the stability analysis in milling Consider a relatively simple dynamic milling model with six finite elements as above in Fig. 4, whose parameters are
listed in Table 1. It is assumed that the motion of this model involves primarily the relative displacement between the cutter and the workpiece. A milling cutter with Nt ¼ 8 straight teeth is used, and the start and exit immersion angles are yst ¼ 0 8C and yex ¼ 90 8C, respectively. Muller’s optimization algorithm was used to predict the stability characteristics of this model. In order to assess the gyroscopic effects of rotating spindle on the stability characteristics of the milling system, a comparison of the eigenvalues between the different milling models with or without gyroscopic effect is made and the results are shown in Fig. 5. From this figure, it should be noted that although the gyroscopic effect of the rotating spindle does not change the instability regions of the milling system, it does change the real parts of the eigenvalues up to about 30% for this model. It should also be noted from this figure that the gyroscopic effect of the rotating spindle would change the resonant frequencies of the system dramatically.
5. Conclusions Fig. 4. The finite element model of the rotating spindle with the milling cutter.
In this paper, a new dynamic milling model of a rotating spindle is developed and the gyroscopic effect of the spindle on the stability characteristics of the milling system is investigated for the first time. Two important conclusions are summarized as follows: (1) Although the gyroscopic effect of the rotating spindle does not change the instability regions in milling, it increases the real parts of the eigenvalues of the system or reduces the critical axial depth of cut. In other words, it makes the stability prediction less conservative. (2) The gyroscopic effect of the rotating spindle usually reduces the chatter frequency of the system, especially as the rotating speed increases.
Appendix A. Elemental matrices of a finite rotating shaft element Fig. 5. A comparison of the stability results from the different models (NGE: without gyroscopic effect; WGE: with gyroscopic effect).
1. Translational mass matrix. ½ M ¼ ½ M0 þ F½ M1 þ F2 ½ M2 , where F is the transverse shear effect
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
(F ¼ 12EI=kAGL2 ): 2 6 6 6 6 6 6 rAl 6 ½M0 ¼ 26 420ð1 þ FÞ 6 6 6 6 6 4 2
3
156 0
156
0 22l
22l 0
4l2 0
4l2
54 0
0 54
0 13l
13l 0
156 0
156
0
13l
3l
2
0
0
22l
4l2
13l
0
0
3l2
22l
0
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM
4l2 3
294
6 0 294 6 6 6 0 38:5l 7l2 6 6 38:5l 0 0 rAl 6 ½M1 ¼ 26 6 0 0 420ð1 þ FÞ 6 126 6 0 126 31:5l 6 6 4 0 31:5l 7l2 31:5l 0 0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 7l2 31:5l
294
0
0
294
0 7l2
0 38:5l
38:5l 0
7l2 0
7l2 3
2
140 6 0 6 6 6 0 6 6 17:5l rAl 6 ½M2 ¼ 26 70 420ð1 þ FÞ 6 6 6 0 6 6 4 0
140 17:5l
3:5l
0 0
0 0
3:5l2 17:5l
140
70 17:5l
17:5l 3:5l2
0 0
0 0
140 17:5l
3:5l2
0
0
3:5l2
17:5l
0
0
17:5l
2. Rotational mass matrix. ½N ¼ ½N0 þ F½N1 þ F2 ½N2 2 36 6 0 36 SYM 6 6 2 6 0 3l 4l 6 6 0 0 4l2 rI 6 3l ½N0 ¼ 26 36 0 0 3l 36 30lð1 þ FÞ 6 6 6 0 36 3l 0 0 36 6 6 2 4 0 3l l 0 0 3l 3l 2
0
l2
0
3l
0
3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
4l2 4l2
0
3
0
6 0 6 6 6 0 6 6 15l rI 6 ½N1 ¼ 26 0 30lð1 þ FÞ 6 6 6 0 6 6 4 0 15l
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 2
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM
15l 0
5l2 0
5l2
0 0
0 15l
15l 0
0 0
0
2
15l
5l
0
0
15l
5l2
0
0
5l2
15l
0
0
5l2
3:5l2
383
384
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
½N2 ¼
rI
4. Stiffness matrix. ½K ¼ ½K0 þ F½K1
30lð1 þ FÞ2 2 0 60 0 6 6 6 0 0 10l2 6 60 0 0 6 6 60 0 0 6 60 0 0 6 6 4 0 0 5l2 0
0
3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 10l
0
2
0 0
0 0
0
0 5l2
0 0
0 0
10l2 0
10l2
3. Gyroscopic matrix. ½G ¼ ½G0 þ F½G1 þ F2 ½G2 rIO ½G0 ¼ 15l 2
3
0 6 36 0 Skew SYM 6 6 6 3l 0 0 6 6 0 3l 4l2 0 6 6 6 0 36 3l 0 0 6 6 36 0 0 3l 36 0 6 6 2 4 3l 0 0 l 3l 0 0 3l l2
0 ½G1 ¼
rIO 15l 2
0
EI l3 ð1 þ FÞ 2
3
12
6 0 12 6 6 6 0 6l 4l2 6 6 6l 0 0 4l2 6 6 6 12 0 0 6l 6 6 0 12 6l 0 6 6 4 0 6l 2l2 0 6l 0 0 2l2 2 0 60 0 6 6 6 0 0 l2 6 60 0 0 EI 6 ½K1 ¼ 3 6 l ð1 þ FÞ 6 0 0 0 6 60 0 0 6 6 4 0 0 l2 0 0
0
SYM
12 0
12
0 6l
6l 0
4l2 0 4l2
7 7 7 7 7 7 7 7 7 7 7 7 7 5 3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM l2 0 0
0 0
0
0
0
0
l2
l2 0
0
0
l2
4l2 10l2
3l
References 3
0
6 0 0 6 6 6 15l 0 6 6 0 15l 6 6 6 0 0 6 6 0 0 6 6 4 15l 0 rIO 15l 2 0 60 6 6 60 6 60 6 6 60 6 60 6 6 40 0
Skew SYM 0 5l2
0
15l 0
0 15l
0
5l
15l 5l2
0 ½G2 ¼
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
½K0 ¼
2
0 0
0
15l
0
0
0
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
15l 5l2 0
3 0 0
Skew SYM 0
0 10l2
0
0 0
0 0
0 0
0 0 0 5l2
5l 0
2
0 0
0
0 0
0 0
0 10l2 0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
[1] S.A. Tobias, W. Fishwick, A Theory of Regenerative Chatter, The Engineer, London, 1958. [2] J. Tlusty, M. Polacek, The stability of machine tools against selfexcited vibrations in machining, ASME Int. Res. Prod. Eng. (1963) 465. [3] H.E. Meritt, Theory of self-excited machine tool chatter, ASME J. Eng. Ind. 87 (1965) 447. [4] I. Minis, T. Yanushevsky, A new theoretical approach for the prediction of machine tool chatter, ASME J. Eng. Ind. 115 (1993) 1. [5] Y. Altintas, E. Budak, Analytical prediction of stability lobes in milling, Ann. CIRP 44 (1995) 357. [6] J. Tlusty, P. Macnell, Dynamics of cutting forces in end milling, Ann. CIRP 24 (1975) 21. [7] J.W. Sutherland, R.E. Devor, An improved method for cutting force and surface error prediction in flexible end milling systems, ASME J. Eng. Ind. 108 (1986) 269. [8] E.J.A. Armarego, R.C. Whitefield, Computer based modeling of popular machining operations for force and power predictions, Ann. CIRP 34 (1985) 65. [9] E. Budak, Y. Altintas, Prediction of milling force coefficients from orthogonal cutting data, in: Proceedings of the ASME Winter Annual Meeting, New Orleans, PED 64 (1993) 453. [10] H.D. Nelson, A finite rotating shaft element using Timoshenko beam theory, J. Mech. Des. 102 (1980) 793. [11] H.D. Nelson, J.M. Macvaugh, The dynamics of rotor-bearing systems using finite elements, ASME J. Eng. Ind. 98 (1976) 593. [12] M. Lalanne, G. Ferraris, Rotor Dynamics Prediction in Engineering, Wiley, UK, 1990.
Study of the gyroscopic effect of the spindle on the stability characteristics of the milling system G.L. Xiong*, J.M. Yi, C. Zeng, H.K. Guo, L.X. Li School of Mechanical Engineering, East China Jiaotong University, Nanchang 330013, PR China
Abstract A new dynamic milling model of a rotating spindle is developed and the gyroscopic effect of the spindle on the stability characteristics of the milling system is investigated for the first time. The results show that although the gyroscopic effect of the rotating spindle does not change the instability regions in milling, it increases the real parts of the eigenvalues of the system or reduces the critical axial depth of cut. In other words, it makes the stability prediction less conservative. # 2003 Elsevier Science B.V. All rights reserved. Keywords: Milling chatter; Rotating spindle; Gyroscopic effect
1. Introduction Machine tool chatter is a self-excited vibration produced by the interaction between the machine tool cutter and the workpiece. In the early research, almost all the researchers voluntarily or involuntarily employed Nyquist stability criterion and the concept of feedback control based on Tobias, Tlusty and Merritt’s work [1–3] to conduct a stability analysis in the frequency domain. According to the authors’ knowledge the gyroscopic effects of the rotating spindle on the chatter frequencies, the mode shapes and the stability characteristics of the milling system have not been considered and incorporated in the milling model in the previous literature. The stability analysis of the milling cutter-rotating spindle system becomes more complicated than the case of stationary milling cutter. The gyroscopic effects of the spindle will make the natural frequency responses become rotating speed-dependent and make one bending vibration mode split into two modes (i.e. the forward-wave and the backward-wave modes). In this case, the frequency response function (FRF) or the system matrices of the milling system become speed-dependent. Therefore, the chatter frequencies and the real parts of the eigenvalues (i.e. stability features) must be solved simultaneously at a specific rotating speed. Although the various methods [4,5] of solving the characteristic equation of the milling system have been well established for the stationary milling models, they cannot be used in the milling model with gyroscopic effects involved. * Corresponding author. E-mail address: [email protected] (G.L. Xiong).
In this paper, a new dynamic milling model including a rotating spindle is developed for the first time. The gyroscopic effect of the rotating spindle is investigated to explore how it relates to the stability characteristics of the milling system.
2. Dynamic milling forces The milling force model for the cutter with a zero helix angle has been well investigated for many decades [6–9]. The milling force model used in this study follows the one in Ref. [5], which is simple but capable of describing the primary characteristics of dynamic milling forces. Fig. 1 shows a cross-section of an end mill cutter, which interacts with a workpiece. The tangential and radial cutting forces applied on the jth tooth in the coordinates R–T rotating with the spindle are given by [5]: Ftj ðtÞ ¼ Kt ðyj Þahðyj Þ;
Frj ðtÞ ¼ Kr ðyj ÞFtj ðtÞ
(1)
where a is the axial depth of cut; Kt(yj) and Kr(yj) the cutting coefficients determined by experiment, which are assumed to be constant for an isotropic material workpiece; and h(yj) the chip thickness that can be expressed in the following form: hðyj Þ ¼ gðyj Þ½Vf sin yj þ ðRj0 ðyj Þ Rj ðyj ÞÞ
(2)
where Vf is the feed per tooth; Rj0 and Rj the radial displacement of the cutter at the previous and present cuts at the position yj, respectively; g(yj) a switch function which determine whether the jth tooth is in or out of cut, namely
0924-0136/03/$ – see front matter # 2003 Elsevier Science B.V. All rights reserved. doi:10.1016/S0924-0136(03)00102-X
380
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
avu ¼ a
Nt X 1 j¼0
2
gðyj Þ½Kt ðyj Þð1 cos 2yj Þ Kt ðyj ÞKr ðyj Þ sin 2yj Þ (7)
auv ¼ a
Nt X
1 gðyj Þ½Kt ðyj Þð1 þ cos2yj Þ Kt ðyj ÞKr ðyj Þsin2yj 2 j¼0 (8)
auu ¼ a
Nt X 1 gðyj Þ½Kt ðyj Þsin2yj Kt ðyj ÞKr ðyj Þð1 cos2yj Þ 2 j¼0
(9) TD
Fig. 1. Dynamic milling force model.
gðyj Þ ¼ 1 when yex yj yst (where yst and yex are the start and exit angular immersions of the cutter, respectively). Otherwise, gðyj Þ ¼ 0. In Eq. (2), the term gðyj ÞVf sin yj represents the static component of the chip thickness. It can be shown easily that the summation of the cutting forces due to these static components is periodic with period 1/ft, where ft is the tooth passing frequency. Therefore, the response to forced vibration may be amplified when ft is near to one of the natural frequencies of the system, which also plays an important role in the initialization of self-excited vibration. However, the periodic cutting force produced by these static components is independent of the response of the system, which does not give any contribution to the regenerative instability mechanism. Thus, it can be removed from Eq. (2). The total dynamic cutting force Fc(t) applied on the cutter teeth can be resolved into two components in the Y and Z directions in the stationary coordinates (XYZ): Fv ðtÞ ¼
N j 1 X j¼0
Fu ðtÞ ¼
N j 1 X j¼0
Fvj ¼
Nj1 X
ðFt j sin yj Fr j cos yj Þ
(3)
j¼0 Nj1 X Fu j ¼ ðFt j cos yj Fr j sin yj Þ
(4)
j¼0
where Nt is the total number of the teeth. Substituting Eqs. (1) and (2) into Eqs. (3) and (4) with Rj ðyj Þ ¼ u sin yj v cos yj and regrouping, gives v v0 Fv ðtÞ avv avu ¼ Fu ðtÞ auv auu u u0 v ¼ ð1 eTD Þ½aðtÞ (5) u
is a time delay operator (i.e. eTD vðtÞ ¼ vðt TÞ ¼ v0 ). e v and u represent the displacements of the spindle in the Y and Z directions, respectively, which are observed from the stationary coordinates. [a(t)] is referred to as the directional dynamic milling coefficient matrix which is periodic with the tooth passing frequency ft.
3. Dynamic finite element model of rotating spindle A finite rotating shaft element using Timoshenko beam theory with gyroscopic terms is shown in Fig. 2. There are two coordinate reference systems used to describe the motion of the element: a fixed coordinates (XYZ) and a coordinate (xyz) rotating with the shaft. The element is considered to be initially straight and is modeled as an eight DOF element. At a given node the element has four degrees of freedom: two displacements v and u, and two slopes about the Z and Y axes which are c and j, respectively. Shear, rotational inertia and gyroscopic effects are included in this element. The elemental kinetic energy consists of both translational and rotational terms, and the potential energy of the element includes elastic bending, shear and axial contributions. The total elemental potential and kinetic energies are developed, then the elemental translational and rotational mass matrices, elemental stiffness matrix and gyroscopic matrix are obtained by applying the extended Hamilton’s principle, which are expressed in the following forms [10,11] (the details of [M0], [M1], [M2], [N0], [N1], [N2], [G0], [G1], [G2], [K0] and [K1] matrices are given in Appendix A).
where avv ¼ a
Nt X 1 j¼0
2
gðyj Þ½Kt ðyj Þ sin 2yj Kt ðyj ÞKr ðyj Þð1 þ cos 2yj Þ (6)
Fig. 2. A typical rotating shaft element.
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
381
(1) Translational mass matrix. ½M ¼ ½M0 þ F½M1 þ F2 ½M2
(10)
where F is the transverse shear effect (F ¼ 12EI=kAGL2 ), EI the bending stiffness per unit of curvature, G the shear modulus, A the cross-sectional area of the shaft, k a constant determined by the shape of the cross-section of the shaft. (2) Rotational mass matrix. ½N ¼ ½N0 þ F½N1 þ F2 ½N2
(11)
Therefore, the total elemental mass matrix is ½M þ ½N. (3) Gyroscopic matrix. ½G ¼ ½G0 þ F½G1 þ F2 ½G2
(12)
(4) Stiffness matrix. ½K ¼ ½K0 þ F½K1
(13)
(5) Disc mass element. In the case of small displacements, the expression for the kinetic energy of a symmetric disc is given by [12]: TD ¼ 12 ðmD ðu2 þ v2 Þ þ IDz ðc2 þ j2 Þ þ IDx ðO2 þ 2OcjÞÞ
(14)
where mD is the mass of the disc, IDx and IDz the moments of inertia of the disc around the X and Z (or Y) axes, respectively. The application of LaGrange’s equations to Eq. (14) gives the mass and gyroscopic matrices in the space-fixed coordinates [12]: 2 38 9 mD 0 0 0 > v > > > > = 0 7 d @T @T 6 6 0 mD 0 7 ¼6 7 > dt @ d_ @d 4 0 0 IDz 0 5> > > >c> ; : j 0 0 0 IDz 2 38 9 0 0 0 0 > v_ > > > > = < > 60 0 0 0 7 6 7 u_ þ O6 7 _ 40 0 0 IDx 5> c> > > > ; : > j_ 0 0 IDx 0 (15)
Fig. 3. A dynamic finite model with a rotating spindle in milling.
displacements at the ith node of the spindle in the Y and Z directions, respectively; vb and ub the corresponding displacements of the bearing kyy, kyz, kzy; and kzz (or cyy, cyz, czy and czz) represent the stiffness (or viscous damping) coefficients of the bearing. The dynamic model shown in Fig. 3 is based on the acceptance of the assumption that the milling machine structure holding the spindle and workpiece clamping system are sufficiently rigid (or they will move together if they move) so that only the relative motion between the tool and the workpiece needs to be considered. Although a continuous spindle is modeled here, it is assumed that the whole axial depth of cut zone deflects the same amount because the stability characteristics of the cutting system are not changed by the very small slope of the spindle tip. Assembling the elemental rotating shaft matrices into the global matrices of the spindle and combining the spindle, the bearing and the milling cutting force models generates the whole dynamic milling model which can be expressed in the following form: s g þ ½GfU_ s g þ ½KfUs g þ ð1 eTD Þ½AðtÞfUs g ¼ f0g ½MfU (17) ^ s g=fU s gg represents the whole diswhere fUs g ¼ ffU ^ sg placement vector of the spindle including the vector fU associated with the degrees of freedom of the spindle tip sg which interacts with the workpiece and the vector fU associated with the remaining degrees of freedom of the spindle. ½M ¼ ½Ms
(18)
where [Ms] is the global mass matrix of the spindle.
T
where d ¼ ðv; u; c; jÞ . (6) Bearing stiffness and damping model. The elastic and viscous damping forces provided by the bearings supporting the spindle can be expressed in the form: kyy kyz FBYi vi vb ¼ kzy kzz FBZi ui u b cyy cyz v_ i v_ b (16) czy czz u_ i u_ b where FBYi and FBZi are the elastic and damping forces of the bearing support at the ith node of the spindle in the Y and Z directions, respectively; vi and ui the
½K ¼ ½Ks þ ½Kb
(19)
where [Ks] is the global stiffness matrix of the spindle, [Kb] the global bearing stiffness matrix associated with the stiffness matrix in Eq. (16). ½G ¼ ½Gs þ ½Cs þ ½Cb
(20)
where [Gs] is the global gyroscopic matrix of the spindle, [Cs] the global structural damping matrix of the spindle, which will be assumed in the form of ½Cs ¼ a½Ms or ½Cs ¼ a½Ms þ b½Ks (a and b are damping ratios), [Cb] the global bearing-damping matrix associated with the damping matrix in Eq. (16).
382
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
Table 1 The parameters of the finite element model of the spindlea Element no.
Length (m) Outer diameter (m) Inner diameter (m) E (Pa) a
1
2
3
4
5
6
0.0450 0.0190 0.0000 4.5Eþ11
0.0400 0.0745 0.0000 2.1Eþ11
0.0500 0.0400 0.0000 2.1Eþ11
0.0200 0.0600 0.0000 2.1Eþ11
0.2000 0.0600 0.0400 2.1Eþ11
0.2000 0.0600 0.0400 2.1Eþ11
All the stiffness coefficients of the bearings located at nodes 5 and 7 are kyy ¼ kzz ¼ 5:0E þ 8 (N/m) and kyz ¼ kzy ¼ 0 (N/m), respectively.
[A(t)] is a time-varying matrix associated with the cutting forces, which can be written in the form: ½aðtÞ 0 ½ AðtÞ ¼ (21) 0 0 where [a(t)] is given by Eqs. (6)–(9). It can be verified easily that [A(t)] is periodic with the tooth passing period T. eTD is a time delay operator.
4. Application of the stability analysis in milling Consider a relatively simple dynamic milling model with six finite elements as above in Fig. 4, whose parameters are
listed in Table 1. It is assumed that the motion of this model involves primarily the relative displacement between the cutter and the workpiece. A milling cutter with Nt ¼ 8 straight teeth is used, and the start and exit immersion angles are yst ¼ 0 8C and yex ¼ 90 8C, respectively. Muller’s optimization algorithm was used to predict the stability characteristics of this model. In order to assess the gyroscopic effects of rotating spindle on the stability characteristics of the milling system, a comparison of the eigenvalues between the different milling models with or without gyroscopic effect is made and the results are shown in Fig. 5. From this figure, it should be noted that although the gyroscopic effect of the rotating spindle does not change the instability regions of the milling system, it does change the real parts of the eigenvalues up to about 30% for this model. It should also be noted from this figure that the gyroscopic effect of the rotating spindle would change the resonant frequencies of the system dramatically.
5. Conclusions Fig. 4. The finite element model of the rotating spindle with the milling cutter.
In this paper, a new dynamic milling model of a rotating spindle is developed and the gyroscopic effect of the spindle on the stability characteristics of the milling system is investigated for the first time. Two important conclusions are summarized as follows: (1) Although the gyroscopic effect of the rotating spindle does not change the instability regions in milling, it increases the real parts of the eigenvalues of the system or reduces the critical axial depth of cut. In other words, it makes the stability prediction less conservative. (2) The gyroscopic effect of the rotating spindle usually reduces the chatter frequency of the system, especially as the rotating speed increases.
Appendix A. Elemental matrices of a finite rotating shaft element Fig. 5. A comparison of the stability results from the different models (NGE: without gyroscopic effect; WGE: with gyroscopic effect).
1. Translational mass matrix. ½ M ¼ ½ M0 þ F½ M1 þ F2 ½ M2 , where F is the transverse shear effect
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
(F ¼ 12EI=kAGL2 ): 2 6 6 6 6 6 6 rAl 6 ½M0 ¼ 26 420ð1 þ FÞ 6 6 6 6 6 4 2
3
156 0
156
0 22l
22l 0
4l2 0
4l2
54 0
0 54
0 13l
13l 0
156 0
156
0
13l
3l
2
0
0
22l
4l2
13l
0
0
3l2
22l
0
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM
4l2 3
294
6 0 294 6 6 6 0 38:5l 7l2 6 6 38:5l 0 0 rAl 6 ½M1 ¼ 26 6 0 0 420ð1 þ FÞ 6 126 6 0 126 31:5l 6 6 4 0 31:5l 7l2 31:5l 0 0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 7l2 31:5l
294
0
0
294
0 7l2
0 38:5l
38:5l 0
7l2 0
7l2 3
2
140 6 0 6 6 6 0 6 6 17:5l rAl 6 ½M2 ¼ 26 70 420ð1 þ FÞ 6 6 6 0 6 6 4 0
140 17:5l
3:5l
0 0
0 0
3:5l2 17:5l
140
70 17:5l
17:5l 3:5l2
0 0
0 0
140 17:5l
3:5l2
0
0
3:5l2
17:5l
0
0
17:5l
2. Rotational mass matrix. ½N ¼ ½N0 þ F½N1 þ F2 ½N2 2 36 6 0 36 SYM 6 6 2 6 0 3l 4l 6 6 0 0 4l2 rI 6 3l ½N0 ¼ 26 36 0 0 3l 36 30lð1 þ FÞ 6 6 6 0 36 3l 0 0 36 6 6 2 4 0 3l l 0 0 3l 3l 2
0
l2
0
3l
0
3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
4l2 4l2
0
3
0
6 0 6 6 6 0 6 6 15l rI 6 ½N1 ¼ 26 0 30lð1 þ FÞ 6 6 6 0 6 6 4 0 15l
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 2
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM
15l 0
5l2 0
5l2
0 0
0 15l
15l 0
0 0
0
2
15l
5l
0
0
15l
5l2
0
0
5l2
15l
0
0
5l2
3:5l2
383
384
G.L. Xiong et al. / Journal of Materials Processing Technology 138 (2003) 379–384
½N2 ¼
rI
4. Stiffness matrix. ½K ¼ ½K0 þ F½K1
30lð1 þ FÞ2 2 0 60 0 6 6 6 0 0 10l2 6 60 0 0 6 6 60 0 0 6 60 0 0 6 6 4 0 0 5l2 0
0
3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM 10l
0
2
0 0
0 0
0
0 5l2
0 0
0 0
10l2 0
10l2
3. Gyroscopic matrix. ½G ¼ ½G0 þ F½G1 þ F2 ½G2 rIO ½G0 ¼ 15l 2
3
0 6 36 0 Skew SYM 6 6 6 3l 0 0 6 6 0 3l 4l2 0 6 6 6 0 36 3l 0 0 6 6 36 0 0 3l 36 0 6 6 2 4 3l 0 0 l 3l 0 0 3l l2
0 ½G1 ¼
rIO 15l 2
0
EI l3 ð1 þ FÞ 2
3
12
6 0 12 6 6 6 0 6l 4l2 6 6 6l 0 0 4l2 6 6 6 12 0 0 6l 6 6 0 12 6l 0 6 6 4 0 6l 2l2 0 6l 0 0 2l2 2 0 60 0 6 6 6 0 0 l2 6 60 0 0 EI 6 ½K1 ¼ 3 6 l ð1 þ FÞ 6 0 0 0 6 60 0 0 6 6 4 0 0 l2 0 0
0
SYM
12 0
12
0 6l
6l 0
4l2 0 4l2
7 7 7 7 7 7 7 7 7 7 7 7 7 5 3 7 7 7 7 7 7 7 7 7 7 7 7 7 5
SYM l2 0 0
0 0
0
0
0
0
l2
l2 0
0
0
l2
4l2 10l2
3l
References 3
0
6 0 0 6 6 6 15l 0 6 6 0 15l 6 6 6 0 0 6 6 0 0 6 6 4 15l 0 rIO 15l 2 0 60 6 6 60 6 60 6 6 60 6 60 6 6 40 0
Skew SYM 0 5l2
0
15l 0
0 15l
0
5l
15l 5l2
0 ½G2 ¼
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
½K0 ¼
2
0 0
0
15l
0
0
0
0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
15l 5l2 0
3 0 0
Skew SYM 0
0 10l2
0
0 0
0 0
0 0
0 0 0 5l2
5l 0
2
0 0
0
0 0
0 0
0 10l2 0
7 7 7 7 7 7 7 7 7 7 7 7 7 5
[1] S.A. Tobias, W. Fishwick, A Theory of Regenerative Chatter, The Engineer, London, 1958. [2] J. Tlusty, M. Polacek, The stability of machine tools against selfexcited vibrations in machining, ASME Int. Res. Prod. Eng. (1963) 465. [3] H.E. Meritt, Theory of self-excited machine tool chatter, ASME J. Eng. Ind. 87 (1965) 447. [4] I. Minis, T. Yanushevsky, A new theoretical approach for the prediction of machine tool chatter, ASME J. Eng. Ind. 115 (1993) 1. [5] Y. Altintas, E. Budak, Analytical prediction of stability lobes in milling, Ann. CIRP 44 (1995) 357. [6] J. Tlusty, P. Macnell, Dynamics of cutting forces in end milling, Ann. CIRP 24 (1975) 21. [7] J.W. Sutherland, R.E. Devor, An improved method for cutting force and surface error prediction in flexible end milling systems, ASME J. Eng. Ind. 108 (1986) 269. [8] E.J.A. Armarego, R.C. Whitefield, Computer based modeling of popular machining operations for force and power predictions, Ann. CIRP 34 (1985) 65. [9] E. Budak, Y. Altintas, Prediction of milling force coefficients from orthogonal cutting data, in: Proceedings of the ASME Winter Annual Meeting, New Orleans, PED 64 (1993) 453. [10] H.D. Nelson, A finite rotating shaft element using Timoshenko beam theory, J. Mech. Des. 102 (1980) 793. [11] H.D. Nelson, J.M. Macvaugh, The dynamics of rotor-bearing systems using finite elements, ASME J. Eng. Ind. 98 (1976) 593. [12] M. Lalanne, G. Ferraris, Rotor Dynamics Prediction in Engineering, Wiley, UK, 1990.