Abstract

This work focuses on the identification of optimal model parameters related to Abrasive Waterjet Milling (AWJM) process. The evenly movement as well as variations of the jet feed speed was taken into account and studied in terms of 3D time dependent AWJM model. This gives us the opportunity to predict the shape of the milled trench surfaces. The required trench profile could be obtained with high precision in lack of knowledge about the model parameters and based only on the experimental measurements. We use the adjoint approach to identify the AWJM model parameters. The complexity of inverse problem paired with significant amount of unknowns makes it reasonable to use automatic differentiation software to obtain the adjoint statement. The interest in investigating this problem is caused by needs of industrial milling applications to predict the behavior of the process. This study proposes the possibility of identifying the AWJM model parameters with sufficiently high accuracy and predicting the shapes formation relying on self-generated data or on experimental measurements for both evenly jets movement and arbitrary changes of feed speed. We provide the results acceptable in the production and estimate the suitable parameters taking into account different types of model and measurement errors.

1. Introduction

The abrasive waterjet (AWJ) machining is a nonconventional low-cost process [1] that was developed to give the opportunity to manufacture complex shapes with difficulty in cutting and milling materials regardless of their properties [2, 3]. This machining technique embraces low cutting forces [4] on the target workpiece (reducing the possible risks of damaging the sample) and does not have heat affected zone. The Abrasive Waterjet Milling (AWJM) process involves the high-speed waterjet produced by the water pump with a small nozzle and abrasive particles included in the jet flow. This forms the circular high-energy jet plume of water and abrasives that impacts on the target surface and erodes the material. The behavior of the material removing and the intensity of the etching rate can be changed and controlled by different machine parameters such as the feed speed of jet movement, pump pressure, and mass flow rate of the abrasive particles. All these parameters can be mathematically characterized by etching rate function which is inaccessible from the experiments and plays key role in the modeling and prediction of the trench surface.

One of the most challenging and crucial questions among the industrial and manufacturing problems which can be interpreted with partial differentiation equations (PDEs) is the identification of the optimal model parameters. The goal is to reproduce the required shapes and processes relying only on the available experimental measurements related to the real systems. In the conditions of complexity and nonlinearity of the problem, the determination process becomes one of the critical questions and leads to involvement of various techniques and approaches.

There were several reported studies in consideration of the direct problem when it is necessary to predict the trench surface with a given model and its set of parameters. Some well-known methods based on the statistical approaches [5] and finite element methods [68] have been used and reported. For solving the direct problem, which particularly involves a nonlinear PDE, some information about the model parameters such as coefficients or energy sources is required, but often most of them are unavailable or unknown in advance and need to be identified.

Even if the direct problems are linear under some considerations, however, the inverse problems of the model parameters identification are usually ill-posed [9, 10] and measurement noise and model errors impose regularization needs. The regularization techniques [1113] in the identification process can assist in performing and coping with such aspects. The posed minimization problem for a cost function (i.e., a mismatch between experimental measurements and modeled estimation) underlies the identification problem. Some common approaches and techniques for various general and particular problems have been used and reported in [1417].

In this paper, we extend the work presented in [18] about the mathematical method to identify required unknown parameters of the generic Abrasive Waterjet Milling (AWJM) model. This model was previously developed and reported in [1922] according to the industrial needs for microwaterjet footprints prediction.

The inverse problem consists of the identification of AWJM model parameters from the experimental observations. These results have to be further used to simulate the required surface profile. Recent research of linear AWJM inverse problems focused on the identification of the beam path has been previously reported in [23].

In our work, the gradient vector of the cost function, which is required for all the family of the gradient descent algorithms, is found numerically by use of the automatic differentiation software TAPENADE [24]. Numerical optimization method based on the limited memory Broyden–Fletcher–Goldfarb–Shanno (BFGS) algorithm [2527], realized in the minimization package N2QN1 from the INRIA MODULOPT library [28], is used to implement the minimization problem. In terms of the adjoint problem, the minimization problem can be represented as an optimal control problem based on the Lagrangian multiplier of the model equations [29, 30].

The parameter identification problem significantly depends on input measurements and is very unstable, which is shown by inclusion of noise in the generated data. In addition, we indicate the influence of using Tikhonov regularization on accuracy of the surface prediction and improvement of the AWJM model parameters identification in case of noisy data.

The paper is organized as follows. Section 2 consists of representation of the mathematical model describing the AWJM process and explanation of the adjoint approach, which is used to obtain the gradient of the cost function. Further in this section we give a short description of gradient descent algorithms, which are used for the minimization problem. In Section 3 we present the model parameters identification and actual numerical results for the moving waterjet with fixed and varying feed speeds. Results based on artificial and experimental observations are provided. Section 4 presents the sensitivity study of the given AWJM model, where the influence of the various measurement noises is studied and several approaches to improving the accuracy in the surface prediction are demonstrated. This paper finalizes by Section 5 with some conclusions and outcomes.

2. AWJM Direct and Adjoint Model

2.1. Proposed AWJM Model

The milling process perpetrated by abrasive waterjet machine is represented as a nonlinear partial differential equation with initial and boundary conditions. This model characterizes the process of the trench surface formation by the jet impact on the workpiece and is suitable for various jet feed speeds independently of the target material properties. To define the problem, we suppose the time interval of the continuous milling process and denote by a bounded domain of where the process takes place.

The proposed Abrasive Waterjet Milling model, coming from previous works [1921] and already partially studied in [18], is presented aswith initial and boundary conditions:

The given AWJM model in (1) describes the trench creation by impact of an abrasive waterjet of radius with the forces caused by etching rate function on the primarily flat surface.

The final form of the trench, which is described here as a function , depends on the different physical parameters such as pump pressure, abrasives mass flow, velocity, and waterjet nozzle diameter. According to the announced mathematical model, these machine settings could be described and represented as set of model parameters , defining the intensity of the jet impact. The vertical position of the jet is fixed during the process relative to the zero level of the workspace, and the intensity of the jet impact continuously depends on the trench depth.

The schematic representation of the problem is illustrated in Figure 1.

2.2. Cost Function and Adjoint Model

In order to identify the optimal AWJM model parameters that are needed for the reconstruction of the required surface, we formulate a minimization problem. Thus, we pose a problem to find such thatwhere the cost function consists of the difference between the model solution and experimental observations and additional regularization termunder the constraint that is the solution of the direct problem, obtained with input set of parameters .

In expression (4), are the experimental measurements, by we denote an a priori estimation of the set of AWJM model parameters , and is the Tikhonov regularization coefficient.

To obtain the solution of minimization problem (3), we consider it an optimal control problem, when the Lagrangian multipliers based approach can be used [29, 30]. Hence, we find the solution of the optimal control problem by examining the critical “points” of the Lagrangian functional associated with the constrained minimization [18], with Lagrangian multiplier under the constraint that is a solution of (1).

To solve the minimization problem and find an optimal solution , the gradient of the cost function is needed for the iterative gradient descent minimization algorithms from the quasi-Newton family.

It is necessary to consider not the continuous but the discrete system to figure out numerically the optimal problem and to find its solution. Actually we have to minimize the discrete cost function which requires the gradient of the discrete cost function. We involve the automatic differentiation software (i.e., TAPENADE) to obtain the gradient of the discretized cost function. Once the gradient is computed we can solve minimization problem (3) using it.

2.3. Minimization Process

Mainly most minimization techniques, which are used to evaluate approximate gradients for constrained problems and to find local minimum of cost functions, are iterative gradient descent algorithms. The quasi-Newton type methods to compute the approximate gradient and descent step of the minimization process have been applied due to complexity and high costs to compute the Hessian on each iteration. To perform the minimization process, we use the N2QN1 minimization package for constrained optimization problems from “MODULOPT” library [28] which considers the L-BFGS update algorithms [18].

The L-BFGS [26] algorithm, named for limited BFGS quasi-Newton approximation, uses an approximation of the inverse Hessian. This method truncates the standard BFGS update to store and use the last values of and , thus reducing the required memory. The update formulas for L-BFGS algorithm can be found in [25].

3. Identification of AWJM Model Parameters

3.1. Evenly Moving Waterjet
3.1.1. Standing with Self-Generated Data

In case of evenly moving abrasive waterjet, we assume the proposed AWJM model in (1) with initial and boundary conditions which describes the formation of the trench surface by the impact of the waterjet beam in the ideal conditions when no measurement or model errors are considered:

In this subsection we study two different cases, when the input data are considered self-generated surfaces or averaged experimental observations.

For the numerical implementation we define a domain , where always depends on the actual experimental parameters or measurements, due to the changes of the jets radius, and . In fact, we restrain the minimization problem to a squared part in the middle of the trench . From the experiments done in collaboration with STEEP project partners, we set  mm. A regular grid of points with the steps  mm is chosen for discretizing . For both cases, the time interval is taken as unit with . The general centred difference approximation for several variables in space and forward difference in time is used for the numerical implementation.

In order to obtain a smooth solution ( instead of ), we change the regularization term to the one with the gradient of the etching rate function:

The value of Tikhonov regularization multiplier suitable for this concrete problem was obtained by L-curve method [31, 32].

The “pseudoexperimental” surface was generated with arbitrary values of model parameters and and etching rate function defined on . To overcome the possible increase of the computing costs and instability in behavior of the etching rate function, we estimate it as the circularly symmetrical projection of the centred vector by the cubic spline interpolation [33].

Based on the demonstrated correctness of the identification process reported in the previous work [18], the determined etching rate function and comparison of the corresponding trench cross-section with the input data are presented in Figure 2. The minimization process is performed using N2QN1 minimizer from the “MODULOPT library” [28].

From the results shown in Figure 2 we can notice that due to the ill-posedness of the inverse problem the final obtained shape of the etching rate function differs from the original but is suitable for AWJM model and leads to very high accuracy in reconstruction of the trench surface (Figure 3) in the ideal conditions under the consideration that no measurement noise or model errors are included.

Further, we base our search of the AWJM model unknowns on the trench surface obtained from the real experimental measurements by extending the average cross-section in the direction of the jet movement (Figure 4). This input data has some peculiarity caused to some extent by specific of the milling process and corresponds to milling process with a jet feed speed of 2000 mm/min and nozzle diameter of AWJ machine of 0.5 mm. This value was used as a background to estimate the model parameter . All the numerical settings stay the same as in the previous experiment except the value of the Tikhonov regularization factor , which is now equal to .

Unlike the previous occasion, we have no estimation of the etching rate function that was used in the produced experiment; thus we start the determination from the zero assumption .

Here, we provide the results of the identification of AWJM model parameters, inaccessible from the experiments, which should be used in the direct simulation to reproduce the required workpiece shape. The identified etching rate function (Figure 5(a)), which describes the formation of the waterjet energy beam, takes acceptable uniform shape, where the highest effect is focused in the centre of the beam. The use of these results in direct simulation gives us the next match between cross-section of the numerical solution and experimental measurements (Figure 5(b)). Thus, the required surface can be numerically reconstructed with an accuracy in terms of norm smaller than 6% by the following expression:

One can observe the mismatch on the edges of the slopes of the trench, but this aspect was not considered and modeled in the used AWJM model in (1). In practice, these effects are explained by the redeposition of the target material and appear as the result of high power of the waterjet impact.

3.1.2. Experimental Measurement Based Identification

According to the considered mathematical model, we assume the constant movement of the jet in straight direction and in this subsection we base our determination process on the original measurements of the real experiments done with waterjet machining tool (Figure 6) but not on the average of the trench profiles.

Usually parameter identification problems induce various difficulties caused by model errors and rather measurement noise. To be able to search unknowns from rough and noisy initial measurements, we include in the AWJM model the error term , which represents the measurement errors as random variables with a Gaussian probability density function and a zero mean.

For the direct simulations, we use the following model:where is the factor corresponding to the percentage of the applied calibrated uncorrelated noise, .

Here the available experimental measurements differ from the previous ones and correspond to waterjet milling process with the jet feed speed of 3000 mm/min. Due to provided data, we define the squared domain by setting  mm with the steps  mm related to the chosen part of the milled trench, where we minimize the cost function.

Results of the determination of the etching rate function and comparison of the reproduced surface with original profile are given in Figures 7(a) and 7(b). We can notice that, respectively, to the decrease of the density of the input measurements we have the small loss in the smoothness of the final function shape which causes not significant loss in the surface restoring precision. Despite this, we still have very high accuracy in the surface prediction using presented identification technique which is expressed as level of error (less than 5%) in terms of norm.

3.2. Waterjet Feed Speed Variations

To extend the possibility of the application of demonstrated identification mechanism in the manufacturing, we include another particular case belonging to the variations of the waterjet feed speed during the milling process. To meet the practical capabilities of waterjet machine, we assume that it accelerates constantly during the movement from the initial position to the final one. For the numerical implementation, it can be described as a change of the time spent by jet beam on each position of the workpiece where we examine the problem.

We define the domain by choosing  mm to capture the most effective area where the trench profile changes due to the jet speed according to available experimental measurements. Previously, the regular grid with the steps  mm was used to discretize the domain to satisfy numerical simulations with the experimental measurements. The existing input data (Figure 8) correspond to the abrasive milling process done with AWJ machine with permanent acceleration from the feed speed of 600 mm/min to 2000 mm/min. In order to level and smooth the noisy input data, a noise filter was applied based on the averaging of measurements of the several trenches done with identical machining parameters. The real feed speed of the waterjet was recorded and correspondingly adjusted to the numerical simulations.

A Tikhonov regularization factor was chosen with L-curve method for this specific problem. This allows us to regularize the solution and to accelerate the minimization process.

Results of the identification of etching rate function and comparisons of the simulated central trench profile with the experimnetal input are shown in Figures 9(a) and 9(b), respectively. These results confirm the possibility of identifying the unknown etching rate function suitable for the AWJ machining to predict the trench surface in case of varied waterjet feed speed. Due to the increase of complexity of the problem and changeability of the surface profile, one can observe some decrease of the accuracy in the surface prediction in comparison with the results related to the evenly moving waterjet. Even so we still have the opportunity to predict the shape formation with acceptable level of error around 7% in terms of norm.

4. Sensitivity Study in Case of Uniform AWJ Movement

4.1. Variety of Input Measurements

Sensitivity study plays one of the key roles in the plenty of the parameters identification problems. By this, it is possible to deeply understand the behavior of the model and improve the correctness of the identification process. Observation of the influences and sensitivity of the AWJM model on measurement or model errors provides us with the opportunity to see what the possibilities of reconstructing the required shape of the trench are regardless of the input data.

In this subsection we demonstrate and compare the numerical results of the proposed approach to identify the etching rate function under different levels of measurement errors in the simulated data. Using the given AWJM model (9), we generated various trench surfaces with predetermined parameters , , and , respectively, to different levels of uncorrelated noise up to 40%. Further, these noisy trenches are used separately or in superposition as the only input of the suggested identification method to find out the unknown AWJM model parameters suitable for the surface prediction requirements. For the numerical implementation, the same as in Section 3.1.2, initial assumptions and parameters of the discretization are used except  mm which are used in order to increase the density of the grid and improve the accuracy. We also assume that initial etching rate function is circularly symmetric and obtained by the cubic interpolation from the centred vector.

As explained above, the simulated input data is obtained by adding a Gaussian white noise of various levels of intensity to the initial trench surface (e.g., Figure 10(a) demonstrates the case of 15% noise), generated by use of the etching rate function (Figure 10(b)).

The general purpose is to identify the unknown parameter . The use of this value in direct AWJM model (5) will produce the closest trench to the initial one. The middle cross-sections of the experimental trenches are named “Target” in the figures in this section. In order to find the smooth solution which is acceptable in manufacturing, we base the minimization process on the cost function in (7) using Tikhonov regularization term on the gradient, which will bring the essential absence of high oscillations. Regularization coefficient has to be reestimated through L-curve method because of the modification of the grid size.

Results of the identification of the etching rate function based on a single trench measurement in case of applied noise with levels of 5%, 15%, and 30% and comparison of the middle cross-sections of noisy, original, and rebuilt trenches are given in Figures 11(a)11(c). Here let us notice that it is still possible to find out the unknown AWJM model parameter even from the measurements with very high level of the noise. The use of these obtained values leads to high enough accuracy in the modeling and prediction of the trench surface despite the not ideal matching of the identified function to the initial one.

Further we assume that there are several different available experimental measurements of exactly the same trench that can be used to identify the unknown AWJM model parameters and to model the required surface. We generate them identically with the same parameters, but the distribution of the noise is always random, so the difference between them is only the random noise applied to the initial surface. To diversify the study we assume two different cases when there are three and ten available measurement inputs, which are shown in Figures 12(a) and 12(b).

Here the identification is based on the minimization of the cost function, which measures the difference between numerical solution and each of the experimental observations. In both cases, our cost function transforms towhere is the number of available trench measurements which are taken as the input for the identification process.

The use of several independent trench measurements leads to the following results for the identification of the unknown function and reproducing the required trench surface (Figure 13) in case of 15% of the noise in the input data, based on ten independent measurements.

The given numerical results of the surface prediction (Figure 13) demonstrate very similar accuracy of matching the modeled surface to the original trench regardless of the use of more inputs and do not bring very high change in the surface reconstruction precision. From the other side, it means that our identification approach allows determining the AWJM model unknowns fairly truly even with only one available input, despite the high oscillations in the data caused by the measurement noise. But, looking more attentively at the identified etching rate function , smoothness of the solution has been increased. This effect may be essential in the further realization of the microwaterjet milling process in the real manufacturing, where the form and behavior of the beam energy parameters need to be strictly determined and suitable for machine parameters.

Considering several independent trench measurements can be interpreted as the averaging of the surfaces in some sense, but to clarify this aspect we demonstrate another case of our identification problem with the use of several superposed trench measurements. In theory, the use of the average of the noisy trenches will provide less rough and noisy input data and will lead to the identification of model parameters more precisely, which in its turn implies the better reproducing of the required surface. Based on that proposition, we introduce the superposition of the experimental observation, taken from the previous test and introduced in the cost function as follows:

A comparison of the various configurations of the cost function is represented in Figure 14.

Holding the acceptable level of accuracy (less than 10%) in the surface reconstruction in comparison with the experimental measurements, using several inputs or either of their averages in (10) and (11) provides better opportunity to match the required trench profile. The use of the average of the input measurements demands in its turn to adapt the regularization coefficient due to the change of the behavior of the input.

The difference between using one and several measurements is not very impressive due to random nature of the noise applied to the input and could be strongly increased by involving hundreds of experimental observations to reduce the influence of the errors and by adjusting the regularization coefficient according to the averaging of the input. The given overview of the identification based on 1, 3, or 10 trenches shows us that, from the other side, we can identify unknowns with reasonable accuracy even with only one trench measurement.

Cost functions (11) and (10) for are not identical, but theoretically they have the same gradient. Nevertheless, numerical implementation of these two cases shows the difference in obtained results and gives the flexibility to find more suitable way for each particular problem. More detailed and precise results of the surface prediction are given in Table 1.

The analysis of the obtained results induces thinking about the particular random distribution of the noise, applied with high level to the original input, which has very high influence on the identification process. The trench surface simulated with the use of the identified etching rate function is quite close to the input (noisy or average of several trenches) in all the cases, which were not aligned and fitted to the original surface due to the fairly random distribution of the noise. It engenders the conclusion that the use of much larger number of measurements can negotiate the noises or make them more uniform, calibrate and fit by this the inputs to the original data, and improve the accuracy of the surface prediction.

Mostly, the use of several trenches (and their average) instead of only one can essentially improve the accuracy in the parameters identification, conducting to reduction of the errors in the surface prediction up to 20% in cases of adverse available inputs. Certainly, it should be noted that sometimes only one measurement is available, and it might be enough to obtain the model parameters required to reconstruct the profile.

4.2. Noise as Model Parameter
4.2.1. Identification of the Measurement Noise

The identification of the unknown AWJM model parameters in Section 4.1 demonstrates a high dependency on the accuracy of the level of the noise. One of the possibilities of improving the quality and accuracy of the identification process is to take into account and identify the measurement errors by considering them as unknown model parameters. We now assume that and we consider the following new cost function:where is the Tikhonov regularization coefficient corresponding to the measurement noise.

To ensure first the correctness and possibility of the identification of the existing noise in the input data, we first use the “true” values of as the initial estimation of the etching rate function and focus on the identification of noise, which is represented as normally distributed random values among all the working domain . Numerically, it leads to the growth of the amount of unknowns and slows down the minimization process. Here and in all the next numerical tests in this section we use the same parameters and assumptions as in Section 4.1.

Firstly, it is necessary to check if the minimizer is able to faithfully determine the randomly distributed values of the model parameter. The initial noise applied to the trench with the level of 5% is shown in Figure 15(a), and the result of the identification is given in Figure 15(b). The main achievement of this experiment is to show the availability of the minimizer to identify quite acceptable values of the simulated measurement errors regardless of error intensity.

Usually there is no information about the behavior and type of the etching rate function that should be used to predict the right form of the surface. Moreover, most of the available measurements are noisy and unclear, so we now need to identify the unknown model parameter and measurement errors at the same time, in order to improve the quality of the surface reconstruction.

We start with some assumptions about the form of the etching rate function , which we have already obtained in Section 4.1, to simplify the continuation of identification process.

The results of the identification of the measurement errors for the cases of 5% and 30% are presented in Figure 16. The comparison of the cross-sections of initial and identified etching rate functions and trench prediction results for the 5% case are given in Figure 17.

Considering the identification result for all the range of noise levels, we can note that, with the decrease of the measurement errors, their influence on the surface formation decreases as well and becomes less significant. It leads to the modification of the form of identified noise, which can be distinctly seen in Figure 16(a), where the noise takes in some sense the shape of the trench surface. This kind of results can be very useful to understand how the final identified form of (Figure 17(a)) should be actualized to reconcile the form of the noise and to improve the accuracy in the surface prediction via direct simulation with the use of optimal AWJM model parameters. From the other side, when the measurement errors are much higher (e.g., 30% in Figure 16(b)), its the values and influence greatly rise up in comparison with etching rate function , and it turns into several problems for the minimizer to correct the exact unknowns in the right ways and inability to adjust model parameters more precisely. But even with actual results the level of errors in the surface prediction is less than 4% in terms of norm.

4.2.2. Removing the Measurement Noise

One more interesting aspect of this work is the ability to improve the surface reconstruction by the improvement of the input data. Assume now that we identified quite acceptable and useful values of the measurement errors (e.g., Figure 15(b)), which affect most strongly the identification of the AWJM model parameter and the accuracy of surface prediction. We can further use this information to reduce the measurement errors, even if it is not ideally determined, from the initial input and obtain much smooth and clear trench (e.g., Figure 18(a) in case of 20%).

After the use of such manipulation, we can perform again the identification of the unknown AWJM model parameter as demonstrated in Section 4.1. By such decision we are able to reduce the influence of the measurement and model errors and really enhance the accuracy of the surface prediction; see Figure 18(b). The complete description of the identification results for all the span of the error levels is given in Table 2. One could see that almost for all the cases our identification approach gives the very high accuracy in the surface prediction regardless of the type of input record and let us predict the surface in direct simulation with level of error less than 3%. In cases of 15% and 20% of noise, we have the increase of the precision in more than two times, while for 30% and 40% of noise we improve the accuracy of the surface prediction in more than four times.

5. Conclusion

Parameter identification is a highly challenging problem in AWJM problem, particularly from the noisy experimental measurements and in case of uneven movement of the waterjet with varied feed speed in 3D case. In this paper we demonstrated the possibility of using the application of inverse problems theory, based on minimization problems, in the real manufacturing problems to estimate the process behavior and forecast the trench shape formation. The general high precision of the AWJM model parameters identification provides good opportunity to predict and simulate the milled trench surfaces regardless of the quality and density of available experimental observations. We showed the capability of the proposed method to cope with different cases independently of type and size of the input data, depth of the milled trench, microwaterjet feed speed, kind of the jet movement, and level of the measurement noise. We gave an overview of how an even minor and insensitive level of noise can affect the accuracy of the results and occasionally leads to considerable errors in surface reconstruction.

We presented the comparison of different approaches of the cost function formulation in accordance with various number of available data, which leads to several particular improvements in cases of high measurement errors. Moreover, we demonstrated the importance of the regularization terms, which have to be considered and carefully adjusted to obtain a more real and precise surface shape. In order to control the surface prediction under noisy conditions, we introduced and implemented a technique to identify the measurement noise independently of the other model parameters and to remove it from the input data of the minimization problem, thereby increasing the quality of the identification. Also, the combined identification of all the possible and not fixed model parameters was presented in this work, to explain how widely this approach can be used.

Competing Interests

The authors declare that there are no competing interests regarding the publication of this paper.

Acknowledgments

The authors would like to acknowledge the funding support of the EU-FP7-ITN (Grant no. 316560) for the works presented as a part of the STEEP ITN project. The authors thank Mr. Pablo Lozano Torrubia from the University of Nottingham for his help in the experimental work.