A Numerical Investigation of a Single-Shot in a DEM-FEM Approach to Shot Peening Simulation

Shot peening (SP) is a controlled and systematic process of surface treatment that has a large number of controllable process parameters that make its application highly challenging. It involves the shooting of small and hard metallic balls at a targeted surface, with the aim of enhancing the fatigue strength of the workpiece under unfavorable service conditions. The compressive residual stress (CRS) induced by this application is expensive to evaluate experimentally. This paper presents a numerical model of the impact of a single-shot on a metallic surface, with the aim to set the stage for a realistic multiple shots peening simulation. The approach proposed herein is a sequential Discrete Element-Finite Element (DE-FE) coupled simulation, based on the use of different types of coefficients of restitution (CoRs) with emphasis on the energetic CoR. The energetic CoR relates the shot/target contact forces to the fractional strain energy needed for localized plastic deformation of the near-surface layer in the workpiece. The generated results of the induced compressive residual stresses (CRS) and equivalent plastic strain (PEEQ) from single-shot simulations are validated with similar results from the literature. Our study clarifies the strain energy aspects of a single-shot impact responsible for the desired effects of CRS and PEEQ, thereby laying the groundwork for accurate and realistic modeling of the SP process via the DEM-FEM approach.


Introduction
Shot peening (SP) is one of the oldest mechanical surface enhancements of components for fatigue applications. Though it involves projecting a stream of energetic shots unto the surface of the component, a single-shot application to shot peening analysis has been the norm [1][2][3][4][5][6][7][8][9] from the start of analytical and numerical approaches. Single-shot simulation plays a very important role in developing a model for the proposed Discrete Element Method -Finite Element Method (DEM-FEM) approach to a realistic simulation of the SP process, in the current research work.
The SP process involves the shooting of small and hard metallic balls at the targeted surface, with the aim of enhancing the fatigue strength of the workpiece under unfavorable service conditions. Therefore, a single shot which is just an impact between two bodies, cannot represent the actual multiple shots as in shot peening but serves the purpose of giving insight into the formation of Compressive Residual Stress (CRS) on the surface of a workpiece subjected to this treatment.
A typical residual stress profile due to shot peening, is presented in Figure 1, been characterized by a surface compressive stress (σ sur ), a near-surface maximum compressive stress (σ max ) at a beneficial beneficial depth (db), a total depth (dT) of compressive stress, and a tensile stress regime in the remaining section of the depth from the treated surface. Typical profiles of the evolving trend of the residual stress as influenced by the thickness of the workpiece, the hardness of the shot, and its velocity was demonstrated in the work of Al-Obaid [10], where a workpiece with a thin section experiences an increase (a kick) in the tensile stress just after the total depth of compressive stress (plastic layer) as in Figure 1. This could aid the external tensile load applied in service, thereby initiating failure of the workpiece from within. That work also indicated that harder shots induce higher compressive stress in the near-surface regime of the workpiece, as well as shots with higher velocity capable of inducing both higher and deeper compressive stress, for enhanced fatigue life.
Al-Obaid [10] obtained a theoretical relationship between the plastic zone and some parameters of the shot and workpiece, which is useful in predicting the residual stress distribution along with the depth of a shot-peened workpiece, without having to resort to experimental measurements. However, there is a crucial shortcoming of the theoretical approach, which produces high tensile stresses at the transition from compressive to tensile stresses within the workpiece, due to huge assumptions made.
James et al. [11] related the fatigue performance of a component to residual stresses arising from the manufacturing and fabrication processes. They posit that residual stresses must be measured over a reasonable depth of several millimeters below the surface, for fatigue life prediction purposes.
However, affordable residual stress measurement practices such as the laboratory X-ray diffraction technique and incremental hole-drilling technique are not generally useful in the light of the above discussed. They embarked on a deeper depth measurement campaign, using synchrotron and neutron diffraction techniques, which are very expensive but capable of assessing the depth of the components, revealing the entire plastic zone and the transition pattern of the residual stress.
The importance of the issues highlighted above, and the limitations associated with the available option of an expensive experiment measurement technique such as synchrotron diffraction method, calls for the development of numerical models that incorporate both residual stress evaluation and fatigue life assessment due to an applied external load. This paper is, therefore, the first step of that endeavor for a realistic shot peening simulation with industrial relevance while maintaining the plasticity concept of the process.
The relevance of the impact created with a single-shot can be seen in the work of Gangaraj et al. [12] in which the estimation of multiple impact coverage was achieved by means of a numerical single-shot indentation, which is also applied in the work of Miao et al. [13].
There are several reports in the literature of single-shot numerical investigations such as Meguid et.al. [14], where a single-shot model was utilized in the study of the effect of shot velocity, size, Typical profiles of the evolving trend of the residual stress as influenced by the thickness of the workpiece, the hardness of the shot, and its velocity was demonstrated in the work of Al-Obaid [10], where a workpiece with a thin section experiences an increase (a kick) in the tensile stress just after the total depth of compressive stress (plastic layer) as in Figure 1. This could aid the external tensile load applied in service, thereby initiating failure of the workpiece from within. That work also indicated that harder shots induce higher compressive stress in the near-surface regime of the workpiece, as well as shots with higher velocity capable of inducing both higher and deeper compressive stress, for enhanced fatigue life.
Al-Obaid [10] obtained a theoretical relationship between the plastic zone and some parameters of the shot and workpiece, which is useful in predicting the residual stress distribution along with the depth of a shot-peened workpiece, without having to resort to experimental measurements. However, there is a crucial shortcoming of the theoretical approach, which produces high tensile stresses at the transition from compressive to tensile stresses within the workpiece, due to huge assumptions made.
James et al. [11] related the fatigue performance of a component to residual stresses arising from the manufacturing and fabrication processes. They posit that residual stresses must be measured over a reasonable depth of several millimeters below the surface, for fatigue life prediction purposes.
However, affordable residual stress measurement practices such as the laboratory X-ray diffraction technique and incremental hole-drilling technique are not generally useful in the light of the above discussed. They embarked on a deeper depth measurement campaign, using synchrotron and neutron diffraction techniques, which are very expensive but capable of assessing the depth of the components, revealing the entire plastic zone and the transition pattern of the residual stress.
The importance of the issues highlighted above, and the limitations associated with the available option of an expensive experiment measurement technique such as synchrotron diffraction method, calls for the development of numerical models that incorporate both residual stress evaluation and fatigue life assessment due to an applied external load. This paper is, therefore, the first step of that endeavor for a realistic shot peening simulation with industrial relevance while maintaining the plasticity concept of the process.
The relevance of the impact created with a single-shot can be seen in the work of Gangaraj et al. [12] in which the estimation of multiple impact coverage was achieved by means of a numerical single-shot indentation, which is also applied in the work of Miao et al. [13].
There are several reports in the literature of single-shot numerical investigations such as Meguid et.al. [14], where a single-shot model was utilized in the study of the effect of shot velocity, size, shape, and the strain hardening coefficient on the plastic zone and residual stress field patterns created due to the indention of the metallic surface. In their work, it is indicated that the surface CRS and the maxima which occurs under the surface together with the compressed layer are influenced by the shot parameters such as velocity, shape, and size. Hong et al. [15] did a 3D finite element dynamic analysis of the single-shot impact on a component showing the effect of parameters such as impact velocity and angle, shot diameter, and initial yield stress of the workpiece on the residual stress distribution.
Al-Hassani [16,17] carried out both static and dynamic single-shot steel ball indentations, respectively, giving only the depth of the plastic zone as a response of the target impinged upon by a single sphere. The residual stress distribution along the depth and the axial direction from the point of impact on the surface were reported by the same author for a single-shot and a repeated single-shot impact at the same location on the surface. He was able to show the importance of strain rate in the distribution of stresses in the target, which is evident in the repeated impact.
With regard to the effects of peening parameters, Mann et al. [18] used a calibrated numerical single-shot via a nano-indentation experiment to study the near-surface hardness and residual stress, which resulted in the prediction of shot peening-induced residual stresses for an aluminum alloy. Hong et al. [15,19] used a single-shot impact on a metallic component to investigate the effects of the various peening parameters on the residual stress distribution. Among other findings, they reported that the shot diameter has a negligible effect on the magnitude of residual stress but that the depth increases linearly with increasing shot diameter. Hirai et al. [20] indicated numerically that the dent formed by a single-shot is influenced by the strain hardening parameter and the yield strength of both the shot and the peened material.
The importance of the materials involved in the SP process can be observed in the work of Sanjurjo et al. [8], who utilized a single shot experimentally and numerically in the analysis of commonly used constitutive material models and their response to monotonic loads. The study indicated that the Johnson-Cook material model gave a dimple diameter closer to that obtained experimentally for duplex stainless steel when compared with the Cowper-Symonds model and modified versions of the Johnson-Cook or Cowper-Symonds models.
All the simulation works highlighted so far are based on the finite element method, and none of these properly predicted the CRS profile obtained experimentally by Mori et al. [21] and Kobayashi et al. [22] for a single-shot. They investigated a single-shot impact on a metallic surface experimentally and reported a difference in the static and dynamic indentation with regard to the residual stress, which is tensile in a single-ball dynamic test at the center of the dent, while that of a static indentation is compressive or near zero at the center. However, as more dents are created close to an isolated one due to multiple shooting, the tensile residual stress at the center of a single dent becomes compressive. Shivpuri et al. [5] and Zion et al. [23] confirmed these findings with numerical simulations supported by experimental results.
Xiao et al. [24] used results from single and repeated impacts to establish a model to calculate the random peening stresses in peen forming operations. They showed that the depth of a repeated impact fluctuates, revealing the loading and rebounding of the impact regime. A relationship was eventually established between the maximum depth of each impact and the diameter of the dimple created.
It should be noted that even in cases where the reported residual stress profiles for a single-shot do not conform to experimental reports, they do, however, give insight into the relationship between the SP parameters and CRS. The use of FEM in simulating the complex SP process for a realistic outcome is currently not practicable due to the huge number of shots (over a million) involved. This is further compounded by the complex velocity distributions of the shots influenced by the nozzle design and settings, coupled with the shot-shot collisions, which reduces their initial velocity at impact.
Current computer resources required to use Finite Element Method (FEM) alone for realistic shot-peening investigations are still not feasible. This paved the way for the use of the discrete element method (DEM) more recently, either alone or coupled with FEM concurrently by Han et al. [25], Jebahi et al. [26], or sequentially by Murgaratnam et al. [27], Bhuvaraghan et al. [28], and Fubin-Tu et al. [29], to overcome the challenges associated with the FEM method.
In their study, Han et al. [25] adopted a soft model for the shot, since both shot and workpiece are physically deformable, instead of the assumed rigid shots model commonly used by many researchers in the literature. Their results conform to experimental profiles obtained for a single-shot, where tensile residual stress is obtained directly under a single-shot indentation. However, they referred to this occurrence as an abnormality associated with the use of a 2-D FE model. This is actually not the case, considering the work of Shivpuri et al. [5] and supported by that of Lewis et al. [23], where similar results were obtained with a 3D FE analysis.
Soady [30] highlighted the issue of the overestimation of CRS in a literature review drawing the attention of the SP research community to the experimental works of Wang et al. [31] and that of Gao et al. [32]. These authors used empirical formulations to relate CRS to the yield stress of the treated workpiece and the SP process parameters, showing that the CRS is a fraction of the yield stress or ultimate tensile stress of the workpiece and not a multiple of those values. The overestimation of CRS is considered to be a pitfall in the strain rate effect included in the material model mostly by the application of the Johnson-Cook constitutive model. Gangaraj et al. [12] did a numerical investigation on the effect of the mesh size on the residual stress convergence and obtained a realistic CRS profile for a single-shot simulation. They established that the ratio of mesh size to dimple size created by the shot should not be larger than 1:10 for the shot and 1:20 for the workpiece, especially in the impact zone.
In light of the above discussions, this paper is aimed at enhancing the DEM-FEM approach previously developed by Edward et al. [33] for a simplified numerical simulation of a realistic SP process. This is needed in order to set the stage for developing a numerical framework that captures all the relevant shot parameters that are of industrial relevance while maintaining the principles of impact mechanics in the SP process simulation.
To the authors' knowledge, the proposed DEM-FEM approach, which is sequential in nature, is one of its kind, where the direct application of imported peening forces from DEM is used in FEM for predicting the workpiece response to a realistic SP process. The work is patterned after the work of Bhuvaraghan et al. [28], where the pressure equivalent of the force was applied to the component using a fixed dent size of a single normal impact as the area of contact via an approximated mesh size. That approach, however, has considerable limitations with respect to a realistic SP process, as highlighted by Fubin-Tu et al. [29], as well as with respect to its accuracy, since dents formation in a realistic SP cannot be uniform. The radius of the dent formed depends on the angle of impact, which is seldom normal to the workpiece surface, more so considering the shot-shot collision effects, which reduce the likelihood of a normal impact.
The main objective of this paper is to establish the basis for a proposed DEM-FEM framework [33] needed to simulate a realistic SP process, by means of single-shot impact analysis from an applied mechanics perspective. This approach utilizes the plastic strain energy component of a single-shot impact, as opposed to the reactive contact force component utilized by previous work [28].
The reactive contact force is responsible for the rebound of the shot and cannot produce the same plastic strain resulting from the impact. In the present work, we convert the plastic strain energy of a single-shot impact into contact forces by use of the energetic CoR instead of the kinematic CoR of the collision. A DEM model is subsequently coupled with combined implicit/explicit FE models, which enhance the direct application of the exported force, thereby overcoming the limitations of the DEM-FEM approach of Bhuvaraghan et al. [28].
The remaining part of this paper is organized as follows: a theoretical background for the numerical formulation of the DEM-FEM model; simulation of a single-shot using the newly proposed model, results, and validation, as well as a discussion of the results and conclusions.

Theoretical Background
In formulating a DEM-FEM model for a realistic SP simulation, the authors consider the mechanics of the process with respect to the impact of a single-shot on the workpiece surface as shown in Figure 2. The kinetic energy of the shot at impact W i , and that at rebound W r , corresponds with the impact velocity v i and the rebound velocity v r at angles α and α' respectively. During impact, contact forces are generated briefly between the bodies at the region of overlap, with deceleration and acceleration of the shot accompanied by rapid energy dissipation. The impact is partially plastic involving energy loss of a fraction of the impact energy of the shot at collision, with plastic deformation of the workpiece. The shot can either be considered rigid, with no plastic deformation or deformable, and then it will experience flattening. shown in Figure 2. The kinetic energy of the shot at impact Wi, and that at rebound Wr, corresponds with the impact velocity vi and the rebound velocity vr at angles α and α' respectively. During impact, contact forces are generated briefly between the bodies at the region of overlap, with deceleration and acceleration of the shot accompanied by rapid energy dissipation. The impact is partially plastic involving energy loss of a fraction of the impact energy of the shot at collision, with plastic deformation of the workpiece. The shot can either be considered rigid, with no plastic deformation or deformable, and then it will experience flattening. Generally, from the first principles, energy balance can be represented as where Wi is the kinetic energy at impact, Wr is the restituted energy, Wp is the energy used for plastic deformation in the workpiece, and Wd is the energy dissipated through elastic stress wave propagation, heat due to friction etc., which is not recoverable. The total kinetic energy loss (WL) by the shot of mass m due to impact is, therefore: It should be noted that is identical for both DEM and FEM simulations for a singleshot, since the impact or initial velocity vi is the same in both cases. The total kinetic energy lost by the ball due to impact with the workpiece can be determined by the CoR obtained from an FE analysis of the impact. To date, researchers have used experiments and finite element analysis (FEA) to obtain the value of CoRs for specific materials. However, in order to focus on the main objective of this paper, we simply refer to a review work by Ahmad et al. [34] where different CoR models, with their pros and cons, were highlighted.
Ahmad et al. [34] highlighted three definitions of CoR, which are the kinematic CoR, kinetic CoR, and the energetic CoR, stating that only the energetic CoR can provide the correct solutions in all cases. This includes impact involving friction and a change in the direction of slip during the collision. For a specific material impacted by a shot, the CoR is calculated based on the normal component of the relative velocities, material properties and geometry of both bodies involved in the collision, as well as the sliding friction. From the review, it was found that an increase in the number of impacts at a specific spot, also increases the reactive contact force and CoR due to surface hardening by plastic deformation. This gives an indication that the maximum reactive contact force exacted by the target (workpiece) on the shot is responsible for the rebound phase of the shot.

CoR Calculations for Single-Shot Impact in DEM.
The coefficient of restitution (CoR) is simply the ratio of the shot velocity (vr) just after a rebound Generally, from the first principles, energy balance can be represented as where W i is the kinetic energy at impact, W r is the restituted energy, W p is the energy used for plastic deformation in the workpiece, and W d is the energy dissipated through elastic stress wave propagation, heat due to friction etc., which is not recoverable. The total kinetic energy loss (W L ) by the shot of mass m due to impact is, therefore: It should be noted that W i = 0.5mv 2 i is identical for both DEM and FEM simulations for a single-shot, since the impact or initial velocity v i is the same in both cases. The total kinetic energy lost by the ball due to impact with the workpiece can be determined by the CoR obtained from an FE analysis of the impact. To date, researchers have used experiments and finite element analysis (FEA) to obtain the value of CoRs for specific materials. However, in order to focus on the main objective of this paper, we simply refer to a review work by Ahmad et al. [34] where different CoR models, with their pros and cons, were highlighted.
Ahmad et al. [34] highlighted three definitions of CoR, which are the kinematic CoR, kinetic CoR, and the energetic CoR, stating that only the energetic CoR can provide the correct solutions in all cases. This includes impact involving friction and a change in the direction of slip during the collision. For a specific material impacted by a shot, the CoR is calculated based on the normal component of the relative velocities, material properties and geometry of both bodies involved in the collision, as well as the sliding friction. From the review, it was found that an increase in the number of impacts at a specific spot, also increases the reactive contact force and CoR due to surface hardening by plastic deformation. This gives an indication that the maximum reactive contact force exacted by the target (workpiece) on the shot is responsible for the rebound phase of the shot.

CoR Calculations for Single-Shot Impact in DEM
The coefficient of restitution (CoR) is simply the ratio of the shot velocity (v r ) just after a rebound to its velocity (v i ) before impact. In an FEM, simulation that includes the dominant physical effects (e.g., plastic dissipation, wave propagation, contact, and friction), the CoR can be computed and can subsequently be used to determine the energy loss. This classical definition of the CoR is known as the kinematic CoR. In a DEM simulation of numerous particles impacting on the workpiece, a CoR must also be defined. However, if the DEM simulation only simulates elastic contact (as is frequently the case), then the CoR obtained from the FEM simulation (e FEM ) would not be the appropriate value to use in the DEM simulation.
In the present work, the CoR to be used in DEM (e DEM ) is computed from the energetic CoR, denoted (E r ), and E r is defined as the ratio of the shot kinetic energy after impact to that before impact.
E r in Equation (3) can be used to obtain the contact force to be exported from the DEM model to the FEM model. The absolute value of (v r /v i ) is the kinematic CoR obtained from a single-shot simulation in FEM, and m is the mass of the shot. It should be noted that E r is the fraction of the impact energy of the shot restituted on its rebound due to the reactive contact force by the workpiece. Therefore: Equation (4) is a very important relationship between the energetic CoR (E r ) and the kinematic CoR of the impact either in DEM (e DEM ) or in FEM (e FEM ). The total energy loss experienced by the shot can be expressed as: It has also been reported that between 80 percent and 90 percent [35,36] of the kinetic energy loss by the shot to the workpiece during a mechanical impact is used in doing plastic work in the workpiece, assuming the shot to be rigid. The assumption of a rigid shot is appropriate for a shot peening process since the shot is much harder compared to the workpiece in most cases.
A detailed analysis of the energy balance for a single-shot impacting a workpiece is obtained in the works of Thornton [37] and Fathallah et al. [36]. Considering the mechanics of a single-shot impact, the kinetic energy of the shot just before impact constitutes its rebound kinetic energy and the total energy lost. Therefore, the energy relation in Equation (4) can be used to obtain the needed contact force to be exported from the DEM model to the FEM model.
A minimum of 80 percent of the total energy loss (W L ) for a single-shot in elastic-plastic impact is used in the plastic deformation of the workpiece. Therefore, equivalent energy via contact forces is required to be exported from DEM to FEM sequentially in order to maintain a realistic simulation of the SP process using the proposed DEM-FEM approach. An appropriate CoR to be implemented in DEM (e DEM ) has to be determined using Equation (7), with the minimum fraction (0.8) of the total energy loss (W L ) chosen instead of 0.9, for a conservative estimate. The use of 0.8 in Equation (7) is justified considering the fractional kinetic energy loss λ = 1 − k, where k is the efficiency coefficient related to the elastic and thermal dissipation in a purely elastic impact. This efficiency coefficient is from the findings of J. Reed, which were reported and referenced by Sherafatnia et al. [38]. Reed stated that k is given a constant value of 0.8 in most works reported in the literature, and also observed that a portion of the energy induced into the workpiece by the shot in purely elastic impact is dissipated as radiation of elastic waves and cannot be recovered. This fraction (λ) of the total kinetic energy loss, therefore, determines the maximum possible value of the CoR for any impact.
From numerical investigations carried out by the authors of the present work, it is evident that the recorded contact forces to be exported to the FEM domain from the DEM domain are responsible for the rebound of the shot during the offloading phase in DE analysis. It should be noted that these reactive contact forces in the DEM domain are not the same in magnitude as the contact forces in the FEM domain, even when the kinematics of the impact is maintained. It is instead required that the exported contact forces from the DEM model must produce the same plastic strain in the FEM model, considering a single-shot analysis with FEM only.
Though the impact velocities v i in both FEM and DEM are the same, the CoR to be implemented in DEM (e DEM ) has to be determined using Equations (4) and (8). Thus, equating the energy restituted in the DEM simulation for a single-shot with the same impact velocity as in FEM domain, to the energy lost due to plastic work in FEM, we have: From Equation (11), a FEM simulation of a single-shot that results in a CoR (e FEM ) of 0.4 will require a CoR (e DEM ) of 0.82 if a collision is implemented in EDEM (a discrete element software) using the Hertz-Mindlin (HM) no-slip contact model [39], in order to obtain the same contact force as in the FEM. This is due to the DEM analysis being purely elastic, but uses the kinematic CoR to mimic the elasto-plastic collision in the FEM. Thus, the DEM serves as a very potent tool in achieving a realistic multiple shots simulation, of a complex shot peening process when coupled with the FEM.

The Use of Contact Laws for a Single-Shot Impact in DEM
Considering a typical single-shot impact shown in Figure 2, for a DEM domain using the EDEM commercial software, the issue of concern in simulating the process in this research is the contact forces generated due to the interaction laws governing the impact. For a detailed analysis of existing contact models that predicts the contact force due to the impact of a particle with a flat surface (workpiece) please refer to the work of Han et al. [25].
The contact models available in EDEM, which are of interest to this work is the Hertz-Mindlin (no-slip) model and the hysteretic spring model (HS). The Hertz-Mindlin (HM) model can calculate contact forces efficiently, but the shot is considered a rigid body. The HS model could be used to accommodate a real plastic deformation effect in the particle. It should, however, be noted that only one of these two models could be implemented at a time when performing a simulation using the EDEM software.
The HM model uses the CoR as indicated in Equation (4) of Section 2.1, in evaluating the fraction of the impact energy, which is restituted, thus describing the impact kinetics. In contrast, the HS model uses the set CoR in a completely different way with the normal and tangential force calculations based on the Walton-Braun theory as detailed in the EDEM user guide [39]. While the CoR of the collision in the HM model takes the value implemented irrespective of the velocity of the shot, the CoR of the collision in HS model changes with the velocity of the shot, irrespective of the fixed CoR implemented.
It should be noted that the reactive contact forces obtained by the use of the HS model could be calibrated to obtain the same plastic strain in FEM. This is made possible by assigning predetermined yield stress for the deformable shot, which is evaluated numerically at the different shot velocities.
The application of the HM model by Bhuvaraghan et al. [28] showed a 27% difference between the maximum contact force obtained for the same single-shot simulation in DEM (87.8N) and FEM (111.6N). This could result in huge discrepancy in the simulated CRS when applied to multiple shots in the SP analysis. However, these authors did introduce a correction factor to calibrate the DEM force to the reactive contact force in FEM, via a look-up table that normalizes the values obtained. Notwithstanding, the force obtained by their approach could not produce the same plastic strain in the FEM, and the simulation time will be extended by the use of a look-up table.
Section 5. of the present work shows that the maximum contact force needed to be generated in the DEM that will produce the same plastic strain in FEM cannot be obtained with the HM model, for a predetermined shot velocity. This suggests a need to properly evaluate the CoR to be used in DEM via the HM contact model in a bid to determine the required force and justifies the preferred use of the HS contact model in the proposed DEM-FEM approach to a realistic SP simulation.
The HS model available in EDEM is calibrated numerically by several trials runs at the different velocities, to generate the required force that will produce the same plastic strain in FEM. A calibrated, yield stress of 2.3 GPa, 2.785 GPa, and 3.2 GPa is obtained for the 0.3556 mm, 0.7 mm, and 1.0 mm shots, respectively.

DE and FE Material Models Characterization
In line with the purpose of this study, the mechanical properties of the shots and workpieces shown in Table 1, were adopted from the literature for all the simulations, both in DEM and FEM. 1240.0 1864 n/a n/a A Johnsons-Cook (JC) constitutive model [41] is selected to capture isotropic hardening and strain rate hardening effects for the workpiece material with thermal softening. The JC model expressed as: is available in the Autodyn Explicit dynamics code in ANSYS. In Equation (12), σ y is the equivalent yield stress, A is initial yield stress, B is the hardening constant, C is the strain rate constant, n is the hardening exponent, and m is the thermal softening exponent.
is a dimensionless strain rate, with · ε 0 a reference strain rate at which A was obtained. is the homologous temperature, T R is the reference temperature, T M is the melting point of the material, and T is the effective temperature. Table 2 summarizes the Johnson-Cook material parameters obtained experimentally and reported in the literature by various authors. These parameters were used together with the corresponding mechanical properties in Table 1. The strain rate effect was not considered for the shot since the focus here is on the response of the workpiece to the impact. More so, the shot in this analysis is considered to be rigid and does not result in any significant difference with strain rate considered. Thus, a purely elastic material model, in accordance with Table 1, was used for the shot throughout this work in the FEM impact analysis. Cast steel shots of 0.7 mm and 1.0 mm sizes are used for the peening of the workpieces made of AISI 4340 material in this work, while the cut-wire shot (CW-14) of size 0.3556 mm is used for the peening of the workpiece made of Inco718 material.

A Process Map for DEM-FEM Coupled Analysis
In the present work, the contact force evaluated from a single-shot impact in the DEM domain using the EDEM commercial software is exported to the FEM domain using ANSYS for the purpose of pre-stressing the same workpiece in Explicit Dynamics (ANSYS Autodyne software) via a transient structural scheme. Figure 3 shows a process flow map of the proposed DEM-FEM analysis considering a single shot.

A Process Map for DEM-FEM Coupled Analysis.
In the present work, the contact force evaluated from a single-shot impact in the DEM domain using the EDEM commercial software is exported to the FEM domain using ANSYS for the purpose of pre-stressing the same workpiece in Explicit Dynamics (ANSYS Autodyne software) via a transient structural scheme. Figure 3 shows a process flow map of the proposed DEM-FEM analysis considering a single shot.
A single-shot impact in EDEM using the calculated CoR obtained from the same single-shot analysis in FEM is carried out to obtain the required contact force. The contact force is subsequently exported to the ANSYS platform as a force load, for CRS and PEEQ evaluation. Equivalent (von Mises) stresses initiated by the imported force load via the ANSYS implicit solver (transient structural) are then applied to the same workpiece in ANSYS Autodyne for explicit analysis. The final results that consist of dent size, CRS, and the Equivalent Plastic Strain (PEEQ) are obtained in Step 5 of the process, but only the CRS and PEEQ are validated with available experimental and numerical results from the literature.

FEM Dynamic Simulation for CoR Calculation
For brevity, an analysis of a single shot simulation of 0.3556 mm CW-14 shot on the surface of an Inco718 workpiece is presented here. However, the same approach applies to other shot types and A single-shot impact in EDEM using the calculated CoR obtained from the same single-shot analysis in FEM is carried out to obtain the required contact force. The contact force is subsequently exported to the ANSYS platform as a force load, for CRS and PEEQ evaluation. Equivalent (von Mises) stresses initiated by the imported force load via the ANSYS implicit solver (transient structural) are then applied to the same workpiece in ANSYS Autodyne for explicit analysis. The final results that consist of dent size, CRS, and the Equivalent Plastic Strain (PEEQ) are obtained in Step 5 of the process, but only the CRS and PEEQ are validated with available experimental and numerical results from the literature.

FEM Dynamic Simulation for CoR Calculation
For brevity, an analysis of a single shot simulation of 0.3556 mm CW-14 shot on the surface of an Inco718 workpiece is presented here. However, the same approach applies to other shot types and targets. A unit cell is modeled using a 0.0072 mm mesh size at the contact surfaces of both shot and workpiece, and more coarse meshes up to 0.03 mm size in other parts of the workpiece. The shot is discretized with 231,065 463,704 elements of 4 noded tetrahedral meshes (element type Tet4) while the workpiece has 463,704 elements of 8 noded hexahedral meshes (element type Hex8), with edge sizing of 20 divisions and a bias factor of 1.1 increments from the surface to the bottom along the four vertical edges. Symmetric boundary conditions are applied to the four vertical faces with the bottom face completely fixed in Figure 4. However, for the bigger particles (0.7 mm and 1.0 mm diameter shots), a workpiece of size 3.6 × 3.6 × 1.29 mm is modelled with a mesh size of 0.017 mm and 0.02 mm, respectively, at the contact surfaces.  Figure 4. However, for the bigger particles (0.7 mm and 1.0 mm diameter shots), a workpiece of size 3.6 × 3.6 × 1.29 mm is modelled with a mesh size of 0.017 mm and 0.02 mm, respectively, at the contact surfaces.
Step 1 of the flow process shown in Figure 3 was executed with a rigid shot of 0.3556 mm, targeted at the midpoint of the top face of workpiece, in the normal direction, as shown in Figure 4. The single shot explicit dynamics run in Ansys mechanical took a total of 218.623 min on four slaves, in an Intel i7 CPU@2.80 GHz of 16 GB Random-access memory (RAM) Notebook. The CoR obtained is used in calculating the required CoR to be implemented in Step 2 of the process. Therefore, the CoR to be implemented in DEM must be properly determined with an accurate FE model, so as to obtain the required contact force that will be exported to FEM, for structural analysis. It should be noted that in determining the needed force at this phase of the schematics in Figure 3, the same single shot simulations are carried out in both the FEM (step1) and DEM (step2) without coupling. This resulted in the data presented in Tables 3 and 4.

Numerical Verification of FE Model
A numerical investigation was carried out to verify the accuracy of the FEM model used in Step 1 to obtain the CoR, to be implemented in the DEM in Step 2.
An elastic impact simulation was first conducted in FEM, where the workpiece was assigned purely elastic properties. This simulation quantifies the fraction of impact energy lost without any plastic deformation of the workpiece. Then the same simulation was repeated but with an elastoplastic workpiece, so as to obtain the fraction of impact energy lost when plasticity is considered. The fraction of impact kinetic energy lost due to induced plastic strain (Ep) could be obtained by subtracting the fraction of energy restituted Erp in the simulation with elastoplastic workpiece from the fraction restituted Ere in a similar simulation with an elastic workpiece. This is expressed as: Step 1 of the flow process shown in Figure 3 was executed with a rigid shot of 0.3556 mm, targeted at the midpoint of the top face of workpiece, in the normal direction, as shown in Figure 4. The single shot explicit dynamics run in Ansys mechanical took a total of 218.623 min on four slaves, in an Intel i7 CPU@2.80 GHz of 16 GB Random-access memory (RAM) Notebook. The CoR obtained is used in calculating the required CoR to be implemented in Step 2 of the process. Therefore, the CoR to be implemented in DEM must be properly determined with an accurate FE model, so as to obtain the required contact force that will be exported to FEM, for structural analysis. It should be noted that in determining the needed force at this phase of the schematics in Figure 3, the same single shot simulations are carried out in both the FEM (step1) and DEM (step2) without coupling. This resulted in the data presented in Tables 3 and 4.

Numerical Verification of FE Model
A numerical investigation was carried out to verify the accuracy of the FEM model used in Step 1 to obtain the CoR, to be implemented in the DEM in Step 2.
An elastic impact simulation was first conducted in FEM, where the workpiece was assigned purely elastic properties. This simulation quantifies the fraction of impact energy lost without any plastic deformation of the workpiece. Then the same simulation was repeated but with an elastoplastic workpiece, so as to obtain the fraction of impact energy lost when plasticity is considered. The fraction of impact kinetic energy lost due to induced plastic strain (E p ) could be obtained by subtracting the fraction of energy restituted E rp in the simulation with elastoplastic workpiece from the fraction restituted E re in a similar simulation with an elastic workpiece. This is expressed as:

Discrete Element-Finite Element Coupled Simulation of Single-Shot
A single-shot simulation was carried out in Step 2 in the flowchart of Figure 3, with the same parameters as in Step 1, making use of the FEM calculated CoR, and facet sizes of 0.02 mm for all simulations in DEM. Facet sizes of 0.04 mm, 0.06 mm, and 0.08 mm are employed in the FEM at steps 4 and 5, with a unit cell target of thickness 1.29 mm, for imported force load of shots 0.3556 mm, 0.7 mm, and 1.0 mm diameter, respectively.
The respective mesh sizes chosen for the FEM analysis in the DEM-FEM combination is obtained by means of a single-shot impact which gives the same dent profile when carried out with FEM only. Contact forces were analyzed and exported per time step from EDEM to ANSYS for implicit analysis (Step 4). The equivalent stress trace is initialized in Step 5 for the pre-stressing of the workpiece, and finally, an explicit analysis is carried out to obtain the workpiece response to the imported force load. The total run time for the multiple implicit/explicit FEM analysis when coupled to DEM is about 10 min for a single-shot force load. Table 3 displays the results obtained from a single-shot impact at different velocities with a workpiece assumed to be purely elastic. This results in different kinematic CoR (e re ), from which the fractional energy restituted (E re ) are obtained by making use of Equation (4). The maximum contact force in FEM (F FEM ) and the corresponding contact force in DEM for the same e re simulation produces a force correction factor at each impact velocity. Table 4 is a replica of Table 3, but for an elastoplastic workpiece material. Tables 3 and 4 resulted in Table 5, where the fractional plastic energy is obtained, and the equivalent CoR (e r ) is calculated.  Tables 3 and 4 for E p and equivalent contact force (F FEM ). From the results displayed in Tables 3-5, the values of fractional plastic energy obtained for single-shot simulations with ANSYS/Explicit dynamics at various shot velocities, were used to verify the accuracy of the model. As a check for the model accuracy, the expression:

Numerical Results and Discussion
E p /E L is satisfied for each simulation at the various shot velocities in Table 5, except for velocities 75 m/s and 100 m/s, where the deviations from the lower limit are within acceptable value of 5 percent reported and referenced by Sherafatnia et al. [38], and it has been highlighted in previous sections of this paper that 80 percent to 90 percent of energy lost by the shot is used for plastic deformation. Figure 5 is generated from these tables, with an indication of a reduction in the kinematic CoR as the impact velocity increases. However, the calculated CoR (e rc ) does not follow the same trend as the kinematic CoR, but increases with velocity. This indicates that the kinematic CoR does not represent the increase in energy transmitted to the component at impact, with increase in the shot velocity. From the results displayed in Tables 3-5, the values of fractional plastic energy obtained for single-shot simulations with ANSYS/Explicit dynamics at various shot velocities, were used to verify the accuracy of the model. As a check for the model accuracy, the expression: Ep/EL is satisfied for each simulation at the various shot velocities in Table 5, except for velocities 75 m/s and 100 m/s, where the deviations from the lower limit are within acceptable value of 5 percent reported and referenced by Sherafatnia et al. [38], and it has been highlighted in previous sections of this paper that 80 percent to 90 percent of energy lost by the shot is used for plastic deformation. Figure 5 is generated from these tables, with an indication of a reduction in the kinematic CoR as the impact velocity increases. However, the calculated CoR (erc) does not follow the same trend as the kinematic CoR, but increases with velocity. This indicates that the kinematic CoR does not represent the increase in energy transmitted to the component at impact, with increase in the shot velocity.  Figure 6 shows the relationship between the contact forces in FEM and the equivalent force from DEM needed to obtain the same plastic deformation. It is observed that the difference in these two plotted quantities increases with the velocity of the shot. This will enhance the ability of shots with higher velocities to induced reasonable plastic strain, which is responsible for residual CRS. This could be the reason why former work by past researchers [28] could not be validated with available experimental results, as observed by [29], even at a velocity of 100 m/s, which is on a high side, for conventional shot peening.  Figure 6 shows the relationship between the contact forces in FEM and the equivalent force from DEM needed to obtain the same plastic deformation. It is observed that the difference in these Metals 2019, 9,1183 13 of 17 two plotted quantities increases with the velocity of the shot. This will enhance the ability of shots with higher velocities to induced reasonable plastic strain, which is responsible for residual CRS. This could be the reason why former work by past researchers [28] could not be validated with available experimental results, as observed by [29], even at a velocity of 100 m/s, which is on a high side, for conventional shot peening. It should be noted that this approach, which involves multiple system (implicit/transient structural and explicit/Autodyn) analysis in ANSYS, is a unique technique of indirectly applying a force load in ANSYS (explicit dynamics), which gives the CRS and PEEQ shown in

Residual stress [MPa]
Residual stress [MPa] It should be noted that this approach, which involves multiple system (implicit/transient structural and explicit/Autodyn) analysis in ANSYS, is a unique technique of indirectly applying a force load in ANSYS (explicit dynamics), which gives the CRS and PEEQ shown in It should be noted that this approach, which involves multiple system (implicit/transient structural and explicit/Autodyn) analysis in ANSYS, is a unique technique of indirectly applying a force load in ANSYS (explicit dynamics), which gives the CRS and PEEQ shown in Figures 7-10. The results obtained are validated with results from the literature.      The resulting CRS and PEEQ profiles in Figures 7 -10 obtained with the present approach are in very good agreement with the numerical results for same single-shot simulations from recent works reported in the literature. The surface CRS, the maximum CRS, and its beneficial depth underneath the surface, the total depth of CRS, are in agreement with results from the literature. These numerical