Tao Li

Sébastien Seguy
email: sebastien.seguy@insa-toulouse.fr

Alain Berlioz
email: alain.berlioz@univ-tlse3.fr

Optimization mechanism of targeted energy transfer with vibro-impact energy sink under periodic and transient excitation

Keywords: Optimal design, Targeted energy transfer, Vibro impact, Nonlinear energy sink, Impact damper

This paper is dedicated to exploit the same optimization mechanism of targeted energy transfer under different types of excitation. Specifically, a linear oscillator (LO) coupled with a vibro-impact (VI) nonlinear energy sink (NES) is analytically studied with an asymptotical method. The optimization mechanism under periodic excitation with a single frequency and under transient excitation is numerically obtained and experimentally validated for the first time. For periodic excitation, the boundary between the regime with two impacts per cycle and that of strongly modulated response (SMR) is proved to be optimal rather than SMR. The chaotic SMR is experimentally observed from the viewpoint of displacement of LO. The above observed mechanism is further applied to explain the optimization mechanism under transient excitation and that under periodic excitation with a range of frequency. It is experimentally verified that the optimization of the latter can be simplified to the optimization under an excitation with a single resonance frequency. For transient excitation, the efficiency of different transient response regimes is experimentally compared, which agrees with the periodic results. Moreover, the efficiency comparison of different lengths of cavity is also experimentally validated. In short, the close relation of optimization under different excitations is clearly demonstrated.

Introduction

In the past ten years, targeted energy transfer (TET) has been extensively studied in the domain of energy control [START_REF] Vakakis | Targeted Energy Transfer in Mechanical and Structural Systems[END_REF][START_REF] Lee | Passive non-linear targeted energy transfer and its applications to vibration absorption: a review[END_REF]. Nonlinearity is used to couple energy absorber and main system instead of linearity like linear spring for tuned mass damper. This kind of energy absorber device is called nonlinear energy sink (NES). Under the context of TET, undesirable energy can be spatially transferred from a main structure to an attached NES in an irreversible way and be efficiently dissipated. The preceding researches have paved the road for the final goal of optimization design for engineering applications. Around this ultimate objective, the researches that explain the phenomena of TET will be briefly introduced and the focus will be laid on the aspect of optimization.

As the original researches about TET [START_REF] Gendelman | Transition of energy to a nonlinear localized mode in a highly asymmetric system of two oscillators[END_REF][START_REF] Vakakis | Energy pumping in nonlinear mechanical oscillators: part ii: resonance capture[END_REF], the subsequent researches are concentrated on NES with cubic nonlinearity, which is termed as cubic NES. Therefore, the optimization study about cubic NES is more mature, though other types of NES also have been studied such as the two recent papers in the domain of sound and magnet [START_REF] Mattei | Nonlinear targeted energy transfer of two coupled cantilever beams coupled to a bistable light attachment[END_REF][START_REF] Benacchio | Design of a magnetic vibration absorber with tunable stiffnesses[END_REF]. In order to optimize the system with cubic NES, different efforts have been dedicated for transient and periodic excitation with a range of frequency respectively. For transient excitation, it is demonstrated [START_REF] Lee | Passive non-linear targeted energy transfer and its applications to vibration absorption: a review[END_REF] that there are three TET mechanisms: fundamental TET, subharmonic TET and TET initiated by nonlinear beating. The last is proved to be most efficient and is closely related to some special orbits in the frequency energy plot (FEP) from the study of the underlying Hamiltonian system [START_REF] Kerschen | Irreversible passive energy transfer in coupled oscillators with essential nonlinearity[END_REF]. The initial energy of system should be high enough to ensure the activation of these special orbits. To insure this, a design criterion to choose the corresponding stiffness is proposed [START_REF] Nguyen | Design criteria for optimally tuned nonlinear energy sinkspart 1: transient regime[END_REF]. By this way, two requirements for sufficient energy dissipation can be guaranteed: the activation of the third TET mechanism and its dura-tion as long as possible. Although the above two aspects are satisfied and design parameters can be obtained by analytical calculation according to this criterion, the requirement of optimum cannot be insured and the underlying mechanism should be explained further.

In return, optimal parameters under transient excitation can be analytically and precisely obtained by a method proposed in [START_REF] Manevitch | Parameters optimization for energy pumping in strongly nonhomogeneous 2 dof system[END_REF][START_REF] Manevitch | Towards the design of an optimal energetic sink in a strongly inhomogeneous two-degree-of-freedom system[END_REF]. For periodic excitation, the response regimes, especially strongly modulated response (SMR), have been deeply studied [START_REF] Gendelman | Attractors of harmonically forced linear oscillator with attached nonlinear energy sink i: description of response regimes[END_REF]. For the purpose of optimization, a tuning procedure to provide the best total system energy suppression is proposed for cubic NES in a range of frequency around the natural frequency of main structure [START_REF] Starosvetsky | Strongly modulated response in forced 2dof oscillatory system with essential mass and potential asymmetry[END_REF] and it is concluded that the quasiperiodic beating regimes can be more efficient than the steady state responses. As a complement, an optimization design procedure is proposed in [START_REF] Gourc | Experimental investigation and design optimization of targeted energy transfer under periodic forcing[END_REF] and its experimental results have validated the aforementioned conclusion by Starosvetsky. Its optimization objective is to avoid the occurrence of detached resonance curve with high amplitude and enlarge the zone of SMR as wide as possible.

From the above analysis, the periodic excitation with one specific frequency is not specially studied. In addition, the optimization mechanism behind different types of excitation and their relation should be studied as a whole. The abovementioned optimization problems also exist for vibro-impact (VI) NES.

VI NES, also named as impact damper that is considerably studied over sixty years [START_REF] Ibrahim | Vibro-Impact Dynamics: Modeling, Mapping and Applications[END_REF], is referred as a ball which is attached to a main structure and can freely move inside its clearance. Its excessive energy can be irreversible transferred and dissipated by consecutive impacts , which is studied like cubic NES in the context of TET [START_REF] Nucera | Targeted energy transfers in vibro-impact oscillators for seismic mitigation[END_REF][START_REF] Nucera | Application of broadband nonlinear targeted energy transfers for seismic mitigation of a shear frame: Experimental results[END_REF][START_REF] Lee | Periodic orbits, damped transitions and targeted energy transfers in oscillators with vibro-impact attachments[END_REF][START_REF] Gendelman | Analytic treatment of a system with a vibroimpact nonlinear energy sink[END_REF][START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF][START_REF] Gendelman | Dynamics of forced system with vibro-impact energy sink[END_REF]. Its optimization study under transient and periodic excitation is also developed respectively.

It is observed that the response regime with two impacts per cycle of main structure is the most efficient [START_REF] Ibrahim | Vibro-Impact Dynamics: Modeling, Mapping and Applications[END_REF]. For transient excitation, it is showed that the free response starting with two impacts per cycle or regimes closed to this regime is more efficient compared to other cases [START_REF] Brown | Vibro-impact attachments as shock absorbers[END_REF]. Similar conclusion is obtained [START_REF] Lee | Periodic orbits, damped transitions and targeted energy transfers in oscillators with vibro-impact attachments[END_REF]22] that TET is the most efficient with a medium clearance when special orbits are activated as the case of cubic NES. During the whole process of TET, transient resonance capture plays the role of energy reduction. For periodic excitation with a range of frequency, charts are developed for engineering applications within the resonance frequency range [START_REF] Bapat | Single unit impact damper in free and forced vibration[END_REF], in which a medium optimal gap (clearance) is obtained for a specific frequency and it is not necessarily optimal at other frequencies. However, the underlying mechanism to find this optimal gap cannot be explained because of the limitation of the used method and is also not experimentally validated. In this same paper, it is found that the optimal value of design parameters for free and forced vibrations differs. This problem also exists for cubic NES. A design procedure [START_REF] Popplewell | A simple design procedure for optimum impact dampers[END_REF] is proposed and the suggested optimum regime is actually corresponding to the minimal point in the slow invariant manifold (SIM) [START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF]. There are some discrepancy between analytical calculations and numerical results because that the proposed analytical solution may not be stable [START_REF] Li | On the dynamics around targeted energy transfer for vibro-impact nonlinear energy sink[END_REF]. The edge of chaos, i.e. the limit between two impacts per cycle and chaos by intermittency which is called chaotic SMR [START_REF] Gendelman | Dynamics of forced system with vibro-impact energy sink[END_REF], is found to be the optimal response regime for self-excited vibrations with impact damper [START_REF] Yoshitake | Quenching of self-excited vibrations by impact damper[END_REF]. Moreover, periodic response regimes and SMR are experimentally observed [START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF][START_REF] Pennisi | Vibro-Impact NES: A Correlation Between Experimental Investigation and Analytical Description[END_REF]. TET through SMR in the vicinity of resonance frequency is also obtained and it is proposed to increase the area of SMR for maximal energy reduction. Therefore, the same difficulty exists for both cubic NES and VI NES.

To exploit the possible same mechanism of optimization for different types of excitation, the same underlying characteristic needs to be found. A method with FEP [START_REF] Kerschen | Irreversible passive energy transfer in coupled oscillators with essential nonlinearity[END_REF][START_REF] Lee | Periodic orbits, damped transitions and targeted energy transfers in oscillators with vibro-impact attachments[END_REF] is used to study the underlying Hamiltonian system, and the results are efficient for week damping. Considering damping, the multiple scales method can be directly used and is proved effective. Based on the already developed asymptotic method [START_REF] Gendelman | Analytic treatment of a system with a vibroimpact nonlinear energy sink[END_REF][START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF][START_REF] Gendelman | Dynamics of forced system with vibro-impact energy sink[END_REF], a general equation governing the variation of motion, e.g., the resonance part of SMR, is developed for both cubic and VI NES [START_REF] Li | Dynamics of cubic and vibro-impact nonlinear energy sink: Analytical, numerical, and experimental analysis[END_REF], in which the interaction force between NES and main structure plays the central role in explaining the variation of different response regimes. As a further study of VI NES, the study of the influence of parameters on dynamics is improved [START_REF] Li | On the dynamics around targeted energy transfer for vibro-impact nonlinear energy sink[END_REF] from the previous study of different parameters in a separate way to the study of different parameters as a whole. By this way, their mutual influence and the same underlying mechanism is clearer.

From the above analysis of optimization researches about cubic NES and VI NES, their optimization studies have been extensively carried out under different excitations respectively. Different design criteria are proposed for specific systems and are relatively mature. Therefore, the main objective of this paper is to study the efficiency of different response regimes in the first place and then to establish the relation of optimization mechanism for different types of excitation.

In this paper, a VI NES coupled to a linear oscillator (LO) is studied under periodic and transient excitations. In Sec. 2, the equation of motion is analytically studied with the multiple scales method. Then, numerical results are demonstrated and the relation between different excitations is analyzed. In Sec. 4, the optimization mechanism is experimentally verified under different excitation types. Finally, the conclusion is addressed. The analytical development of a system comprising a VI NES and a LO as presented in Fig. 1 is carried out here for both transient and periodic excitation. Since the transient case is a special case of periodic excitation (i.e. zero amplitude of outside force), only the latter is analyzed. The objective here is to apply the asymptotic method used in a series of papers [START_REF] Gendelman | Analytic treatment of a system with a vibroimpact nonlinear energy sink[END_REF][START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF][START_REF] Gendelman | Dynamics of forced system with vibro-impact energy sink[END_REF][START_REF] Li | On the dynamics around targeted energy transfer for vibro-impact nonlinear energy sink[END_REF][START_REF] Li | Dynamics of cubic and vibro-impact nonlinear energy sink: Analytical, numerical, and experimental analysis[END_REF] to obtain the possible optimization design criterion of the system with damping rather than its underlying Hamiltonian system. Its motion between impacts is described by the following equation:

ẍ + ελ 1 ẋ + x = εG sin Ω τ + ε 2 λ 1 GΩ cos Ω τ ε ÿ = 0 ∀ |x -y| < b (1)
The corresponding physical parameters are expressed as follows:

ε = m 2 m 1 , ω 0 2 = k 1 m 1 , f 0 = ω 0 2π , τ = ω 0 t, λ 1 = c 1 m 2 ω 0 , Ω = ω ω 0 , G = F ε
where x, m 1 , c 1 and k 1 are the displacement, mass, damping and stiffness of the LO respectively. y and m 2 are displacement and mass of VI NES. The dots denote differentiation with respect to dimensionless time τ. b represents the clearance. x e (t) is the displacement imposed on the base by shaker. εG sin Ω τ and ε 2 λ 1 GΩ cos Ω τ represent the force contribution of displacement and that of velocity of excitation respectively imposed on LO. The latter term is conserved here to demonstrate its physical meaning and is ne-glected during the following analysis because of its small magnitude.

When |x -y| = b, impacts occur. The relation between after and before impact is obtained under the hypothesis of the simplified shock theory and the condition of total momentum conservation:

x + = x -, y + = y - ẋ+ + ε ẏ+ = ẋ-+ ε ẏ-, ẋ+ -ẏ+ = -r ẋ--ẏ-, for |x -y| = b (2)
where r is the restitution coefficient and the superscripts + andindicate time immediately after and before impact. New variables representing the displacement of the center of mass and the internal displacement of the VI NES are introduced as follows:

v = x + εy, w = x -y (3)
Substituting Eq. (3) into Eqs. ( 1) and ( 2), the equation between impacts in barycentric coordinate are given as:

v + ελ 1 v + ε ẇ 1 + ε + v + εw 1 + ε = εG sin Ω τ + ε 2 λ 1 GΩ cos Ω τ ẅ + ελ 1 v + ε ẇ 1 + ε + v + εw 1 + ε = εG sin Ω τ + ε 2 λ 1 GΩ cos Ω τ ∀|w| < b (4)
and the impact condition (2) is rewritten as:

v + = v -, w + = w -, v+ = v-, ẇ+ = -r ẇ-, for |w| = b (5)
Multiple scales are introduced in the following form:

v(τ; ε) = v 0 (τ 0 , τ 1 , . . .) + εv 1 (τ 0 , τ 1 , . . .) + . . . w(τ; ε) = w 0 (τ 0 , τ 1 , . . .) + εw 1 (τ 0 , τ 1 , . . .) + . . . τ k = ε k τ, k = 0, 1, . . . (6)
A detuning parameter (σ ) representing the nearness of the forcing frequency Ω to the simplified natural frequency of the LO is introduced:

Ω = 1 + εσ (7)
Substituting Eqs. ( 6) and [START_REF] Kerschen | Irreversible passive energy transfer in coupled oscillators with essential nonlinearity[END_REF] into Eqs. ( 4) and ( 5), equating coefficients of like power of ε and only conserving the first two orders:

Order

ε 0 : z Π(z) M(z)
Fig. 2 Representation of the non-smooth functions Π (z) and M(z)

D 2 0 v 0 + v 0 = 0 D 2 0 w 0 + v 0 = 0, ∀|w 0 | < b (8) v + 0 = v - 0 , w + 0 = w - 0 , D 0 v + 0 = D 0 v - 0 , D 0 w + 0 = -rD 0 w - 0 , for |w 0 | = b (9)
Order ε 1 :

D 2 0 v 1 + v 1 = -2D 0 D 1 v 0 -λ 1 D 0 v 0 -w 0 + v 0 +G sin (τ 0 + σ τ 1 ) (10)
For order ε 1 , only the term related to LO is conserved and will be used later. Combining the first and second order, the SIM and fixed points could be obtained. From the analysis of order ε 0 , v 0 represents an ideal undamped harmonic oscillator expressed as follow:

v 0 = C(τ 1 ) sin (τ 0 + θ (τ 1 )) (11)
From the standpoint of w 0 , Eq. 8 and Eq. 9 represent a harmonically forced impact oscillator with symmetric barrier. For the regime of 1:1 resonance with two impacts per cycle, its solution can be searched in the following form:

w 0 = C(τ 1 ) sin (τ 0 + θ (τ 1 )) + 2 π B(τ 1 )Π (τ 0 + η(τ 1 )) (12)
Where Π (z) is a non-smooth sawtooth function [START_REF] Pilipchuk | Closed-form solutions for oscillators with inelastic impacts[END_REF]. This folded function and its derivative are depicted in Fig. 2 and are expressed as follows: According to Eq. ( 12) and ( 13), impact occurs at T 0 = π/2-η + jπ with j = 0, 1, 2, . . .. The impact condition |w 0 | = b is rewritten with Eq. ( 12) as:

Π (z) = arcsin(sin z), M(z) = dΠ dz = sgn(cos z) (13)
C cos (η -θ ) = b -B (14)
Rewriting now the inelastic impact condition (9) yields:

C (1 + r) sin (η -θ ) = 2 π B (1 -r) (15)
Combining Eqs. ( 14) and [START_REF] Nucera | Targeted energy transfers in vibro-impact oscillators for seismic mitigation[END_REF], a relation between B and C is obtained as follows:

C 2 = 1 + 4 (1 -r) 2 π 2 (1 + r) 2 B 2 -2bB + b 2 (16)
An example of SIM with b = 1 and r = 0.6 is presented in Fig. 3. The stability of SIM is analyzed by the approach used during the asymptotic analysis in [START_REF] Gourc | Targeted energy transfer under harmonic forcing with a vibro-impact nonlinear energy sink: analytical and experimental developments[END_REF] that is originally used in [START_REF] Masri | On the stability of the impact damper[END_REF][START_REF] Bapat | Stable periodic motions of an impact-pair[END_REF]. By this stable analysis method, the stable branch is defined by the condition that the modulus of all the eigenvalues of a certain matrix relating two consecutive impacts is less than unity. This stability analysis can also be accomplished by direct numerical integration of Eqs. ( 8) and [START_REF] Manevitch | Parameters optimization for energy pumping in strongly nonhomogeneous 2 dof system[END_REF].

In order to obtain the fixed points or study the nonstationary motion evolution of system on the SIM, i.e. the analysis of the transient resonance capture under transient excitation and the transient resonance transition of intermittent chaos (i.e. SMR), Eq. ( 10) at the next order of approximation is analyzed. To identify terms that produce secular terms, the function of w 0 is expanded in Fourier series in the following way:

w 0 =C(τ 1 ) sin (τ 0 + θ (τ 1 )) + E (τ 1 ) sin (τ 0 + ζ (T 1 )) + RFC (17)
where RFC represents the rest frequency components except the simplified natural frequency of LO. The component C (τ 1 ) is decided by the motion of LO and E (τ 1 ) is totally related to the motion of VI NES. E (τ 1 ) is also related to the periodic impact force. The Eq. ( 17) is a more general and relaxed analytical description with respect to the motion of VI NES compared to the former studies of 1:1 resonance with two symmetric impacts per cycle.

Substituting Eqs. ( 11), ( 12) and ( 17) into Eq. ( 10) and eliminating terms that produce secular terms give:

D 1 C = - 1 2 λ 1 C - 1 2 E sin (Θ ) + 1 2 G sin (η) D 1 η = 1 2 G cos (η) /C - 1 2 E cos (Θ ) /C + σ (18)
where

Θ = ζ -θ η = σ τ 1 -θ (19)
Θ represents the phase difference related to LO and VI NES. η represents the phase difference related to LO and outside excitation.

The fixed points can be obtained by equating the left side of Eq. ( 18) to zero and then combining it with Eq. [START_REF] Nucera | Targeted energy transfers in vibro-impact oscillators for seismic mitigation[END_REF]. In this way, the fixed points (number and position) can be obtained. Compared to the classic non-asymptotic method, its functionality is twofold. Firstly, its position can be used to judge the type of response regime. Secondly, its value can be precisely calculated for 1:1 resonance, which may be related to optimal response regime. Meanwhile, the corresponding equivalent force E between LO and VI NES can be used to analyze the performance. The above analysis will be studied further in the next section by numerical simulation.

Optimal TET and response regimes by numerical simulations

For engineering application, the length of cavity is more controllable than other parameters. Therefore, it will be chosen as the design parameter. About the possible response regime, it has been studied and is quantized by the number of impacts per cycle (i.e. z) as demonstrated in Fig. 4. The displacement of LO (i.e. x), the relative displacement (i.e. w) and the projection of motion of system into SIM are demonstrated for different response regimes in the first, second and third column respectively. For a specific system described by Fig. 1 under periodic excitation, the following response regime will occur consecutively with the increase of the length of cavity (i.e. b): (a) chaos with infinite number of impact per cycle, i.e. the regime with more than two impacts per cycle with z > 2, (b) two asymmetric impacts per cycle with z = 2, (c) two symmetric impacts per cycle with z = 2, (d) chaos by intermittency, i.e. SMR with z < 2, (e) chaos with no duration of two impacts per cycle. The characteristic is clearly showed in the second column in Fig. 4. Of course, the above regimes are not exhaustive. For the case of transient excitation, the above regime will occur one after another for a sufficient high initial energy except SMR. With the knowledge of variation mechanism of response regimes, which response regime is most optimal for periodic excitation and which response regime should be set as the starting point for transient excitation will be studied in this section. Therefore, the system under periodic and transient excitation will be studied consecutively to see the optimal efficiency of TET and its corresponding response regimes by varying the length of cavity b. Meanwhile, the important role of SIM will be demonstrated.

For all simulation, the following parameters are fixed except that it is specially pointed out: λ 1 = 1.43, ε = 0.84, r = 0.6.

Periodic excitation with one single frequency

Parameters and initial conditions for periodic excitation are fixed as follows: G = 0.02, σ = 0, x 0 = 0, ẋ0 = 0, y 0 = 0, ẏ0 = 0.

For resonance excitation, the corresponding amplitude of displacement of LO possessing the same energy during stable periods is calculated for both steady and not stationary response. The ratio of this amplitude with and without VI NES is a central and comparable index. This index is denoted as Ae. Meanwhile the ratio between the maximal amplitude of LO with VI NES and that of LO without VI NES is termed as Am. The difference between the above two indices can be distinguished from the displacement of LO for SMR, in which the amplitude is no longer stable.

With the increase of the length of cavity b, different response regimes have been observed and classed according to their efficiency as showed in Fig. 5: area A with more than two impacts per cycle, area B with two impacts per cycle (asymmetric or symmetric), area C as SMR with partial duration of two impacts per cycle, area D as chaos without any duration with two impacts per cycle. It is found that only the response regime with two impacts per cycle is steady with the same Am as Ae. It is observed that the optimal value lies at the boundary between resonance regime with two impacts per cycle and SMR. The following conclusion can be made: the response regime with two impacts per cycle is more efficiency than the case with more than two impacts per cycle. However, the result of efficiency comparison of SMR and the two impacts per cycle depends on the specific conditions. For the same case with two impacts per cycle, the symmetric case is more efficient than the asymmetric case.

To better understand the above numerical results, the analytical results will be demonstrated here. Specifically, the fixed points obtained from Eq. ( 16) and Eq. ( 18) (the stable state of this equation) and its variation in the SIM will be used to judge the response regime and the precise calculation of the value for system design.

With the analytical results in the second section, two fixed points (circle for stable point and cross for unstable point) can be obtained. According to its position, three different areas are observed as showed in Fig. 6: with more than two impacts per cycle for fixed point in area A, with two impacts per cycle (full duration or partial duration) in area B and no duration of two impacts per cycle for area D.

For small value of b, the response regime with more than two impacts per cycle occurs and the localization of fixed points is showed in Fig. 6(a), actually these two points only theoretically exists and cannot precisely describe the actual case, in which the amplitude usually is not steady.

With the increase of b, relative locations between fixed points and SIM will change. In the case showed in Fig. 6(b), SIM will move up and results in the move of one fixed point from area A to B and another fixed point in the red unstable branch. In this case, the response regime with two asymmetric and symmetric impacts per cycle will consecutively occurs until one fixed point meet the minimal point in the SIM. Meanwhile, the efficiency of energy reduction is increased during this process as showed by the preceding numerical results. One important thing is that three important points can be analytically obtained: the boundary point between area A and area B, the boundary point between two symmetric and asymmetric impacts per cycle and the minimal point at SIM.

When one fixed point passes the minimal value of SIM, two cases can exist and denoted by circle for case 1 and ellipse for case 2 as showed in Fig. 6(c) and SMR will occur. For case 1, the stable point passes from the thick stable branch to the thin unstable branch. In this case, the optimization point is localized in the minimal point of SIM. In case 2, the unstable point in the thin unstable branch will move into the thick stable branch and finally meet with the other fixed point. In this case, SMR will be more efficient than the case with two impacts per cycle. Anyway, the above two cases can be analytically calculated. In [START_REF] Li | On the dynamics around targeted energy transfer for vibro-impact nonlinear energy sink[END_REF], the influence of frequency of excitation on the existence of these two cases is explained and it is showed that the case 1 applies for the excitation with the natural frequency of LO or even high frequency. A further increase of b, these two fixed points will meet and then disappear, there will be no fixed points as showed in Fig. 6(d) and then chaos occurs.

From the above analysis, the variation of E can be conversely inferred from the variation of C with the variation of b. For the optimal case with minimal C, E related to interaction force is the largest. Moreover, its variation mechanism is inverse to that of C. For SMR, the transient interaction force E can be inferred from the corresponding steady value and can be used to explain the variation mechanism of SMR as has been done in [START_REF] Li | Dynamics of cubic and vibro-impact nonlinear energy sink: Analytical, numerical, and experimental analysis[END_REF]. In this paper, this value will be used to establish the bridge between transient and periodic excitation.

In summary, the efficiency of different response regimes is numerically obtained and analytically explained. The optimal response regime is the limit between that with two impacts per cycle and SMR for the excitation with the natural frequency of LO.

Transient excitation

In order to compare with the periodic excitation, the transient excitation for the same system is studied. Eq. ( 1 With the variation of b, one response type with b = 0.01 is chosen to demonstrate all possible stages during the whole process. Different tools are used to explain the dynamics from different points of view as showed in Figs. 7 and8. The displacement of LO with and without VI NES, displacement of VI NES and relative displacement are demonstrated respetively from Fig. 7(a-c). In Fig. 7(a), the efficiency of VI NES is clearly observed by the difference between blue and red curve. In Figs. 7(b) and (c), the transient resonance between LO and VI NES is clearly showed. In Fig. 7(d), the motion is projected to SIM and the ability of SIM to govern the transient resonance is showed.

Evidently, the response undergoes different stages during the whole process, which is showed by the relative displacement during the whole process as showed in Fig. 8(a). At the starting place, the regime with more than two impacts per cycle (i.e. z > 2) is excited as showed in Fig. 8(b). Then the two asymmetric and symmetric impacts per cycle consecutively appears in the stage 2 as demonstrated in Fig. 8(c). In Fig. 8(d), the stage with sparse impacts occurs when the energy of system is low. The above three parts are similar to stage A, B and D in the case of periodic excitation just without SMR and the change direction of regimes is fixed. The instantaneous energy distribution is showed in Fig. 8(e). From about τ = 100 to τ = 200, there exits strongly energy exchange. The accumulated energy dissipated by LO and VI NES is showed in Fig. 8(f). The above mentioned period of strong energy exchange can be observed from different perspectives. The energy decay inclination of LO is almost linear in Fig. 7(a), the projection of motion is controlled by SIM at these moments demonstrated in Fig. 7(d), the response regime is two impacts per cycle showed by the relative displacement in Fig. 8(c). In Fig. 8(f), the accumulated energy dissipated by the damp- It is observed that the most efficient case is obtained when the response starts from the stage 2, i.e. two impacts per cycle, which is well demonstrated by Fig. 9(b) with b = 0.05, the decay inclination is linear that can be seen from the variation of displacement denoted by the blue curve in its middle subfigure. It should be pointed out that the value of b is so intentionally chosen that these limited points can almost represent the typical starting point in the SIM. It is specially noticed that the starting point of the above optimal response for transient excitation is different from that of the optimal response regime for periodic excitation. The conclusion is that the optimal starting point should be located a little away from the lowest point in SIM and the specific distance remains unknown and will be discussed later.

Discussion

In section 2, the general Eq. ( 17) governing the variation of motion is developed to further explain the transient transition process of response. This kind of process exists for SMR and response under transient excitation. In section 3.1 and 3.2, the optimization mechanism of TET for periodic excitation with one frequency and transient excitation is obtained. The frequency of periodic excitation is chosen the same as the natural frequency of LO in order that the optimization mechanism under these two different excitations is comparable. In this section, the relation between transient and periodic excitation about the mechanism of optimization is established at the first place. Then the obtained results are used to explain the existing problems in the previous studies or to provide possible design criteria for the future study.

Transient and periodic excitation

For periodic excitation, the steady state response regime is identified by two impacts per cycle. The value of E in Eq. ( 17) can be analytically calculated and is related to one stable point in the blue stable branch in SIM. Of course, every point in the stable branch relates to a value of E and can be analytically obtained, though it may not actually exist. From the variation law of efficiency for different lengths of cavity under fixed periodic excitation (i.e. constant G), the value of C becomes smaller as the stable fixed point moves from up to down in the right branch. The sum of right hand side of Eq. ( 17) should be zero for the points in the branch. Therefore the value of E increases until the maximal value corresponding to the minimal point of SIM. A further increase of length of cavity b from the case of two impacts per cycle, the SMR showed in Fig. 4(d) appears and the corresponding actual value of E varies around an analytical calculated value during the whole process.

For transient excitation, there does not exist SMR. Response regime continuously transits from one type to another type with the decrease of C. However, the transient response can be approximately described by the steady state type, i.e. the above interaction force can be applied. Its effectiveness can be seen from the precise description of SIM for the resonance parts of the motion projection of SMR and transient response. We can imagine an outside force for any moment during the whole process as if the outside excitation is abruptly added to hold this transient regime. The equivalent value of E can be obtained by the analytical results in the section 2. However, this value of the imagined force will not be constant. If the original starting place is set to the minimal point of SIM, the value of E will be the largest with a corresponding largest outside force. Yet, its duration will be very short. In contrast, if the starting point in SIM is designed to a very high point in the right red branch, its efficiency will also be very low considering the low efficiency at the start. Therefore, the design criterion is that the efficiency of TET should be not only high at the beginning but also last as long as possible. Theoretically, this optimal point can be obtained, but its analytical calculation should be very complex. Although it is still not analytically obtained in this paper, its mechanism is further explained. For the former study about optimization mechanism of TET, most work is concentrated in the activation of TET. However, only the condition of activation is not enough for optimal objective because it should also possess a long enough duration. In [START_REF] Nguyen | Design criteria for optimally tuned nonlinear energy sinkspart 1: transient regime[END_REF], one practical design criterion for cubic NES is proposed with applicability though may not optimal. But the underlying mechanism is not developed in detail. The design contradiction claimed in [START_REF] Bapat | Single unit impact damper in free and forced vibration[END_REF], i.e. different values for the same parameter under transient and periodic excitation, can also

Periodic excitation with one frequency and with a range of frequency

To study the optimization mechanism for excitations with a range of frequency in the vicinity of the natural frequency of a LO, the optimization mechanism under one excitation with a specific frequency should be studied at the first place. When the frequency response function of LO is obtained during a range of frequency, the amplitude at different frequency points is different. The value at the point of natural frequency is maximal. If an optimal value of b is chosen according to the above analysis for the resonance point, the problem is that the value of b will not be still optimal for the other points and will result in intermittency chaos (i.e. SMR) or chaos without duration of two impacts per cycle. If this value is optimally designed for other frequencies, the same dilemma occurs. Therefore, there does not exist one value of b optimal for all frequency. Although the value of b can be chosen to make the SMR appear for all frequency, it does not mean that this value is optimal. If the objective is to control the resonance value of LO, the problem is simplified to the case with one frequency. This design objective is reasonable and practical.

Experimental validation and analysis

The objective of this part is to experimentally verify the obtained optimization mechanism in the preceding section. For periodic excitation, the focus is to validate the periodic case with one single frequency at the first place and then to verify the results under a range of frequency. As for transient excitation, the study with experimental validation is already partly done by Bapat [START_REF] Bapat | Single unit impact damper in free and forced vibration[END_REF], the experimental results here will provide more details about the efficiency of different transient response regimes and will be compared to the former study results.

Periodic excitation

Experimental setup

The global experimental configuration is showed in Fig. 11(a). A ball is put inside a cavity of LO and can freely move inside. The whole system is embedded on 10 kN electrody- The cylinder and the other side cover are made of hardened steel. The physical parameters of the system have been identified by performing modal analysis and are summarized in Table 1.

Single frequency excitation

The frequency of excitation slowly varies from 7.82 Hz to 7.84 Hz, which can be considered almost fixed to the value 7.83 Hz. This value is closed to the natural frequency f 0 = 7.86 Hz. Its acceleration is fixed to 0.06 g. The response without VI NES is recorded as reference and the results with VI NES is obtained by varying the value of b from 5-50 mm. It should be pointed out here the unit of b expressed in the experimental results is mm which is different from m for numerical results. The results are showed in Fig. 12. The time history of displacement of LO is showed, but only the stable period regardless of the starting and ending parts is used for calculation and comparison. Specifically, the time period 20-70 s is considered. In Fig. 12 This kind of chaotic SMR [START_REF] Gendelman | Dynamics of forced system with vibro-impact energy sink[END_REF] is experimentally validated from the point view of VI NES by the existence of intermittent strong impacts [START_REF] Li | Dynamics of cubic and vibro-impact nonlinear energy sink: Analytical, numerical, and experimental analysis[END_REF]. Of course, it can be described directly by the strongly modulated envelope of x as showed here. To compare the efficiency of different b, the data during 20-70 s denoted by blue dotted rectangle is used to calculate the equivalent amplitude (Ae) possessing the same energy, which is represented by the black dotted line. For not steady state response and specially SMR, the maximal amplitude (Am) during this process is recorded which is the vertical value of point C in this case. The calculated results are showed in Fig. 12(d) and the variation law observed in Fig. 5 is verified. One difference is that the numerically observed area D is not experimentally observed because of the limit of the length of cavity. The efficiency of the regime with two impacts per cycle (b = 30 mm) and that of the SMR (b = 32.5 mm) are almost the optimal value and the optimal value b can be inferred to be located between them. The obtained value Ae = 40.16 % is very interesting from the point of energy reduction.

Frequency band excitation

To study the influence of b around the natural frequency of LO and the possibility of application of the above observed mechanism, the system periodically excited under a range of frequency is studied. The starting frequency f s and ending frequency f e during sweep is showed in Table 1 and the acceleration is still fixed to 0.06 g.

The displacement of LO is recorded for different b and is showed in Fig. 13. The shift of resonance peak between different experimental results comes from the difference of the starting record time. The maximal value of every sweep is indicated by different values of Y as showed in Fig. 13(a) and (b). For the frequency around resonance, the same mechanism described in the former section is well demonstrated by the value of Y. Its performance is also optimal for a medium value of b. The transition process of response for different value of b close to the optimal response is showed in Figs. 13(c-e). The optimal value locates between b = 27.5 mm (i.e. two impacts per cycle) and b = 30 mm (i.e. SMR). For b > 30 mm, the efficiency will decrease for all frequency with the transition of response regime from SMR with long duration of two impacts per cycle to that with short duration and until totally no two impacts per cycle. In return, the efficiency for the resonance frequency will decrease with the decrease of b, but the response regime in other excitation frequencies will shift from SMR to two impacts per cycle and even to chaos (i.e. many impacts per cycle and different from the above mentioned one), which means that their efficiency will increase at first and then decrease. In this way, the optimal mechanism for a range of frequency is clear. If the optimization objective around resonance frequency is to control the maximal amplitude, the system can be optimized at the resonance frequency according to the above optimization mechanism. In this optimal case, SMR will occur for other frequency close to resonance frequency and its duration of two impacts per cycle will increase with the decrease of the distance of this frequency from resonance frequency. Therefore, the optimal design for a range of frequency is directly related to that of resonance frequency. Although it is optimally designed for a range of frequency, VI NES will not work at its optimal state for other frequency except resonance frequency.

Transient excitation

Experimental setup

The same LO coupled with VI NES in Fig. 11 is used for the transient case as showed in Fig. 14 , so the parameters are almost the same as showed in Table 1 except a little addition of mass owing to the ring bolt. The LO is installed to a cast iron bench. One laser sensor and one accelerometer are used to measure the displacement and acceleration of LO respectively, meanwhile impact moments can be detected. For all tests, the initial displacement of LO is drawn to around 20 mm which is realized by a stretching device displayed in Fig. 14(b). The location of the ball is random inside the cavity and it is supposed that it will not influence the stable response. The initial velocity of ball and LO is zero. Only the length of cavity is varied to see the transition of response regime for one specific length of cavity and further to compare the efficiency of different lengths of cavity.

Efficiency of different transient response regimes

The response of LO with b = 15 mm is demonstrated in Fig. 15. The time history of displacement is displayed in Fig. 15(a) and it is observed that the displacement is not symmetric at the initial stage until the time around t1 = 0.9254 s. However, the acceleration in this period is symmetric showed in Fig. 15(b). The pulse value of acceleration denotes an impact moment. According to the number of impacts during one cycle of LO, different response regimes can be detected. An enlarged view of acceleration is showed in Fig. 15(c). The transition of impact number (i.e. z) can be clearly identified. The area 1,2 and 3 are enlarged in Fig. 15(d-f) respectively. In area 1, the regime with three impacts per cycle disappears. In area 2, the regime with two impacts per cycle escapes into one impact during some cycles and will come back to the regime with two impacts per cycle. In area 5, the regime with two impacts per cycle escapes and cannot come back again. The different marked points correspond to different areas respectively. From the regime with three impacts per cycle to the regime with no impacts per cycle, different response regimes exist as represented in the SIM. To compare their efficiency, the decay rate of displacement amplitude during the period [t1,t2] displayed in Fig. 15(a) is used.

An enlarged view of Fig. 15(a) is showed in Fig. 16(a). The local minimal displacements are marked out and represented by a i and plotted in Fig. 16(b), which also represent its envelop. The variation of absolute value of difference (i.e. da = a i+1a i ) is demonstrated in Fig. 16(c). It is observed that the value of da varies during the whole process and cannot be approximated by a constant value. This means that the inclination rate is more complex than expected and cannot be described as linear, though it is truly different to the case without VI NES and is closer to a linear decay rate. Therefore, the index of linear decay rate in the former study of impact damper may cause misunderstanding and it also should be noticed for other kinds of NES. Then the relative decay rate δ is used to compare different response regime. Specifically, δ is defined and calculated as follows:

δ = da * 2 a i+1 + a i (20)
The results are showed in Fig. 16(d). In zone 1, regime with three impacts per cycle exists. In zone 2 and 4, regime with two impacts per cycle exists. In zone 3, occasional short escape of two impacts per cycle appears. In zone 5, almost no impacts occur. The blue broken polyline can approximately describe the variation tendency of the efficiency of different response regimes. Except the phenomenon of occasional out of the excitation of VI NES, the obtained efficiency law agrees well with the former conclusion from the study of periodic excitation. The occasional out of excitation is probably related to the friction between ball and LO and it is not considered in our study. It means that the experimental results can just cover and validate the results in the zone A and B marked out in Fig. 10 because of the limitation of b.

Conclusion

A LO coupled with a VI NES is analytically studied with multiple scales method in Sec. 2. A general equation is developed to explain SMR and transient response. In Sec. 3, the optimization mechanism is studied under periodic and transient excitation respectively by numerical simulation. The optimal response regime is obtained for periodic excitation with one single frequency and the starting response regime is obtained for transient excitation. The obtained results are used to explain the optimization relation between different excitation types. In the next section, the obtained optimization mechanism under periodic with one single and a range of frequency and under transient excitation are experimentally validated.

The interaction force related to impacts is analytically developed and is used as the bridge to connect the steady state response and nonstationary response under periodic and transient excitation. The response regime at the boundary between the regime with two impacts per cycle and SMR is numerically and experimentally proved to be the most efficient response regime for the periodic excitation, which is contrary to the traditional impression that SMR is the more efficient. The conclusion also applies for other frequency since the variation mechanism of response regime should be the same. Moreover, this optimal value can be analytically calculated and is practical for engineering applications. The optimization mechanism and optimal design under transient excitation is clarified by that mechanism under periodic excitation. Its requirement is that design parameters should be chosen to make not only the response regime with high TET exist at the beginning but also last as long as possible. The optimization objective with the classic method to test the efficiency of VI NES under periodic excitation with a range of frequency is clarified again and the reason why there is no optimal parameters (e.g., b) for all frequency is also explained.

By experimental validation, if VI NES is designed optimal for LO under the excitation with resonance frequency, the performance will be optimal as a whole for a range frequency around resonance frequency. From experimental results under transient excitation coupled with VI NES, it is observed that the traditional concept of the linear decay rate of displacement amplitude of main structure does not apply and the actual situation is more complex than a linear description, which should be noted for other types of NES. The efficiency of different transient response regimes is also compared from a new point of view and the conclusion agrees with that of periodic excitation.

Although the theoretical, numerical and experimental studies in this paper are not quantitatively connected, the results are enough to accomplish its objectives. Considering that VI NES is studied under the context of TET, the above conclusion may be generalized to other types of NES.

Fig. 1

 1 Fig. 1 Schema of a LO coupled with a VI NES under periodic excitation

B C 2 Fig. 3

 23 Fig. 3 SIM of VI NES: one stable branch in blue thick line and two unstable branches in red thin line

Fig. 4

 4 Fig. 4 Displacement of LO in the first column, relative displacement in the second column and motion of system projected into SIM in the third column for different response regimes [25]: (a) z > 2; (b) z = 2 and asymmetric; (c) z = 2 and symmetric; (d) z = 2 and z < 2 for SMR; (e) z < 2

Fig. 5 Fig. 6 Fig. 7 A

 567 Fig. 5 Efficiency of TET for different response regimes under periodic excitation: blue asterisk represents the average amplitude of LO and red upside-down triangle represents the maximal amplitude of LO

Fig. 8 A

 8 Fig. 8 A typical response with b = 0.01: (a) the whole duration of the relative displacement; (b) stage 1 with z > 2; (c) stage 2 with z = 2; (d) stage 3 with z < 2 ; (e) the instantaneous energy distribution between LO and VI NES; (f) the energy dissipated by LO and VI NES

Fig. 9

 9 Fig. 9 Responses starting from different response regimes by the variation of b: (a) start from stage 1 with b = 0.01; (b) start from stage 2 with b = 0.05; (c) start from stage 3 with b = 0.1

Fig. 10

 10 Fig. 10 Efficiency comparison of different response regimes with the variation of b: percent of energy dissipated by VI NES represented by red cross and time needed for the 95% energy reduction of LO denoted by blue top-down triangle

Fig. 11

 11 Fig. 11 Experimental setup: (a) the global configuration; (b) a detailed view of VI NES

Fig. 12

 12 Fig. 12 Displacement of LO and TET efficiency of VI NES under periodic excitation with a fixed frequency and different b: (a) small b with z > 2 or z = 2; (b) large b with z < 2 (c) b = 35mm with SMR; (d) TET efficiency

Fig. 13

 13 Fig. 13 The displacement history of the LO and the TET efficiency of VI NES under the periodic excitation with a range of frequency around f 0 for different b: (a) a small b with z > 2 or z = 2 at the point of resonance; (b) a large b with z < 2 at the point of resonance; (c) a detailed view for b = 25 mm around the point of resonance; (d) a detailed view for b = 27.5 mm around the point of resonance; (e) a detailed view for b = 30 mm around the point of resonance

f s -f e 6 . 5 - 9

 659 Hz shaker acceleration 0.06 g namic shaker. The displacement of the LO as well as the imposed displacement of the shaker are measured using contactless laser displacement sensors. The acceleration of LO and shaker is measured by an accelerometer. A detailed view of the VI NES is presented in Fig. 11 (b). It simply consists of a closed cavity of length d + 2b, where d is the diameter of the ball and b can be adjusted by a cylinder in the cavity.

Fig. 14

 14 Fig. 14 Experimental setup: (a) the global configuration; (b) the displacement regulation device; (c) a detailed view of VI NES

Fig. 15

 15 Fig. 15 Response of LO with b = 15mm : (a) time history of the displacement; (b) time history of the acceleration; (c) an enlarged view of acceleration; (d) an enlarged view of area1; (e) an enlarged view of area2; (f) an enlarged view of area3

4. 2 . 3

 23 Efficiency of different lengths of cavity Except b = 15 mm, other values of b have been chosen to find the optimal value of b and the corresponding starting response regime. The time history of displacement for b = 20 mm is used to display the envelop extraction of each b and it is showed in Fig. 17

Fig. 16

 16 Fig. 16 Efficiency of different response regimes with b = 15mm : (a) enlarged view of time history of displacement; (b) time history of local minimal displacement; (c) decay rate; (d) relative decay rate

Fig. 17

 17 Fig. 17Efficiency compare with different b: (a) envelop of b = 20 mm; (b) imposed envelops of displacement of LO

  Fig. 17Efficiency compare with different b: (a) envelop of b = 20 mm; (b) imposed envelops of displacement of LO

Table 1

 1 Parameters of the experiment

		Physical Parameters
	m 1	4.7 kg	c 1	3.02 Ns/m
	k 1	11.47 * 10 3 N/m m 2	32 g
	b	0 -50 mm	r	0.6
		Reduced Parameters
	ε	0.76 %	λ 1	1.91
	f 0	7.86 Hz
		Single Frequency Test
		f r	7.82 -7.84 Hz
	shaker acceleration	0.06 g
		Frequency Band Test

Acknowledgements The authors acknowledge the French Ministry of Science and the Chinese Scholarship Council under grant No. 20130449 0063 for their financial support.


Oleksii O Shtehin

Vincent Wagner
email: vincent.wagner@enit.fr

Sébastien Seguy
email: s.seguysebastien.seguy@insa-toulouse.fr

Yann Landon

Gilles Dessein
email: g.desseingilles.dessein@enit.fr

Michel Mousseigne
email: m.mousseignemichel.mousseigne@univ-tlse3.fr

Stability of ball-end milling on warped surface: semi-analytical and experimental analysis

Keywords: chatter, milling, stability lobes, warped surface. Stability of ball end milling on warped surface: semi-analytical and experimental analysis

published or not. The documents may come   L'archive ouverte pluridisciplinaire

Introduction

High-performance machining of crooked spatial surfaces is a complex task, which solution are supremely important for various type of industry like casting, molding, dies production and monolithic components in the aeronautical industry. In most cases, finishing of such surfaces is carried out with ball end mills on CNC machines. Ball end milling of inclined surfaces can be considered as particular case of crooked spatial surfaces machining. This process is widely used, which is caused by competitive prices of machines and accessories and by the wide spectrum of the problems that can be solved with this machining process.

However, the main limitation of this machining method consists in difficulty of tilt and lead angles control. Dynamical properties of ball end milling process are specific and its dynamic behaviour prediction is a knotty question. Particularly, it concerns ensuring process stability.

Chatter vibrations are relevant for all types of cutting. Their genesis consists in regenerative effect. One of the first studies about chatter vibration was applied in the case of orthogonal turning [START_REF] Tobias | Theory of regenerative machine tool chatter[END_REF]. Stability analysis of the corresponding dynamic system leads to the well-known stability lobes. This approach was adapted for milling system where value and orientation of cutting force are changing [START_REF] Altintas | Analytical prediction of stability lobes in milling[END_REF]. On this approach only the mean cutting force is analyzed, this method is called zero order stability solution [START_REF] Budak | Analytical models for high performance milling -Part I: cutting forces, structural deformations and tolerance integrity[END_REF][START_REF] Budak | Analytical models for high performance milling -Part II: process dynamics and stability[END_REF]. However, this analytical method is widely used, faster and easy to adjust in a real process. Thus, this method was applied for ball end milling in 3 axis milling [START_REF] Altintas | Analytical prediction of stability lobes in ball end milling[END_REF]. In this study, the chip thickness variation in both radial and axial directions was considered, and the variations of the cutting force coefficients for each point on the cutting edge were taken into account. This approach was widespread in the third axial direction [START_REF] Altintas | Analytical prediction of three dimensional chatter stability in milling[END_REF]. The obtained three dimensional chatter stability model was extended to 5-axis ball end milling [START_REF] Ozturk | Investigation of lead and tilt angle effects in 5-axis ball-end milling processes[END_REF] and to positioned 5-axis ball end milling [START_REF] Mousseigne | Predicting the dynamic behaviour of torus milling tools when climb milling using the stability lobes theory[END_REF] by adding the effect of lead and tilt angles on the process. An automatic adjustment of tool axis orientations was developed to avoid chatter along the tool path in 5-axis milling [START_REF] Chao | Chatter free tool orientations in 5-axis ball-end milling[END_REF].

On the same way, improvements were made by analysing the exact force expression occurring in milling [START_REF] Davies | Stability prediction for low radial immersion milling[END_REF][START_REF] Bayly | Stability of interrupted cutting by temporal finite element analysis[END_REF][START_REF] Insperger | Stability of up-milling and down-milling, part 1: alternative analytical methods[END_REF][START_REF] Insperger | Updated semi-discretization method for periodic delay-differential equations with discrete delay[END_REF][START_REF] Merdol | Multi frequency solution of chatter stability for low immersion milling[END_REF]. These multi-frequency stability solutions are based on improved mathematical methods for differential equations with delayed terms to model milling. In peripheral milling, with small radial depth of cut, interrupting cutting generates a period doubling instability (flip) [START_REF] Mann | Stability of up-milling and down-milling, part 2: experimental verification[END_REF]. This special cutting condition is observable only at high speed on the first lobes of stability lobes. Flip lobes thus become negligible for high orders. These approaches are also highly effective for fine modelling of the stability dynamic with the cutting tool run-out [START_REF] Insperger | On the chatter frequencies of milling processes with runout[END_REF][START_REF] Zhang | An efficient approach for milling dynamics modeling and analysis with varying time delay and cutter runout effect[END_REF], the helix angle [START_REF] Zatarain | Analysis of the influence of mill helix angle on chatter stability[END_REF], the geometrical defects [START_REF] Seguy | Chatter in interrupted turning with geometrical defects: An industial case study[END_REF], the multiple chatter frequency [START_REF] Insperger | Multiple chatter frequencies in milling processes[END_REF], the non-linear cutting laws [START_REF] Stépán | Nonlinear dynamics of high-speed milling subjected to regenerative effect[END_REF][START_REF] Stépán | Delay, parametric excitation, and the nonlinear dynamics of cutting processes[END_REF], the spindle speed variation [START_REF] Seguy | Suppression of period doubling chatter in high-speed milling by spindle speed variation[END_REF].

The accuracy of this multi-frequency method has been again improved recently [START_REF] Insperger | Full-discretization and semi-discretization for milling stability prediction: Some comments[END_REF][START_REF] Ozoegwu | Hyper-third order full-discretization methods in milling stability prediction[END_REF][START_REF] Tang | Accurate and efficient prediction of milling stability with updated full-discretization method[END_REF]. In 5axis milling this improved full discretization method was applied to construct posture stability graphs to guide the selection of cutter postures during tool-path generation [START_REF] Wang | Chatter-free cutter postures in five-axis machining[END_REF].

Other approach for chatter modelling is the numerical simulation of the motion equation, called time domain simulation [START_REF] Campomanes | An improved time domain simulation for dynamic milling at small radial immersions[END_REF][START_REF] Liu | Modelling the machining of peripheral milling[END_REF][START_REF] Surmann | Simulation of milling tool vibration trajectories along changing engagement conditions[END_REF][START_REF] Arnaud | Simulation of low rigidity part machining applied to thin walled structures[END_REF]. In this case, the fine modelling at the scale of chip formation give accurate cutting force for general case and the machined surface only for system with a very simple dynamic [START_REF] Surmann | Simulation of milling tool vibration trajectories along changing engagement conditions[END_REF][START_REF] Arnaud | Simulation of low rigidity part machining applied to thin walled structures[END_REF]. Time simulation is able to model rapid phenomenon that are not periodical and highly complex to be modelled with frequency based approaches, like the ploughing effect [START_REF] Ahmadi | Analytical stability lobes including nonlinear process damping effect on machining chatter[END_REF]. The dynamics for 5-axis in flank milling was modelled with hemispherical tools for surface machining [START_REF] Ahmadi | Machining chatter in flank milling[END_REF]. This highly comprehensive time-based modelling allows process geometry, cutting force and stability to be obtained for large axial depth of cut in flank milling.

Although significant research has been reported in the stability of 5 axis milling [7-9, 27, 33] there is a little effort reported for the accurate stability prediction of this process. This paper is devoted to study of warped surfaces machining on milling machines with ball end mills. The dynamic model presented is based on the classical stability lobe theory, but new improvements are incorporated in order to improve the accuracy. The non-linear effect of radial allowance on the contact angle was modelled. A more complete experimental analysis was conducted in order to study the effect of machining parameters on the stability of climb milling.

All the aspects of the dynamical model and the original averaging method are presented in Section 2. An innovative experimental approach, to study stability during inclined surface milling, is described in Section 3. The experimental results are analysed to estimate the accuracy of the model in Section 4. Finally, concluding remarks are addressed.

Modelling

Tool immersion: contact angle

The contact patch between end mill and workpiece is described by start angle 𝜙 𝑠𝑡 and exit angle 𝜙 𝑒𝑥 . Angle between exit and start angle is a contact angle 𝜙 𝑐 = 𝜙 𝑒𝑥 -𝜙 𝑠𝑡 (Fig. 1). For machining with flat end mills contact angle is constant for variable axial depth of cut 𝑎 𝑝 value and it is determined only by radial depth of cut 𝑎 𝑒 . Furthermore, for down-milling exit angle always equals π (Fig. 1a) and for conventional cutting (up-milling) start angle always equals 0 (Fig. 1d) [START_REF] Budak | Analytical models for high performance milling -Part I: cutting forces, structural deformations and tolerance integrity[END_REF]. For curved surfaces machining with ball end mills this dependency is usually respected for contour milling and is not respected for copy milling. Hence, it is necessary to determine precisely start and exit angles as a function of current values of tool radius 𝑟 [mm], lead angle 𝜑, radial allowance 𝑎 𝑝,𝑟 [mm] and radial depth of cut 𝑎 𝑒 [mm]. As is shown in a Fig. 1, contact angle decreases as lead angle increases for constant radial allowance. It is significant that down milling and conventional milling in ball end machining are not pureas opposed to flat end machining, because loading and unloading of cutting edge occurs more graduallywith a time stretched impact.

For upward copy milling start and exit angles may be defined [START_REF] Shtehin | Vyznachennya kutiv vrizannya ta vykhodu pry obrobtsi pokhylykh poverkhon sferychnymy kintsevimi frezamy [Definition of start and exit angles in ball end milling of inclined surfaces[END_REF] with Eq. (1) for up-milling, Eq. ( 2) for down-milling and Eq. (3) for slot-milling, correspondingly:

{ 𝜙 𝑠𝑡 = 𝜋 2 -arctg √1 -𝐾 2 𝐾 • sin𝜑 𝜙 𝑒𝑥 = 𝜋 2 + arctg ( 𝑎 𝑒 𝑟 • sin𝜑 √ 4 - 𝑎 𝑒 2 𝑟 2 ) (1) a) b) c) d) e) f) { 𝜙 𝑠𝑡 = 𝜋 2 -arctg ( 𝑎 𝑒 𝑟 • sin𝜑 √ 4 - 𝑎 𝑒 2 𝑟 2 ) 𝜙 𝑒𝑥 = 𝜋 2 + arctg √1 -𝐾 2 𝐾 • sin𝜑 (2)
{ 𝜙 𝑠𝑡 = 𝜋 2 -arctg √1 -𝐾 2 𝐾 • sin𝜑 𝜙 𝑒𝑥 = 𝜋 2 + arctg √1 -𝐾 2 𝐾 • sin𝜑 (3)
In these equations coefficient 𝐾 = 1 -𝑎 𝑝,𝑟

𝑟

. Mutual values of radial allowance and lead angle, wherein cutter axis is excluded from cutting process, can be calculated [START_REF] Cosma | Horizontal path strategy for 3D-CAD analysis of chip area in 3-axes ball nose end milling[END_REF]:

𝑎 𝑝,𝑟 𝑐𝑟𝑖𝑡 = 𝑟(1 -cos𝜑) (4)
As far as possible, it is necessary to machine with radial allowance values less than 𝑎 𝑝,𝑟 𝑐𝑟𝑖𝑡 (Fig. 3).

Cutting forces and dynamic motion

Regenerative effect in milling consists in influence of previous cutting tooth passing on chip formation during current cutting tooth passing. Inasmuch as machine-tool-workpiece system oscillates with some frequency and magnitude, instantaneous chip thickness includes static (Fig. 2a) and dynamic (Fig. 2b) components. Static component is a function of spindle speed 𝑛, feed per tooth 𝑓 𝑧 and depths of cut variables. Dynamic component ℎ 𝑑 (𝜙 𝑗 ) of instantaneous rotational angle 𝜙 𝑗 of 𝑗 𝑡ℎ tooth depends upon the difference between previous and current positions of tool and workpiece with a period 𝜏.

ℎ 𝑑 (𝜙 𝑗 ) = [([𝑥(𝑡) -𝑥(𝑡 -𝜏)]cos𝜙 𝑗 + [𝑦(𝑡) -𝑦(𝑡 -𝜏)]sin𝜙 𝑗 )]𝑔(𝜙 𝑗 ) (5)
Here 𝑔(𝜙 𝑗 ) is a so-called switch function, which determines whether the tooth is cutting or not. This step function depends on contact angle:

{ 𝜙 𝑗 ∈ [𝜙 𝑠𝑡 ; 𝜙 𝑒𝑥 ] → 𝑔(𝜙 𝑗 ) = 1 𝜙 𝑗 ∉ [𝜙 𝑠𝑡 ; 𝜙 𝑒𝑥 ] → 𝑔(𝜙 𝑗 ) = 0 (6)
For the Single Degree Of Freedom (SDOF) dynamic model (Fig. 2) only 𝑦-axis component must be taken into account. So, [𝑥(𝑡) -𝑥(𝑡 -𝜏)] component equals 0. 𝑘 𝑦 is the modal stiffness and 𝜁 𝑦 is the modal damping of the system in 𝑦-direction. Therefore, cutting forces also have static and dynamic components. Static component can be neglected in dynamic analysis of cutting process [START_REF] Stépán | Delay, parametric excitation, and the nonlinear dynamics of cutting processes[END_REF]. Dynamic component describes 𝜏 periodic perturbations of the system. Hence, cutting process can be described by the next delay differential equation in generalized coordinates.

𝑞̈(𝑡) + 2𝜁𝜔 𝑛 𝑞̇(𝑡) + 𝜔 𝑛 2 𝑞(𝑡) = 𝑔(𝜙 𝑗 ) 1 𝑚 𝐹 𝑞 (𝑞(𝑡), 𝑞(𝑡 -𝜏)) (7)
Here 𝑚 is the modal mass, 𝜁 is the modal damping, 𝜔 𝑛 is the natural frequency of the mode and 𝐹 𝑞 is the actual force.

For the multi-dimensional case Eq. ( 7) becomes a matrix equation:

[𝐌]𝑞̈(𝑡) + [𝐂]𝑞̇(𝑡) + [𝐊]𝑞(𝑡) = [𝐅] (8)
Here freedom of the system (usually 1 or 2) and its symmetry. These values may be defined by modal analysis or by an experiment.

[𝐅] is matrix of dynamic (linear or non-linear) milling forces.

Cutting force models during ball end milling is more complicated than during flat end milling.

It is caused by several reasons. Effective cutting velocity is changing if radial allowance and lead angle are changing. As a result magnitude and vector direction of cutting force components are changing too. Thus, due to increase of lead angle, the radial component also increases. Normally axial components of cutting force are determined experimentally for different initial cutting conditions. After that it is possible to calculate specific cutting force coefficients, radial and tangential cutting force components.

For linear SDOF system 𝑦-axis component of the linear cutting force 𝐹 can be calculated in classical representation as:

𝐹 𝑦 = 𝐾 𝑦 𝑎 𝑝 ℎ 𝑑 (𝜙 𝑗 ) (9)
Here 𝑎 𝑝 is the axial depth of cut [mm] and 𝐾 𝑦 is an experimentally identified specific cutting force coefficient [N/mm 2 ]. In such manner tangential and radial components of the cutting force for 𝑗 th tooth are:

{ 𝐹 𝑡𝑗 = 𝐾 𝑡 𝑎 𝑝 ℎ 𝑑 (𝜙 𝑗 ) 𝐹 𝑟𝑗 = 𝐾 𝑟 ′𝑎 𝑝 ℎ 𝑑 (𝜙 𝑗 ) (10)
Here 𝐾 𝑡 and 𝐾 𝑟 ′ are the specific tangential and radial cutting force coefficients [N/mm 2 ],

respectively. A dimensionless radial cutting force coefficient 𝐾 𝑟 will be used from this point on, which is:

𝐾 𝑟 = 𝐹 𝑟 𝐹 𝑡 (11)
With reference to tangential and radial components 𝑦-axial cutting force component in matrix form is:

{𝐹 𝑦𝑗 } = [-cos𝜙 𝑗 -sin𝜙 𝑗 ] { 𝐹 𝑡𝑗 𝐹 𝑟𝑗 } (12)
After substitution Eq. ( 10) into Eq. ( 12):

𝐹 𝑦𝑗 = 1 2 𝑎 𝑝 𝐾 𝑡 𝐴 0 Δ𝑦 (13)
Where 𝐴 0 is directional dynamic milling force coefficient, which depends on the contact patch, with a period 𝑁 2𝜋

(𝑁 is number of teeth) and ayy the dynamic coefficient in the direction y:

𝐴 0 = 𝑁 2𝜋 𝑎 𝑦𝑦 = 𝑁 2𝜋 1 2 [cos2𝜙 -2𝐾 𝑟 𝜙 + 𝐾 𝑟 sin2𝜙] 𝜙 𝑠𝑡 𝜙 𝑒𝑥 (14)
For non-linear SDOF system 𝑦 component of the non-linear cutting force 𝐹 can be calculated in accordance with the experimentally verified three-quarter rule. This non-linearity takes into account the friction between cutting edge and workpiece during the cutting process.

𝐹 𝑦 = 𝐾 𝑦 𝑎 𝑝 (𝑓 𝑧 -𝑦(𝑡) + 𝑦(𝑡 -𝜏))

3/4 (15)
Here 𝑓 𝑧 is the feed per tooth [mm]. Using the first few terms of a Taylor series:

Δ𝐹 𝑦 ≈ 𝑘 1 (𝑦(𝑡 -𝜏) -𝑦(𝑡)) - 1 8𝑓 𝑧 𝑘 1 (𝑦(𝑡 -𝜏) -𝑦(𝑡)) 2 + 5 96𝑓 𝑧 𝑘 1 (𝑦(𝑡 -𝜏) -𝑦(𝑡)) 3 (16
)

𝑘 1 = 3 4 𝐾 𝑦 𝑎 𝑝 √𝑓 𝑧 4 (17)
Hence, the equation of motion for SDOF system is:

𝑦̈(𝑡) + 2𝜁𝜔 𝑛 𝑦̇(𝑡) + 𝜔 𝑛 2 𝑦(𝑡) = 1 𝑚 Δ𝐹 𝑦 ( 18
)

Stability analysis

During contact between cutting tooth and workpiece, system endures forced oscillations.

Then, when tooth leaves workpiece, system endures damped free vibrations. In machining classical stability analysis of linear dynamic system is based on Lyapunov stability [START_REF] Altintas | Analytical prediction of stability lobes in milling[END_REF]. For this the equation of motion must be written in this characteristic form first. All eigenvalues 𝜆 of the transfer matrix must be negative or equal zero for Lyapunov stability:

ℝ𝑒(𝜆 𝑖 ) ≤ 0 (𝑖 ∈ ℕ * ) (19)
The determinant of the resulting characteristic equation is:

𝑑𝑒𝑡[[𝐼] + 𝜆[𝐺 0 (𝑖𝜔 𝑐 )]] = 0 (20)
Oriented transfer function 𝐺 0 (𝑖𝜔 𝑐 ) of SDOF system describes with directional dynamic milling force coefficient 𝐴 0 and direct transfer function 𝐺 𝑦𝑦 (𝑖𝜔 𝑐 ). Direct transfer function 𝐺 𝑦𝑦 (𝑖𝜔 𝑐 ) characterizes dynamic compliance of the system as a response to periodic force with chatter frequency 𝜔 𝑐 .

𝐺 0 (𝑖𝜔 𝑐 ) = 𝐺 𝑦𝑦 (𝑖𝜔 𝑐 )𝐴 0

For linear SDOF system chatter free limited axial depth of cut is:

𝑎 𝑝,𝑙𝑖𝑚 = 2𝜋 𝑁𝐾 𝑡 𝑎 𝑦𝑦 ℝ𝑒[𝐺 𝑦𝑦 (𝑖𝜔 𝑐 )] (22)
Radial allowance as a function of axial depth of cut can be determined (Fig. 3):

𝑎 𝑝,𝑟 = 𝑟 [1 -𝑐𝑜𝑠 (𝑎𝑟𝑐𝑐𝑜𝑠 (𝑐𝑜𝑠𝜑 - 𝑎 𝑝 𝑟 ) -𝜑)] (23)
In the case of horizontal surface (𝜑 = 0) radial allowance 𝑎 𝑝,𝑟 , determined with Eq. ( 23), will be equal to axial depth of cut 𝑎 𝑝 .

Fig. 3 Difference between radial allowance and axial depth of cut in a case of inclined surface ball end milling Direct transfer function 𝐺 𝑦𝑦 , or flexibility in the 𝑦 direction may be determined as:

𝐺 𝑦𝑦 = 𝜔 𝑛 2 𝑘 𝜔 𝑛 2 -𝜔 𝑐 2 + 𝑖2𝜁𝜔 𝑛 𝜔 𝑐 (24)
In ball end milling of curved surfaces increase of radial allowance leads to non-linear increase of contact angle, which follows from dimensionless solution of Eq. (1-3), e.g. for slot-milling of surfaces with lead angle 30 ○ and 45 ○ , respectively (Fig. 4). In milling, directional dynamic milling force coefficient also changes. , must be vanished.

6.

Steps 2-5 should be repeated for the rest of specified range with a definite stage of iteration. Precision of this method depends on the iteration step. In some cases, non-linear special aspects of ball end milling have profound effects on cutting process. This is particularly so important when high-speed machining is coupled with lowimmersion milling. In this case, critically short time of cut (equal oscillation period or even less) may cause self-interruption of cutting process and thus may cause period-doubling of oscillation. Besides, non-linearity of cutting force inherent to ball end milling is a result of curvature of cutting edge. This non-linearity can be taken into account in a tight range using linearization.

Experimental tests

Cutting force coefficient

For approval of the proposed iterative method and analysis of non-linear effects, a series of 46° and 134° (lead angle is 30°), and 55° and 125° (lead angle is 45°), respectively. Time of cut equals 0.004 s (lead angle is 30°) and 0.003 s (lead angle is 45°). Test results for up-milling with lead angles 30° and 45° are presented in Fig. 13 and Fig. 14.

Quasi-periodic stable oscillations are typical for up-milling. This can be explained by additional impact load, which is inherent in conventional milling, when tooth of the cutting tool leaves workpiece. Thus, magnitude of free flight oscillation grows up instantaneous after impact excitation (tooth leaving of workpiece). As previously noted, during ball-end milling

Discussion

More than 80 cutting tests have been conducted to emphasize the behaviour of inclined surface. A quite good agreement is observed between simulations and machining tests. The discrepancies between the experiment and the modelling are observed for some spindle speed, in particular for up-milling (lead angle 30°) at 4600 rpm, see Fig. 13. It is presently not clear what causes these discrepancies. A possible cause could be due to the non-linear variation of damping of the system. Its variation in turn may be based on a host of factors: variations of friction force, elastic properties, flank wear growth, temperature increase in cutting zone, and so on and so forth.

Conclusions

This paper presents a study of warped surfaces machining on milling machines with ball end mills. A specific model based on the classical stability lobe theory was use and improved with more important aspects. The non-linear effect of radial allowance on the contact angle was integrated by an original averaging method. The cutting coefficients are updated during changing of effective radius and cutting velocity. This modelling allows predicting the stability lobe on inclined surface more accurately but easy to adjust in a real industrial process.

An original experimental procedure was developed in order to compute the cutting coefficient for various inclined surfaces. In real context, the cutting forces are measured by a standard dynamometer with various inclinations. The cutting coefficients are non-constant, and large amplitude of variation is observed as function of the tilt surface.

More complete experimental analysis was conducted in order to study the effect of machining parameters on the stability of climb milling. The lead angle has a direct impact on cutting stability, the greater the angle decreases more the stability increases.

The comparison between experiment and simulation show good correlation for the prediction of the depth of cut without chatter. Therefore, the improvement in the simulation of inclined surface inclined surfaces machining on milling machines with ball end mills needs the development non-linear damping models.

List of Figures Captions

List of Tables

Fig. 1

 1 Fig. 1 Start and exit angles in end mill machining, a) down-milling with flat end mill, b) down-milling with ball end mill (lead angle 30 ○ ), c) down-milling with ball end mill (lead angle 45 ○ ), d) up-milling with flat end mill, e) up-milling with ball end mill (lead angle 30 ○ ), f) up-milling with ball end mill (lead angle 45 ○ )

Fig. 2

 2 Fig. 2 Single Degree Of Freedom model of cutting with ball end mill, a) static representation, b) dynamic representation

  [𝐌] is the modal mass matrix, [𝐂] is the modal damping matrix and [𝐊] is the modal stiffness matrix. Matrix parameters have solutions in accordance with number of degrees of a) b)

Fig. 4

 4 Fig. 4 Non-linear effect of radial allowance on the contact angle

Fig. 5 Fig. 6

 56 Fig. 5 Iterative plotting of stability lobes diagrams for ball end milling

  tests were carried out. Milling pre-tests of workpiece (steel C35) with different lead angles (from 0° up to 60°) were carried out on a 5-axis milling machine Deckel Maho DMU 50 eVolution. During these pre-tests spindle speed was 3900 rpm and 4485 rpm; radial depth of cut was 0.5 mm and 1.0 mm. Pre-tests were carried out with Sandvik Coromant 2-flute ball end mill CoroMill 216 R216-10A16-050 with inserts R216-10 02 E-M 1010; mill diameter is Ø10 mm and free length is 49 mm (modal stiffness k = 7519 N/mm). Cutting scheme was slotting, climb and conventional milling. In the process of cutting, force components F x , F y and F z were measured with Kistler 9257BA dynamometer as shown in a Fig.7. Feeding was in the y-direction. Tangential, radial and axial cutting force components were obtained as: helix angle of the cutting tool is ψ hl = -15 °.

Fig. 7

 7 Fig. 7 Cutting force experimental setup

Fig. 8

 8 Fig. 8 Tangential (blue line), radial (green line) and axial (red line) cutting force components

Fig. 9

 9 Fig. 9 Cutting force components for various lead angles (15° 30° 45° 60°) in slot-milling

Fig. 10 showsFig. 10 2 ( 26 )Fig. 11

 101022611 Fig.10shows the experimental setup. During machining, the vibration of the main system (1), which includes ball end mill (2) and workpiece (3) was measured by the velocimeter laser (4), the signal was collected by a conditioner[START_REF] Altintas | Analytical prediction of stability lobes in ball end milling[END_REF] and stored in a personal computer[START_REF] Altintas | Analytical prediction of three dimensional chatter stability in milling[END_REF].

Fig. 12

 12 Fig. 12 Experimental stability for slot-milling (lead angle 45°), tests results marked with green ○, limited stable tests results are marked with blue ▽, unstable tests results are marked with red ✕

1 2 3

 12 there is time stretched impact excitation. This kind of stability is according to supercritical Hopf bifurcation.

Fig. 13 Fig. 14

 1314 Fig. 13 Experimental stability for up-milling (lead angle 30°), stable tests results are marked with green ○; limited stable tests results are marked with blue ▽; unstable tests results are marked with red ✕

Fig. 1 Fig. 2 Fig. 3 Fig. 5 Fig. 6 Fig. 7 Fig. 9 Fig. 14

 123567914 Fig. 1 Start and exit angles in end mill machining, a) down-milling with flat end mill, b) down-milling with ball end mill (lead angle 30°), c) down-milling with ball end mill (lead angle 45°), d) up-milling with flat end mill, e) up-milling with ball end mill (lead angle 30°), f) up-milling with ball end mill (lead angle 45°) Fig. 2 Single Degree Of Freedom model of cutting with ball end mill, a) static representation, b) dynamic representation Fig. 3 Difference between radial allowance and axial depth of cut in a case of inclined surface ball end milling Fig. 4 Non-linear effect of radial allowance on the contact angle Fig. 5 Iterative plotting of stability lobes diagrams for ball end milling Fig. 6 Comparison of classical approach and proposed method for stability lobes in slot-milling with 30° and 45° Fig. 7 Cutting force experimental setup Fig. 8 Tangential (blue line), radial (green line) and axial (red line) cutting force components Fig. 9 Cutting force components for various lead angles (15° 30° 45° 60°) in slot-milling Fig. 10 Experimental setup: a) photo of setup; b) principal scheme of setup Fig. 11 Experimental stability for slot-milling (lead angle 30°), stable tests results are marked with green ○; limited stable tests results are marked with blue ▽; unstable tests results are marked with red ✕

  Determine stability boundary for the current radial allowance. For this 𝑎 𝑝,𝑙𝑖𝑚 a corresponding 𝑎 𝑝,𝑟,𝑙𝑖𝑚 must be calculated with Eq. (22) and Eq.(23) respectively. Compare every 𝑎 𝑝,𝑟,𝑙𝑖𝑚 value with initial radial allowance value 𝑎 𝑝,𝑟 (𝑖). All of them, which do not satisfy the condition 𝑎 𝑝,𝑟 (𝑖) -

	𝑎 𝑝,𝑟 | 𝑚𝑖𝑛 𝑚𝑎𝑥 .
	2. Choose stage of iteration ∆𝑎 𝑝,𝑟 for determining of every value of the range as
	𝑎 𝑝,𝑟 (𝑖) = 𝑎 𝑝,𝑟 (𝑚𝑖𝑛) + 𝑖∆𝑎 𝑝,𝑟	where	𝑖 = {0,1, … , 𝑁}	and	∆𝑎 𝑝,𝑟 =
			∆𝑎 𝑝,𝑟 2	≤ 𝑎 𝑝,𝑟,𝑙𝑖𝑚 ≤ 𝑎 𝑝,𝑟 (𝑖) +
	∆𝑎 𝑝,𝑟
	2

1. Firstly, it is necessary to determine range covered of radial allowance values 𝑎 𝑝,𝑟 (max)-𝑎 𝑝,𝑟 (min) 𝑁 3. For the first value of the range must be determined start and exit angles 𝜙 𝑠𝑡 (𝑎 𝑝,𝑟 (𝑖)), 𝜙 𝑒𝑥 (𝑎 𝑝,𝑟 (𝑖)) with Eq. (1-3).

4. Respective directional dynamic milling force coefficients must be obtained. 5.

  axial component is less than radial component. Hence, it can be expected that cutting stability will go down if lead angle is increased.

	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥
	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥	𝜙 𝑠𝑡	𝜙 𝑒𝑥

Table 1

 1 in Table 1. Modal parameters identified

		k	fn
	1.2 %	2815 N/mm	864 Hz

Table 1

 1 Modal parameters identified

		k	fn
	1.2 %	2815 N/mm	864 Hz