On Multi-Objective Based Constitutive Modelling Methodology and Numerical Validation in Small-Hole Drilling of Al6063/SiCp Composites

Discrepancies in capturing material behavior of some materials, such as Particulate Reinforced Metal Matrix Composites, by using conventional ad hoc strategy make the applicability of Johnson-Cook constitutive model challenged. Despites applicable efforts, its extended formalism with more fitting parameters would increase the difficulty in identifying constitutive parameters. A weighted multi-objective strategy for identifying any constitutive formalism is developed to predict mechanical behavior in static and dynamic loading conditions equally well. These varying weighting is based on the Gaussian-distributed noise evaluation of experimentally obtained stress-strain data in quasi-static or dynamic mode. This universal method can be used to determine fast and directly whether the constitutive formalism is suitable to describe the material constitutive behavior by measuring goodness-of-fit. A quantitative comparison of different fitting strategies on identifying Al6063/SiCp’s material parameters is made in terms of performance evaluation including noise elimination, correlation, and reliability. Eventually, a three-dimensional (3D) FE model in small-hole drilling of Al6063/SiCp composites, using multi-objective identified constitutive formalism, is developed. Comparison with the experimental observations in thrust force, torque, and chip morphology provides valid evidence on the applicability of the developed multi-objective identification strategy in identifying constitutive parameters.


Introduction
Particulate Reinforced Metal Matrix Composites (PRMMCs) present low-temp properties, high strength-to-weight ratio, good wear resistance, and low thermal expansion coefficient, but poor machinability [1][2][3]. Due to the high-cost and time-consuming experiments in R&D of PRMMC components for achieving desired machining quality, FE simulations have been employed to investigate the mechanics of cutting, optimization of cutting processes, and redesign of cutting tool [4]. To produce the realistic simulation results, it is vital to develop a proper constitutive model of MMCs with capability to capture the mechanical behavior in cutting. A key focus for the establishment and application of a proper constitutive model is the correlation of experimental data, since it should be capable of representing the material mechanical behavior in the range of strains, rates of strain, and temperature in machining [5], and so it is necessitated to derive the material constitutive models and identify the parameters from a numbers of experimental observations in views of engineering applications. The interdependence and intercoupling amongst the constitutive variables, such as deformation history, strain rate, and temperature lead to the difficulties in quantifying material constants. It is considerably questionable to improve the classical phenomenological material models only based on some available experiment data, regardless of the reliability and measurement errors [6][7][8]. These improvements to the classic constitutive model in the original scientific hypothesis were determined simply by experimental data statistics. These improved models lacked in significant difference in nature and physical significance in constitutive equation form, and so rarely adopted in finite element modelling afterwards. The solutions to the constitutive identification for a certain material may consist in the fulfillment of the applicability for classical material model [9]. This is necessitated to make scientific assumption and the reasonable identification of constitutive parameters with suitable numerical tools.
Most engineering materials exhibit varying deformation behavior at low (quasi-static) and high (dynamic) deformation rates or temperatures. This brings great challenges to the establishment and identification of the material constitutive model that can describe the mechanical response under different loading conditions. The Johnson-Cook (J-C) constitutive model is one of the most commonly used semi-empirical phenomenological ones for describing the plastic deformation behaviors at high strain, high strain rate, and high temperate, especially suitable for the simulation of machining processes [10]. Its distinct advantages in application over other models rely on its exactness at macro-scale and simplicity since the constitute model, characterized only using a few parameters to be identified, can describe the material behavior rather well.
In the J-C constitutive equation, the flow stress σ is assumed to be a multiplicative function of equivalent plastic strain ε p , normalized strain rate ε 0 with . ε 0 and . ε being the reference strain rate and strain rate and temperature terms. The three terms, respectively, represent respectively the elasto-plastic, viscosity, and thermal softening effects of the material plastic behavior, as described by Johnson and Cook [11].
The material constants (A, B, C, n, m) are identified by fitting stress-strain data from quasi-static compression and split-Hokinson pressure bar (SHPB) experiments. T * is a normalizing temperature term expressed by where T is the material temperature, T melt the melting point, and T room the room or reference temperature. One of the formalism characteristics for J-C model in Equation (1) lies in identifying the decoupling relationship between strain rate and deformation temperature.
In determining the empirical constants in J-C constitutive formalism, different strategies can be found in published literatures [12][13][14][15]. Amongst them, one strategy widely used in identifying material model for metals and alloys is an ad hoc strategy, with the following steps [15].
Identify and fix the initial yield stress A from the quasi-static test at the reference strain rate Identify the strain rate hardening coefficient C by fitting linearly the data of σ/ A + Bε n p against ln . ε * under the same strain at room/reference temperature for SHPB data from: Average a set of m values calculated under different strains or rates of strain, and find an estimation of m.
One difficulty in determining the material model, however, lies in the identification of initial yield stress since in some cases the yield point of most materials is not explicitly identified. This difficulty brings a challenge to the determination of yield point, particularly for metal matrix composites [16,17]. So, here it is more appropriate to think of the coefficient A as a variable to be optimized. Furthermore, during the process of identification using the proposed methodology in Ref. [15], another difficulty in calculating strain rate-hardening coefficient C is that at different strain rates the averaged C could not adapt to each other across the strain range [18], so it is infeasible to find a suitable average value with small variance. This may be partially attributed to the randomness from the test data, or the application and exactness of material models in different material types, time, and length scales. The identifying in thermal-softening parameter m employs the same solving procedure as determining the coefficient C. To accomplish this determination using ad hoc strategy, the experimental data at different temperatures, but identical strains rates are adopted. Again, the same difficulty may come into being: the coefficient being of different orders of magnitude or large variance arrived at from different strains or strain rates data. The extended formalism of J-C constitutive model with more fitting parameters, despite applicable efforts, would increase the difficulty in identifying constitutive parameters. The major difficulties in widespread application of constitutive model for industrial simulation consist in a larger number of mechanical tests for classical identification methodology of constitutive model, and the ensuing hard identification job of material parameters involving in constitutive characterization. Additionally, large fluctuation in the stress-strain curves, especially at high strain-rate loading, may influence the reliability and accuracy of fitting constitutive model. Thus, based on experimental data from quasi-static and dynamic tests, this work is focused on the identification strategy of J-C constitutive model parameters, and provides evidence on the applicability of the developed multi-objective optimization strategy in identifying model parameters. The feasibility of constitutive model obtained using weighted multi-objective strategy is verified through the comparison of experimental and finite element results of cutting forces and chip morphologies.

Test Procedure
Al6063 matrix composites reinforced with 65% volume fraction of SiC particulate (Al6063/SiC p /65p) are investigated to focus on the mechanical responses in quasi-static compression and SHPB tests. Figure 1 shows the microstructure of Al6063/SiC p /65p composites, with uniformly distributed SiC particles in the absence of obvious particulate clustering and preferred orientation. At low strain rates, the temperature rise caused by plastic strain at room temperature is considered to be negligible, while at high strain rates, there exists coupling between strain rate and deformation temperature since a great amount of heat is not timely transferred. To deduce the uncoupled relationship between temperature and strain rate, the isothermal mechanical behavior at elevated temperature is to be considered in the dynamic tests at high strain rates [19]. The quasi-static compression tests were carried out at a universal compression machine CMT4305 (MTS Systems Corporation, Eden Prairie, MN, USA). The specimens were machined into the cylinders of Φ6 mm × 9 mm in size. During quasi-static compression tests, the deformation rates were maintained at 10 −4 , 10 −2 , and 10 0 /s, respectively. The dynamic compression tests under high strain rates were conducted on cylinder specimens of Φ2 mm × 2 mm in size by means of SHPB. Both the quasi-static and dynamic compression tests are repeated at least three times. The SHPB testing setup at high temperatures was provided by the School of Aeronautics at Northwestern Polytechnical University (Xi'an, China). The detailed descriptions of high-temperature SHPB system and its workflow have been given in [20]. In SHPB tests, the stress-strain measure at the strain rates ranging from 10 2 to 10 4 was performed at room temperature. Also, the mechanical responses at loading temperatures of 300 °C and 500 °C at high strain rates were taken into account. For the hightemp SHPB tests, the specimens were first preheated to the desired test temperatures through the heating unit shown in Figure 2, and then tested in accordance with the conventional SHPB procedure [21].  The quasi-static compression tests were carried out at a universal compression machine CMT4305 (MTS Systems Corporation, Eden Prairie, MN, USA). The specimens were machined into the cylinders of Φ6 mm × 9 mm in size. During quasi-static compression tests, the deformation rates were maintained at 10 −4 , 10 −2 , and 10 0 /s, respectively. The dynamic compression tests under high strain rates were conducted on cylinder specimens of Φ2 mm × 2 mm in size by means of SHPB. Both the quasi-static and dynamic compression tests are repeated at least three times. The SHPB testing setup at high temperatures was provided by the School of Aeronautics at Northwestern Polytechnical University (Xi'an, China). The detailed descriptions of high-temperature SHPB system and its workflow have been given in [20]. In SHPB tests, the stress-strain measure at the strain rates ranging from 10 2 to 10 4 was performed at room temperature. Also, the mechanical responses at loading temperatures of 300 • C and 500 • C at high strain rates were taken into account. For the high-temp SHPB tests, the specimens were first preheated to the desired test temperatures through the heating unit shown in Figure 2, and then tested in accordance with the conventional SHPB procedure [21]. The quasi-static compression tests were carried out at a universal compression machine CMT4305 (MTS Systems Corporation, Eden Prairie, MN, USA). The specimens were machined into the cylinders of Φ6 mm × 9 mm in size. During quasi-static compression tests, the deformation rates were maintained at 10 −4 , 10 −2 , and 10 0 /s, respectively. The dynamic compression tests under high strain rates were conducted on cylinder specimens of Φ2 mm × 2 mm in size by means of SHPB. Both the quasi-static and dynamic compression tests are repeated at least three times. The SHPB testing setup at high temperatures was provided by the School of Aeronautics at Northwestern Polytechnical University (Xi'an, China). The detailed descriptions of high-temperature SHPB system and its workflow have been given in [20]. In SHPB tests, the stress-strain measure at the strain rates ranging from 10 2 to 10 4 was performed at room temperature. Also, the mechanical responses at loading temperatures of 300 °C and 500 °C at high strain rates were taken into account. For the hightemp SHPB tests, the specimens were first preheated to the desired test temperatures through the heating unit shown in Figure 2, and then tested in accordance with the conventional SHPB procedure [21].   Figure 3a,b illustrates the stress-strain curves of Al6063/SiC p /65p composites in quasi-static and dynamic SHPB compression tests, respectively. It is to note that in quasi-static mode, Al6063/SiC p /65p composites show rate dependence to stress response, which differs from pure Al6063 alloy with no rate-dependent mechanical behavior in quasi-static compression [22]. This is mainly due to the strengthening effect caused by the impedance of the Al6063 matrix dislocation motion by SiC particulates and the transition of load bearing role from Al6063 matrix to SiC reinforcement in high volume fraction SiC reinforced composites [23]. The strain-rate hardening effect for the deformation rates ranging from 10 −4 to 10 −2 is not significant, so the reference strain rate in fitting constitutive model for Al6063/SiC p /65p composites is identified as 10 −2 . As illustrated in Figure 3b dynamic mode, the mechanical behavior of Al6063/SiC p /65p composites is shown to be of significant rateand temperature-dependence. The higher the deformation rate is, the greater the data fluctuation of the measured stress-strain curves is; the higher the deformation temperature is, the less the fluctuation is. Therefore, the measurement noise at high strain rates and low temperatures become more considerable in dynamic mode, which influences the reliability and accuracy of fitted data in constitutive identification.  Figure 3a,b illustrates the stress-strain curves of Al6063/SiCp/65p composites in quasi-static and dynamic SHPB compression tests, respectively. It is to note that in quasi-static mode, Al6063/SiCp/65p composites show rate dependence to stress response, which differs from pure Al6063 alloy with no rate-dependent mechanical behavior in quasi-static compression [22]. This is mainly due to the strengthening effect caused by the impedance of the Al6063 matrix dislocation motion by SiC particulates and the transition of load bearing role from Al6063 matrix to SiC reinforcement in high volume fraction SiC reinforced composites [23]. The strain-rate hardening effect for the deformation rates ranging from 10 −4 to 10 −2 is not significant, so the reference strain rate in fitting constitutive model for Al6063/SiCp/65p composites is identified as 10 −2 . As illustrated in Figure 3b dynamic mode, the mechanical behavior of Al6063/SiCp/65p composites is shown to be of significant rate-and temperature-dependence. The higher the deformation rate is, the greater the data fluctuation of the measured stress-strain curves is; the higher the deformation temperature is, the less the fluctuation is. Therefore, the measurement noise at high strain rates and low temperatures become more considerable in dynamic mode, which influences the reliability and accuracy of fitted data in constitutive identification.

Multi-Objective Strategy Considering Weighted Measurement Errors
The J-C model in Equation (1) involves a set of five parameters ( = ( , , , n, m )) to be determined by matching numerous sets of measured data at quasi-static and dynamic tests. In this paper, the above parameters set identification can be turned into minimising the weighted summation of squared of aggregate errors or residuals function between the parameterized constitutive equation as a function of three independent variables ( = ( , , )) and the sets of experimental data. That is to solve nonlinear least squares problems over all the data points. Here, a Chi-squared error criterion is used as a measure of parameter optimization.

Multi-Objective Strategy Considering Weighted Measurement Errors
The J-C model in Equation (1) involves a set of five parameters (P = (A, B, C, n, m)) to be determined by matching numerous sets of measured data at quasi-static and dynamic tests. In this paper, the above parameters set identification can be turned into minimising the weighted summation of squared of aggregate errors or residuals function between the parameterized constitutive equation as a function of three independent variables (X = (ε p , . ε, T)) and the sets of experimental data. That is to solve nonlinear least squares problems over all the data points. Here, a Chi-squared error criterion is used as a measure of parameter optimization.
Materials 2018, 11, 97 where σ Exp is the experimental observation data of flow stress, σ Model is the constitutive model function of an independent variable vector X and five material parameters vector P to be identified (The relative error of A against the offset yield point σ(ε p = 0.002) is within a range of 10%, C and n are in the range of 0-1, B and m are non-negative), N is the total number of experimental data points, and W is the weighting matrix with the diagonal weighting factor ω i corresponding to the observation point i. The weighting matrix should be set to the measure errors of any observations, rather than simply to uniform weighting factor that is preferred in [16]. This is due to the fact that the measure noise under different combination of strain rates and temperature differs considerably in the total observations. Since J-C material model does take into account different mechanical responses under the quasi-static and dynamic loadings, the Chi-square error can be subdivided into the quasi-static term quas−static and dynamic term χ 2 dynamic as Thus, a seemingly bi-objective nonlinear least squares equation involving the quasi-static and dynamic terms. According to Equations (9) and (11), an overall objective function is represented as additive relationship with different weighting factors by Equation (12).
The fluctuation in some flow stress curves corresponding to high strain-rate loading conditions observed in Figure 3b influences the reliability and accuracy of fitting constitutive model, and the same phenomena have been reported in some references on particulate reinforced metal matrix composites [24][25][26]. The major reason should be because bad signal-to-noise ratios give rise to noisy experimental measurement, such that the experimental data at high rates are difficult to represent the material behavior accurately [26]. Through residuals analysis, it is found that the measure errors of experimental points under each loading condition are almost of the same order of magnitude, but not under other conditions. Additionally, the measure errors in each condition are shown to behave as a normal distribution. These facts hold essential implication that the measure errors under each combined condition of deformation rates and forming temperatures can be assumed to be identical, and so the experimental deviation between measured data and "true" value are treated as random measured noise that satisfies Gaussian distribution N 0, θ 2 .
ε,T . The distribution of measured noise is of varying standard deviation θ . ε,T with strain rates and temperatures. The relationship between the experimental observations σ Exp and fitted model data σ Model is satisfied, as follows.
Since the measure error covariance θ 2 .
ε,T in each condition before optimization analysis are unknown, the covariance analysis of measured error is to be performed, under the assumption of identical measure error in each condition. Similarly, So, based on the definition of weighting matrix, this seemingly bi-objective optimization actually in practice is a multi-objective minimization problem relying on the number of loading condition combinations. Thus, the following aim is to minimize the bi-objective reduced Chi-squared error summation for the flow stress data experimentally determined from different loading conditions. An updated Levenberg-Marquardt formula based on the gradient scaling is employed in solving Equation (12) [27].
where λ Norm is a normalized damping factor to the leading diagonal element in J T WJ + λI from standard Levenberg-Marquardt formula. The component J ij at the ith observation with respect to the jth parameter in the J matrix for J-C model can be defined as Materials 2018, 11, 97 8 of 22 Initialization and update of damping factor λ, Jacobian matrix J and step length h as well as convergence and error criteria are described detailedly in Appendix A. Figure 4 show the flowchart of multi-objective optimization strategy for identifying J-C parameters.

Start
Initialize all variables, including, σ Model (P), τ , X, σ Exp , , P int , P min , P max , Con et al. Set J to 0, and the counter of iteration N it and function evaluation N fe to 0.

No
Yes Figure 4. Flowchart of multi-objective optimization strategy for identifying J-C parameters.

Identification Results and Discussion
Based on the relationship of ln(σ − A) against ln ε in quasi-static mode, σ/(A + Bε n ) − 1 against ln . ε * and ln 1 − σ/(A + Bε n ) 1 + C ln . ε * against ln T * in dynamic mode, the ad hoc strategy in Equations (3)-(5) is employed for identifying material parameters A, B, n, C, m, as shown in Figure 5. A quantitative comparison of different fitted strategies and solving algorithms for Al6063/SiC p /65p material model identification are summarized in Table 1. When compared to ad hoc optimization strategy, multi-objective strategy applied in J-C model optimization leads to a higher predictability with a high value of R squared (R 2 ). This indicates that it is not always necessary to modify such a classic constitutive model under scientific hypothesis that the parameter coefficient A is a variable to be optimized rather than offset yield point. Besides, the overall fit standard error of 21.88 MPa indicates the model fitting errors are of the same order as the measured ones. That is to say, the multi-objective optimization model can be relatively capable of fitting experimental noise.

Identification Results and Discussion
Based on the relationship of ln( − ) against lnε in quasi-static mode , σ/( + ) − 1 against ln * and ln 1 − σ/( + )(1 + ln * ) against ln * in dynamic mode, the ad hoc strategy in Equations (3)- (5) is employed for identifying material parameters , , n, , m, as shown in Figure 5. A quantitative comparison of different fitted strategies and solving algorithms for Al6063/SiCp/65p material model identification are summarized in Table 1. When compared to ad hoc optimization strategy, multi-objective strategy applied in J-C model optimization leads to a higher predictability with a high value of R squared (R 2 ). This indicates that it is not always necessary to modify such a classic constitutive model under scientific hypothesis that the parameter coefficient A is a variable to be optimized rather than offset yield point. Besides, the overall fit standard error of 21.88 MPa indicates the model fitting errors are of the same order as the measured ones. That is to say, the multi-objective optimization model can be relatively capable of fitting experimental noise.  Figure 6 shows the experimental observations and model prediction using J-C model parameters that were identified by multi-objective optimization strategy. It can be seen from the Figure 6a that multi-objective model prediction agrees equally well with the experiment observations in quasi-static and dynamic deformation modes. The higher visual difference between the observations and model prediction at high deformation rates such as 1800 s −1 and 2200 s −1 is mainly attributed to the experimental measurement resulting from poor signal-to-noise ratio, which easily occurs in the SHPB tests for particulate reinforced metal matrix composites. Hence, the important significance of introducing the varying weighting factors under different loading conditions into formulating multiobjective function is manifested in this aspect. When compared to the ad hoc strategy, the multiobjective optimization one enables the identified J-C model to fit the observations under high temperatures and high strain rates more accurately, in Figure 6b.    Figure 6 shows the experimental observations and model prediction using J-C model parameters that were identified by multi-objective optimization strategy. It can be seen from the Figure 6a that multi-objective model prediction agrees equally well with the experiment observations in quasi-static and dynamic deformation modes. The higher visual difference between the observations and model prediction at high deformation rates such as 1800 s −1 and 2200 s −1 is mainly attributed to the experimental measurement resulting from poor signal-to-noise ratio, which easily occurs in the SHPB tests for particulate reinforced metal matrix composites. Hence, the important significance of introducing the varying weighting factors under different loading conditions into formulating multi-objective function is manifested in this aspect. When compared to the ad hoc strategy, the multi-objective optimization one enables the identified J-C model to fit the observations under high temperatures and high strain rates more accurately, in Figure 6b.

Materials and Methods
An accurate and reliable set of J-C material parameters under various combination of strain rates and temperatures could give rise to similar chip morphologies and cutting forces to experimental observations in cutting [28][29][30]. The workpiece parts in contact with the drill from the lip and axis centre tend to experience a wide range of strains, strain rates, and temperatures. So, the feasibility of constitutive model obtained using weighted multi-objective strategy is verified through the comparison of experimental and finite element results of drilling.

Materials
The J-C material model for Al6063/SiCp/65p composites is given in Table 1. For damage initiation and evolution, a combination of Johnson-Cook and shear failure criteria is used for representing ductile fracture caused by crack nucleation, growth, and voids coalescence in Al matrix and shear facture induced by local shear band in machining. A Johnson-Cook criterion for damage onset is modified for the formula put forward by Johnson [31], as , ̅ = + exp(− ) 1 + ln( * ) 1 + * where is the equivalent plastic strain at the onset of damage, and − are the material damage parameters, is the stress triaxiality with = / wherein is the hydrostatic pressure, and is the shear strength. Here the Johnson-Cook damage criterion is used to describe the ductile damage of Al6063 matrix tearing.
In terms of macro-and microscopic strain-stress relationships under homogeneous strain boundary conditions shown in Figure 7, the effective strain of representative volume element (RVE) for two-phase constituent composites can be defined by [32].
where 〈 〉 and 〈 〉 is the averaged strains of Al6063 matrix and SiC phase. 〈•〉 denotes the volume average of physical and mechanical properties. For notational simplicity, let = (1 − ).

Materials and Methods
An accurate and reliable set of J-C material parameters under various combination of strain rates and temperatures could give rise to similar chip morphologies and cutting forces to experimental observations in cutting [28][29][30]. The workpiece parts in contact with the drill from the lip and axis centre tend to experience a wide range of strains, strain rates, and temperatures. So, the feasibility of constitutive model obtained using weighted multi-objective strategy is verified through the comparison of experimental and finite element results of drilling.

Materials
The J-C material model for Al6063/SiC p /65p composites is given in Table 1. For damage initiation and evolution, a combination of Johnson-Cook and shear failure criteria is used for representing ductile fracture caused by crack nucleation, growth, and voids coalescence in Al matrix and shear facture induced by local shear band in machining. A Johnson-Cook criterion for damage onset is modified for the formula put forward by Johnson [31], as where ε D p is the equivalent plastic strain at the onset of damage, and d 1 − d 5 are the material damage parameters, η is the stress triaxiality with η = p/q wherein p is the hydrostatic pressure, and τ max is the shear strength. Here the Johnson-Cook damage criterion is used to describe the ductile damage of Al6063 matrix tearing.
In terms of macro-and microscopic strain-stress relationships under homogeneous strain boundary conditions shown in Figure 7, the effective strainε of representative volume element (RVE) for two-phase constituent composites can be defined by [32].
where ε Al and ε SiC is the averaged strains of Al6063 matrix and SiC phase. · denotes the volume average of physical and mechanical properties. For notational simplicity, let α t = (1 − V SiC ).
Materials 2018, 11, 97 11 of 22 Figure 7. Localized relationship between representative volume element (RVE) for multi-particle inclusion and Unit Cell for single-particle inclusion.
In macro scale, due to the polycrystalline aggregates characteristic and randomly oriented distribution of SiC particle in as-cast Al6063 composites, the mechanical responses of Al6063/SiCp/65p composites are roughly considered macroscopically isotropic. The plastic strain part on both sides for Equation (22) should be equal, and SiC particulate is elastically deformed, and therefore the failure criterion for Al6063 matrix composites can be approximately deduced using Equations (21) and (22) Equation (23) is normalized to obtain the failure parameters.
In macro scale, due to the polycrystalline aggregates characteristic and randomly oriented distribution of SiC particle in as-cast Al6063 composites, the mechanical responses of Al6063/SiC p /65p composites are roughly considered macroscopically isotropic. The plastic strain part on both sides for Equation (22) should be equal, and SiC particulate is elastically deformed, and therefore the failure criterion for Al6063 matrix composites can be approximately deduced using Equations (21) and (22) as Equation (23) is normalized to obtain the failure parameters.
The Johnson-Cook damage parameters d Al 1 − d Al 5 for 6063 aluminum alloy were given in [33]. The Johnson-Cook failure parameters for Al6063/SiC p /65p composites are deduced by Equation (25), as shown in Table 2. The shear criterion is a phenomenological one for describing shear band localization, with the following form [34] ε S p θ s , .ε p (26) and θ s = (q + k s p)/τ max (27) where ε S p is the equivalent plastic strain at damage initiation, q the effective Mises stress, and k s material parameter.
A scalar-valued damage factor ω is a measure of damage initiation. For J-C criterion in Equation (21), have For shear criterion in Equation (26), have When the damage factor ω in Equation (28) or Equation (29) reaches 1, the onset of damage is deemed to initiate.
Provided that material damage is initiated, the stress-strain law originally based on J-C material model no longer describes the material deformation behavior accurately [35]. To avoid energy-dissipated effect arising from refining mesh, a Hillerborg's fracture energy proposal is introduced to lower mesh dependency by employing a stress-displacement law after material damage occurs. The Hillerborg's concept is defined through fracture energy equation G f .
where the element characteristic length L is introduced to represent the effective plastic displacement u p after damage initiation, σ y is yield stress, ε d p is the plastic strain at the onset of damage, ε f p and u f p are the plastic strain and displacement at failure, respectively. The overall damage variable D is a comprehensive measure of all active damage criteria. The damage evolution law in a maximum sense is employed for predicting damage, which is calculated as the maximum of the individual damage variable d j (j = 1 or 2, corresponding to Johnson-Cook damage and shear damage).
A linear softening law of flow stress with equivalent plastic displacement u p is assumed as and in terms of Equation (30), have After damage initiation, the effective flow stressσ is given as Therefore, when the value of D reaches 1, the material failure occurs, accompanied by element deletion.
For Al6063/SiC p /65p composites involving α-SiC hard phase and Al6063 matrix, the specific heat capacity can be estimated according to the rule of mixtures by Equation (35).
where C SiC p and C Al6063 p are the specific heat capacity of α-SiC and Al6063. The physical and mechanical properties of Al6063/SiC p /65p composites for material modeling are listed in Table 3.

Finite Element Modelling of Drilling
A three-dimensional (3D) drilling model based on FE-based appoach has been built for simulating small-hole manufacturing process of Al6063/SiC p /65p composites. For workpiece modelling, the workpiece part neighbored to the tool tip are accounted for in drilling model, and a cone-like concave machined surface is prefabricated on the workpiece surface, such that stable drilling can be arrived at as soon as possible, as shown in Figure 8a. Modelling the workpiece in such a way makes it convenient to generate the structured mesh in uncut chip along the spiral cutting path so as to facilitate the formation of chip, and to enable the improvement of computational efficiency. Meanwhile, on account of the limited workpiece's model size and the possible action of the reflected stress waves from the workpiece's cylindrical boundary surface where more extra material exists beyond this boundary in the actual drilling infinite elements are incorporated into the boundary domain. To facilitate the implementation of infinite element modelling in which is limited to hexahedral meshes and sweep algorithm in ABAQUS, the workpiece is modelled as a cylinder with 4 mm in diameter and 1.5 mm in height. To equilibrate the computation cost and accuracy, the modeling for drill bits are only focused on the tool tip part involving the realistic cutting of the workpiece, and the sharp edges is taken into account. The PCD tool is deemed as a rigid body. The boundary conditions for the axial feed and rotation of tool are imposed on a reference point on the tool center axis, as shown in Figure 8b. The geometric simplification of drill bit and the assembly configuration of drilling model are illustrated in Figure 8b. An 8-node linear brick C3D8R with reduced integration and hourglass control is adopted for mesh division of drilling model, with the minimum mesh size of 10 µm. where is the critical shear yield stress, generally √3 times lower than the tensile yield stress. The cutting forces and chip morphology were obtained numerically by running FE model of drilling Al6063/SiCp/65p composites based on the constitutive formalism determined using multiobjective optimization strategy. The thrust force is the applied resultant force on cutting tool along the z direction.   It is worth noting that at least five-layer locally meshed elements along the feed direction are allocated for feed per revolution, only by such meshing strategy can the chip be formed, as shown  Figure 8c. This is mainly ascribed to material damage and element deletion involved in drilling simulation, as well as edge radius if the blunted cutting tool is to be accounted for.
The stick/slip contact model across tool-chip interface is adopted, as suggested by Zorev [36], for tool-chip contact definition. The contact formation for sticking/sliping along local tangent directions is defined as in terms of the division of slipping (τ ≤ τ crit ) and sticking regions, where τ crit is the critical shear yield stress, generally √ 3 times lower than the tensile yield stress. The cutting forces and chip morphology were obtained numerically by running FE model of drilling Al6063/SiC p /65p composites based on the constitutive formalism determined using multi-objective optimization strategy. The thrust force is the applied resultant force on cutting tool along the z direction.

Experimental Set-Up and Signal Processing
To validate the built constitutive equation and drilling model for Al6063/SiC p /65p composites, the drilling experiments of Al6063/SiC p /65p composites were performed on a CNC machining center DMU 80 monoBLOCK (EAMTM, Brussels, Belgium), as shown in Figure 9. The measurement of the thrust forces and torque during drilling was accomplished by a rotating 4-Component Dynamometer (RCD) Kistler 9123 (Kistler Instrument China Ltd., Shanghai, China). Due to semi-enclosed characteristics and the deformed complexity of drilling operation, the contact between the newly-formed chips by element deletion and machined surface in 3D drilling is difficult to define, so the chips flow away from the workpiece side wall, and do not impede the tool operation. This is in contradiction to the fact that in realistic drilling the newly-formed chips is inconvenient to remove, thus increasing cutting forces. To better approximate the simulation conditions, a peck drilling operation (with 0.05 mm step size) that facilitates the chip removal is employed to reduce the interference from the chips. Based on special consideration of drilling modelling, only the force signals at stable cutting are extracted for comparison. To eliminate the interference when the forcing frequencies approach natural frequencies of piezoelectric dynamometer (Figure 10a), the Type I low-pass Chebyshev filter was used for processing the experimentally acquired signals of thrust force and torque, with the cutoff frequency, as suggested by Szymon [37].
where the cut-off frequency is associated with the cutting tool kinematics represented by tooth passing frequency . The already filtered signals were dealt with shift compensation for avoiding zero drift. After signal processing, the acquired thrust force signals at rotational speed of 1500 rpm and feed of 50 mm/min are shown in Figure 10b. Due to semi-enclosed characteristics and the deformed complexity of drilling operation, the contact between the newly-formed chips by element deletion and machined surface in 3D drilling is difficult to define, so the chips flow away from the workpiece side wall, and do not impede the tool operation. This is in contradiction to the fact that in realistic drilling the newly-formed chips is inconvenient to remove, thus increasing cutting forces. To better approximate the simulation conditions, a peck drilling operation (with 0.05 mm step size) that facilitates the chip removal is employed to reduce the interference from the chips.
Based on special consideration of drilling modelling, only the force signals at stable cutting are extracted for comparison. To eliminate the interference when the forcing frequencies approach natural frequencies of piezoelectric dynamometer (Figure 10a), the Type I low-pass Chebyshev filter was used for processing the experimentally acquired signals of thrust force and torque, with the cut-off frequency, as suggested by Szymon [37]. (37) where the cut-off frequency f c is associated with the cutting tool kinematics represented by tooth passing frequency f zo . The already filtered signals were dealt with shift compensation for avoiding zero drift. After signal processing, the acquired thrust force signals at rotational speed of 1500 rpm and feed of 50 mm/min are shown in Figure 10b.

Validation of Thrust Force and Torque
Based on our drilling modelling strategy, the simulation results during stable cutting is take into account for comparison with the experimental data. When compared to the experimental observations in Figure 11, the maximum relative errors for the simulated thrust force and torque at rotational speed of 1500 rpm and feed of 50 mm/min are 11.00% and 22.27%, within the acceptable range. The errors are expected to decline as the hole diameter increases since the effects of tool geometric errors induced by sharpening and the dynamometer's measurement error on forces disturbance are lessened. The great fluctuation of simulated force illustrated by error bar may be partly due to uneven occurrence of element failure and partly due to the tool-chip contact instability

Validation of Thrust Force and Torque
Based on our drilling modelling strategy, the simulation results during stable cutting is take into account for comparison with the experimental data. When compared to the experimental observations in Figure 11, the maximum relative errors for the simulated thrust force and torque at rotational speed of 1500 rpm and feed of 50 mm/min are 11.00% and 22.27%, within the acceptable range. The errors are expected to decline as the hole diameter increases since the effects of tool geometric errors induced by sharpening and the dynamometer's measurement error on forces disturbance are lessened. The great fluctuation of simulated force illustrated by error bar may be partly due to uneven occurrence of element failure and partly due to the tool-chip contact instability caused by the discontinuous chips formation and rake angle's variation during cutting. The numerically and experimentally obtained thrust force and torque under other cutting conditions are compared in Figure 11. Figure 11a, which demonstrates the thrust force dependence on the drilling conditions, shows that there are the better machining conditions at rotational speed of 2000 rpm and feed of 75 mm/min. The developed FE model allows for presenting this dependence in the much more visual and usable form. The comparative results show both the simulation results of thrust force and torque are within the acceptable error ranges, and the constitutive formalism determined using multi-objective strategy can perform well in force prediction for small hole drilling.

Validation of Chip Morphology
The mechanical behaviors of the workpiece play a significant role in mechanics of chip formation and evolution where the workpiece material is subjected to plastic deformation and shear along the primary shear zone by dint of cutting tool [38]. The detailed chip formation and evolution mechanisms are investigated in Figure 12. Under the shearing and extruding actions from the lips and chisel edge accompanied by element failure, two approximatively centrosymmetric pieces of chip are gradually formed along with spiral downward cutting motion of drill bit. The chip is formed under the combined effects of the lips and chisel edges, and flows out along the rake face. The curling of chip is partly due to the velocity and deformation difference between the upper surface layer and the lower surface layer, and partly due to the compression from the rake face. Because the cutting lips of PCD tool is composed of two straight edges, uniform and straight morphologies of the chip root that is formed by the cutting lips are concluded to be triggered by shear failure. Less amount of chip formed by the chisel edge is mainly attributed to plastic extrusion of chisel edge with large negative rake angle, and resultant element failure is suspected of being triggered by Johnson-Cook damage criterion. Since the intersection between the cutting lips and chisel edge divides the chip into two segments, the segment formed by cutting lips is the main chips for comparison in Figure 13. The maximum equivalent plastic strain occurs right round the intersection of the margin and lip for drill bit. This is relevant to the highest cutting speed and low hydrostatic pressure here. Eventually, discontinuous chips under two damage criteria are formed, as is shown in Figure 13a. Details of the chip morphologies from simulation and experiment observation are captured in Figure 13b-e. From Figure 13b,c, the lamella morphologies from simulation and experiment observation are clearly observed in the free surface. Since the lamella morphologies result from the workpiece material' nonuniform deformation sheared by drill bit, the interval for such microstructure is in proportion to

Validation of Chip Morphology
The mechanical behaviors of the workpiece play a significant role in mechanics of chip formation and evolution where the workpiece material is subjected to plastic deformation and shear along the primary shear zone by dint of cutting tool [38]. The detailed chip formation and evolution mechanisms are investigated in Figure 12. Under the shearing and extruding actions from the lips and chisel edge accompanied by element failure, two approximatively centrosymmetric pieces of chip are gradually formed along with spiral downward cutting motion of drill bit. The chip is formed under the combined effects of the lips and chisel edges, and flows out along the rake face. The curling of chip is partly due to the velocity and deformation difference between the upper surface layer and the lower surface layer, and partly due to the compression from the rake face. Because the cutting lips of PCD tool is composed of two straight edges, uniform and straight morphologies of the chip root that is formed by the cutting lips are concluded to be triggered by shear failure. Less amount of chip formed by the chisel edge is mainly attributed to plastic extrusion of chisel edge with large negative rake angle, and resultant element failure is suspected of being triggered by Johnson-Cook damage criterion. Since the intersection between the cutting lips and chisel edge divides the chip into two segments, the segment formed by cutting lips is the main chips for comparison in Figure 13. The maximum equivalent plastic strain occurs right round the intersection of the margin and lip for drill bit. This is relevant to the highest cutting speed and low hydrostatic pressure here. Eventually, discontinuous chips under two damage criteria are formed, as is shown in Figure 13a. Details of the chip morphologies from simulation and experiment observation are captured in Figure 13b-e. From Figure 13b,c, the lamella morphologies from simulation and experiment observation are clearly observed in the free surface. Since the lamella morphologies result from the workpiece material' nonuniform deformation sheared by drill bit, the interval for such microstructure is in proportion to feed per revolution during drilling. In the chip back surface, the similar morphological characteristics of the simulated chips to realistic ones again verify the validity of the drilling model and the feasibility of the developed multi-objective constitutive identification strategy with varying weighted consideration to each combination of deformation rates and temperatures.

Conclusions
In this paper, the mechanical behaviors of Al6063/SiCp/65p composites in quasi-static compression and SHPB tests are investigated to obtain the stress-strain responses for identifying J-C constitutive constants. Aimed at hard identification of material parameters using conventional ad hoc method, a weighted multi-objective optimization strategy for identifying J-C constitutive parameters is proposed, coupled with updated Levenberg-Marquardt implementation algorithm. When compared to ad hoc optimization strategy, multi-objective strategy applied for J-C parameters optimization leads to a higher predictability of experimental observations and better capability of fitting experimental noise. It is found that it is not always necessary to modify such a classic J-C constitutive model under scientific hypothesis that parameter the coefficient A is a variable to be optimized rather than offset yield point.

Conclusions
In this paper, the mechanical behaviors of Al6063/SiCp/65p composites in quasi-static compression and SHPB tests are investigated to obtain the stress-strain responses for identifying J-C constitutive constants. Aimed at hard identification of material parameters using conventional ad hoc method, a weighted multi-objective optimization strategy for identifying J-C constitutive parameters is proposed, coupled with updated Levenberg-Marquardt implementation algorithm. When compared to ad hoc optimization strategy, multi-objective strategy applied for J-C parameters optimization leads to a higher predictability of experimental observations and better capability of fitting experimental noise. It is found that it is not always necessary to modify such a classic J-C constitutive model under scientific hypothesis that parameter the coefficient A is a variable to be

Conclusions
In this paper, the mechanical behaviors of Al6063/SiC p /65p composites in quasi-static compression and SHPB tests are investigated to obtain the stress-strain responses for identifying J-C constitutive constants. Aimed at hard identification of material parameters using conventional ad hoc method, a weighted multi-objective optimization strategy for identifying J-C constitutive parameters is proposed, coupled with updated Levenberg-Marquardt implementation algorithm. When compared to ad hoc optimization strategy, multi-objective strategy applied for J-C parameters optimization leads to a higher predictability of experimental observations and better capability of fitting experimental noise. It is found that it is not always necessary to modify such a classic J-C constitutive model under scientific hypothesis that parameter the coefficient A is a variable to be optimized rather than offset yield point.
A 3D FE simulation in small hole drilling of Al6063/SiC p /65p composites was implemented, based on the constitutive model determined using weighted multi-objective strategy and failure criteria in a maximum sense. The maximum prediction error of no more than 11.00% for thrust forces indicates the reliability and accuracy of the proposed drilling model. It can be found from numerical and experimental results that there are the better machining conditions at rotational speed of 2000 rpm and feed of 75 mm/min. The developed FE model allows for presenting the thrust force dependence on the drilling conditions in the much more visual and usable form. Accurate prediction of thrust force is of great importance for investigating drilling induced edge defects at hole exit and entrance in PRMMCs or delamination defects in laminated composites. Hence, this proposed FE model can be used to optimize cutting processes to improve drilling performance and hole quality, based on the cutting forces effect on drilling induced defect. Accurate prediction of chip formation contributes to the disclosure of material removal mechanism. Comparative results of thrust force, torque, and chip morphology also verify the feasibility of the developed multi-objective constitutive identification strategy with varying weighted consideration to each combination of deformation rates and temperatures.
However, while applying Brayden updating algorithm in the J-C model optimization, numerical unstability and divergence problems arise in the solving process. And the reason for nonconvergence, according to Broyden [32], is that in the 1st and 2M Para iterative updating of Jacobian matrix, a bad approximation results from χ 2 (P) > χ 2 (P + h). So in these iterations, the Brayden rank-1 updating algorithm is replaced by finite differences, and therefore the function evaluation is required to judge M Para or 2M Para [40].

A.3. Update of Step Length h
For the implementation of the algorithms above, the update of the step length h can be determined by comparing the quality of χ 2 (P + h) against χ 2 (P). The gain ratio Q(h), as a measure of suitable or unsuitable h, is a ratio between the actual and the expected improvement of an L-M update in the cost function [41]: for Equation (19), have The values for h i is identified by judging whether Q(h) > 3 or not, wherein 3 is a specified threshold that determines the step acceptance in L-M algorithm.

A.4. Convergence Criteria
The following convergence criteria are employed for iterative termination. 1 Gradient criterion max J T W σ Exp − σ Model < 1 (A7) 2 Step length increment criterion max|h i /P i | < 2 (A8) where 2 and 3 are the specified thresholds for iterative termination.

A.5. Error Criteria
The coefficient of determination, R 2 and overall fit standard error OFSE are used to measure the overall quality of fit and unbiased statistical error of model predictability and reliability, respectively [42].
where N total is the total experimental data point, k the number of material parameters to be fit.