Strain-rate-dependent model for the elastoplastic dynamic contact of sphere and plate

A force-indentation contact model is presented for the dynamic contact loading of elastoplastic particle and plate to incorporate the material’s strain-rate-dependent plasticity, built theoretically from the well-known Hertz contact law and Hill’s solution for elastic and elastoplastic quasi-static contacts. A theoretical relation of the relative impact velocity and plastic strain rate is introduced to solve the model’s parameters. A Johnson–Cook strain rate dependence is included into the model to consider dynamic effects. We validate the model using finite element analysis and show that the model can accurately simulate the force-indentation relation. The impact responses of plate simulated by applying the model combined with a substructure technique are validated using finite element analysis and laboratory test. With the aid of the model, a significant decrease in contact pressure during fully plastic indentation and the independence of dynamic contact-loading path upon loading rate are observed.


Introduction
The study of the dynamic contact law between sphere and plate is a fundamental problem in mechanics, involved in several physical phenomena and engineering applications, such as impact damage of spacecraft and satellite panels [1,2], ship collisions [3], lubrication [4] and shot peening [5]. Hertz solved analytically the forceindentation relation that is well understood for the static contact between elastic spheres [6]. For the application of contact law to the sphere-plate contact, the indentation is replaced by the relative displacement between the sphere and plate [7,8].
Earlier work on contact law for static elastoplastic contact focused on modeling the force-deformation relation between an elastoplastic sphere and a rigid plane, or between a rigid sphere and an elastoplastic half space [9][10][11][12][13]. The models have been used to inform experiments in a limited range of material properties, geometries, or loading cases [14][15][16]. Static contact laws can be applied to solve impact-contact problems by assuming quasi-static contact process [17][18][19]. However, due to some of their restrictive assumptions, the elastoplastic contact models are unsuitable for dynamic contact between materials with strain rate dependence, for which few models have been presented. A strain-rate-dependent model has been built for dynamic compression of elastoplastic spheres [20], and the parameters of an empirical force-indentation relation were determined by finite element (FE) analysis. Alfredsson and Nordin [21] extended an elastoplastic model by fitting model parameters from FE data, to predict the impact of shot-peening materials with strain rate dependence. A theoretical model with less empirical parameters is still required to predict more general characteristics of dynamic contact between strain rate dependent materials.
The objective of this study is to propose a theoretical model for the dynamic contact between elastoplastic sphere and plate. The model derives the quasi-static contact law from Hertz contact theory and Hill's solution [22]. The strain rate dependence is incorporated by introducing dynamic yield stress into the model. We validate Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. the model with FE analysis. Some contact phenomena are predicted and explained by the model. The dynamic responses of plate to low and moderate velocity impacts are simulated by the present model with a substructure technique, and are validated using the FE analysis and laboratory test.

Model description
A dynamic contact is considered to be quasi-static if the stress wave propagation has more than one reflection within the contact duration [23]. A quasi-static elastoplastic contact of sphere and plate has four distinct stages: (I) elastic, (II) elastic-plastic and (III) fully plastic indentations, and (IV) unloading stages [9,24]. In stage I, the contact is elastic and follows the Hertz contact law [6]. Stage II starts at the inception of contact yield with an indentation d Y [11]. The spherical plastic region [25,26] is encased within the surrounding elastic materials, with small observable permanent indentation [11]. The evolution of plastic region has been solved analytically by Hill [11,22]. Stage III starts at the inception of fully plastic yield with an indentation d p [9]. The plastic deformation is no longer contained beneath the contact surface. The contact pressure over the contact surface is approximately uniform, and the mean contact pressure is well known to be constant for static contact [27]. Stage IV is considered to be elastic [11], and obeys the Hertz elastic contact law [28,29].
In Stage I, the contact pressure is distributed parabolically over the contact area with contact radius a. The contact force P is estimated by indentation d according to the Hertz contact law [6]: The effective Young's modulus is The elastic indentation ends at ( ) when the contact plasticity first occurs according to the von Mises criterion. d Y is given by and s Y is the static yield strength of the softer material [30].
In Stage II, the contact pressure distribution flattens gradually as the plastic region evolves. The material's strain rate dependence begins to influence the contact behavior by raising the material yield stress. The material dynamic yield stress s dY is updated instantaneously with plastic strain rate. To derive the force-indentation relation in Stage II, Hill's solution [22] for the expansion of the elastic-perfectly plastic spherical shell is applied by replacing s Y with s . dY Similar to Stronge's derivations [9], the contact force is derived as the following:

* *
/ are the dynamic yield indentation and dynamic yield contact force, respectively. In the derivation of equation (2), the Stronge's geometry relation [9] without considering the pile-up deformations is applied, and is re-written as d d = + a a 0.5 1 dY dY * respectively. The application of Hill's solution implies a condition that the static yield stress uniformly distributes throughout the plastic region. However, s dY distributes non-uniformly, due to the complex spatial distribution of the strain rate. Therefore, we assume that s dY is the average value of the plastic region at an instant indentation d. Unlike those parameters for static contact, the parameters d , dY a dY and P dY change with the evolving of plastic region.
Stage III begins at d d = , p and the contact pressure is considered to distribute uniformly. The ratio of average contact pressure p to s dY is considered to be the constant 2.8 with increase of d, as selected for static contact [9,31]. The constant measured in [32] is 2.87. d p is selected to be d 84 Y as selected by Stronge [9] for static contact, and is considered to be independent upon the loading rate. Other different selections of d p can be found in [33][34][35]. Then, the contact force is derived as the following:  The unloading stage (Stage IV) is considered to be purely elastic, and follows the Hertz elastic contact law. Here, we take the formulation presented by Jin et al [24]. The contact force is In equation (4), d r is the permanent indentation, and R e * is the effective radius of unloaded curvature. By keeping continuity in contact force at the transition between the loading and unloading stages [36], where d u and P u are the indentation and contact force just before unloading.
Different types of strain-rate dependence can be used [21,37,38] to determine s . dY Here, the Johnson-Cook type strain-rate dependence is applied,  is the equivalent plastic strain rate, e 0  =10 −3 s −1 is the quasi-static strain rate at which the static yield stress s Y is measured, and C is an empirical parameter determined experimentally [20,39]. The Johnson-Cook model can capture the behavior of moderate-speed sphere impact experiments [20,40].
In the present contact law, only s dY remains undetermined. To determine s d in equation (5), e p  must be defined for sphere-plate contact. Some definition of strain rate has been made for sphere-sphere contact [20]. Since all knowledge of local strain rates is lost in the analytical description, the average plastic strain rate e p  in local contact plastic region is defined approximately based on a uniaxial compression model (see figure 1), instead of using the real equivalent plastic strain rate in equation (5).
Considering a rigid sphere compressed into a thin layer supported by a rigid flat and neglecting shear effect between adjacent elastic-plastic elements, the vertical displacement is Then, the average plastic strain rate e p  is defined as Hence, where v s is the velocity of the impactor centroid, v 1 is the velocity of impacted point at the plate bottom, and h is the plate thickness. For the contact between the two elastoplastic materials, the average plastic strain rate is calculated by replacing the indentation of the rigid sphere with the total indentation of the two materials. The expression of e , p  equation (7), for the sphere-plate impact is different from the representative strain rate of the sphere-sphere impact [20] and sphere-plane impact [21]. As e p  is obtained, s dY in equation (5) is modified at the present step by the rate at which the sphere and plate are moving together. Then, the contact model parameters are updated as the simulation progresses. Thus, the contact force at a given indentation can be obtained by the present contact model. For an impact event, an iterative equilibrium can be achieved at each time step, and the model parameters can be identified step by step.
For impact velocity much smaller than wave speed, according to the examination of energy balance during the elastic [41,42] and elastoplastic [43][44][45] impacts between sphere and massive target, the energy dissipated by elastic waves is only a few per cent, carried mainly by the Rayleigh wave. At least 90% of the initial kinetic energy is dissipated in plastic deformation [44,45]. Wu et al [23] found that the energy dissipation due to stress wave propagation is less than 3% of the total impact energy, if there is more than one reflection at boundary during the contact. Hence, the dynamic contact of plate can satisfy the quasi-static condition more easily.
The selection of invariant d p to plastic strain rate avoids the rollback of indentation stage.d p is a geometric parameter that represents the inception of plastic region on the contact surface.

Simulations
To study the dependences of force-indentation response and wave propagation on material strain rate sensitivity, an FE model is developed in LS-DYNA. The materials are modelled as elastic-perfectly plastic with the Johnson-Cook strain-rate dependence. The von Mises yield criterion is used for the yield condition. An isotropic rectangular plate impacted eccentrically by a sphere without friction is simulated for different impact velocities v 0 up to 20 m s −1 . The plate is simply supported along two opposite edges. The material for sphere and plate is SUJ2 steel [37,46] with remarkable strain rate sensitivity. = C 0.0884 is determined from the experimental data [46] and e = -10 0 3  s −1 . The geometrical sizes and material properties are listed in table 1. The eccentric impact point is located at the position that is 40 mm away from the plate center in the length direction (X) and 90 mm away from the plate center in the width direction (Y). The output interval is 5×10 −7 s in the simulation. The sphere and plate are meshed by 146 000 and 235 012 Solid 164 elements, respectively, with nodes spaced at 0.1 mm in contact area. The meshed model for highly accurate simulations was explained in detail in [24].
Due to the wave propagation and high deformation concentration in local contact region, a fine mesh and a small time step are required for the accurate FE analysis. Thus, the FE analysis becomes computationally expensive. To extend the simulation of sphere-plate impact dynamics, the dynamic substructure method (DSM) [24] is used and implemented in Visual C++ 6.0. DSM was derived from the plate finite element theory by a combination of static condensation and mode superposition, to reduce the computational degrees of freedom. The plate is considered to be globally elastic, and connects the impactor by force-indentation relation governed by an elastic-plastic contact model, which avoids the dense meshes in contact region. The differential equations couple the plate dynamics based on the Kirchhoff plate theory and sphere motion according to the Newton's second law, which are solved by the Newmark timeintegration algorithm. DSM captures the local deformation relying on the contact model for computational efficiency, sacrificing the knowledge regarding the intricacies of interaction. Here, DSM incorporating the present strain-rate-dependent contact model is called DSM-VP. DSM in [24] incorporating a rate-independent contact model is called DSM-EP, and is applied for comparison in this study.

Experiment
To validate the strain-rate-dependent contact model experimentally, an experimental setup was designed to measure the temporal signals of displacement, stress and contact duration for sphere-plate impact. A schematic diagram of the experimental setup is shown in figure 2(a). The eccentric impacts of a narrow plate and a hard sphere were tested at four positions with =   Figure 4 shows the comparison of the force-indentation curves obtained by the FEM, DSM-VP and DSM-EP, respectively. The geometry and material properties are the same for all simulation methods. In figure 4, the maximum indentations for seven impact velocities range from Stage II to Stage III. To obtain the DSM-VP results, equations (1)-(4) provide the explicit relation of the contact force and material dynamic yield stress. The model parameters are updated step by step using equations (5) and (7) based on the sphere-plate approach velocity. As the relation of contact force and velocity is implicit and nonlinear, the sphere-plate approach velocity is solved iteratively according to the solution procedure given in [24]. There is no modification for DSM-EP including the rate-independent contact model [24]. The DSM-EP results are used as reference to evaluate the degree of strain rate dependence. Figure 4 shows that the force-indentation relation exhibits a large strain rate dependence, even for low velocity impact. The material strain rate dependence significantly influences the contact law, and causes an obvious increase of contact stiffness and contact force, and a decrease of indentation and permanent indentation. The curves obtained by FEM and DSM-VP are in good agreement. By considering the dependence of contact law on the   relative approach velocity between the particle and plate, it is seen that DSM-VP can effectively simulate the local contact behavior of strain-rate-dependent material during impacts.
Observing the curves obtained by FEM and DSM-VP in figure 4, it is seen that there is a common forceindentation relation (loading path) in loading stage obeyed by different loading rates (or different impact velocities), similar to that for static contact predicted by DSM-EP. The same phenomenon is also observed in the experimental result [40]. A common loading path also exists among the contact radius-indentation curves in figure 5(a). It seems unusual for the dynamic contact of materials with strain rate dependence, because the plastic strain rate varies significantly from impact velocities, as seen in figure 5(b). In figure 5(b), the FEM curves are calculated according to equations (5) and (7), where the approach velocities are extracted from the database of LS-DYNA.
The independence of loading path on loading rate is interpreted by examining the relation of P and e . It can be concluded that the contact loading law is dependent upon the loading rate for a given impact system. The contact stiffness of the loading law might be dependent upon the dynamic constitutive behavior, such as Johnson-Cook material model, which is dependent upon the parameter C. Therefore, there might be a direct relation between the contact stiffness and parameter C. It suggests that the contact stiffness can be used as a measurement for the rate-dependent parameter of a material model, and the result of force-indentation is insensitive to the selection of loading rates. Figure 7 shows the FE simulation results of contact pressure at impact velocities of 3 m s −1 , 10 m s −1 and 20 m s −1 , respectively. The contact pressure distributions along the contact radius at different indentations are plotted. The contact pressure is found to be independent upon the impact velocity. The contact pressure decreases generally with indentation (or contact radius) in Stage III. Comparing the pressures at indentations of d 3.92 p and d 2.2 , p about 15.2% decrease is measured. It is quite different from the results of the static contact [48], for which the average contact pressure p remains approximately constant in Stage III.
The theoretical basis of the present model offers an interpreting of the occurrence of decreasing contact pressure in Stage III. In figure 5(a), DSM-VP predicts that e p  has a significant decreasing phase during contact loading process, resulting in a simultaneous decrease of s dY with indentation. The present contact model assumes that the ratio of p to s dY is constant in Stage III. To keep a constant value of s p , dȲ / it is expected that p should decrease obviously with the decrease of s . dY In addition, if p is defined as dynamic hardness [28], an obvious measurement error might happen due to the significant variation of contact pressure in Stage III.
To validate the present model further,   Although the contour patterns are complicated, the two contours match perfectly. At = t 170 μs, the stress disturbance initiated from the impact zone has propagated to the entire plate. The high stress is still concentrated in the local contact region around the impact position. Outside the local contact zone, the maximum von Mises stress is 181.4 MPa, which is below the static yield stress, 345 MPa. The global response of the plate behaves elastic. By including material strain rate dependence in the contact model, it can be seen that DSM-VP can effectively simulate the impact response of plate.
The above validation is for the impact of the sphere on the plate with the same materials, and the impact is called Case 1. To further validate the present model, two more impacts called Case 2 and Case 3 are simulated. Case 2 is the impact of a hard sphere with s Y =1380 MPa on a soft plate with s Y =345 MPa. Case 3 is the impact of a soft sphere with s Y =345 MPa on a hard plate with s Y =1380 MPa. For the impacts of the two cases, the other material and structural parameters are the same as those for impact of Case 1. For the proposed model in the present study, the yield strength is selected to be that of the softer material, thus the model provides the same predicting results for all the three cases. Figure 10 shows the comparisons of the maximum contact force, maximum indentation, and permanent indentation obtained by the FEM and DSM-VP for the three cases under different impact velocities. The relative error between the DSM-VP and FEM results is listed in table 3. Figure 10(a) and table 3 show that the maximum contact force P m predicted by DSM-VP is in good agreement with the FEM data for the three cases. The maximum relative errors are merely 5.3%, 8.6% and 6.5% for the impacts of Case 1, Case 2 and Case 3, respectively. Figure 10(b) and table 3 show that the maximum indentation d m predicted by DSM-VP is in good agreement with the FEM data for Case 1 and Case 3 with impact velocity v 0 20 m s −1 , and for Case 2 with impact velocity v 0 15 m s −1 . The maximum relative errors are 11.2%, 12.5% and 9.8% for the impacts of Case 1, Case 2 and Case 3, respectively. For Case 2 with v 0 =20 m s −1 , the error is 23.3%. Figure 10(c) and table 3 show that DSM-VP overestimates the permanent indentation. For the impacts of Case 1 with v 0 15 m s −1 , the permanent indentation predicted by DSM-VP matches the FEM data well, and the maximum relative error is 11.8% as v 0 15 m s −1 , and is 16.1% as v 0 =20 m s −1 . For the impacts of Case 2 and Case 3, the maximum relative errors reach 38.8% and 20.3%, respectively.
According to the above validations, it can be seen that the present model shows high accuracy in simulating the low-velocity and moderate-velocity impacts between the identical and dissimilar materials. The present model can predict the maximum contact force and maximum indentation with high accuracy. However, for the impact of a hard sphere on a soft plate and the impact of a soft sphere on a hard plate, the present model predicts the permanent indentation with relatively lower accuracy.

Experimental validation
Using DSM-VP and DSM-EP, we simulate the impact tests with v 0 = 3.13 m s −1 at four impacted points, as illustrated in figure 2(b). The geometrical size and material property are listed in table 2, and the Johnson-Cook  Figure 11 shows the force-indentation curves obtained by DSM-VP and DSM-EP. The maximum indentation reaches contact stage II. The impacts at different positions exhibit the same contact behavior. The material strain rate dependence shows a significant influence even for low velocity impact. It makes a smooth transition at the unloading point, as plotted in the inset I. It is seen from the inset I that the contact force reaches to maximum faster than the indentation. This delay can be attributed to the material strain rate dependence [46], and it does not appear in the case of rate-independent model, as shown in the inset II. Figures 12 and 13 plot the stress and displacement responses obtained by the experiments and DSM-VP, respectively, for four impact points. The stress is recorded by five strain gauges mounted on the plate bottom surface. For the impact at point D, different waves that sweep over strain gauges are denoted as S2, S3, S4, S5, and S12, respectively, according to the theoretical identification method [24,49]. Figure 12 illustrates that wave    propagations have been captured experimentally and numerically, and the DSM-VP results are in good agreement with the test data. Figure 13 shows that the displacement responses predicted by DSM-VP match the experimental results perfectly. It illustrates that DSM-VP can effectively simulate the low velocity impact response of plate.

Conclusion
We have presented a model for the dynamic contact of elastoplastic sphere and plate. The model parameters are determined theoretically. A Johnson-Cook strain-rate dependence is used to modify the effective yield stress based on the relative velocity of the sphere and plate. An analytical description of average plastic strain rate is provided based on a simplified model of plate contacted by a rigid sphere. The present model is validated by the experimental data and FE simulation results of the impacts between the identical and dissimilar materials. The predictions of the model are shown to agree with the results of FE simulations and experimental data. The model predicts the maximum contact force and maximum indentation with high accuracy. For the prediction of permanent indentation, the present model shows high accuracy in the case where the sphere and plate have the same hardness, but shows relatively lower accuracy in the case where the sphere is more rigid than the plate and in the contrary case. The model captures well the impact response of sphere and plate for different materials and for a range of loading conditions. With an aid of the present model, a significant decrease in contact pressure during fully plastic indentation, and the independence of dynamic contact-loading path upon loading rate are observed.