Metamodel-based multi-objective optimization of a turning process by using finite element simulation

ABSTRACT This study investigates the advantages and potentials of the metamodel-based multi-objective optimization (MOO) of a turning operation through the application of finite element simulations and evolutionary algorithms to a metal cutting process. The objectives are minimizing the interface temperature and tool wear depth obtained from FE simulations using DEFORM-2D software, and maximizing the material removal rate. Tool geometry and process parameters are considered as the input variables. Seven metamodelling methods are employed and evaluated, based on accuracy and suitability. Radial basis functions with a priori bias and Kriging are chosen to model tool–chip interface temperature and tool wear depth, respectively. The non-dominated solutions are found using the strength Pareto evolutionary algorithm SPEA2 and compared with the non-dominated front obtained from pure simulation-based MOO. The metamodel-based MOO method is not only advantageous in terms of reducing the computational time by 70%, but is also able to discover 31 new non-dominated solutions over simulation-based MOO.


Introduction
The optimization of metal cutting processes in turning operations has been studied extensively in the literature. However, optimization studies including two or more objectives, i.e. multi-objective optimizations (MOOs), are limited. The approaches towards solving a multi-objective optimization problem (MOOP) were divided into two groups by Deb (2001). The first group of methods use the a priori approach. In this approach, the MOOP is transformed into a single-objective optimization by adopting a preference vector, generated by some higher-level information. The a priori approach has been employed in multiple studies (Aggarwal et al. 2008;Bhushan 2013;Bouacha et al. 2014;Chabbi et al. 2017;Khamel, Ouelaa, and Bouacha 2012;Tebassi et al. 2016), where a desirability function transforms multi-response turning problems into a single-objective optimization problem. Alternatively, several researchers adopted the theories of grey relational analysis (GRA) to solve multiple performance characteristics optimization problems based on the a priori technique (Bouzid et al. 2014;Gok 2015;Kazancoglu et al. 2011;Pawade and Joshi 2011).
The second group of methods use the a posteriori approach. This approach involves obtaining a set of solutions in the form of a trade-off front, where the desired solution is selected according to some higher-level information concerning the problem. In this approach, the decision maker can gain a better understanding of the decision variables and objectives, and the relations between the two. In addition, it provides the freedom to analyse the results before selecting the optimal solution. Evolutionary algorithms (EAs), owing to their characteristic of using a population of solutions that evolve in each generation, are well suited for the a posteriori approach in solving MOOPs.
To the best of the present authors' knowledge, there are only a limited number of studies on the MOO of turning operations. Most of the research is focused on single-objective optimization, or where an MOO is converted to a single objective by following the a priori method, as mentioned earlier. Researchers often perform a limited number of turning experiments designed by methods such as factorial, Taguchi, etc. The cutting parameters are optimized with statistical tools e.g. ANOVA, based on the experiments. Moreover, in some literature the experiments are approximated by using basic metamodelling methods, such as the response surface method (RSM), and the optimization is carried out on the metamodels instead.
Recent applications (2007)(2008)(2009)(2010)(2011) of evolutionary optimization techniques in the optimization of machining parameters have been reviewed by Yusup, Zain, and Hashim (2012). Out of all the articles reviewed in Yusup, Zain, and Hashim (2012), only seven studies focus on the MOO of turning or multi-pass turning operations using evolutionary algorithms Baskar 2010, 2011;Cus, Balic, and Zuperl 2009;Durán, Barrientos, and Consalter 2007;Srinivas, Giri, and Yang 2009;Sultana and Dhar 2010;Xi and Liao 2009). They all incorporate mathematical models of the turning operation or an approximate model from a limited number of experiments and only two studies optimize the objectives simultaneously based on the a posteriori approach. In one study by Durán, Barrientos, and Consalter (2007), the production rate and production cost of turning operations were optimized by using genetic algorithms, and a Pareto-optimal front was obtained. Sultana and Dhar (2010) minimized cutting temperature and cutting force in turning AISI-4320 steel by using an MOO algorithm based on genetic algorithms (GAs), subject to keeping the surface roughness less than a constant value. An experimental study along with predictive models were implemented using the response surface method. The cutting variables considered were cutting speed, feed rate, pressure and the flow rate of high pressure coolant. In a study by Pytlak (2010), which was not included in the aforementioned review article, unit production cost, time per unit, and the resultant cutting force of hard finish turning operations of hardened 18HGT steel were optimized. The modified distance method (MDM) is based on evolutionary computations generated the Pareto-optimal front from the objective function equations. The best solution from the non-dominated set was selected by using the hierarchical optimization method.
Reviewing the literature from 2012, a number of studies considering the MOO of turning operations are found. Metamodel-based optimization of surface roughness and tool wear rate in the turning of titanium metal matrix composites was conducted in Aramesh et al. (2013). Experiments using a three level factorial method having three variables were incorporated to build a metamodel using Kriging. The strength Pareto evolutionary algorithm was used for optimization. A non-dominated sorted genetic algorithm (NSGA-II) was employed in Ganesan and Mohankumar (2013) to find the optimal cutting parameters in CNC turning by optimizing three objective functions, minimum operating time, minimum production cost and minimum tool wear, expressed with mathematical equations. In a study by Thepsonthi and Özel (2014), the cutting force and tool wear of a high performance micro-milling process were generated by running 2D finite element simulations (in DEFORM-2D software). They were used as inputs for optimizing tool path and process parameters along with burr formation and surface roughness data extracted from experiments. The multi-objective particle swarm optimization (MOPSO) technique was employed for performing optimization. Santos et al. (2015) simultaneously optimized the machining force and chip thickness ratio when turning aluminium alloys by building second-order polynomial approximation models from a central composite experimental design. The genetic algorithm available in MATLAB ő software was used for optimization. Satyanarayana et al. (2015) used Taguchi's full factorial to design 27 experiments for turning a super alloy Inconel 718. The experiments were adopted to find two mathematical expressions for material removal rate and surface roughness, where NSGA-II obtained the optimal solutions in terms of cutting parameters. A later study (Satyanarayana, Reddy, and Ruthvik Nitin 2017) optimized tool wear and tool temperature by conducting a three level full factorial (27 experiments). Using regression models, two analytical expressions were formed as the objective functions, and a genetic algorithm was used for optimization. Abbas et al. (2016) performed an MOO study of cutting parameters in a turning operation. Surface quality and material removal rate were the objectives, while speed, feed rate and depth of cut were the variables. The Pareto-optimal front was obtained from the non-dominated solutions of a five level factorial experiment (125 experiments). In addition, the results were compared with a metamodel-based MOO (EGO) with a total of 80 experiments (64 initial and 16 validation experiments). In Klancnik et al. (2016), optimal surface roughness, cutting forces and tool resistance time were obtained by optimizing three mathematical models constructed with a gravitational search algorithm from 15 experiments.
In all the aforementioned sources, mathematical models, experimental work or metamodels that were built on experiments formed the basis of the study. On the other hand, research work in which FE simulations are coupled to MOO algorithms are very rare. Umer et al. (2014) minimized cutting force and tool-chip interface temperature by using three different surrogate models-RSM, radial basis functions (RBFs) and neural networks-and a multi-objective optimization genetic algorithm (MOGA2) implemented in MODEFRONTIER. The objectives were obtained by modelling the oblique turning process in ABAQUS FEM software. Recently, Sadeghifar et al. (2018) performed an MOO on radial turning operation by using the a priori approach (a weighted sum method). A hybrid optimization, coupling genetic algorithm and sequential quadratic programming, optimized a combined function of tool temperature, cutting force and residual stresses. The FE model was approximated by RSM and used in the optimization process.
The present authors' review of the literature shows that multi-objective optimization of turning operations based on FE simulations, built on the a posteriori approach by employing EAs, is a relatively new research topic. This is more noticeable when the optimization is based on metamodels. In almost all existing studies, the metamodels represent the experiments and the Pareto-optimal fronts are obtained from a metamodel-based MOO. The first drawback is that there is no proof that the metamodels employed in the previous studies represent the experiments accurately. Yet, the experiments are complex and time consuming. Hence it is not possible to perform a posteriori EA-based MOO with the actual experiments and find the true Pareto-optimal front. Therefore the second drawback, or a gap in the existing literature, is that the accuracy of the non-dominated solutions obtained from the metamodel-based MOO cannot be easily assessed, since the true Pareto-optimal front is unknown.
In the present research work, the intention is to fill this gap by analysing the performance and advantages of a metamodel-based MOO of a turning operation in comparison with an FE simulationbased MOO of the same operation. A recent study by Amouzgar et al. (2018a) developed a new framework for the automated MOO of a machining process based on FE simulations. The framework was demonstrated by optimizing a metal cutting process in turning AISI-1045, using an uncoated K10 tungsten carbide tool. The aim of the MOO was to minimize tool-chip interface temperature and tool wear depth, which are obtained from FE simulations, while maximizing the material removal rate. The optimization was based on five variables, including three geometrical variables (clearance angle, rake angle and tool cutting edge radius) and two process variables (feed rate and cutting speed). The strength Pareto evolutionary algorithm SPEA2 was employed as the multi-objective optimizer and the Pareto-optimal front was found after 17 generations. Without the use of the metamodel proposed in the current work, the optimization process took about two weeks.
In this work, the MOO of the same turning operation (cutting AISI-1045 with an uncoated carbide tool) is performed by approximating the FE simulations with metamodels. For this purpose, the Latin hypercube sampling method is utilized to create 100 designs of experiments (DoEs) with the same five variables: clearance angle, rake angle, tool cutting edge radius, feed rate and cutting speed. FE simulations of the turning operation in DEFORM-2D software are run to find the responses for the DoEs. Thereafter, seven different well-known metamodelling methods are constructed for two of the objectives (tool-chip interface temperature and tool wear depth). The best metamodelling method for each objective is selected by using two performance metrics. The MOO is performed on the selected metamodels with the same algorithm (SPEA2) and the trade-off front is obtained. The potential of metamodel-based MOO is discussed by analysing and comparing the non-dominated solutions on the front with the trade-off front of the simulation-based MOO. Thus, this work differs from the existing literature in the following ways: (i) multi-objective optimization without combining the objectives, (ii) multiple metamodelling methods with parameter tuning for finding optimum parameters, (iii) objective functions that are obtained by FE simulation rather than empirical formulae, (iv) validation of optimal solutions by simulation, and (v) a real-world turning operation simulated by state-of-the-art FE software.
In the next section, a brief overview of FE simulation of the turning operation in DEFORM-2D is presented. Then, Section 3 is dedicated to defining the seven metamodelling methods and formulating the performance metrics that are used to select the best method. In Section 4, the SPEA2 and the parameters in the algorithm are presented and the overall MOOP is formulated. In Section 5, parameter settings for the metamodelling methods are described and the results are presented and discussed. Finally, the concluding remarks are presented in Section 6.

FE simulation of the turning process
The FE model is based on the Lagrangian method. Figure 1 illustrates a schematic diagram of a 2D cutting operation. The variables considered in the study are the tool geometry parameters including clearance angle γ , rake angle α and tool cutting edge radius r, and process parameters including cutting speed ν and feed rate f , which are indicated in Figure 1. Three objectives are considered, which are extracted from the simulation after 7 mm of cut length, i.e. minimizing tool wear depth, minimizing the tool-chip interface maximum temperature and maximizing the material removal rate (MRR). MRR is calculated by where ν (mm/sec) is the cutting speed, f (mm/rev) is the feed rate and d (mm) is the depth of cut. The depth of cut is in the third dimension; thus, it is kept constant. The variables and objectives are shown in Table 1, and Table 2 shows the upper and lower limits of the variables.

Boundary conditions
The tool is fixed and the cutting is done by the movement of the workpiece towards the tool. The vertical speed of the workpiece is set to zero (fixed in the y-direction), and horizontal speed (in the  Variable  (Jaspers 1999). x-direction) is assigned to the nodes on the bottom edge of the workpiece. The two edges of the workpiece that are in contact with the tool have heat exchange with the environment and the temperature of the nodes in the two other edges are kept constant at 20 • C. In the same way, heat exchange for the tool is defined for the two edges in contact with the workpiece. The temperature of the nodes on the other two edges is set to room temperature.

Material properties
The turning process simulates the cutting of a plain carbon steel workpiece (AISI-1045) by an uncoated tungsten carbide (K10) cutting tool. The tool is assumed to be rigid and the material properties of the workpiece are shown in Table 3. The JC constitutive model is used to simulate the workpiece material behaviour within the range of strain rate, strain and temperature during the process by where σ is the flow stress, is the true strain,˙ is the true strain rate,˙ 0 is the reference true strain rate, T is the workpiece temperature, T 0 is the ambient temperature, T m is the workpiece material melting temperature and A, B, C, n and m are the model constants. The material model data for AISI-1045 are taken from Jaspers (1999), wherein the calibration of the JC constitutive model was carried out by using the SHPB high strain rate test. Table 4 shows the JC parameters for the workpiece incorporated in this study.

Thermal properties
The thermal properties of AISI-1045 taken from Halim and Maria (2008) and K10 acquired from Klocke (2011) are shown in Table 5. In order to reach the thermal steady state condition during the 7 mm of cut length used in this study, an ideal contact condition is assumed. This is accomplished by setting the heat transfer coefficient in DEFORM-2D to a high value (h int = 100, 000 kW/m 2• C). Table 5. Thermal properties of the workpiece and the tool material.

Contact and friction model
The contact between the chip and workpiece is assumed to be governed by the sliding model, i.e. τ = μσ n , where μ = 0.8 is the Coulomb friction coefficient (sliding coefficient) and σ n is the interface normal pressure. The tool and workpiece are assumed to undergo a hybrid contact consisting of two friction models, namely the sticking-sliding model defined by where L st and L sl are the sticking and sliding contact lengths, respectively. The constant sticking and sliding coefficients assigned for this study are m = 1.0 and μ = 0.6, respectively.

Tool wear model
Usui's wear rate model (Usui, Shirakashi, and Kitagawa 1984) is adopted in this study. The model is defined byω where σ n is the normal stress, T is the temperature and V s is the sliding velocity of the predicted nodal data of the tool contact surface. The wear constants D 1 and D 2 are given in Maekawa et al. (1989) for plain carbon steels and uncoated tungsten tools: D 1 = 7.8 × 10 −9 , D 2 = 5.302 × 10 3 .
Tool wear depth and maximum nodal temperature is extracted from simulations after more than half of the workpiece length has been cut by the tool, e.g. 7 mm of cut length. The FE step length, number of steps and total simulation time are defined accordingly for each simulation. The third objective, e.g. MRR, is calculated by using Equation (1).

Mesh convergence
In FE modelling, the size of the mesh affects the accuracy considerably. Therefore, based on the result of a mesh convergence study with six different mesh sizes and numbers of elements, a suitable mesh is selected and shown in Table 6. The simulation time for this mesh size varies between 250 and 600 min. The large variation of computational time is caused by different variable combinations for each simulation. Thus, by considering the capability of running several simulations at the same time within the developed framework in Amouzgar et al. (2018a), the computation time for obtaining the responses of the 100 DoEs, used for constructing the metamodel, was less than 48 hours.

Metamodelling
In this study, seven different methods are employed to construct metamodels of total tool wear depth and tool-chip interface temperature. The metamodels are evaluated by two performance metrics and the best method is chosen as the solver for MOO.

Metamodelling methods
The recently developed radial basis functions with a priori bias (RBF pri ) along with the conventional RBF (RBF pos ), the response surface method (RSM), Kriging (KG), support vector regression (SVR), neural networks (NNs) and multivariate adaptive regression (MARS) are used. A summary of each method is described in the Supplemental Data for this article, which can be accessed at https://doi.org/10.1080/0305215X.2019.1639050

Performance metrics
The two error metrics used to evaluate the metamodelling accuracy are as follows.
(1) Normalized root mean squared error (NRMSE) is given by where . . , f n s ), f andf are the true and predicted function values, n s and n t are the numbers of sample and test points. Root mean squared error (RMSE) indicates the deviation of the metamodel output from the actual function response and provides an overall indication of the global accuracy. RMSE is typically of the same order as the actual function values and indicates the relative performance of the methods across different functions independently. Therefore, to compare the methods over different test functions, the normalized value NRMSE is used. (2) Rank error (RE) is defined in Joachims (2005) and given by where SwappedPairs = |{(i, j) : When metamodels are used in evolutionary optimization algorithms, the ranking among pairs of solutions is crucial, rather than the actual function values. Therefore, rank error will demonstrate the performance of the metamodelling methods in this regard.
The two performance measures are compared by adopting the k-fold (in this study k = 10) crossvalidation method. The sample points are divided into k sets, and surrogate models are trained multiple times using samples in k − 1 sets, while leaving out the samples in one set as test points. This process is repeated k times for each metamodel and the median of each set of k performance metrics is used for comparison study.

Multi-objective optimization
MOO requires several evaluations of multiple objectives. There are several methods for solving MOOPs; however, evolutionary algorithms (EAs) are one of the most popular methods and have been applied in different fields (Das et al. 2011; Reddy, Abhyankar, and Bijwe 2011; Reddy, Bijwe, and Abhyankar 2014). EA methods start with a population that evolves through iterations by employing evolutionary operators, namely selection, crossover and mutation, that mimic the evolutionary mechanisms observed in nature. In this study, an in-house implementation is employed of the wellknown strength Pareto EA (SPEA2) (Zitzler, Laumanns, and Thiele 2001) as the MOO solver with the following parameter values: • Initial population size: 300 • External set population size: 300 • Crossover probability: [0.7, 0.95] • SBX crossover distribution index (Deb and Kumar 1995): [2,15] To cope with the stochastic nature of evolutionary algorithms like SPEA2, each optimization run was replicated 10 times, each replicated run conducted with randomly chosen parameter values within the ranges mentioned above. Next, all 10 sets of non-dominated solutions were combined and the best set of non-dominated solutions, with a population of 300 solutions, was found and reported as the Pareto-optimal set of the problem. The SPEA2 method has been shown to have some advantages compared with other existing techniques for problems with more than two objectives (Zitzler, Laumanns, and Thiele 2001). The method has previously been used for optimizing realworld engineering problems (Amouzgar, Cenanovic, and Salomonsson 2015;Amouzgar, Rashid, and Strömberg 2013).
The overall formulation of the optimization problem is Min ω(x), where x = [γ , α, r, ν, f ] is the vector of design variables and x L and x U are the lower and upper bounds of the design variables as defined in Table 2. The first objective function value is calculated by using Equation (1) and the second and third objectives are obtained through FE simulations of the turning process as described in Section 2.

Comparison of metamodelling methods
The iterative Latin hypercube sampling method is used to create a total of 100 DoEs. The FE simulation of the cutting process is run for all DoEs to find the corresponding tool wear depth and maximum tool-chip interface temperature. Three out of the total 100 simulations aborted before reaching the stopping criterion, which was the length of cut (7 mm). Therefore, 97 DoEs and their responses are employed to evaluate the performance of each metamodel with the aim of finding the best method to integrate with the SPEA2 algorithm. However, to perform an unbiased comparison, a parameter tuning process should be performed for all metamodelling techniques. The parameters that have to be specified by the user for each metamodelling method have a great impact on the accuracy of the models. Therefore, a parameter selection procedure should be applied to all methods. Hold-out validation, which is an effective and simple method, is used for tuning the parameters. The sampling data are grouped into two subsets, e.g. training and testing sample sets. Each metamodelling technique is trained by using the training set for different parameter combinations. The testing set is used to find the optimal parameters, which are the parameter combination that produce the lowest NRMSE for the validation set. The parameter combinations are generated by using  grid search to explore the entire parameter space. This procedure is done for both objectives, i.e. tool wear depth and interface temperature.
The parameters needing to be tuned for each metamodelling technique and their related ranges for the grid search are shown in Table 7. In general, there are two types of parameter requiring to be tuned, namely continuous and categorical parameters. The continuous parameters are taken from the range {10 i |i = −4, −3, −2, −1, 0, 0.3, 0.5, 0.7, 0.9, 1, 2, 3}.
The polynomial and the regression function in RBF pri and KG are selected from polynomials of order zero, one and two, represented respectively by poly0, poly1 and poly2. The most common type of categorical parameter is adopted for the basis functions in RBF pri , the correlation function in Kriging and the activation function in neural networks. The number of hidden nodes in the hidden layer of NN is given by an integer in the range [1,10] After determining the optimal parameter configurations, the performance metrics of all metamodelling methods with their corresponding optimal parameters are calculated by using k-fold cross-validation. The results are shown in Tables 8 and 9 for the tool wear depth and tool-chip interface temperature.
The RBF approaches and Kriging are the best methods with regards to RE for both objectives. In the NRMSE case, the RBFs and Kriging performed better than the other three methods (SVR, NN and MARS). However, for tool wear depth, RBF pri is the most accurate method, and Kriging predicts the interface temperature with the lowest error. Thus, to construct the metamodels employed in the MOO study, RBF pri is selected for tool wear depth and Kriging is selected for interface temperature.
RBF and Kriging are methods that use both integration and interpolation for fitting the metamodel. Based on the LibSVM documentation (Hsu, Chang, and Lin 2003), SVR requires a rigorous two-step parameter tuning process which, for the sake of a fair comparison, is not implemented in the parameter setting procedure of this study. Neural networks are known to be very sensitive to the parameters, especially the number of hidden layers. On the other hand, each additional hidden layer drastically increases the number of free parameters (weights) of the network to be learned through back-propagation. This makes the learning algorithm susceptible to getting stuck at sub-optimal parameters. MARS uses several parameters, and attempting to tune all of them is a rigorous study in itself. Hence, all parameters except for the cost penalty factor (c) are fixed to their recommended values.
The aforementioned results are consistent with the results obtained in a detailed comparison study by Amouzgar et al. (2018b) of the same metamodelling methods on several benchmark functions and a turning operation.  2. Pareto-optimal front in objective space obtained from simulation and metamodel-based MOO.

Comparison of simulation-MOO and metamodel-MOO
A simulation-based MOO on the same problem with identical settings was carried out to find the Pareto-optimal front in a recent study by Amouzgar et al. (2018a). In that study, the optimization algorithm converged after 17 generations. The non-dominated solutions of generation 17 reported as the Pareto-optimal front are shown in Figure 2 with * and called 'Simulation-MOO' in the legend. Each generation of the study was completed in around 20 hours (40 FE simulations in each generation) and the 17 generations, i.e. a total of 680 simulations (17 × 40), were executed in two weeks.
In this study, the metamodel-based MOO process is run with Equation (1) as the expression for MRR (the first objective) and RBF pri and Kriging as metamodels for tool wear depth and interface temperature (the second and third objectives), respectively. The final Pareto-optimal set containing the 300 non-dominated solutions was obtained from the overall 3000 non-dominated solutions generated from the 10 MOO runs, each running for 500 generations. The final 300 non-dominated solutions that create the Pareto-optimal front obtained from metamodelling-based MOO are also shown in Figure 2. By comparing the two fronts, the accuracy of the metamodels in capturing the true Pareto-optimal front (simulation) is evident. However, there are two regions in the front where the metamodels were relatively inaccurate. First, the solutions between MRR values of 500 to 800 and second region is the lower end of the front, where the metamodels predicted negative values of tool wear depth.
To validate the metamodel-based Pareto-optimal front, 100 non-dominated solutions (the validation set) are selected from the total 300 solutions based on a nearest neighbour density estimation technique (Zitzler, Laumanns, and Thiele 2001). The true second and third objective values of the validation set are obtained by running FE simulations, where 97 out of the 100 solutions were valid and they are depicted in Figure 3. The first region of the front mentioned above is still not covered by the validated metamodel-based Pareto-optimal front. However, the second region is now rectified by the updated metamodel-based front, where the solutions with negative ω values have shifted to the  positive side of the tool wear depth axis. Now, the validation set can be used as test samples to calculate the performance metrics for RBF pri and Kriging, which are shown in Table 10. Higher error values in Table 10 compared with the error measures obtained from the k-fold method in Tables 8 and 9 are observed. However, the errors still confirm a relatively good accuracy of the two metamodelling methods. An outcome of a metamodel-based MOO is the ability to discover new non-dominated solutions or find solutions that dominate the results obtained from simulation-based MOO. To investigate this aspect, the non-dominated population from simulation-based MOO (60 solutions) are combined with the validated solutions obtained from the metamodel-based MOO (97 solutions). The fast nondominated sorting algorithm is applied to the set to find the solutions with rank one in the front. The updated Pareto-optimal front along with the simulation-based MOO front is depicted in Figure 4. Figure 5 illustrates the relationships between the variables and objectives of the simulation-based MOO solutions and the newly discovered solutions by using metamodel-based MOO. The results are promising because 31 new non-dominated solutions are discovered of which 9 dominate the solutions obtained from the simulation-based MOO and the other 22 are completely new optimal solutions. This is more intriguing when the noticeable reduction in the computation time is considered. The original study required 680 FE simulations (the population size of 40 × 17 generations) while this study was carried out with 200 FE simulations, reducing the total computation time from two weeks to four days. The saved computational time can be employed in advancing the process, by using the adaptive metamodel-based design optimization as described in Wang and Shan (2007).    set of non-dominated solutions will be discovered. The process can continue with the intention of exploring more and better solutions.
The advantage of utilizing metamodels in MOO can also be assessed by calculating the hypervolume indicator (Zitzler et al. 2003). The hypervolume indicator is a measure to compare the non-dominated sets in each generation of optimization, or between two optimization methods. The hypervolume indicator gives the hyper-volume between the solutions in the Pareto-optimal front and a reference point. Figure 6 shows the hypervolume for the minimization of a two-objective problem. It is obvious that a front with higher hypervolume is closer to the true Pareto-optimal front, when the nadir point (Deb, Miettinen, and Chaudhuri 2010) is selected as the reference point. Comparison of hypervolume in each generation of both simulation and metamodel-based MOO, plotted in Figure 7, shows consistently higher values for solutions obtained by metamodel-based MOO. The graph of simulation-based MOO shows that the run was terminated after 17 generations because there was very little improvement in the hypervolume value. At the time of termination, the simulation-MOO had used 680 simulations. On the other hand, metamodel-based MOO (which used metamodels built using only 200 simulations) continued to improve the hypervolume up to 30 generations. The reduction in runtime due to the use of fewer simulations is far greater than the time required for running 13 (= 30 − 17) additional generations.

Concluding remarks
Metamodel-based MOO of cutting AISI-1045 with an uncoated carbide tool in a turning operation with three objective functions and five variables was performed. Two of the objectives were evaluated by finite element simulations. RBF pri and Kriging among a total of seven metamodelling methods were selected, based on a performance study, to approximate the objective functions. The metamodels were updated by validating a set of selected solutions with FE simulations from the Pareto-optimal front obtained after running the SPEA2 algorithm. Afterwards, the solutions in the combined set consisting of the initial DoE and the validation set were ranked and the non-dominated solutions obtained. The comparison of these solutions with the Pareto-optimal front created in the recent simulation-based MOO disclosed the benefits of using metamodels in multi-objective optimization studies of computationally expensive and complex simulations of turning operations. The optimization method based on metamodelling was able to find more non-dominated solutions in a fraction of the time when compared with optimization methods that are based entirely on FE simulations.
In conclusion, employing metamodels in MOO studies of turning operations with FE simulations and even experimental work is efficient (giving time reduction) and effective (giving better solutions). Combining physical experiments, FE simulations and metamodels in performing MOO of more complex machining operations is an area of interest for future work. Furthermore, metamodelbased MOO of other manufacturing processes such as drilling and metal casting is in the present authors' future plans.