Ensemble Kalman filtering for force model identification in milling

Available online at www.sciencedirect.com

Available online atonline www.sciencedirect.com Available at www.sciencedirect.com Available online at www.sciencedirect.com

ScienceDirect ScienceDirect

Procedia CIRP 00 (2019) 000–000 Procedia CIRP 00 (2017) 000–000 Procedia CIRP 82 (2019) 296–301 Procedia CIRP 00 (2019) 000–000

www.elsevier.com/locate/procedia

www.elsevier.com/locate/procedia www.elsevier.com/locate/procedia

17th CIRP Conference on Modelling of Machining Operations

17th CIRPfiltering Conferencefor on Modelling of Machining Operations in milling Ensemble Kalman force model identification Ensemble Kalman filtering for force model a,∗ b Max Schwenzer Stemmler , Muzaffer Aybidentification , Thomas Bergsa , in Dirkmilling Abelb 28th, Sebastian CIRP Design Conference, May 2018, Nantes, France a b for Machine Toolsa,∗ and Production Engineering (WZL), bRWTH Aachen University, 30, 52074 Aachen, Germany Max Schwenzer , Sebastian Stemmler , Muzaffer Ayb , Campus-Boulevard Thomas Bergs , Dirk Abel b

a Laboratory

Institute of Automatic Control, RWTH Aachen University, Campus-Boulevard 30, 52074 Aachen, Germany A new methodology to analyze and physical30, 52074 architecture Laboratory for Machine Tools and Production Engineering the (WZL), functional RWTH Aachen University, Campus-Boulevard Aachen, Germany of Institute of Automatic Control, RWTH Aachen University, Campus-Boulevard 30, 52074 Aachen, Germany existing products for an assembly oriented product family identification a

b

Abstract

Paul Stief *, Jean-Yves Dantan, Alain Etienne, Ali Siadat Mechanistic Abstract force models are popular to describe the force in cutting technology. Process simulation, process optimization, and process control rely on the accuracy of these models. Standard identification techniques are not capable of identifying a mechanistic force model on-line and École d’Arts et Métiers, et Métiers ParisTech, LCFCProcess EA 4495, 4 Rue Augustin Fresnel, Metz 57078,and France Mechanistic forceNationale models Supérieure areit popular to describe theArts force in cutting technology. optimization, process control in hard real-time. However, is necessary to adjust the model to increasing tool wear, e.g.simulation, in a modelprocess predictive controller for force control in rely on the of these models. Standard identification techniques are not capable of identifying a mechanistic model on-line and milling. Thisaccuracy work introduces the ensemble Kalman filter to the field of force model identification in cutting technologyforce – enabling for the first hard real-time. However, it 3is87necessary adjust themodels. model to increasing wear, e.g. in a model controller for force in *in Corresponding author. Tel.: +33 37 mechanistic 54 30; to E-mail address: [email protected] time a continuous parametrization of force The approachtool shows high accuracy and predictive fast convergence in spite of thecontrol presence milling. This work introduces theapproach ensembleisKalman to the field of 1000 forcerandom model identification in cutting technology – enabling theideal, first of measurement noise. The novel validatedfilter statistically using initial distributions and different ensemble sizes.for The time a continuous parametrization of with mechanistic models. The approach shows highofaccuracy and fast convergence spiteconverges of the presence simulated force signals is augmented differentforce levels of noise (signal-to-noise ratios 50, 15, five). Nevertheless, theinfilter within of measurement The novelofapproach is validated statistically using 1000 random initial distributions and different ensemble sizes. The ideal, three, eight, and noise. 32 revolutions the tool respectively. Abstract simulated force signals is augmented with different levels of noise (signal-to-noise ratios of 50, 15, five). Nevertheless, the filter converges within three, eight, 32 revolutions tool respectively. c 2019  The and Authors. Publishedofbythe Elsevier B.V. In today’s business environment, the trend towards more product and customization unbroken.ofDue to this development, of Peer-review under responsibility of the scientific of Thevariety 17th CIRP Conference onisModelling Machining Operations, inthe theneed person © 2019 The Authors. Published by Elsevier B.V. committee agile and reconfigurable production systems emerged to Dr cope with various andMsaoubi. product Toofdesign and optimize production c the  2019 The Authors. Published byOzturk Elsevier B.V. of Conference Dr Erdem and Co-chairs Tom andproducts Dr Rachid Peer-review underChair responsibility of the scientific committee ofMcleay The 17th CIRP Conference on families. Modelling Machining Operations systems as well as responsibility to choose the of optimal productcommittee matches, of product analysis are on needed. Indeed, most of theOperations, known methods to Peer-review under the scientific The 17th CIRPmethods Conference Modelling of Machining in the aim person Keywords: Ensemble Kalman filter;family force model; identification; mechanistic force model; analyze a product orChair one product on and themodel physical level. Different product however, of the Conference Dr Erdem Ozturk Co-chairs Dr Tom Mcleay and Drfamilies, Rachidmilling Msaoubi. may differ largely in terms of the number and nature of components. This fact impedes an efficient comparison and choice of appropriate product family combinations for the production Ensemble Kalman filter; force model; model identification; mechanistic force model;functional milling and physical architecture. The aim is to cluster Keywords: system. A new methodology is proposed to analyze existing products in view of their pling step, namely ensemble filter (EnKF). Its unique these products in new assembly oriented product families for the optimization of existing assembly lines and Kalman the creation of future reconfigurable assembly systems. Based on Datum Flow Chain, the physical structure of the products 1. Introduction features areis analyzed. Functional subassemblies are identified, and pling architecture step, namely ensemble Kalman (EnKF). Its unique a functional analysis is performed. Moreover, a hybrid functional and physical graph (HyFPAG) is thefilter output which depicts the 1. Introduction features are similarity between product families design support to both, production system identification, planners and product designers. An illustrative • Continuous In cutting technology, force is by theproviding key to understand the proexample of a nail-clipper used tocombine explain the proposed methodology. An industrial case study on two productcapability, families of steering columns of • On-line and hard real-time cess. Mechanistic forceismodels a physically-motivated thyssenkrupp Presta France is then carried out to give a first industrial evaluation of the proposed approach. • Continuous identification, In cutting technology, force is the key to understand the pro• The low amount of data that is required. basis function with empirical coefficients to account for mate©cess. 2017Mechanistic The Authors.force Published by combine Elsevier B.V. • On-line and hard real-time capability, a physically-motivated rial properties. They tunemodels a given model structure to the material Peer-review under responsibility the scientific of the 28th CIRP Design 2018. • TheConference low amount of data that is required. basis function with empiricalofcoefficients tocommittee account for mate-

Those features are in particular important if the model is of the workpiece, of the tool, and to the tool geometry. Mecharial properties. Theyare tune a given modelthe structure to the material used for process control [17]. An interesting approach in this nistic force models used to model cutting process. They Keywords: Assembly; Design method; Family identification Those features are incontrol particular important if the themaximum model is of the workpiece, of the tool, and to the toolmaterial geometry. Mechafield is model predictive (MPC) to control enable a purely geometric simulation of the removal to used for process control [17]. An interesting approach in this nistic force models are used to model the cutting process. They of the active force in milling (misleadingly often also referred estimate the underlying process force and the tool load. field is model predictive control (MPC) to control the maximum enable a purely geometric simulation of the material removal to to as “cutting force”). Stemmler et al. [16] use a mechanistic New developments use mechanistic models in order to of themodel active to force in and milling (misleadingly often also referred themilling underlying process force and theeliminate tool load.[16] the 1.estimate Introduction of the product range characteristics manufactured and/or force optimize the feed of a machining center in orcontrol the process to reduce [2] or to as “cutting force”). Stemmler et al. [16] use a mechanistic New developments use mechanistic models in order to assembled in athis system. this context, the main challenge in der to track given forceInreference as close as possible. The overshoot of the active force when the engagement condition force to optimize isthenow feednot of only a machining center instep orcontrol milling process to reduce [2]inoraeliminate [16]state the Due the to the fast development the domain of modelling and analysis to copein with MPC model solves a constrained optimization problem everysingle changes. The model coefficients represent certain tool der to track a active givenproduct force reference as close product as possible. The overshoot of adapted theand active the engagement condition communication an force ongoing ofBut digitization and products, a the limited range or existing families, describing force for a finite horizon. and must be when thewhen tooltrend wears. identifying the MPC solves a constrained optimization problem in every step changes. The model coefficients represent a certain tool state digitalization, manufacturing enterprises important but This also to be able to analyze and toofcompare products define work handles the problem how to update thetomechaspecific coefficients is expensive in both, are timefacing and effort. Most describing active force finite horizon. and must bein adapted when thereal-time tool wears. But identifying the challenges today’s market environments: a capable. continuing new families. It canfor beaobserved that classical existing nisticproduct forcethe model continuously. Nevertheless, the authors conapproaches are far from being or on-line This work handles the problem of howbroader, to clients update thefeatures. mechaspecific coefficients is expensive in both, time and effort. tendency towards reduction of product development timesMost and product areapplication regrouped in of or sider thefamilies scope of asfunction much since this ennistic force modelcalibrate continuously. Nevertheless, thehardly authors conapproaches areintroduces farlifecycles. from being orthere on-line ables to directly mechanistic force model for the first shortened product Inreal-time addition, is ancapable. increasing However, assembly oriented product families are to find. This paper a novel approach to recursively estisider application much broader, sincemonitorthis entime. The work of may contribute process modeling, demand customization, being at the timeatin a global On the the scope product family level,asto products differ mainly in two mate theofcoefficients of mechanistic forcesame models every samables to directly calibrate mechanistic force model for the first This paper introduces a novel approach to recursively estiing, optimization and(i) control. competition with competitors all over the world. This trend, main characteristics: the number of components and (ii) the ∗ Corresponding time.ofThe work may(e.g. contribute to process modeling, monitormate the of mechanistic force models at every samauthor.the Tel.: +49-241-80-28021; fax: +49-241-80-22293. which is coefficients inducing development from macro to micro type components mechanical, electrical, electronical). ing, optimization and control. E-mail address: [email protected] (Dirk Abel). markets, results in diminished lot sizes due to augmenting Classical methodologies considering mainly single products ∗ Corresponding author. Tel.: +49-241-80-28021; fax: +49-241-80-22293. product varieties (high-volume to low-volume production) [1]. or solitary, already existing product families analyze the E-mail address: [email protected] c 2019 The 2212-8271  Authors. Published by Elsevier (Dirk B.V. Abel). To cope with this augmenting variety as well as to be able to product structure on a physical level (components level) which Peer-review under responsibility of the scientific committee of The 17th CIRP Conference on Modelling of Machining Operations, in the person of the Conference identify possible potentials inDrthe existing c 2019 2212-8271  The optimization Authors. Published by Mcleay Elsevierand B.V. Chair Dr Erdem Ozturk and Co-chairs Dr Tom Rachid Msaoubi. causes difficulties regarding an efficient definition and production system, it is important to have a precise knowledge comparison of different product families. Addressing this Peer-review under responsibility of the scientific committee of The 17th CIRP Conference on Modelling of Machining Operations, in the person of the Conference Chair Dr Erdem Ozturk and Co-chairs Dr Tom Mcleay and Dr Rachid Msaoubi. 2212-8271 © 2019 The Authors. Published by Elsevier B.V. Peer-review under responsibility of the by scientific of The 17th CIRP Conference on Modelling of Machining Operations 2212-8271 © 2017 The Authors. Published Elseviercommittee B.V. 10.1016/j.procir.2019.04.028 Peer-review under responsibility of the scientific committee of the 28th CIRP Design Conference 2018.



Max Schwenzer et al. / Procedia CIRP 82 (2019) 296–301 Author name / Procedia CIRP 00 (2019) 000–000

2. State of the art

h. This is the most criticized drawback of this identification method, e.g. by [10, 9, 24, 15].

2.1. Mechanistic force model In general, two different model structures are distinguished in mechanistic force modeling: lumped-mechanism models merge all effects in a single coefficient, whereas dualmechanism models use two coefficients to separately describe the effect of shearing and of ploughing. Dual-mechanism models tend to be more accurate [21, 1, 23]. The most popular representatives of the dual-mechanism models are the linear model according to Altintas and Lee [3] and the exponential formulation according to Kienzle [11]. Fi,Altintas = b · (Kie + Kic · h) Fi,Kienzle = ki · b · h

297 2

1−mi

i ∈ t, r, p,

i ∈ t, r, p.

(1) (2)

Where h denotes the uncut chip thickness, b the uncut chip width, and the coefficients ki and mi , or Kie and Kic respectively, represent material-specific tuning. The same functional relationship is assumed for all components of the force. For the rest of this paper, the index i denotes the tangential (t), radial (r), and passive (p) component of the force and the coefficients. Although most of the work on mechanistic force models is conducted assuming a linear model structure (Eq. 1), the coefficients themselves are often assumed to be a function of the uncut chip thickness h [9, 20, 18, 23, 22]. This turns the former linear model into a non-linear model. Therefore, this paper uses directly the non-linear force model Eq. 2. To account for the helical cutting edge of an end-mill, the cutter is divided along the z-axis into disk elements, for which a straight cutting edge is assumed. Each disk element is twisted towards the underlying element in the direction of the helix angle. This approximates the helix like a spiral staircase. 2.2. Force model identification Two different approaches exist to identify mechanistic force models. The methods of • average forces, and • instantaneous uncut chip thickness. The method of average forces originates from turning where the chip geometry does not change, which results in a static force vector. To apply this approach to the intermittent cut of milling, the force signals need to be averaged over a multiple of revolutions – because the arithmetic mean of a cutter revolution corresponds to the static force component. Averaging the force per revolution eliminates the influence of the radial deviation of the tool [20]. All models assume the force to scale proportionally with the uncut chip width b. However, several experiments have to be conducted at different feed rates f to depict the influence of the uncut chip thickness

The method of instantaneous uncut chip thickness was developed to overcome this central drawback making more use out of the measured data. It even shall improve the overall accuracy of the identified models. Essentially, it performs a curve fit on the measured force signals. The coefficients are adjusted iteratively to obtain an optimal least-squares fit. Usually the method uses global optimization algorithms, such as evolutionary algorithms [9], or particle swarm optimization [24] – independently whether it fits a linear model [9], or a nonlinear model [24]. Local optimization algorithms have only rarely been examined in an unconstrained [14] or constrained version [15]. Using local constrained algorithms for curve fitting have the advantage of allowing an identification at runtime [15]. Gonzalo et al. [8] compare the method of average forces with the method of instantaneous forces. As a reference for the coefficients, they identify them in turning and define results as the “single point of truth”. They conclude that the method of instantaneous forces shows a better conjunction to the “true” coefficients obtained from turning. Nevertheless, the improvement is small but they argue that the method of instantaneous forces has a better physical credibility due to the correspondence to the turning coefficients. Adem et al. [1] support the finding that an optimization-based curve fit is more accurate than the average forces approach. An interesting combination of both methods is to average the force measures at the same instantaneous chip thickness over every tooth [20, 4, 19]. The approach can be intuitively exemplified assuming a constant tooth pitch  Nz −1   1  360 ◦ F h (ϕ) = F h ϕ+ ·i , Nz i=0 Nz

(3)

where Nz is the number of teeth and ϕ is the cutter rotation position. To average the force for all values of the instantaneous uncut chip thickness, an 1/Nz -revolution of the cutter rotation ◦ must be considered ϕ(t) ∈ [0, 360 Nz ]. This method combines the insensitivity to the radial deviation, which is obtained through averaging the forces per revolution, and the high utilization rate of the force signal (resulting in a low number of required experiments) through using the instantaneous uncut chip thickness. Nevertheless, none of the methods allows an on-line identification in hard real-time.

3. Approach In the following section, the function and principle of the EnKF and its implementation to a mechanistic force model are described. Afterwards the set-up of the simulative experiments and the evaluation criteria are discussed.

Max Schwenzer et al. / Procedia CIRP 82 (2019) 296–301 Author name / Procedia CIRP 00 (2019) 000–000

298

3.1. Ensemble Kalman filter

Table 1. Force coefficients for the Kienzle force model

The EnKF is a fully non-linear estimator, which is introduced by Evensen [5]. It propagates an ensemble of individual state vectors in time and takes its mean as the best-guess. Essentially, it is a particle filter with an Markov chain Monte Carlo simulation as its backbone to represent the evolution in time of the probability density of the states [6]. Using the mean value as the best-guess is the standard in literature [5, 6, 7, 13]. Nevertheless, Evensen [5] himself noted that for large ensembles, the use of the median is also conceivable. An aspect of particular charm is that the filter uses the variance within the ensemble as covariance matrices and; therefore, do not have to be tuned manually. Covariance matrices represent how much to trust the model or the measurements. Assuming an ensemble matrix Eq. 4 with n individual state vectors x – indicated by the superscript. The index (·)k|k−1 denotes a value at time step k based on the information of the previous time step k − 1. First, in the prediction phase, the EnKF uses the nonlinear model function to propagate every member of the ensemble forward in time using the identical input uk for every ensemble, Eq. 5, and determines its covariance matrix Pk|k−1 , Eq. 5.

Xk|k−1 =



1 xk|k−1 ,



3

..., 

n xk|k−1



Xk|k−1 = f Xk−1|k−1 , uk   T  Pk|k−1 = Xk|k−1 − X k|k−1 · Xk|k−1 − X k|k−1

(4) (5) (6)

In this way, the spreading of the ensemble to the mean X k|k−1 represents the error variance. This removes the need to choose the state covariance matrix manually. Second, in the update phase, the EnKF updates the predicted state vectors through the measurements. Distorting the measurements with artificial noise εk creates an individual measurej . It weights the ment vector for every ensemble member xk|k−1 influence of the measurements by the Kalman gain Kk −1  Kk = Pk|k−1 · HkT · Hk · Pk|k−1 · HkT + Rk ,    j j j = xk|k−1 + Kk · zk + εkj − Hk · xk|k−1 . xk|k

(7) (8)

Where Hk denotes the measurement or output matrix. The covariance matrix of the measurement noise Rk consists simply of the artificial noise that is added to the measurements:

with Rk = εk · εTk ,   1 εk = εk , . . . , εnk .

(9) (10)

upper bound lower bound assumed values

kt

mt

kr

mr

1800 800 1200

0.6 0.05 0.3

1200 600 900

0.3 0.01 0.2

3.2. Practical remarks on EnKF implementation Initial ensemble. The choice of the initial ensemble affects the time in which the filter converges. It is recommended to create an ensemble by applying perturbations to a best-guess value [6]. In accordance to [15], we define a technologically reasonable range of the coefficients (Tab. 1) and initialize the ensemble by uniformly distributing them within half of this range around its center. Note, that a uniform random distribution does not reflect potential characteristics or correlations between the states. The initial values are determined once beforehand and kept fixed for every ensemble size and combination of measurement noise. This ensures consistency. To statistically evaluate the performance of the EnKF, Ninit = 1000 distinct initial ensembles were created. With this, a general conclusion of the convergence behavior of this novel identification method can be derived. Artificial measurement noise. At every sample point the measurements are individualized for every ensemble member by adding artificial noise. This distortion increases the robustness of the filter, which relies on a properly distributed ensemble. This paper adds white Gaussian noise to the measurement value and normalize the covariance matrix. The power of the added noise is 30 % of the power of the ideal input signal, see section 3.3. This translates to 12 dB and 10.5 dB for the tangential and radial force respectively. Since white noise (i.e. equal power over time and frequencies) is used, it creates high errors at low values of the uncut chip thickness h where the force measurement is close to zero. Ploughing and stability. Ploughing is the part of the cut where a defined cut fails due to the low value of the uncut chip thickness h ≈ 0. However, in the cyclic intermittent cut of the milling process, ploughing inevitably occurs. Many works on the identification of force models explicitly exclude small uncut chip thicknesses, e.g. h ≤ 0.06 mm [8] or h < 0.01 mm [19]. In turning and with the method of average forces, K¨onig et al. [12] examine only uncut chip thicknesses h ≥ 0.1 mm. This is because the models themselves are only valid for a defined cut. Therefore, the EnKF is not updated when the sum of the uncut chip N hl < 0.05 mm. thickness of all N disk elements is small l=1 By this, we essentially exclude samples where ploughing occurs along the whole cutting edge. Though, individual sections (represented by disk elements) of the cutting edge still can be affected by ploughing. It is assumed that this improves the robustness, stability, and convergence time of the filter since the artificial noise that the EnKF adds to its measurement input lets the filter measure a force value even if the uncut chip thickness is zero.



Max Schwenzer et al. / Procedia CIRP 82 (2019) 296–301 Author name / Procedia CIRP 00 (2019) 000–000

Table 2. Overview of the tool geometry parameters and the process parameters

Tool geometry

Process parameter

Diameter Number of teeth Helix angle

D Nz

10 mm 2

β

45 ◦

Rake angle

γ

20 ◦

Feed Cutting velocity Depth of cut Width of cut

f vc

0.12 mm 2.44 m s−1

ap

2 mm

ae

2 mm

299 4

Using a sample frequency of f s = 10 kHz, a revolution translates roughly to 129 samples. To analyze robustness, the mean of the absolute error per revolution ∆Fi,abs,rev for all Ninit = 1000 initial ensembles is presented. To examine the convergence time of the identification, the instantaneous weighted error is defined as Eq. 12. The absolute error is weighted by the instantaneous uncut chip thickness hk . The convergence time kcnv of an ensemble j should be the time step or number of samples from which on the weighted error, Eq. 12, stays below 3 % of the maximum force per revolution

lim

  hk j · F(ki,1.1 , mij )k − Fi,k . maxk hk 3 · max Fi,k , < ∀k ≥ kcnv . k 100

j = ∆Fi,h,k

(12)

j ∆Fi,h,k

(13)

k>kcnv →∞

4. Results and discussion S/N: 50,

S/N: 15,

S/N: 5

Fig. 1. One exemplary cutter revolution of the simulated tangential force including measurement noise of different signal-to-noise ratios

3.3. Simulated experiments To test the filter, this paper simulates 42 revolutions of a typical orthogonal milling operation with an flat end-mill, Tab. 2. The engaged part of the tool axis was divided into N = 23 disk elements to account for the helix angle β. No radial deviation was assumed. The exponential force model according to Kienzle, Eq. 2, was used with the coefficients shown in Tab. 1. Artificial noise was added to the ideal simulation signal. The different signal-to-noise ratios (S/N) were assumed to be 50, 15, and five respectively, Fig. 1. All calculations were performed with MATLAB R2017a from The MathWorks on an AMD Ryzen7-2700 computer running Windows 10. For the sake of conciseness, we limit the results to the tangential component of the force and the coefficients (index i = t). Nevertheless, the approach equally applies to the radial and passive force component and let derive identical results.

3.4. Evaluation criteria This paper analyzed different ensemble sizes at different signal-to-noise ratios of the simulated force signals regarding accuracy, robustness, and convergence. For accuracy, we examine the mean absolute error per revolution for a certain ensemble j

j = ∆Fi,abs,rev

1 · 360 ◦



360 ◦ 0◦

Fϕ (ki,1.1 , mi ) − Fi,ϕ dϕ.

(11)

Fig. 2 depicts the evolution of the mean error per revolution ∆Fi,abs,rev of the Ninit = 1000 initial ensembles. The lines denote the mean error of different ensemble sizes and the circles represent the variation of the error within an ensemble size (it is the scaled standard deviation × σ). It indicates how sensitive a particular ensemble size is to its initial distribution. A small circle means that the mean error ∆Fi,abs,rev is a good measure to evaluate how fast and accurate the particular ensemble size can identify the model coefficients. Even at low signal-to-noise ratios, the EnKF is able to accurately describe the tangential force component. The black circles are particularly striking: they represent the standard deviation of the smallest ensemble size (ten members). They show a very high variance and an unsteady convergence. Zooming into the first 15 revolutions, Fig. 3, shows that for ideal signals, small ensemble sizes are favorable in terms of convergence. However, for distorted signals, small ensemble sizes converge worse showing high variations. The results indicate that a larger ensemble size results in a higher “inertia” of the ensemble, holding it from converging fast. However, assuming measurement noise, a larger ensemble size results in a more stable convergence. The higher inertia holds the EnKF from falling for distorted measurement signals. The convergence time kcnv confirms these observations, Fig. 4. The top diagram emphasizes that the ensemble size acts similar to an inertia: the larger it is, the slower it converges at ideal measurements. The 42 simulated revolutions translate to 5400 samples. If the EnKF is unable to identify proper coefficients that cause the error to converge (Eq. 13) within this time, it is considered to diverge. One can see that only from ensemble sizes larger than 50, the error converges reliably for all initial ensembles. The higher the signal quality, the smaller can the ensemble size be chosen. The ensemble size of ten is an exception as it suffers from convergence at the majority of all initial ensemble distributions

Max Schwenzer et al. / Procedia CIRP 82 (2019) 296–301 Author name / Procedia CIRP 00 (2019) 000–000

S/N: 50

S/N: 15

Ensemble sizes:

# revolutions / 10, 30,

×60

60,

S/N: 15

S/N: 5

Ensemble sizes:

# revolutions / 10, 30,

60,

200 100

1,500

×10

1,000 500 0 4,000 2,000

100

Fig. 2. Evolution of the mean error per revolution for different signal-to-noise ratios (top: S/N 50, middle: S/N 15, bottom: S/N 5) and selected ensemble sizes. Size of the circles represent the scaled standard deviation

S/N: 50

300

0 ×30

kcnv / -

S/N: 5

5

400 kcnv / -

Zoom, Fig. 3

kcnv / -

300

0

10

30

50

70

90

Fig. 4. Convergence time of all Ninit initial ensembles represented as boxplots (whiskers: 5 %, 95 %; box: 25 %, 75 %, · maximum/minimum) for different signal-to-noise ratios (top: S/N 50, middle: S/N 15, bottom: S/N 5)

×60

×30

×10

Mean,

100

Fig. 3. Zoom into the evolution of the mean error per revolution for different signal-to-noise ratios (top: S/N 50, middle: S/N 15, bottom: S/N 5). Size of the circles represent the scaled standard deviation

(high maximum at S/N 50 and a large whisker at S/N 15). Assuming that a signal-to-noise ratio of five seldom occurs with a properly configured measurement chain, an ensemble size of 50 is a good trade-off between convergence and stability. Fig. 5 illustrates the evolution of the coefficients for an ensembles with 50 members. It depicts the evolution of the mean (black curve) and the maximum value (red curve) and the minimum value (blue curve) of all initial distributions (Ninit = 1000). Using a simulation instead of measurement data enables us to evaluate whether the filter converges to the correct values (gray lines). The model coefficients kt and mt quickly converge

Minimum,

Maximum

Fig. 5. Evolution of the states that represent the tangential force coefficients for an ensemble size of 50. (top: S/N 50, middle: S/N 15, bottom: S/N 5)

for a signal-to-noise ratio of 50 and only a little slower for significantly worse signal-to-noise ratios. It is astonishing that every initial distribution converges to the correct pair of values. This suggests that the problem is well-posed, i.e. that there exist only one unique pair of values to describe a certain force signal. 5. Conclusion This paper introduces ensemble Kalman filtering to identify mechanistic force models. The methodology allows for a continuous estimation for model coefficients. The EnKF demonstrates a unique robustness to measurement noise and



Max Schwenzer et al. / Procedia CIRP 82 (2019) 296–301 Author name / Procedia CIRP 00 (2019) 000–000

a fast convergence behavior in spite of drastic measurement noise. With an ensemble size of 50, the EnKF converges in 95 % of the initial ensemble distributions within 298, 561, and 2038 samples for a signal-to-noise ratio of 50, 15, and five respectively. In general, the larger the ensemble size, the slower the convergence time but the higher the robustness against measurement noise. An ensemble size of 50 seems to be a good compromise between robustness and convergence time. The problem of force model identification is of well-posed nature. All initial distributions converge to the ideal model coefficients. Only the simulation of artificial force signals made it possible to draw this essential conclusion. Though the approach was analyzed statistically on a sound number of initial ensembles (1000), several aspects remain for further examination. The current implementation uses a fixed power of the artificial measurement noise within the EnKF. A formulation relative to the uncut chip thickness is desirable to limit the amount of introduced perturbation at low forces. Furthermore, a validation on real measurements is pending. Finally, for an on-line identification within the process, varying engagement conditions must be tested: when the ensemble converges, can it still react to changes of the coefficients? The system will soon be tested in conjunction with a force controller for milling. Acknowledgements The authors would like to thank the German Research Foundation DFG for the kind support within the Cluster of Excellence “Internet of Production” (Project ID: 390621612). References [1] Adem, K.A.M., Fales, R., El-Gizawy, A.S., 2015. Identification of cutting force coefficients for the linear and nonlinear force models in end milling process using average forces and optimization technique methods. The International Journal of Advanced Manufacturing Technology 79, 1671– 1687. doi:10.1007/s00170-015-6935-3. [2] Altintas, Y., Aslan, D., 2017. Integration of virtual and on-line machining process control and monitoring. CIRP Annals - Manufacturing Technology 66, 349–352. doi:10.1016/j.cirp.2017.04.047. [3] Altintas, Y., Lee, P., 1996. A general mechanics and dynamics model for helical end mills. CIRP Annals - Manufacturing Technology 45, 59–64. doi:10.1016/S0007-8506(07)63017-0. [4] Dotcheva, M., Millward, H., Lewis, A., 2008. The evaluation of cuttingforce coefficients using surface error measurements. Journal of Materials Processing Technology 196, 42–51. doi:10.1016/j.jmatprotec.2007.04.136. [5] Evensen, G., 1994. Sequential data assimilation with a nonlinear quasigeostrophic model using monte carlo methods to forecast error statistics. Journal of Geophysical Research 99, 10143. doi:10.1029/94JC00572. [6] Evensen, G., 2003. The ensemble kalman filter: theoretical formulation and practical implementation. Ocean Dynamics 53, 343–367. doi:10.1007/s10236-003-0036-9. [7] Gillijns, S., Mendoza, O., Chandrasekar, J., De Moor, B., Bernstein, D., Ridley, A., 2006. What is the ensemble kalman filter and how well does it work?, in: 2006 American Control Conference, IEEE. p. 6 pp. doi:10.1109/ACC.2006.1657419. [8] Gonzalo, O., Beristain, J., Jauregi, H., Sanz, C., 2010. A method for the identification of the specific force coefficients for mechanistic milling simulation. International Journal of Machine Tools and Manufacture 50, 765– 774. doi:10.1016/j.ijmachtools.2010.05.009.

301 6

[9] Grossi, N., 2017. Accurate and fast measurement of specific cutting force coefficients changing with spindle speed. International Journal of Precision Engineering and Manufacturing 18, 1173–1180. doi:10.1007/s12541-0170137-x. [10] Jayaram, S., Kapoor, S., DeVor, R., 2001. Estimation of the specific cutting pressures for mechanistic cutting force models. International Journal of Machine Tools and Manufacture 41, 265–281. doi:10.1016/S08906955(00)00076-6. [11] Kienzle, O., 1952. Die Bestimmung von Kr¨aften und Leistungen an spanenden Werkzeugmaschinen, in: VDI-Z, VDI-Verlag. pp. 299–305. [12] K¨onig, W., Essel, K., Witte, L., Verein Deutscher Eisenh¨uttenleute, 1982. Spezifische Schnittkraftwerte f¨ur die Zerspanung metallischer Werkstoffe. Verlag Stahleisen. OCLC: 64323901. [13] Mandel, J., 2009. A brief tutorial on the ensemble kalman filter. arXiv:0901.3725 [physics] arXiv:0901.3725. [14] Rubeo, M.A., Schmitz, T.L., 2016. Mechanistic force model coefficients: A comparison of linear regression and nonlinear optimization. Precision Engineering 45, 311–321. doi:10.1016/j.precisioneng.2016.03.008. [15] Schwenzer, M., Auerbach, T., D¨obbeler, B., Bergs, T., 2018. Comparative study on optimization algorithms for online identification of an instantaneous force model in milling. The International Journal of Advanced Manufacturing Technology URL: http://link.springer.com/10.1007/s00170-018-3109-0, doi:10.1007/s00170-018-3109-0. [16] Stemmler, S., Abel, D., Schwenzer, M., Adams, O., Klocke, F., 2017. Model predictive control for force control in milling. IFAC-PapersOnLine 50, 15871–15876. doi:10.1016/j.ifacol.2017.08.2336. [17] Thombansen, U., Buchholz, G., Frank, D., Heinisch, J., Kemper, M., Pullen, T., Reimer, V., Rotshteyn, G., Schwenzer, M., Stemmler, S., Abel, D., Gries, T., Hopmann, C., Klocke, F., Poprawe, R., Reisgen, U., Schmitt, R., 2018. Design framework for model-based self-optimizing manufacturing systems. The International Journal of Advanced Manufacturing Technology 97, 519–528. URL: http://link.springer.com/10.1007/s00170-018-1951-8, doi:10.1007/s00170-018-1951-8. [18] Wan, M., Zhang, W.H., Dang, J.W., Yang, Y., 2009. New procedures for calibration of instantaneous cutting force coefficients and cutter runout parameters in peripheral milling. International Journal of Machine Tools and Manufacture 49, 1144–1151. doi:10.1016/j.ijmachtools.2009.08.005. [19] Wan, M., Zhang, W.H., Qin, G.H., Wang, Z.P., 2008. Consistency study on three cutting force modelling methods for peripheral milling. Proceedings of the Institution of Mechanical Engineers, Part B: Journal of Engineering Manufacture 222, 665–676. doi:10.1243/09544054JEM1085. [20] Wan, M., Zhang, W.H., Tan, G., Qin, G.H., 2007. New algorithm for calibration of instantaneous cutting-force coefficients and radial run-out parameters in flat end milling. Proceedings of the Institution of Mechanical Engineers, Part B: Journal of Engineering Manufacture 221, 1007–1019. doi:10.1243/09544054JEM515. [21] Wang, J.J.J., Chang, H.C., 2004. Extracting cutting constants via harmonic force components for a general helical end mill. The International Journal of Advanced Manufacturing Technology 24, 415–424. doi:10.1007/s00170-003-1664-4. [22] Wang, L., Si, H., Guan, L., Liu, Z., 2018. Comparison of different polynomial functions for predicting cutting coefficients in the milling process. The International Journal of Advanced Manufacturing Technology 94, 2961– 2972. doi:10.1007/s00170-017-1086-3. [23] Wei, Z.C., Guo, M.L., Wang, M.J., Li, S.Q., Liu, S.X., 2018. Prediction of cutting force in five-axis flat-end milling. The International Journal of Advanced Manufacturing Technology doi:10.1007/s00170-017-1380-0. [24] Zhang, Z., Li, H., Meng, G., Ren, S., Zhou, J., 2017. A new procedure for the prediction of the cutting forces in peripheral milling. The International Journal of Advanced Manufacturing Technology 89, 1709–1715. doi:10.1007/s00170-016-9186-z.