Modeling, identification and application of process damping in milling of titanium alloy

2.1. Formation mechanism of process damping

The 2 DOF model shown in Fig. 1(a), is always used to describe the dynamics of the milling system [3, 4]. ${\mathrm{\Phi }}_j$ is the angular position of the $j$th cutting tool, $\mathrm{\Omega }$ is the spindle speed. The vibration displacements of tool center $o^{\text{'}}$ determined in the workpiece coordinates (i.e., $xoy$) can be recorded as $\text{(}{x}_{tool}\text{,}{y}_{tool}\text{)}$, the wave caused by vibration leave on the workpiece surface in $r$ direction (i.e., chip thickness or radial direction), $t$ direction represent the cutting speed direction. The interference between flank face and wave will occur when the chatter amplitude grows larger. Flank face intend into the wave left behind in radial direction as shown in Fig. 1. As a result, indentation forces (i.e., $F_v^d$ and $F_f^d$), arise in normal, $v$, and frictional, $f$, directions on the tool flank face as shown in Fig. 1(b), which create a damping effects. The damping forces in $v$ and $f$ directions acting on each tooth are oriented and summed up as a resultant damping force in $x$ and $y$ directions.

Fig. 1End milling system dynamics

a)

b)

Fig. 2Formation mechanism of indentation area

a)

b)

Process damping is proportional to indentation area. In Fig. 2(a), $\alpha$ is the clearance angle of tool, angle $\gamma$ represent the slope of wave left in radial direction, $u_r^{wave}$ is the radial coordinate values of wave, which represent the radial vibration displacement. The instantaneous slope of wave, represent the radial vibration speed (i.e., ${\dot{u}}_r^{wave}$). It is assumed that the tool does not indent the previous wave as shown in Fig. 2(a). Therefore, the interference between flank face and workpiece occurs only as the tool slide down the wave, where ${\dot{u}}_r^{wave}<0$. The occurrence of interference is determined by the relationship between $\alpha$ and $\gamma$. If ${\dot{u}}_r^{wave}$ satisfies the conditions:

The interference will occur. In Eq. (1),$v_c$ is the cutting speed, $s$ represents tangential arc length of tool edge, $t$ is time here. It can be known from Eq. (1) that interference is easy to occur as $v_c$ and $\alpha$ decreases, indentation area and process damping will also increase. The law is also shown in Fig. 2(b).

2.2. Calculation of indentation area

The calculation of indentation area $U\left(t\right)$shown in Fig. 1(a) is rather complicated. As shown in Fig. 3, the cutting tool is divided into a number of axial elements, elemental axial width is $dz$. In axial level $z$, point $A$ represents the current position of the tool edge, which can be also seen in Fig. 4(a), correspondingly, $u_{r,0}^{wave}$ is the current radial coordinate values of wave. As shown in Fig. 4(a), $B_i$ ($i=$1, 2, 3...) represent the point on the surface wave generated by the tool edge in the previous time step, corresponding radial coordinate values of wave calculated in previous time steps is $u_{r,i}^{wave}$. For every $B_i$, $C_i$ is the corresponding point on the flank face, $u_{r,i}^{flank}$ is its radial coordinate values. The distance between $B_i$ and $C_i$ is $d_i$. $d_0$ is zero. The indentation area $U\left(t\right)$ in every axial level is bounded by the tool flank face and the workpiece surface wave. This area is calculated at each simulation time step numerically.

Fig. 3Digitized end mill model

Fig. 4Calculation of indentation area

a)

b)

The calculation of $u_{r,i}^{flank}$ is calculated according to Eq. (2), here, $∆s=v_c∆t$ is the tangential arc length between adjacent two time point. $∆t$ is the simulation time step:

The intersection point of flank face and wave, $D$, needs to be searched in the calculation process. As shown in Fig. 4(b), when $d_i<\text{0}$, the calculation of indentation area can stop. The number of time step front and back intersection can be recorded as $n-1$ and $n$ respectively, here $d_{n-1}>\text{0}$, $d_n<\text{0}$, thus, indentation area $U\left(t\right)$ can be calculated by summing up the discrete areas according to Eq. (3):

2.3. Effects of anti-vibration clearance angle on process damping

It can be seen from Eq. (1) that interference is much easier to occur and process damping is more significant as cutting speed decreases, smaller clearance angle can lead to the same effects. However, in the tool design, severe extrusion is easy to occur and the surface quality will become very poor if the clearance angle is too small. Therefore, transition edge needs to be adopted in design. As shown in Fig. 5, the tool modified have two clearance angle, ${\alpha }_1$ and ${\alpha }_2$, ${\alpha }_1$ is designed to be smaller to suppress vibration, which can be named anti-vibration angle, the length of its corresponding transition edge is $W$.

Fig. 5Effects of anti-vibration clearance angle on indentation area

Fig. 6Calculation of indentation area with anti-vibration clearance angle

The calculation of indentation area follows the methods provided in Section 2.2, the main parameters in Fig. 6 is the same as ones in Fig. 4(a). Minor modification needs to be done because the adoption of anti-vibration clearance angle. The calculation of $u_{r,i}^{flank}$ follows Eq. (4), then the calculation of indentation area also follow the Eq. (3):

3.1. Equations of milling dynamics

Dynamic cutting force model is of great importance in analysis of milling dynamics. At each axial level $z$, the location of the$j$th cutting tooth is determined by its angular position ${\mathrm{\Phi }}_j$ (as shown in Fig. 1):

In Eq. (5), $N$ is the number of teeth, $\beta$ is the helix angle (as shown in Fig. 3), $R$ is the tool radius. $\theta$ is the angular position of tool bottom (axial level $z=\text{0}$), calculated from angular velocity $\mathrm{\Omega }$ and the time $t$: $\theta =\mathrm{\Omega }t$.

The dynamic cutting force is modeled using the linear-edge force model [20], where elemental tangential and radial forces, $F_t$ and $F_r$, at flute, $j$, are expressed in terms of cutting pressure coefficients, $K_c$, cutting edge force coefficients, $K_e$, instantaneous chip thickness, $h$, and elemental axial width, $dz$ (as shown in Fig. 3):

The differential cutting forces are then oriented in $x$ and $y$ directions using the forces in $r$ and $t$ directions as follows:

For the $j$th cutting flute, differential indentation forces acting on the axial level$z$, in normal, ${dF}_v^d$ and in frictional ${dF}_f^d$ directions are given as:

where $U_j(t,z)$ is the indentation area at axial level $z$, which is calculated by methods provided in Section 2.2-2.3. $\mu$ is the coefficient of contact friction between the flank face and the workpiece surface, $K^d$is indentation coefficient. The process damping force in $r$ and $t$ directions are expressed in terms of the forces in $v$ direction and the friction coefficient $\mu$ (as shown in Fig. 1(b)):

The differential process damping forces are then oriented in $x$ and $y$ directions using the process damping forces in $r$ and $t$ directions as follows:

Thus, equations of milling dynamic can be modeled:

where, $S$ is the numbers of axial elements (${Sdz=a}_p$, $a_p$ is axial cutting depth), milling at small radial immersion is adopted in the paper, generally, only one tooth is cutting workpiece in one tooth passing period. ($m_x\text{,}{m}_y$), ($c_x\text{,}{c}_y$), ($k_x\text{,}{k}_y$) are the modal mass, modal damping, modal stiffness respectively (as shown in Fig. 1(a)), which can be obtained from modal tests and analysis. ($x_{tool}\text{,}{y}_{tool}$) is the vibration displacement of tool, which is also shown in Fig. 1(a). $g$ is a unit step function which determines whether the tooth is in or out of the cutting, which will be discussed in Appendix. Besides, the calculation of dynamic chip thickness $h$ in Eq. (6) and $u_{r,i}^{wave}$ in Eq. (2) are rather complicated, which will be also discussed in detail in Appendix.

3.2. Determining chatter stability lobes from time domain simulations

The Eq. (11) is solved using the classical 4th order Runge Kutta method, where time increment $∆t$ is adopted as 10^-5 seconds. Using the following steps, the stability lobes are evaluated by finding the limiting axial cutting depth through a range of spindle speeds:

1. Modal parameters of machine tool structure, cutting tool geometry, feed rate, tool orientation, entry angle and exit angle, cutting force coefficients, and a starting axial cutting depth, $a$, are specified.

2. A range of spindle speeds and a spindle speed step size are specified.

3. For a given spindle speed, a static time domain simulation is run (suppressing any tool and workpiece vibrations) and the maximum static uncut chip thickness, $h_{s,max}$, is stored.

4. A second time domain simulation is run. For given spindle speed and cutting depth, Eq. (11) is solved, vibration displacements and dynamic chip thickness $h$ is obtained. The largest dynamic chip thickness of all cutting points on the tool, $h_{d,max}$ is stored for the last few revolutions of the simulation.

5. The non-dimensional chatter parameter, $\eta$, is evaluated as:

6. If $\eta$ is greater than a pre-determined limit (1.25 is used in this thesis), the process is unstable, otherwise the process is stable. The criterion to judge is provided in [21].

7. If the process is stable, $a_{min}$ is set to the current value of $a$, $a$ is doubled and steps 3-6 are repeated until chatter occurs. Then $a_{max}$ is set to the value of $a$ when chatter occurred.

If the process is unstable, $a_{max}$ is set to the current value of $a$, $a$ is halved and steps 3-6 are repeated until the process is stable. Then $a_{min}$ is set to the stable axial cutting depth, $a$.

8. Once the range $a_{min}$ to $a_{max}$ is found between which the limiting axial cutting depth lies, a bisection search is performed with $a_{min}<a<a_{max}$, repeating steps 3-6 until the limiting axial cutting depth, $a_{lim}$, is found within a given tolerance.

4.1. Analysis of milling dynamic equations based on Routh criterion

If nonlinear effects such as process damping are ignored, linear time-delay model such as zero order solutions are always adopted in generating stability lobes [4], which can predict stability limit accurately at higher speed. According to classical zero order solutions [3], the equations of milling dynamics are shown in Eq. (13). ${\alpha }_{xx}$, ${\alpha }_{xy}$, ${\alpha }_{yx}$, ${\alpha }_{yy}$ shown in Eq. (13) stand for the average directional coefficients, which relate the dynamic displacements to the dynamic cutting forces. $a_p$ is the axial cutting depth. The physical meanings of other parameters in Eq. (13) are the same as ones in Eq. (11), ${K_t=K}_{tc}$, $\left(x,y\right)=\left(x_{tool},y_{tool}\right)$:

The paper presents one new method to analyze the root reason of the occurrence of chatter.

Firstly, the Eq. (13) is converted from time to frequency domain:

Eq. (14) can be written simply:

where, $\mathbf{Z}\left(\omega \right)$is the dynamic stiffness matrix of milling system. The characteristic equations of milling system can be obtained from Eq. (15):

The expansion of Eq. (16) is shown as follows:

Additional damping elements:

Additional stiffness elements:

In Eq. (13), ${\alpha }_{xy}>\text{0,}$${\alpha }_{yx}<\text{0,}$ in general [3]. For regenerative chatter, the phase $\beta ={\omega }_{c}T\in \left(\pi ,2\pi \right)$ [3], thus $\sin \left({\omega }_{c}T\right)<\text{0}$. For additional damping terms, $c_{xy}>\text{0}$, $c_{yx}<\text{0}$, for additional stiffness terms, $k_{xy}<\text{0}$, $k_{yx}>\text{0}$. It is the non-diagonal elements that cause chatter, which can be verified by Routh criterion [22]. In Eq. (17), ${j\omega }_c$ can be replaced by $s$, and the equation can be transformed into Laplace domain:

Eq. (20) can be expanded:

where:

The first row of Routh criterion table is:

According to Routh criterion [22], the system becomes unstable when there are negative elements in Eq. (23). In Eqs. (18)-(19), $k_{xy}<\text{0,}$$c_{yx}<\text{0,}$ thus $k_{xy}{c}_{yx}>\text{0,}$ then $k_{yx}>\text{0,}$$c_{xy}>\text{0}$, then ${c_{xy}k}_{yx}>\text{0}$, as a result, $a_1$ decrease. The coefficient $\mathbf{L}\mathbf{s}\left(\mathrm{4,1}\right)$ decreases obviously as a result of $a_1$ decreasing. Meanwhile, the other elements in Eq. (23) have no obvious trends to decrease. Thus, the limiting axial cutting depth can be predicted when the $\mathbf{L}\mathbf{s}\left(\mathrm{4,1}\right)$ become negative as $a_p$ increases. The prediction method will be verified by milling experiments.

4.2. Identification of process damping

From Section 4.1, it could be concluded that the way to suppress chatter can be attributed to increase diagonal elements$\mathrm{}{c}_i$, $k_i$ ($i=x\text{,}y$) or eliminate non-diagonal elements $c_{ij}$, $k_{ij}$ ($i=x\text{,}y$), these ways can all increase the $a_1$ in Eq. (23) to increase stability. In [14-17], the effects of process damping is considered as an additional damper $c_x^p$, $c_y^p$. Obviously, the process damping can improve limiting axial cutting depth by increasing $c_i$ ($i=x\text{,}y$).

Fig. 7Variations in absolute stability limits with cutting speed

In general, according to linear time-delay model ignoring process damping, there is hardly any stability lobes and stability limits are very small at low speed as shown in Fig. 7, which is due to the many tightly packed vibration waves at each tooth period [3]. As discussed in Section 2, the indentation area and process damping increase as the cutting speed decreases, therefore the absolute stability limit commonly varies with cutting speed as shown in Fig. 7. The stability limit decreases at higher speeds, where the effect diminishes. For given speed $n_s$, $a_{lim}^{hs}$ is the stability limits predicted ignoring process damping, which can be generated by time simulation methods provided in Section 3, where process damping force is ignored in Eq. (11). Whereas $a_{lim}$ is generated when process damping force is considered in Eq. (11). The process damping is identified using the differences between $a_{lim}^{hs}$ and $a_{lim}$.

Considering the symmetry of tool structures, the modal parameters ($m_i$, $c_i$, $k_i$, ${\xi }_i$, $i=x\text{,}y$) in two different directions can be assumed as equal: $m_x=m_y=m$, $c_x=c_y=c$, $k_x=k_y=k$. Structural damping ratio $\xi =c/\left(2\sqrt{mk}\right)$. Thus, for given spindle speed, the process damping can be identified following the steps as below:

1) Generate $a_{lim}$ and $a_{lim}^{hs}$.

2) Entitle value: ${\xi }_1=\xi$.

3) Increase damping ratio: ${\xi }_1={\xi }_1+∆\xi$ ($∆\xi$ is taken as 0.0001).

4) For given ${\xi }_1$, obtain new critical cutting depth $a_{lic}$. In the step, axial cutting depth $a$, increases from $a_{lim}^{hs}$ until $\mathbf{L}\mathbf{s}\left(\mathrm{4,1}\right)$ is negative, then $a_{lic}=a$.

5) If $a_{lic}<a_{lim}$, go to step 3, otherwise, ${\xi }_t={\xi }_1$. Where, ${\xi }_t$ is total damping ratio, which is the sum of structural damping ratio and process damping ratio.

Then, the process damping $c_p$ can be obtained as follows: