In typical metallic contacts, stresses are very high and result in yielding of the material. Therefore, the study of contacts which include simultaneous elastic and plastic deformation is of critical importance. This work reviews the current state-of-the-art in the modeling of single asperity elastic–plastic contact and, in some instances, makes comparisons to original findings of the authors. Several different geometries are considered, including cylindrical, spherical, sinusoidal or wavy, and axisymmetric sinusoidal. As evidenced by the reviewed literature, it is clear that the average pressure during heavily loaded elastic–plastic contact is not governed by the conventional hardness to yield strength ratio of approximately three, but rather varies according to the boundary conditions and deformed geometry. For spherical contact, the differences between flattening and indentation contacts are also reviewed. In addition, this paper summarizes work on tangentially loaded contacts up to the initiation of sliding. As discussed briefly, the single asperity contact models can be incorporated into existing rough surface contact model frameworks. Depending on the size of a contact, the material properties can also effectively change, and this topic is introduced as well. In the concluding discussion, an argument is made for the value of studying hardening and other failure mechanisms, such as fracture as well as the influence of adhesion on elastic–plastic contact.
Contact mechanics is critically important for a wide variety of applications, such as bearings, machine interfaces, mechanical seals, cams, electrical contacts, gears, tire/road contacts, joint structures, impact mechanics, and biomechanical systems. Physical surfaces are rarely perfectly flat where they come into contact. Surface curvature or roughness often causes the contact areas to be extremely small, and therefore, the resulting contact pressures and stresses are usually relatively high. In many cases, this results in failure or yielding in the contact regions. In metallic contacts especially, there is often plastic deformation. Such plastic deformation must be considered along with elastic deformation, and hence, metallic contacts are often labeled as elastic–plastic or elastoplastic.
This paper will summarize the state-of-the-art in the modeling of elastic–plastic contacts over a wide variety of geometries and loading situations. Specifically, it will focus on contact models for single peaks (often referred to as asperities), which can either be applied directly to single peak contact problems such as gears and rolling elements or integrated into statistical and multiscale mathematical frameworks that consider the complex case of rough surface contact. Due to the breadth of research into single asperity contact, cases of rough surface contact are only discussed briefly in Sec. 10.
Analytical solutions of contact problems are limited to very special geometries and constitutive behaviors [1–5]. Therefore, most works on elastic–plastic contact employ the finite element method (FEM) or a combination of elasticity theories and approximations (often referred to as semi-analytical models) to account for elastic–plastic behavior. Typically, a finite element model is created with the appropriate geometry, boundary conditions, mesh, and material properties. Next, a parametric study is performed by varying the material properties and geometrical parameters. Then, the outputs of the analysis, usually contact force, deformation, and contact area, are normalized in a logical way that consolidates them as much as possible. This normalization is often relative to a known elastic solution (such as Hertz contact theory) or to the initiation of plastic deformation. The initiation of plastic deformation can sometimes be found analytically by using a particular yield criterion (such as the von Mises criterion) and an appropriate elastic solution. The results of the finite element analysis (FEA) are next studied for trends and fitted to phenomenological equations. To be acceptable, such curve fits should be within a few percent of agreement with the results of the independent computational analysis. In some cases, they are also validated against experimental measurements. The goal of the phenomenological equations is to predict the behavior of elastic–plastic contact in science and engineering applications. One of the most important applications is to account for the asperity interactions between rough surfaces, which govern friction, plastic deformation, fatigue, and wear between solid surfaces.
Only a few previous reviews have been published on the topic of this paper [6–9]. These papers are over 16 years old, with minimal overlap with the current review, given considerable progress in the field of elastic–plastic contact mechanics has occurred since their publication. To achieve an even more thorough coverage over the area of contact mechanics, one could consider all four of these reviews in addition to the current paper. The review by Barber and Ciavarella  concentrates on elastic contact cases that include friction, wear, and roughness. The contemporaneous review by Adams and Nosonovsky  covers similar topics but concentrates on the forces occurring between contacting surfaces, including adhesion in elastic contact. The earlier work by Bhushan in 1998  reviews the elastic–plastic contact of rough surfaces. The most closely related work to the current paper is the Bhushan  review of single asperity contacts from 1996 .
Bhushan and Peng  also reviewed the narrow area of numerical methods (including the boundary element method and the finite element method) applied to the contact mechanics of multilayered rough surfaces and should be referred to when considering that topic, although a considerable amount of work has been occurred since its publication in 2002 [11–16]. The current review will not focus on the effects of coatings or layered contacts. Recent work in the area includes research into the inception of yielding in spherical shells , layered elastic–plastic spherical contacts [13–16], and the indentation of thin coatings for characterization of their material properties [18,19]. Recently, Chen et al.  completed a work on a universal model of elastic–plastic contact of thin layered or coated spherical surfaces.
This paper is organized logically based first on the type of loading (i.e., normal or tangential) and then on the geometry of the asperity. For each asperity geometry, the regimes of elastic, elastic–plastic, and fully plastic contact may also be examined individually. Much more research has been conducted on spherical asperity geometries than any other shape. Accordingly, these are covered more extensively than other geometries in this paper. First, the results of purely normal contact problems, with no tangential load and shear traction, will be analyzed. These problems are much simpler to solve because there is only one load on the system that is coupled directly to the pressure or normal traction distribution. In contrast, when both normal and tangential tractions are present, their effects are coupled through the contact area, which makes such problems more difficult to solve.
Normal Cylindrical Contact
Surprisingly, little work has been published on the contact of a cylindrical geometry against a flat (also known as line contact), as shown in Fig. 1. In Fig. 1, δ is known as the interference or normal displacement, F is the contact load or force, b is the contact area half width, and L is the length of the cylinder. Perhaps this lack of work is because this contact condition is not seen much in practical situations due to misalignment and geometrical imperfections. Also, it cannot be used to consider asperities on rough surfaces, unless they are highly anisotropic. Nevertheless, such contact applies to many cases, e.g., wheel contacts, gear teeth interaction, and cylindrical roller bearings.
Many have suggested that most practical instances of cylindrical contact fall into the case of plane strain (when the ends of the cylinder are held rigid) such that there is no strain outside of the cross section perpendicular to the cylinder's axis of symmetry [5,20]. Alternatively, one might consider a cylindrical rolling element bearing to be closer to the case of plane stress, since it is unconfined at the ends and therefore should have no stress acting outside of the plane perpendicular to the cylinder's axis of symmetry. Similarly, applications such as gears and cams should be in the plane stress state because they are usually lubricated (low friction) and are unrestrained at their ends. However, even with a significant amount of friction, if two identical cylinders are in contact, there will be no relative motion between them. Therefore, two identical cylinders in contact which are unconfined at their ends will be in the plane stress state. This is the case for the meshing of two identical gears.
In the case of contact between different materials, the contact yield strength is typically assumed to equal the yield strength of the material with the lower value of
where Eq. (3a) is valid as shown for plane strain and axisymmetric problems. Equation (3a) applies to plane stress problems, such as cylindrical contact with unconfined ends and negligible friction, by setting and to zero. Note that Eqs. (3a) and (3b) are also valid for spherical contact and are often used when considering both elastic and elastic–plastic contact.
where is the average pressure in the fully plastic regime.
It is also interesting to note that the average normal pressure (hardness) for a cylindrical indentation contact in plane stress predicted by FEA is much smaller relative to the case of a spherical indenter (around 1.9 to 2.0 times the yield strength of the flat ). This also agrees with Cinar and Sinclair , but the magnitudes are slightly different. For the plane stress case, it is clear that hydrostatic stress cannot form in the contact, and therefore, plastic deformation from “distortion” based stress (i.e., deviatoric stress) occurs much more readily.
Normal Spherical Contact
Single asperity spherical contact (also known as point or Hertzian contact) is arguably the most important problem in contact mechanics and has been the concern of researchers for more than a century [5,8]. Many analytical, experimental, and numerical studies have been performed to simulate and predict various contact properties, such as the real radius of contact, stress distribution, and contact force. Spherical contact models are used by many researchers from different fields such as tribology, mechanical impact [31–35], and electrical contact [11,31,32,36–54]. Even though the spherical geometry is often used to consider asperity contact, it can also be used to consider the flattening of particles between surfaces, such as in anisotropic conductive films  or in the presence of wear particles  and nanoparticles . As with the previously discussed case of cylindrical contact, spherical contact can be divided into three regimes: fully elastic, elastic–plastic, and fully plastic .
In the elastic regime, a closed-form analytical solution is available , although it is still approximate since it assumes a parabolic geometry. For the fully plastic regime, an analytical solution is only possible under simplifying assumptions, such as neglecting the effects of sink-in and pile-up. Although there have been many efforts to provide an analytical solution for the case of elastic–plastic contact, none of the presented models captures the full complexity of this condition. In the elastic–plastic regime, most of the published models simplify the analysis by studying the contact between a sphere and a flat and apply the finite element analysis to model this contact, starting with the works of Sinclair and Follansbee [56,57]. These models can be categorized into three main groups, namely, flattening models, indentation models, and more comprehensive models. In the flattening models, the flat surface is considered rigid, and the whole deformation occurs in the sphere. On the other hand, in the indentation models, the sphere is rigid, and the flat is deformable. Even though Jackson and Kogut  explicitly stated the differences between flattening and indentation, and the need to bridge between them, separate research on flattening and indentation has persisted without considering these differences. Only a few papers have considered general cases that account for deformation of both of the objects [32,41,59,60].
In Secs. 3.1–3.3, a summary of the investigations of the elastic regime, the fully plastic regime, and some of the most well-known elastic–plastic models is presented. The fully plastic regime is considered before the elastic–plastic regime, since it can be considered as the limit of the more complicated elastic–plastic cases. Finally, the results of original FEM simulations of three different contact cases representing flattening, indentation, and the intermediate case (by varying the ratios of yield strengths between the two bodies) are compared with the reviewed models.
For the elastic contact of two spheres, a closed-form analytical solution was provided by Hertz  and aroused considerable interest during the past century. Although the Hertzian theory approximates the spherical surfaces by parabolas, it is still very accurate for spheres. The Hertzian theory  can only be used for elastic contacts where no plasticity occurs. For contact involving metals, yield first occurs beneath the contact interface when a critical load/interference is applied and the stress state exceeds the strength of the material.
Fully Plastic Regime.
Next, the fully plastic regime of normal spherical contact is examined because, along with the elastic case, it is another limit of the elastic–plastic problem. The fully plastic regime is considered reached, when plastic deformation has completely engulfed the contact area at high loads. Many researchers [5,41,52,61–73] have studied the average pressure in the fully plastic regime to match theoretical predictions against the experimental indentation results of Brinell  and Meyer . Tabor's book on the hardness of materials has a very thorough but dated summary of this area . On a side note, a popular (but erroneous) truncation model for the fully plastic regime is often incorrectly attributed to Abbott and Firestone . The authors have obtained an original copy of this reference, and it appears to be on the bearing area curve used to characterize surfaces. Abbott and Firestone did not mention plastic deformation in their work.
As described earlier, flattening describes the case of a deformable sphere pressed against a rigid flat surface, while indentation describes the case of a rigid sphere pressed into a deformable flat surface. In the elastic regime, these can be considered equivalent contact conditions, but in the elastic–plastic and fully plastic regimes they differ significantly. Most work on spherical fully plastic contact has considered indentation as a way to characterize hardness (i.e., the average pressure during fully plastic contact). Hardness has historically been used as a way to define the strength of a material on its surface and actually originally derived from simple scratch tests . However, it is now mostly defined as a measure of a surfaces resistance to normal indentation and varies between specific methodologies.
The hardness is defined by Meyer  as the contact force divided by the projected area of contact in the direction of the indentation, , where is the contact force, and is the radius of contact (i.e., hardness is the average pressure, . Alternatively, Brinell  suggested that hardness is the contact force divided by the total surface area of the contact.
Conventionally, from the experimental results of the Brinell indentation test by Tabor , the popular ratio of (although sometimes this is approximated as 3) was accepted for the relationship of the average pressure or hardness to the yield strength of fully plastic spherical contacts without strain hardening. This average pressure definition is now referred to as the conventional hardness. To reduce the effect of strain hardening, Tabor used work hardened materials in all of his tests. In 1944, Ishlinskii  provided a new way of drawing the slip lines that helped to provide an analytical model for the indentation of a rigid sphere into a fully plastic flat. Tabor's hardness ratios was also verified theoretically using slip-line theory by Ishlinskii , who suggested that the hardness ratio varies between 2.61 and 2.84 and does not significantly depend on the indentation depth. In 1972, Lee et al.  came to a similar conclusion via finite element simulations and experimental measurements of spherical indentation. However, recent studies [11,79–81] have shown the dependence of the average normal pressure on the indentation depth and deformed geometry.
Deviations from the conventional hardness value were first observed experimentally by Chaudhri et al.  when studying the compression or flattening of metal spheres. Mesarovic and Fleck [64,79], Ye and Komvopoulos , and Kogut and Komvopoulos  found the same trend using finite element analysis for both indentation and flattening. The effect of strain hardening was shown to counter the effect of the average pressure deviating from the conventional hardness in some cases. However, it was not theoretically understood or explained until the works of Jackson and Green [52,70,71,82]. Jackson and Green theorized that as a sphere is either flattened or indented into a surface, its geometry approaches that of a compressed column, and therefore, the contact pressure also approaches the yield strength, which is approximately one-third of the pressure obtained from a shallow spherical contact (i.e., the conventional hardness).
Storåkers et al.  later used theoretical similarity methods in an attempt to better describe viscoplastic contact and also included a discussion of spherical indentation. Their analysis accounted for time-dependent strain and could be applied to creep problems. They also included hardening exponents in their equations. The resulting equation for hardness or average pressure does not produce the trend of a decreasing hardness with deformation.
Based on the work of Ye and Komvopoulos , Kogut and Komvopoulos  showed that the average pressure in elastic–plastic spherical indentation often reaches a maximum value that is less than . They studied the dependence of the maximum average pressure on the material properties and proposed a phenomenological equation based on the ratio of the elastic modulus to the yield strength. However, as was shown in several subsequent papers described later, the deformed geometry is what directly affects the average pressure.
Note that Eq. (8a) is from the original work by Jackson and Green (JG) , and Eq. (8b) is from work by Jackson et al. (JGM) [53,84]. The reason that the new equation was created is that theoretically when a/R approaches 1, the geometry approaches that of a compressed column and the average pressure should approach the yield strength. However, Wadwalkar et al.  showed that due to a barreling effect in flattening, R must be modified for large deformations for the average pressure to approach the yield strength when a/R approaches 1.
Their solution is based on the assumption of plane strain rather than axisymmetry (as these are mechanically very similar).
The comparison of Eqs. (6)–(10) with new FEM results for flattening and FEM results presented in Ref.  for indentation is shown in Fig. 2. The new FEM results for flattening employ a mesh of 29,541 elements and are for a material with , to approximate the rigid perfectly plastic case and the fully plastic limit of the average pressure to yield strength ratio. The cases of indentation and flattening are compared separately in Figs. 2(a) and 2(b), respectively. It has to be noted that the FEM results presented in Ref.  do not represent a rigid perfectly plastic material, but are close to the fully plastic case by applying a very small yield strength relative to the elastic modulus. For indentation, Eqs. (9) and (10) appear to be closest to the FEM results. For flattening, all three models appear to be in approximate agreement when , with Eq. (8b) showing the best match. They all deviate for larger values of . This may be due to large deformations of the sphere changing the effective radius, R. Notably, for low values of , the average pressure to yield stress ratio approaches 2.84 , which differs from Eq. (9). Above , the original JG equation (Eq. (8a)) deviates from the FEM results, but the updated version (JGM, Eq. (8b)) shows relatively good agreement. Two versions of the FEM results are also provided. In the first version, uses the original radius of the sphere. The second version adjusts to use the deformed radius of the sphere calculated from the base of the sphere. Both Eqs. (8b) and (9) fall in between these two versions. The JGM model (Eq. (8b)) follows the adjusted model more closely because it assumes that as contact radius, approaches the effective radius of the sphere, the pressure ratio will decrease to 1. The model by Olsson and Larsson (Eq. (9)) agrees better with the unadjusted FEM results because it was fit to other unadjusted FEM results.
During contact, the material usually begins yielding at very small interferences, especially for metals. Therefore, the elastic regime only covers a small range of interferences. On the other hand, the fully plastic regime happens at very large interferences. For the majority of contacts, part of the material is in the elastic regime, and part of it undergoes plastic deformations. This is called the elastic–plastic regime. So far, there are no closed-form solutions for elastic–plastic contacts that are derived from fundamental principles.
A majority of the studies on elastic–plastic contact use finite element analysis to model, simulate, and mathematically describe the contact. In Secs. 3.3.1–3.3.5, some of the most well-known spherical contact models have been reviewed in three groups: flattening, indentation, and the general case of plastic deformation in both the flat and the sphere.
The indentation of a flat by a rigid hemisphere has been studied for over a century by many researchers. As described earlier, the study of indentation contact sought to develop models for testing material hardness. In the indentation models, contact between a rigid sphere and a deformable flat is considered as shown in Fig. 3. Here, is the penetration, is the radius of real contact area, is the radius of truncated contact area which neglects sink-in and pile-up deformations in the flat, and is the radius of the sphere. In most experimental indentation studies, the indenter is considered to be rigid. However, elastic deformation in the indenter can cause prediction errors of up to 12%, as shown by Rodriguez et al. . Since, in reality, there is no ideally rigid body, there is a need for a criterion in which the contact can be approximately defined as indentation. Tabor  suggested that the yield strength of a spherical indenter should equal 2.5 times that of the flat, to prevent plastic deformation of the indenter. Ghaednia et al.  verified that if the yield strength of the sphere is 1.7 times larger than the yield strength of the flat, the contact can be practically designated as indentation.
In 1971, Hardy et al.  pioneered using the FEM to analyze the indentation problem. They investigated the elastic–plastic indentation qualitatively, but did not provide any predictive formulations. Using an improved meshing scheme, Sinclair et al.  and Follansbee and Sinclair  also studied indentation using the finite element method and used this to obtain some approximate phenomenological formulations. They compared their numerical predictions with experiments  and other existing theories, such as Ishlinskii's  and Richmond's  slip-line theories. Per Sinclair and Follansbee, Hertzian theory applies to the elastic regime when , where is the yield strength of the flat. Similar to Johnson , Sinclair and Follansbee assumed the average normal pressure in the fully plastic regime to equal the conventional hardness, , once the deformation reaches . However, as noted before, the assumption of a constant hardness is overly simplistic and makes predictions that can be erroneous.
where a is calculated from Eq. (14). For the elastic–plastic large deformation phase, when , Eq. (15) is used while setting . Komvopolous and Ye also noted a decrease in the average pressure to yield strength ratio at large deformations. Kogut and Komvopolous  reformulated Komvopolous and Ye's model and not only considered loading, but the unloading phases as well.
Bartier et al.  also performed an experimental and numerical verification of analytical models in terms of contact radius in an indentation contact. They showed that the theoretical models, dependent only on the strain hardening exponent, overestimate the dimensionless contact radius. Another observation in their study was the coefficient of friction that had a non-negligible effect on the contact radius for large penetration depths.
By dividing the contact problem into three phases, Brake  recently formulated a semi-analytical model. The three phases are the elastic, elastic–plastic, and fully plastic phases. The elastic phase follows the Hertz theory, and the fully plastic phase follows Johnson's theory , which assumes a constant conventional hardness. Brake considered the elastic–plastic phase using a polynomial interpolation between the two limiting phases. The model was verified by the impact experiments of Kharaz and collaborators [88–90]. In a more recent work , Brake used the same methodology with different transition functions to connect the elastic phase and fully plastic phase. His formulation for the fully plastic phase uses the Meyer hardness exponent. Brake  also considered the frictional contact for the oblique impact case. This model shows an improvement over the previous model presented in Ref. .
The case where a rigid flat is in contact with a deformable sphere is referred to as flattening and is schematically represented in Fig. 4. Note that is the real radius of contact and is the relative normal displacement (i.e., the interference or penetration depth) of the sphere into the flat, if the sphere was not deforming.
As with indentation, there are still no closed-form solutions of elastic–plastic flattening problems, but there are a large number of models. Johnson  experimentally observed the contact of spheres and cylinders with a rigid flat in the fully plastic regime in 1968. The first elastic–plastic contact model was proposed by Chang, Etsion, and Bogy , often known as the CEB model. In the CEB model, the sphere remains elastic up to a critical interference. For larger interferences, the model imposes volume conservation of the deformed sphere during plastic deformation. It also assumes that the average pressure on the sphere is the conventional hardness and does not vary once the fully plastic regime is reached. Additionally, the contact force to interference relation, as well as the contact area to interference relation, contains a discontinuity at the transition point from elastic to plastic deformation. Zhao et al. attempted to correct this discontinuity by using a polynomial template for the transition between the elastic and plastic regimes . Unfortunately, their correction was shown to compare very poorly to finite element predictions .
Thornton and Ning [42,43] formulated theoretical models of elastic–plastic contact with the aim of implementing them for the study of impact. With the development of the finite element method and access to computational resources, modeling capabilities were advanced by many researchers, including Li and Thornton , Wu et al. , Kogut and Etsion , and Jackson and Green . Shankar and Mayuram  studied the transition from the elastic–plastic regime to the fully plastic regime based on the evolution of the plastic core. Lin and Lin  also examined the problem and compared their finite element results to several of the models that will be discussed later.
For the contact of two elastic–plastic bodies, it has been suggested that the yield strength of the weaker material is used for the effective yield strength in Eq. (16a) . A better methodology is to use the lowest of the two materials .
Etsion's research group has also characterized the inception of yielding for coated or layered contacts [13–15]. Their work on layered surfaces predicts that the inclusion of a coating can actually increase the resistance of the surface to plastic deformation past that of a continuous surface composed entirely of either the substrate or coating material. There appears to be an optimal coating thickness that can be utilized by the coating industry to improve their design methodology. For instance, it has already been demonstrated that optimizing the coating thickness on machining tools can significantly prolong their lives .
After is exceeded, the contact is considered to be in the fully plastic regime and the average pressure is assumed to equal to the conventional hardness, i.e., . Wang  improved the formulation of Kogut and Etsion, especially for the unloading phase that was presented in Ref. . The predictions of Eq. (20), now referred to as the Kogut and Etsion model, were later confirmed using an optical measurement of the contact area .
The symbol /Sys in Eq. (21) captures the geometrically varying hardness and is equivalent to given by Eq. (8b). Substitution of Eq. (8b) into Eq. (21) is one of the biggest differences between Jackson and Green's model and previous models, since the former predicts that the hardness or the average pressure decreases as the contact radius increases. With larger deformations, the shape of the sphere approaches that of a cylindrical column loaded against a flat surface. In the fully plastic regime, the average pressure then decreases from the conventional hardness value of approximately as found by Tabor  and toward a value of .
Since the range of loads considered by Jackson and Green  were limited, Eqs. (9) and (21)–(23) are not valid beyond . Expanding on this work, Wadwalkar et al.  considered much larger loads and deformations at which the ratio of the contact radius to sphere radius was almost unity. Later, Sahoo and Chatterjee  also studied the effect of on the contact. Etsion's research group also removed the assumption of frictionless contact by considering full-stick between the surfaces [88,104]. Later, Megalingum and Mayuram  also revisited the spherical contact problem and confirmed that the hardness limit does not hold in the fully plastic regime. They provide a set of lengthy equations in the Appendix of their work for modeling spherical contact. The authors of the current review paper could not successfully include this model in the comparison with the other models.
Several researchers have also investigated the unloading of elastic–plastic spherical contacts [40,53,54,106]. In these works, the load is applied as before, but then it is removed incrementally to observe how the sphere recovers. In many of these works, it is assumed that the sphere recovers elastically and the Hertz model can be used to consider the elastic recovery of the sphere. Still, once the sphere is completely unloaded, a complex state of residual stresses remains in the body of the sphere due to the sustained plastic deformation. This residual stress is similar to that aimed for in the industrial practice of shot peening. Altering tensile and compressive stress regions are effective at stopping the propagation of cracks and, consequently, increasing fatigue life (i.e., cracks do not propagate in a compressive stress region).
Kadin et al.  found that all subsequent loading stages were elastic for both materials, with and without hardening, unless the original load is exceeded. This suggests that, in many repeated surface contacts, the deformation becomes elastic after the first few loading cycles. Kadin et al. also postulated that a more realistic kinematic hardening law might predict continuing hardening and plastic deformation. The unloading curves predicted by the modeling of Etsion's group were later validated using an optical technique to measure contact area . The experimental measurements of the unloaded and deformed radius of the spheres did not agree with the modeling predictions, perhaps due to differences in the boundary conditions at the surface.
The General Case (Combined Flattening and Indentation).
It was previously often assumed that the flattening and indentation models could be used for all contact cases. However, recently Ghaednia et al.  showed that there could be up to a 25% difference between these models, and that there is a transition phase between the flattening and indentation regimes. Jackson and Kogut  also found similar differences but did not provide a general model. Ghaednia et al. defined a new variable, the yield strength ratio, Sy*, defined as the ratio of the yield strength of the sphere to that of the flat. They showed that the yield strength ratio could be used to formulate the transition between the two limiting models, i.e., flattening and indentation. The deformations of both bodies, contact radius, contact force, and hardness were analyzed. For the convenience of the comparison, the authors kept the minimum yield strength constant and changed the yield strength ratio. In this way, previous flattening and indentation models will result in the same predictions for all of the yield strength ratios. Ghaednia et al.  also provided phenomenological equations that were more general than the indentation and flattening models. Due to the length of these equations, they are listed in Appendix A.
A recent paper by Olsson and Larsson  provided a model for the elastic–plastic contact of two spheres, rather than a sphere and flat. The two spheres can also have different radii, which make this model more general. The model is based on a semi-analytical formulation and fitted to the finite element results. It appears to be aimed for use in modeling the compaction of powder. Unfortunately, the model is not in a closed form. However, it may be useful to predict the behavior of the general case. Olsson and Larsson's paper did provide a closed-form prediction of the variation of the average pressure (i.e., hardness) as a function of the deformed geometry in the fully plastic regime. Their equation was compared to others in Sec. 3.2.
Figure 5 shows the FEM predictions of the contact radius for different yield strength ratios compared to several of the discussed models. For this case, the Young's modulus and Poisson's ratio of both the sphere and the flat are GPa and , respectively. The minimum yield strength is set to MPa, assigned to either the flat, the sphere, or both objects. The applied displacement is held constant at . The results show an interesting trend when the yield strength ratio is close to 1. At this point, the contact radius decreases sharply. The reason for such a trend is that the two objects then distribute and share the induced deformation thereby reducing the contact radius. Essentially, when the yield strengths are nearly equal, the amount of total plastic deformation approximately doubles at the same load compared to the ideal flattening or indentation cases.
Figure 6 shows the transition of the average normal pressure from flattening to indentation, as predicted by Ghaednia et al. . For this case, the modulus of elasticity (Young's modulus) and Poisson's ratio of both the sphere and the flat are GPa and , respectively. The minimum yield strength is set to Results show a sudden jump from the flattening regime with to the indentation regime with . Note also that, for both cases, the average pressure never reaches the conventional value of the hardness or . This is again due to relatively large deformations and contact shapes morphing from that of a sphere toward a column geometry. Essentially, the angle between the surfaces at the edge of the contact is increased from very small values toward 90 deg.
Figure 7 shows the Ghaednia et al.  finite element predictions of the contact force for different yield strength ratios compared to several earlier models. For this case, the modulus of elasticity (Young's modulus) and the Poisson's ratio of the sphere and the flat are GPa, and , , respectively. The minimum yield strength is set to MPa. The results show a slight decrease in force near a yield strength ratio of and then a sudden increase in force in the indentation regime. This behavior during the transition is not captured by flattening or indentation contact models. However, the model provided by Ghaednia et al. agrees closely with the FEM predictions and captures the intermediate behavior between indentation and flattening.
Numerical Comparison of Sphere Models.
Few of the well-known contact models for flattening and indentation that have been discussed in Secs. 3.3.1 and 3.3.2 have been compared to each other and to the FEM results. In this work, new numerical results are obtained from finite element modeling. The same model and meshing described in Ref.  have been used, but the cases are different from that presented in Ref. . The contact of a sphere with radius mm has been modeled for two different cases: indentation, in which the flat is the weaker object (with lower yield strength), and flattening, in which the sphere is the weaker material.
Five different contact models have been studied here: the Jackson and Green  and Kogut and Etsion models  from the flattening group, the Komvopolous and Ye  and Kogut and Komvopolous  from the indentation group, and the Ghaednia et al. model , which is the only model that provides a closed-form formulation that considers the effect of deformation of both of the objects. A comparison between these models for impact experiments can also be found in Ref. .
In the case of indentation, the yield strengths of the flat and sphere are MPa and MPa, respectively. This satisfies the criteria in Refs.  and  for an indentation, since the yield strength of the sphere is ten times larger than that of the flat. The modulus of elasticity and Poisson's ratio of both the flat and the sphere are GPa and , respectively. The contact force and radius of contact have been analyzed for deformations from to . Note that the radius of contact here is the radius as projected onto the original flat plane. One should consider that because of the relatively small yield strength of the flat, the deformations are mostly in the elastic–plastic regime.
Figure 8 shows the contact radius as a function of displacement during indentation. The FEM results shown in this figure are obtained for a model developed for this specific comparison rather than from any of the models used in previous publications. Both of the flattening models [52,93] show larger contact radius values than the FEM results, while both of the indentation models predict smaller radii of contact. This is because the pure flattening models underpredict the pressure during fully plastic contact, while pure indentation models overpredict it. The predictions of the model in Ghaednia et al.  show the best match with the FEM results.
Figure 9 shows variations in the contact force due to varying displacement during indentation. The Kogut and Etsion model  shows a large difference with the FEM results, likely because the formulation provided in their model is for very small deformations compared to the deformations considered here. The Jackson and Green model  predicts only slightly lower values than predicted by the FEM analysis, as would be expected from a model developed for flattening contact. Kogut and Komvopoulos  also predicts lower values than predicted by the FEM analysis. Again, the Ghaednia et al.  and Ye and Komvopoulos models  agree most closely, which is as expected because this case can be classified as indentation.
Next, flattening is considered. In this case, the yield strength of the flat and sphere are MPa and MPa, respectively. The Young's modulus and the Poisson's ratio of the flat and sphere are GPa, GPa and respectively. Similar to the previously considered indentation case, the deformations are mostly in the elastic–plastic regime.
Figure 10 shows the predicted changes in the contact radius due to the change in the displacement. Just as with indentation, both of the flattening models [52,93] predict larger values than predicted by the FEM analysis. In contrast, the indentation models [11,80] predict smaller values. Ghaednia et al.  and the Jackson and Green model  show better matches with the FEM results. These differences in prediction between the models are expected, since the finite element analysis used in this paper considers the deformation of both of the objects, while all the models but Ghaednia et al. do not. For example, for the flattening models, even very small deformations on the flat can decrease the contact radius significantly.
Figure 11 depicts the comparison of the predicted contact force between the models and the FEM results during flattening. As with the indentation case, the Kogut and Etsion model shows large differences for large deformations. For this case, the models of Ghaednia et al., Jackson and Green, and Kogut and Komvopoulos all show good agreement with the FEM results. Ye and Komvopoulos also showed a good match for smaller displacements .
It can be concluded from the numerical comparison that for the flattening case, the Jackson and Green  and Ghaednia et al. models  can provide better predictions for the contact force. For the indentation case, Ye and Komvopoulos  and Ghaednia et al.  show better agreement with the FEM. While the comparison here is not extensive, it has been performed for a wider range of material properties than in previous works.
Spherical Contact With Strain Hardening.
One of the major simplifications of almost all of the studies, whether experimental, theoretical, or numerical, is neglecting the material strain hardening. Tabor used work hardened materials  in his tests to exclude the effect of strain hardening on the Brinell hardness test. In the aforementioned theoretical and numerical works, the material is often assumed to be elastic–perfectly plastic.
where HP varies from to for most practical materials. Hence, the upper limit of could be one third of E. Brizmer et al.  analyzed the effect of a tangent modulus (2% of the Young's modulus) on a sphere with bilinear material properties contacting a rigid flat. The tangent modulus was also used by Shankar and Mayuram  to improve the effort in Kogut and Etsion  to account for the strain hardening effect. Interestingly, Kogut and Etsion  stated in their paper that even a high amount of hardening (the tangent modulus was 10% of the Young's modulus) caused only a fairly small deviation from the elastic–perfectly plastic results. Shankar and Mayuram  proposed a new empirical formulation similar to the one developed by Kogut and Komvopoulos  for three different values of the tangent modulus: , E, and . They showed that the interference ratio at which the contact changes to fully plastic also changes based on the ratios and . Later, Sahoo et al.  studied the effect of the tangent modulus on the flattening case of spherical contact. They analyzed the effect of strain hardening on the contact area and contact force and discussed the difficulties and perhaps impossible task of developing a predictive formulation when hardening is important. Sahoo et al. also confirmed that small amounts of hardening () do not influence the results significantly compared to the perfectly plastic case, but for HP = 0.1, the hardening results differ by up to 17% from the perfectly plastic case. For HP = 0.5, the load and area differ by up to 52% and 33%, respectively.
This reduces the hardness, and as a result, the contact force values decrease for the fully plastic regime. Note that there appears to be no fundamental mechanics based derivation for Eq. (26) . Instead, Eq. (26) is based on the similar equation for combining the elastic properties of two contacting bodies, which is derived from the theory of elasticity. However, a similar derivation is not possible for combining the plasticity properties of two contacting bodies. Brakes formulation has been verified with experimental results of indentation tests considering the effect of strain hardening.
Here, n is the strain hardening exponent, and the critical values should be obtained for the full-stick case using Eq. (19). Additional equations were also provided for the permanent deformation of the sphere after unloading.
Elliptical or ellipsoid contacts, i.e., contacts between surfaces with different principal radii or axes of curvature, have also been investigated, although not nearly as thoroughly as spherical contacts. Such elliptical contact could, for example, arise from the contact of two cylinders whose axes are not aligned. In the elastic regime, the case of elliptical contact can be analyzed using the Hertzian theory [6,15]. Horng  performed such analysis with a technique similar to that in the work of Chang et al.  on spherical contacts, which relied on a volume conservation approach for the elastic–plastic contact regime. Jeng and Wang  followed a similar methodology but with a template inspired by Zhao et al. . Jamari and Schipper  also followed a similar approach, and their model was verified with the experimental data of Chaudhri et al. . Lin and Lin  provided one of the first finite element analyses of elastic–plastic elliptical contacts and derived a phenomenological description that appears to be the most accurate formulation available for making predictions on the elastic–plastic elliptical contact. Their equations for predicted contact area, contact force, and deflection are provided in Appendix C.
Two-Dimensional Sinusoidal Contact
where L is the width of the wavy surface. Note that the contact area of the 2D wavy surface is rectangular, similar to cylindrical contact, and thus resembles a “line” contact.
To the best of our knowledge, Eq. (31) has not been presented previously. Gao et al.  provided an alternative to predicting the point of initial yield using a finite element model of elastic–plastic two-dimensional sinusoidal contact. In their work, the elastic limit is predicted by . Although of the same form as Eq. (31), this differs from Eq. (31) by approximately a factor of 2. Of course, the approximations in Eq. (31) are limited to small contact areas.
Gao et al.  noticed some clear differences between sinusoidal and cylindrical contacts. In the former case, at heavy loads, the average pressure rose to approximately 5.8 times the yield strength and was not limited by the conventional hardness of three times the yield strength. These high pressures were required to press the surfaces together into complete contact. Krithivasan and Jackson  also observed this increase in pressure for three-dimensional (3D) sinusoidal surfaces (see Sec. 7) and stated that there is no upper limit to the pressure as the amplitude is increased. They deduced that the absence of such a limit was due to hydrostatic stresses dominating as complete contact is reached. Later, Sun et al.  also noted the same increase in pressure when using dislocation plasticity theory to describe the plastic deformation in flattened two-dimensional sinusoidal surfaces. These observations have important ramifications when considering rough surface contacts. For instance, if a self-affine fractal description of a surface is assumed, this will result in zero contact area and infinite contact stress between contacting elastic–plastic rough surfaces . The experimentally observed phenomenon of asperity persistence  can also be explained by this mechanism.
Axisymmetric Sinusoidal Contact
and is graphically shown in Fig. 12. The elastic–plastic axisymmetric sinusoidal asperity that is described here considers interactions with adjacent asperities and also the effect of the substrate at the base of the asperity. Although the asperity contact described here may not be perfectly periodic due to the applied axisymmetric boundary conditions, it does effectively include interactions with adjacent asperities by having a confined boundary at the outer radius of the base. This is drastically different from single spherical asperity models which do not have this confinement and lateral interaction. Figure 12 shows the schematic and boundary conditions considered by this model.
The model has been developed for a broad range of material properties (yield strength in the range of 25 MPa < Sy < 40,000 MPa, Young's modulus in the range of 50 GPa < E < 400 GPa, and amplitude to wavelength ratios in the range of 0.00005< Δ/λ < 0.0125). Bilinear hardening was also assumed by using a tangent modulus of 0.01E. From the analysis of the finite element results, it is found that the dimensionless contact pressure depends on the dimensionless parameter, . Figure 13 shows the effect of (not ) on the dimensionless contact pressure–area relation for GPa and . It is clear that as increases, the pressure required to gain complete contact (with the entire surface in contact) increases.
This fit agrees with the curve very well and differs from it on average by 2.34% and overall by no more than 5.75%. Also note that can be much greater than conventional hardness for this geometry as well. This same phenomenon was also noted in other works considering the elastic–plastic contact of wavy surfaces [4,116,117,125,126].
This agrees closely with the finite element predictions and differs from these on average by 3.69%.
Song et al.  recently examined the same case of elastic–plastic contact of a sinusoidal surface and included the effect of hardening. They also observed that the pressure was limited to three times the yield strength when no hardening is allowed (in contrast to the results presented in the previous paragraphs). In contrast, when hardening was considered, the pressure observed by Song et al. increased past , which is in contradiction with the results in Ref. . To illustrate the possible increase in pressure with relatively little hardening, Eq. (33) is plotted in Fig. 15 as a function of , while normalizing the pressure by the yield strength. Here, typical properties for a steel alloy were assumed ( GPa and MPa), and the amplitude over the wavelength was varied in the range of . Not only does the critical pressure deviate from three times the yield strength but also may far exceed this value.
Another recent paper on the topic is by Liu and Proudhon , which also includes the effect of strain hardening. Some of the results shown in this work appear to indicate that the average pressure does surpass the yield strength by up to a factor of 6, well past the conventional value of three. However, the authors attribute this to hardening and still state that the pressure reaches a limit of approximately three times the yield strength. They do note the decrease in the pressure with large deformations, as discussed earlier in this paper in Sec. 3.2 on spherical contact.
Three-Dimensional Sinusoidal Contact
The asperities of rough surfaces are typically described by spherical or parabolic geometries in most rough surface contact models. However, the geometries of actual asperities are much different than this, especially at the base of the asperities. This is especially important when an asperity is under a heavy load. Therefore, wavy or sinusoidal surface models are important for the consideration of rough surface contacts used for modeling electrical resistance, friction, and wear .
where is given in Eq. (16b) and f = 1/λ.
The following results build mostly on a finite element simulation of three-dimensional elastic–plastic sinusoidal contact conducted by Krithivasan and Jackson . In their analysis, the range of material properties and geometries was somewhat limited ( and ). The studied range of was limited due to convergence difficulties of the nonlinear FEM model. The current work expands on this range.
The contact of sinusoidal or wavy surfaces also shows some interesting differences when compared to spherical contact. In a spherical contact that is dominated by plastic deformation, the average pressure decreases as the load and the ratio of the contact radius to the radius of curvature increase (see Sec. 3.2). This pressure is approximately 2.84 multiplied by the yield strength for small contact areas and decreases toward a limiting value equal to the yield strength as the load is increased . However, this trend does not occur in sinusoidal contacts. For the fully plastic case, Manners  used slip-line theory to show that as complete contact was approached, the pressure could grow beyond any bounds. This was due to the stress state becoming hydrostatic and therefore not allowing plastic deformation. In actual contacts, the surface pressure will remain bounded because elastic deformation allows for the surfaces to be pressed together. The average pressure does, however, increase relative to the yield strength and appears to have no upper limit. It increases with the amplitude of the sinusoidal surface being flattened.
Combined Normal and Tangential Loading
Sphere Normal and Tangential Loading, Presliding.
The elastic–plastic contact of spheres under combined normal and tangential loading has been studied for a long time, beginning with the classic works of Cattaneo  in 1938, Mindlin in 1949 , and Mindlin and Deresiewicz in 1953 . Interest in this area focused on deducing the underlying mechanics of friction. Mindlin used a predefined friction coefficient between two surfaces. He set an upper limit on the local shear stress, equal to the local normal stress multiplied by a predefined coefficient of friction. Known as the local Coulomb friction law, whenever the computed shear stress exceeds the upper limit, local slip takes place. Therefore, the sliding of the entire surface then occurs when all shear stresses over the contact area reach the upper limit, thereby satisfying the (global) Coulomb friction law. Mindlin also obtained the surface shear stress distribution for the full stick and partial slip conditions. Keer et al.  followed Mindlin's approach, obtaining the region for sliding of elastic bodies in contact. Hamilton  found the yield inception by using Hertz contact pressure and the Mindlin shear stress distribution. Sackfield and Hills  modified the stress distribution by considering the effect of the shear stress on surface displacements for two dissimilar cylinders.
Bowden and Tabor  presented a different approach, which considered the start of surface slip in relation to the mechanical properties rather than a local friction law as in Ref. . They used a failure mechanism related to the material properties to determine the sliding inception. They suggested that the tangential load at sliding inception was equal to the real contact area times the material shear strength. Courtney-Pratt and Eisner  measured the contact area of a metallic sphere pressed against a smooth flat. They observed an increase in the contact area when the tangential load was increased. Tabor  defined this phenomenon as “junction growth,” explaining that a contact area has already yielded plastically under a given preload and must grow when it is subjected to an additional tangential loading. This is because the tangential loading can reduce the mean contact pressure and require additional area to accommodate the additional shear stresses.
At approximately the same time as Kogut and Etsion, Zhang et al.  considered the sliding contact of a sphere against a flat surface using the finite element method. They were able to map the different regimes of elastic, elastic–plastic, and fully plastic contact with the sphere subjected to combined normal and tangential loading. They also examined the effect of the friction coefficient and tangential loading on the stress distribution, contact area, and load.
Later, Li et al.  used a similar methodology as Ref.  to account for larger deformations on the asperity level by incorporating the effects of the Jackson and Green model . Their model considered the statistical rough surface contact in the static friction model.
where is the local Coulomb friction coefficient.
Several experimental investigations of static friction and junction growth under combined normal and tangential loading in spherical contact have been reported in the literature. Static friction and contact area at sliding inception of a deformable sphere loaded against a rigid flat were measured by an experimental system employing load cells and displacement transducers in Ref. . This work confirms that the effective static friction coefficient does indeed decrease with normal load as predicted by the theoretical results.
Later, a novel test rig, which enables real time and direct in situ measurements of static friction, contact area, and relative displacement, was developed in Ref. . This rig allowed for the contact area between a sphere and a flat to be measured optically while under load. The initial results of the rig for normal loading were confirmed by Hertz theory and the Kogut and Etsion model . The rig was also used to confirm the theory of unloading of elastic–plastic spherical contact . Using this test rig, measurements of friction force and contact area at the instant of sliding inception were reported in Ref. . In this work, several different materials and spheres of varying sizes were loaded against a hard sapphire flat. It was observed that the static friction coefficient is strongly dependent on the normal load.
To investigate the static friction occurring in tin-plated electrical connectors, measurements were also conducted on tin interfaces using an automated tilting plane apparatus . The experiment measurement was of a nominally flat on flat interface. The static friction was observed to decrease with normal load. Good agreement between the measurements and the static friction model created by Li et al.  was obtained once the existence of the tin oxide on the surface was accounted for.
The phenomenon of junction growth, during which the contact area grows as tangential load is applied, was also confirmed experimentally using the same rig . Junction growth was observed to increase the contact area by up to 45%. Based on the measurements in Ref. , Etsion  revisited the Cattaneo–Mindlin [134,155] concept of interfacial slip in tangentially loaded compliant bodies and presented an alternative approach that treats sliding inception as a failure mode involving material plastic strain. Etsion argued from recent finite element models and experiments that the original assumptions made by Cattaneo–Mindlin were nonphysical but allowed the analytical solution to be developed.
Normal and Tangential Loading of Sinusoidal Contact.
This section described original results from finite element analysis of the sinusoidal geometry described by Eq. (38) rather than the typically assumed spherical geometry . In these computations, normal loading is applied first with the sinusoidal surface and the flat assumed to be in full stick. Next, the tangential loading increases gradually while the normal preload remains constant. The maximum friction shear stress criterion is used to determine the sliding inception. Local sliding occurs when the frictional shear stress at one element on the contact area reaches . The sliding of the whole asperity occurs when all the elements on the contact area slide.
Figure 17 presents typical results for the instantaneous dimensionless tangential load, , as a function of the dimensionless tangential displacement of the flat, , for several different combinations of material properties, geometrical parameters, normal pressure, and critical shear strength. In each case, as the tangential loading progresses, increases gradually until a constant value of is approached, corresponding to the static friction coefficient. As shown in Fig. 17, all the model parameters have some effect on the static friction coefficient.
Figure 18 shows the isolated effects of different parameters on the static friction coefficient. It can be seen in Figs. 18(a)–18(d) that the static coefficient of friction decreases as E,, and increase and, in contrast, increases with . An interesting finding is that when the geometric ratio, , is very small, it is very difficult to initiate sliding, and the static friction coefficient could be greater than 1, as shown in Fig. 18(d). This would be analogous to the high friction obtained from smooth clean surfaces in contact. For convenience, the nondimensional parameter, , proposed by Gao et al.  may be used to characterize the static coefficient of friction. From the finite element analysis, it was found that the static coefficient of friction decreases as increases.
Figure 18(e) shows the static coefficient of friction decreasing with increasing normal contact pressure. For low- and medium-strength contact pressures, the trend is the same as the spherical contact. However, under high contact pressures, the static coefficient of friction shows a different behavior. Under these pressures, the tangential load can cause complete contact for some cases, and in such a condition, the maximum shear force will not increase, no matter what the normal force is (i.e., the static friction coefficient will decrease continually).
In Ref. , the value of the critical shear stress was set to , representing the shear strength of two identical metal surfaces bonded together. However, often differs from this value as a function of temperature conditions, the presence of contaminants, and lubrication. The effect of on the static friction coefficient is shown in Fig. 18(f).
Similar to the equations for the sliding inception of a sphere against a flat, Eq. (58) only represents the interface of a single asperity contact, not a randomly rough surface. While it cannot be considered the same as the coefficient of friction between two objects with rough surfaces, it may be considered as a qualitative friction coefficient that provides the expected trends of a rough surface contact. For example, Eq. (58) does effectively predict the static friction for a periodic sinusoidal surface with many points of contact loaded tangentially.
Sliding Contact of Spheres and Cylinders.
A handful of papers have examined the sliding contact of two interacting curved surfaces. The case of a curved surface sliding against a flat surface has also been studied extensively, as noted in Sec. 8.2, but may not be representative of asperity interactions between two rough surfaces with similar magnitudes of roughness. This is because when two rough surfaces come into contact, most of the contact will occur between the asperities on both surfaces. Models of sliding of a curved surface against a flat would only be applicable to a rough surface when in contact with a much smoother surface.
For cylindrical line contact, Komvopoulos  studied the contact of a layered flat against a sliding rigid cylinder using FEM. One of the first papers to report using finite elements to study a sliding interaction between spheres is Faulkner and Arnell . Their study neglected the recovery phase of the interaction and therefore could not be related to friction. Jackson et al.  later expanded the investigation to include the recovery of the surfaces after contact. This study included both a semi-analytical and finite element based approach in order to derive phenomenological equations describing the forces during the interaction. Boucly et al.  studied the sliding and rolling interaction of two-dimensional asperities. Vijaywargiya and Green  concurrently studied the sliding interaction and energy loss occurring between cylinders in contact. Mulvihill et al.  continued this work on the interaction of sliding spheres and cylinders, but also included the effect of adhesional shear strength at the interface. This work by Mulvihill et al. is very significant because it includes the effect of both plowing and adhesive friction. Shi et al.  investigated sliding interaction between hemispherical asperities and also included lateral misalignment (i.e., in the previous studies, the centers of the sphere's passed directly above each other). Finally, Zhao et al.  recently included the effect of strain hardening into the sliding interaction of asperities. Others have also used molecular dynamics to analyze these types of interactions at the nanoscale [165–167]. Nonetheless, this is an area that still requires a great deal of work. If scientists and engineers hope to one day predict friction at a fundamental level, it will begin with modeling elastic–plastic asperity contact. In some of the above works, there appears to be great promise as there is already agreement with some experimental observations.
Adhesive Contact in the Elastic–Plastic Regime
In most of the previously discussed works, there are only repulsive tractions acting between the surfaces of the asperities. However, in the actual contact of surfaces, there is also an attraction that acts on the surfaces when they are out of contact. This is due to the nature of the interaction between the atoms of the surfaces, as described by their potentials. Surfaces will be attracted to each other until they are close enough that the atoms on one surface begin to repel atoms on the opposing surface. Such adhesion will result in the contact area being larger than if only repulsive tractions are considered.
There is a large amount of research investigating the adhesion of elastically deforming surfaces, but only a relatively small amount of work investigating elastic–plastic adhesion of contacting surfaces. This may be due to the fact that adhesive contact has a significant effect only on softer materials (e.g., rubber and polymer). In contrast, elastic–plastic contact often occurs between metals which are rigid enough to effectively resist the adhesive stress. Nonetheless, this section will provide a brief summary of research on elastic–plastic adhesive contacts.
As was the case with elastic–plastic spherical contact, probably the first work on elastic–plastic adhesive contact was performed by Chang et al. . In the CEB model, a semi-analytical solution is found by expanding the Hertz elastic contact solution to assume fully plastic volume-conserving contact once plastic deformation has initiated. Adhesion was included in the model by using the well-known Derjaguin, Muller, and Toporov (DMT) model  for elastic contacts. The DMT model simply adds the adhesive pressures to the outside of the contact region of the spherical contact, which then effectively adjusts the contact force. The inclusion of adhesion increases the contact area relative to a model that only considers repulsion.
where is the work of adhesion, is the maximum adhesive tensile stress, is the average elastic–plastic pressure, is the elastic–plastic contact radius, and and are both approximated by the hardness (corresponding to the fully plastic contact case). Generally, as increases, more plasticity is expected, while represents the transition from an unstable pull-off to a malleable one, during which the material between the surfaces yields rather than cracks (S = 1 is the approximate boundary between the two mechanisms).
Kogut and Etsion  followed in the footsteps of Chang et al.  by including adhesion directly in their semi-analytical solution in Eq. (20) (which was fitted to finite element results without adhesion) using the DMT model .
Sahoo and Banerjee  focused on elastic–plastic adhesive contact between rough surfaces, but also provided a methodology for including adhesion in elastic–plastic contact following the well-known Johnson, Kendall, and Roberts (JKR) model . The JKR model includes the adhesive pressure within the contact area but neglects it outside of the contact area. The JKR model is considered more applicable to soft materials as opposed to the DMT model which is more applicable to harder surfaces. Sahoo and Banerjee's work showed that adhesion reduced the stiffness between the surfaces and suggested that roughness can reduce the pull-off force required to separate surfaces.
Du et al.  examined the adhesive contact of an elastic–plastic sphere by including an adhesive stress based on the Lennard-Jones potential in the finite element model and may have been the first to do so. This approach allowed for a clear observation of the bifurcation of the pull-off mechanism. In some cases, the pull-off was sudden and brittle, while in other cases there was a ductile failure and a neck of plastically deformed material formed between the surfaces. Du et al. also suggested that the parameter proposed by Mesarovic and Johnson  could be used to predict which pull-off mechanism would occur. In a follow-up work, Du et al.  examined this transition and elastic–plastic spherical adhesion more closely. They found that the plasticity could be related to the applied interference relative to the critical interference. They also found that the transition from a brittle to ductile separation took place at approximately . Note that their calculations assumed that was equal to the hardness.
Later, Kadin et al.  performed a very similar analysis using FEM and spring elements to consider the adhesion between the surfaces. Here, again, the adhesive stress was based on the Lennard–Jones potential. Stress contours of the contacts during and after unloading were used to conclude that plastic deformation initiates just below the surface near the edge of the contact radius, as opposed to in an adhesionless contact where plastic deformation initiates just below the center of contact. The authors suggested that the yield strength should be substituted for in the expressions for S, and that the well-known Tabor parameter  was an important variable for characterizing elastic–plastic adhesive contact, in addition to S and the normalized deflection. High values of S and normalized deflection and low values of the Tabor parameter appeared to result in the most plastic deformation and hysteresis. Low values of the Tabor parameter appeared to cause more necking during pull-off.
Rough Surface Contact
This section is meant to only briefly introduce rough surface contact and ways that single asperity models are included in rough surface contact frameworks. There are many published papers in this area, and the references mentioned here are not exhaustive of the literature. If one seeks a more in-depth review of rough surface contacts, they can use one of the previous reviews on the topic [6–9,179], although these are somewhat out of date.
All surfaces are rough to some degree. The area in actual contact, labeled as the real contact area, is often much smaller than what it appears to be. It is very difficult to include this roughness in a model due to its complexity of randomness and features over multiple scales. Nonetheless, many researchers have approached this problem, and there are therefore many models that can be used to make predictions. It is very difficult to measure the real area of contact with precision to verify the models [151,180–188]; therefore, all of these models should be used with great care and caution.
After some early turmoil [189,190], in recent years researchers appear to be coming to a consensus of what a valid elastic rough surface contact model looks like and which of those in the literature meet these expectations [191–196] (except perhaps when adhesion is considered ). This type of progress has not yet been made with regards to elastic–plastic contact of rough surfaces. A few works have made some comparison between the various models, but it is extremely difficult to obtain a truly accurate deterministic model of elastic–plastic rough surface contact. Therefore, here the various kinds of models will be briefly introduced along with how they relate to single asperity models.
There are many possible ways to classify rough surface contact models, but the essential types are deterministic models, the hardness model, statistical models, stacked multiscale models, truncation models, and diffusion-based models. Many researchers also use fractal as a category, especially following the seminal work by Majumdar and Bhushan [198–200]. In this work, we do not consider fractals to be a category, because fractals actually characterize the multiscale nature of the surfaces and not the contact mechanics methodology used. Surfaces are now generally accepted to have features of roughness on many scales, and can, in some cases, exhibit a fractal structure. Many subsequent works [86,201–215] have employed fractals to characterize rough surfaces, typically using one the aforementioned models to handle the contact mechanics. Notably, natural surfaces may be rough but still not adhere strictly to a fractal description [216–219] or may be difficult to effectively characterize experimentally due to instrument limitations .
Each of the methodologies listed earlier for rough surface contact is summarized later and some representative works are discussed.
Deterministic models make very few assumptions about the surface geometry and are aimed at solving the “complete” problem. Liu et al.  provide an in-depth summary of some of the various methods used to solve deterministic rough surface contact problems. There have been a few attempts at using finite element analysis to describe elastic–plastic rough surface contacts, but a consensus and thorough comparison to the existing theories has not yet been achieved. For example, when implementing deterministic models that require a mesh, it may be difficult to obtain adequate mesh density to reach mesh convergence and a reliable solution.
Liu et al.  used the FEM with plastic deformation and the simplex algorithm to consider cylindrical and 2D rough surface contact in plane strain. Somewhat different from other works, Jacq et al.  and Sainsot et al.  used a semi-analytical approach to solve the elastic–plastic problem. Their method relied on a modified boundary element method to account for plasticity and is more efficient than finite element models. It should be noted that this semi-analytical approach may be limited to cases where the plastic deformation is localized and not severe. Similarly, Wang et al.  performed a deterministic elastic–plastic rough surface analysis based on a semi-analytical approach that used the conjugate gradient method and the fast-Fourier transform. They validated their results using commercial finite element software. A similar work solved the coupled elastic–plastic contact and frictional heating thermal problem using a fast-Fourier transform solution methodology .
Pei et al.  used finite elements to investigate the elastic–plastic contact of self-affine fractal surfaces. They also provided fitted relationships for pressure and contact area as a function of surface and material properties. Liu et al.  also used finite elements to consider the elastic–plastic contact of rough surfaces in a microelectromechanical system. Sahoo and Ghosh  examined the deterministic contact of elastic–plastic fractal surfaces and made some comparisons to experimental measurements with mixed results. Thompson and Thompson [227,228] discussed the methodologies used to construct and implement rough surface contact finite element models with the aim of predicting thermal contact resistance. Specifically, they showed how to effectively incorporate real rough surfaces into a finite element mesh. Wang et al.  attempted to accelerate the elastic–plastic deterministic solution time by using a multilevel approach . Megalingam and Mayuram  constructed a finite element based deterministic model of rough surface contact and made some comparisons to a hybrid statistical/deterministic model. In the hybrid model, they used the real surface geometry but considered individual asperity contacts as described by elastic–plastic spherical contact models.
The so-called hardness model is one of the first and most basic rough surface contact models. Here, the contact area is predicted by dividing the applied load by the hardness, assumed to equal three times the yield strength. Based on the repeated discussions in this review of how pressure during fully plastic contact varies widely from this conventional hardness, one might conclude that this model has become antiquated. However, it can provide a quick order of magnitude estimate of the contact area between rough metallic surfaces. This model most likely was originally developed by Bowden and Tabor , but its provenance is uncertain.
where is the contact force asperity model. For and , many of the elastic–plastic asperity models discussed in Secs. 3, 4, 6, and 7 of this review can be used. Note that an analytical solution of the original elastic and Gaussian version of the Greenwood and Williamson model is available in Ref. .
Although Greenwood and Williamson did not include an elastic–plastic model, many subsequent works have done this successfully [3,92,96,111,112,236–239]. It does appear that the roughness of the surface tends to average or reduce the influence of the chosen single asperity model (there is sometimes not much variation between rough surface contact models that employ different asperity models). Nonetheless, the accuracy of the chosen elastic–plastic asperity model is important and influences the resulting effectiveness of the rough surface contact model . Although the statistical family of models is the most widely known and used, they are often criticized for considering the asperities to be all the same size and radius of curvature. Indeed, in its original form, it did not consider the multiscale nature of surfaces [241,242]. In addition, most statistical works employ spherical asperity models, which also results in the asperities being effectively isolated and unable to interact with adjacent asperities. However, refinements to the model have been proposed, and it is still considered a very useful model [202,243,244]. For example, the omission of adjacent asperity interaction might be corrected by using a wavy asperity model or a corrected asperity interaction model . Shoulder to shoulder asperity contact when the asperities are misaligned could also be an important mechanism to consider [237,238]. Recently, Beheshti and Khonsari [101,102] considered the effect of an underlying curvature of a larger scale feature, such as a cylinder, on the elastic–plastic contact using a statistical model. Since they consider multiple scales of features, their analysis leads to the next type of models which inherently consider multiple feature scales.
Stacked Multiscale Models.
where Ar is the real area of contact, is the contact area of a single asperity on a certain scale of roughness, η is the real asperity density, An is the nominal contact area, and the subscript i denotes a specific asperity scale level, with imax denoting the smallest scale level considered. Since each scale bears the same load, , the single asperity load at the scale can be related to the total load by . This assumes that load is distributed uniformly across all of the asperities at a given scale, which may be a valid assumption for elastic–plastic contacts since the pressure distribution is often uniform. The spherical and sinusoidal elastic–plastic asperity contact models discussed in Secs. 3, 4, 6, and 7 can easily be employed in the multiscale model, as has been done in several previous works [119,185,247–251]. Gao and Bower  included elastic–plastic deformation in a similar multiscale model while assuming a Weierstrass fractal description of the surfaces. Goedecke et al. [205,206] later incorporated a fractal description of the surfaces into the multiscale model along with elastic–plastic deformation and creep, with the aim of predicting dynamic friction.
One of the earliest models to employ a fractal description to describe rough surfaces also used a simple truncation of the rough surface geometry to determine the contact area [198–200]. In fact, many refer to this truncation model as a fractal model, but it is categorized differently here because other families of models can also employ fractals. The model determines the load on each contact patch using either elastic Hertz contact or, in the elastic–plastic regime, by assuming the pressure to equal the conventional hardness. Contrary to intuition, this resulted in the contact becoming more elastic as the size of the contact patches grew, because they have a larger radius of curvature. Although this model did advance the field at the time, it has been severely criticized for this counterintuitive behavior . Nonetheless, the model did show reasonable agreement with some experimental measurements.
Later, Chung and Lin  modified this fractal truncation model to include the elliptical shape of the asperity contacts. They also made impressive comparisons to experimental measurements and other elastic–plastic models. In most cases, all of the data appeared to be in agreement.
Since surfaces are now considered to be multiscale in nature and perhaps fractal, they can contain many scales or layers of asperities. Since in real surfaces, these scales should possess a continuous rather than discrete spectrum, Persson  uses a diffusion theory to model rough surface contact. The general concept is that the pressure diffuses through the different scales of asperities. Persson's diffusion model has been shown to be very effective for elastic rough surface contact [255,256], but has not been thoroughly evaluated for elastic–plastic contact. Persson also appears to account for plasticity in a very simplistic way that is similar to the hardness model. The pressure appears to be simply truncated at the strength of the surface material.
There have been other attempts at comparing and verifying the accuracy of these various models, in addition to those discussed earlier. A comparison between various statistical models that incorporate several different spherical contact formulations was made by Jackson and Green . Additionally, a comparison between fractal and statistical models was made by Kogut and Jackson . This area clearly requires much more attention.
Equations (65) and (66) have very similar forms even though they are derived from different methods for considering the roughness of the surface. As with the original plasticity index, higher values of suggest that a surface is dominated by plastic deformation. When using Eq. (66), there is a clear transition to the elastic regime at .
In some cases, it is important to consider the scale of contact when considering what material properties to include in a model. As a contact decreases below the microscale, the relevant materials properties, such as elastic modulus and yield strength, can change significantly from the bulk properties. The largest and most influential change usually occurs with the yield strength. Generally, with decreases in scale, the yield strength increases well past its bulk value until it reaches a theoretical limit of approximately a tenth of the elastic modulus . For instance, the experimental results of Greer and Nix  showed how the yield strength of different diameter pillars in compression effectively changed with their diameter. Burek performed a similar study for tin .
Here, is the bulk yield strength of the material, is a characteristic length scale for the material, and is the depth of plastic indentation. To employ strain-gradient plasticity in an asperity model, one needs to iteratively solve one of the elastic–plastic asperity models considered in this review while coupling it to Eq. (68) by the yield strength and plastic indentation depth. Additional information about scale effects and nanomechanics can be also found in a book by Bhushan .
Early research to consider scale-dependent properties in contact mechanics is found in Polonsky and Keer [263,264]. The authors modeled plastic deformation by considering matrix dislocation nucleation and motion, and so their model was inherently scale dependent. A later study by Bhushan and Nosonovsky  was likely the first to employ the strain-gradient method in surface interactions. They used a scale-dependent yield strength model based on strain-gradient plasticity and considered how this can affect contact and slip at small scales. Jackson  followed this by considering scale-dependent strength in an elastic–plastic statistical rough surface contact model. The scale-dependent strength tends to increase with smaller scales, and this work showed that the contact area, consequently, decreases when scale-dependent strength is included .
Almeida et al. [267,268] considered how scale-dependent properties could be important in the contact of microswitches in microelectromechanical systems. Their work used a stacked multiscale model rather than the statistical framework to consider rough surface contact and electrical contact resistance. Similarly, Adams et al.  considered how the friction force changes with the scale of the asperity contact area using a friction dislocation theory proposed by Hurtado and Kim [270,271] and also accounted for the effect of adhesion. A multiscale rough surface contact model was later formulated to predict the thermal contact resistance between surfaces while considering scale-dependent strength [185,272]. Jackson et al.  used a similar methodology to consider the multiscale rough surface electrical contact resistance. There, the scale dependence of both the strength and electrical resistivity was considered. The work suggested that not all rough surface contacts would be influenced by the scale-dependent properties.
As an alternative to considering scale-dependent plasticity as a result of scale-dependent yield strength, others have used FEMs structured into individual crystalline grains in the material. This type of model should inherently account for the size effect of the grain boundaries. Therefore, each grain is modeled as a set of elements with a different set of anisotropic material properties which correspond to their crystal orientation, Most of the research using this technique aims to analyze fretting fatigue. For example, Goh et al.  considered a rigid cylinder loaded normally and tangentially into a polycrystalline elastic–plastic flat surface. The grains were randomly oriented, about 1/200 the size of the cylinder radius, and consisted of 100 elements on average. The model made predictions of elastic shakedown after repeated load cycles and the onset of ratchetting and fretting.
Researchers have also used single crystalline structures with dislocations to consider the scale-dependent behavior of surface contact. Nicola et al.  considered a simplified rendering of a rough surface contact that was a periodic set of evenly spaced flat punches. They found that friction had a negligible effect on the results, and that there was a scale-dependent effect when the contacts were on the scale of microns (which is similar to strain-gradient theory). This scale effect was observed by comparing conventional plasticity predictions to the discrete dislocation predictions and predicts that the strength of the material effectively increased for small-scale contacts.
Sun et al.  used a similar discrete dislocation methodology to consider the flattening of a sinusoidal metallic surface. They also observed a size dependence, and that the pressure required to flatten the surface could nearly reach the theoretical limit of strength (approximately 1/10th of the elastic modulus). This was observed while decreasing the wavelength and holding the amplitude of the features constant, which also effectively increases the amplitude to wavelength aspect ratio of the asperity. Krithivasan and Jackson , Gao et al. , and Saha and Jackson  also noticed that the required pressure increased with this aspect ratio without the inclusion of any scale effects (as noted in Sec. 5). In addition to the scale effects, such an increase could therefore also be partially attributed to the effect noted by Krithivasan and Jackson  and Manners , where to reach complete contact (flattening) the stress becomes hydrostatic due to the periodic boundaries and plasticity reduction. Later, Sun et al.  also noted that the increase was due to the dislocation plasticity being limited by periodicity (i.e., the regular spacing between adjacent asperities). This might be understood as the same mechanism as described in Ref.  using von Mises plasticity, which theorizes that plastic deformation is only caused by deviatoric stresses. Hence, the observed increase in the pressure could be obtained with a combination of these effects. These effects also explain the observed phenomenon of asperity persistence from the 1970s  (where experiments on complete contact required pressures much higher than expected).
Very recently, Ng Wei Siang and Nicola  investigated contact between a 2D wavy surface against a flat, in which the surfaces were again rendered as single crystal structures (referred to as a discrete dislocation plasticity analysis). Therefore, their analysis is an effective model only of a very small-scale contact. Ng Wei Siang and Nicola also investigated the differences between flattening and indentation as discussed previously in Sec. 3. They suggested that under certain conditions, the contact of two elastic–plastically deforming surfaces can be reduced to a case of a deformable surface against a rigid surface, but that this reduction should be performed with great care. They also noted some scale effects in the plastic deformation, but these diminished as the load was increased.
If one is considering the contact of nanoparticles, the scale-dependent strength becomes very important. Nanoparticles have been shown by experiments to possess an elevated strength [278–280] due to dislocations easily traveling to the surface, rendering a strong crystalline lattice in the particle. Many nanoparticles are spherical in shape, and so an elastic–plastic spherical model might be used in unison with scale-dependent properties to consider nanoparticle contact. Ghaednia and Jackson  applied such a model to the contact between nanoparticles dispersed between rough surfaces. They theorized that the particles are strong enough to increase the separation between the surfaces and therefore decrease the real area of contact which can decrease the friction.
This work provides a broad overview of the current state-of-the-art in modeling elastic–plastic asperity contact. Of course, this is a continuously growing area, and the current work is not exhaustive. For instance, it did not discuss in detail the interaction of coated (i.e., layered) surfaces [13–15] or time-dependent creep and relaxation [125,281–284]. Nonetheless, it is clear that elastic–plastic interactions are critically important for many engineering applications and not only for rough surface contacts. Due to the small size of typical surface contacts, it is very likely that plastic deformation will occur, especially in metallic contacts. Since most surfaces are now considered to be rough over many scales, similar to a fractal, it is not always obvious how to implement elastic–plastic models since the theory of elastic superposition no longer applies. Although many current models do exist [50,102,236], experimental testing of these models is one of the most important areas of future research.
This review highlighted the observation that the pressure during plastic contact is not strictly governed by the conventional hardness. Indeed, many works show that the average pressure can range from approximately the yield strength to an order of magnitude or more past it.
Another section examined the case where both surfaces are able to deform elastic–plastically, which differs from most previous analyses that only consider the deformation of either the asperity (flattening) or the flat (indentation). Interesting phenomena not captured by other models occur when the strengths of both materials are similar.
The derivation of the critical amplitude to allow plastic deformation of a sinusoidal contact was also corrected, and the critical values at which initial plastic deformation occurs in two-dimensional sinusoidal contacts were presented for perhaps the first time. However, it appears that the correction does not result in a drastic change in the results, and previous models of elastic–plastic sinusoidal contact are still very accurate.
Additionally, tangentially loaded contacts were summarized for spherical contacts and sinusoidal contacts. The offspring of these models may one day allow engineers to predict the friction based on the fundamental mechanics occurring between the surfaces.
Representative work in elastic–plastic adhesive contacts was given an overview. The work in this area has evolved from semi-analytical approximations to a more deterministic approach that incorporates adhesion in the elastic–plastic finite element models. Researchers have observed that there are different mechanisms of pull-off that depend on the stability and ductility. In the authors' opinion, this area of elastic–plastic adhesion requires the most focus and work in the future.
In another section, the concept of rough surface contact was introduced along with the various available models. This section is meant mostly as an introduction to the subject and to give an idea of how the individual asperity models can be used in a rough surface contact. The plasticity index was also reviewed as a way to qualitatively determine how much plastic deformation will occur in a rough surface contact.
Finally, Sec. 11 reviewed the subject of scale-dependent properties, and specifically scale-dependent yield strength. This is especially important in rough surface contacts where the small-scale asperities can possess strength much higher than bulk values. This will effectively reduce the area of contact between surfaces. Still, for some surfaces, the effect does not appear to be that important.
The modeling of elastic–plastic contact is clearly still a growing area that will continue to be a critical component of research. If the ultimate goal of tribologists is to predict friction, elastic–plastic contact mechanics will play a critical part in it. There are several areas that could benefit greatly from future research, for example, sliding elastic–plastic asperity contact, including adhesion in elastic–plastic contacts, as well as the inclusion of strain hardening in contact models. In addition, the analysis of scale-dependent properties is an area that is very important for rough surface contact, but that has received little attention.
Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation.
Division of Civil, Mechanical and Manufacturing Innovation, National Science Foundation (Grant No. 1362126).
In this section, we review the formulation developed by Ghaednia et al.  for the single asperity elastic–plastic contact of a deformable sphere with a deformable flat. The formulation is divided into five steps: finding deformations on each object, updating the radius of curvature, calculating the contact radius, calculating the hardness, and finding the contact force. The displacement of the center of the sphere is , where and denote the deformations of the sphere and the flat, respectively. These parameters are normalized as , , and . Using the Ghaednia et al. formulation, the deformation ratio of the sphere, , may be expressed in terms of the normalized displacement, , and the yield strength ratio, , as shown here
The contact radius has been divided into phases and subphases based on and respectively, where is found from Eq. (16).
where the fully elastic forces, , , follow Hertzian theory, the fully plastic force, , and and are the weights for the fully elastic and fully plastic regimes, respectively, and are calculated as follows:
Since, in the current case, wavelengths in the x and y directions are equal, β = α = 2πf.
where represents the dimensionless z coordinate of the function's maximum value.