A Review of Elastic–Plastic Contact Mechanics

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 [8] concentrates on elastic contact cases that include friction, wear, and roughness. The contemporaneous review by Adams and Nosonovsky [9] covers similar topics but concentrates on the forces occurring between contacting surfaces, including adhesion in elastic contact. The earlier work by Bhushan in 1998 [7] reviews the elastic–plastic contact of rough surfaces. The most closely related work to the current paper is the Bhushan [6] review of single asperity contacts from 1996 [6].

Bhushan and Peng [10] 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 [17], 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. [16] 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.

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 $C1Sy.$

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 $ν1$ and $ν2$ 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 $p¯$ 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 [29]). This also agrees with Cinar and Sinclair [20], 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.

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 [51] or in the presence of wear particles [55] and nanoparticles [50]. As with the previously discussed case of cylindrical contact, spherical contact can be divided into three regimes: fully elastic, elastic–plastic, and fully plastic [5].

In the elastic regime, a closed-form analytical solution is available [1], 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 [58] 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 [1] 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 [1] 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.

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 [74] and Meyer [75]. Tabor's book on the hardness of materials has a very thorough but dated summary of this area [67]. On a side note, a popular (but erroneous) truncation model for the fully plastic regime is often incorrectly attributed to Abbott and Firestone [76]. 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 [77]. 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 [75] as the contact force divided by the projected area of contact in the direction of the indentation, $F/πa2$⁠, where $F$ is the contact force, and $a$ is the radius of contact (i.e., hardness is the average pressure, $p¯)$⁠. Alternatively, Brinell [74] 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 [67], the popular ratio of $p¯/Sy=2.84$ (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 [27] 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 [27], 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. [78] 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. [81] when studying the compression or flattening of metal spheres. Mesarovic and Fleck [64,79], Ye and Komvopoulos [11], and Kogut and Komvopoulos [80] 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. [83] 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 [11], Kogut and Komvopoulos [80] showed that the average pressure in elastic–plastic spherical indentation often reaches a maximum value that is less than $2.84Sy$⁠. 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) [52], 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. [72] 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. [2] 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 $E/Sy=2×104$⁠, 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. [2] 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 $a/R<0.4$⁠, with Eq. (8b) showing the best match. They all deviate for larger values of $a/R$⁠. This may be due to large deformations of the sphere changing the effective radius, R. Notably, for low values of $a/R$⁠, the average pressure to yield stress ratio approaches 2.84 [27], which differs from Eq. (9). Above $a/R=0.4$⁠, 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, $a/R$ uses the original radius of the sphere. The second version adjusts $a/R$ 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, $a,$ 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, $a$ is the radius of real contact area, $a′$ is the radius of truncated contact area which neglects sink-in and pile-up deformations in the flat, and $R$ 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. [59]. 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 [66] 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. [41] 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. [85] 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. [57] and Follansbee and Sinclair [56] also studied indentation using the finite element method and used this to obtain some approximate phenomenological formulations. They compared their numerical predictions with experiments [67] and other existing theories, such as Ishlinskii's [27] and Richmond's [65] slip-line theories. Per Sinclair and Follansbee, Hertzian theory applies to the elastic regime when $p¯/Syf≤1.1$⁠, where $Syf$ is the yield strength of the flat. Similar to Johnson [5], Sinclair and Follansbee assumed the average normal pressure in the fully plastic regime to equal the conventional hardness, $p¯=2.8Syf$⁠, once the deformation reaches $aE′/RSyf=50$⁠. 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 $Q>400$⁠, Eq. (15) is used while setting $A=A′/0.71$⁠. Komvopolous and Ye also noted a decrease in the average pressure to yield strength ratio at large deformations. Kogut and Komvopolous [80] reformulated Komvopolous and Ye's model and not only considered loading, but the unloading phases as well.

Bartier et al. [87] 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 [47] 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 [5], 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 [48], 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 [48] also considered the frictional contact for the oblique impact case. This model shows an improvement over the previous model presented in Ref. [47].

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 $a$ 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 [91] 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 [92], 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 [73]. Unfortunately, their correction was shown to compare very poorly to finite element predictions [52].

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 [44], Wu et al. [45], Kogut and Etsion [93], and Jackson and Green [52]. Shankar and Mayuram [94] studied the transition from the elastic–plastic regime to the fully plastic regime based on the evolution of the plastic core. Lin and Lin [95] 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) [96]. A better methodology is to use the lowest $C2Sy$ of the two materials [22].

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 [98].

After $δ/δc=110$ 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., $2.8Sys$⁠. Wang [99] improved the formulation of Kogut and Etsion, especially for the unloading phase that was presented in Ref. [54]. 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 [100].

The symbol $HG$/S_ys in Eq. (21) captures the geometrically varying hardness and is equivalent to $p¯/Sy$ 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 $2.84Sys$ as found by Tabor [67] and toward a value of $Sys$⁠.

Since the range of loads considered by Jackson and Green [52] were limited, Eqs. (9) and (21)–(23) are not valid beyond $a/R=0.41$⁠. Expanding on this work, Wadwalkar et al. [72] considered much larger loads and deformations at which the ratio of the contact radius to sphere radius was almost unity. Later, Sahoo and Chatterjee [103] also studied the effect of $E′/Sy$ 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 [105] 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. [107] 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 [106]. 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.

Several unloading models were successfully used to describe the elastic recovery phase of spherical impacts [53]. The resulting models compared well with experimental measurements of the coefficient of restitution [31,32].

It was previously often assumed that the flattening and indentation models could be used for all contact cases. However, recently Ghaednia et al. [41] 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 [58] also found similar differences but did not provide a general model. Ghaednia et al. defined a new variable, the yield strength ratio, S_y*, 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. [41] 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 [60] 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 $Es=Ef=200$ GPa and $νs=νf=0.3$⁠, respectively. The minimum yield strength is set to $Sy,min=503$ MPa, assigned to either the flat, the sphere, or both objects. The applied displacement is held constant at $0.04R$⁠. 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. [41]. For this case, the modulus of elasticity (Young's modulus) and Poisson's ratio of both the sphere and the flat are $Es=Ef=200$ GPa and $νs=νf=0.3$⁠, respectively. The minimum yield strength is set to $Sy,min=200MPa.$ Results show a sudden jump from the flattening regime with $p¯/Sy=2.4$ to the indentation regime with $p¯/Sy=2.7$⁠. Note also that, for both cases, the average pressure never reaches the conventional value of the hardness or $p¯/Sy=2.84$⁠. 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. [41] 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 $Es=71.7$ GPa, $Ef=200GPa$ and $νs=0.33$⁠, $νf=0.3$⁠, respectively. The minimum yield strength is set to $Sy,min=300$ MPa. The results show a slight decrease in force near a yield strength ratio of $Sy*=1$ 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.

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. [41] have been used, but the cases are different from that presented in Ref. [41]. The contact of a sphere with radius $R=1$ 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 [52] and Kogut and Etsion models [93] from the flattening group, the Komvopolous and Ye [11] and Kogut and Komvopolous [80] from the indentation group, and the Ghaednia et al. model [41], 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. [32].

In the case of indentation, the yield strengths of the flat and sphere are $Syf=100$ MPa and $Sys=1000$ MPa, respectively. This satisfies the criteria in Refs. [41] and [66] 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 $E=200$ GPa and $ν=0.3$⁠, respectively. The contact force and radius of contact have been analyzed for deformations from $δ/R=0.0005$ to $δ/R=0.06$⁠. 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. [41] 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 [93] 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 [52] 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 [80] also predicts lower values than predicted by the FEM analysis. Again, the Ghaednia et al. [41] and Ye and Komvopoulos models [11] 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 $Syf=530$ MPa and $Sys=165$ MPa, respectively. The Young's modulus and the Poisson's ratio of the flat and sphere are $Ef=201$ GPa, $Es=200$ GPa and $νf=νs=0.3,$ 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. [41] and the Jackson and Green model [52] 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 $(δ/R≤0.03)$⁠.

It can be concluded from the numerical comparison that for the flattening case, the Jackson and Green [52] and Ghaednia et al. models [41] can provide better predictions for the contact force. For the indentation case, Ye and Komvopoulos [11] and Ghaednia et al. [41] 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.

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 [66] 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.

Here, n is the Meyer hardness exponent [63,66] and is related to the strain hardening exponent, $nϵ$⁠, by $nϵ=n-2$⁠. It should be noted that $n=2$ for elastic–perfectly plastic materials and $n=3$ for elastic contact.

where HP varies from $0$ to $0.5$ for most practical materials. Hence, the upper limit of $Et$ could be one third of E. Brizmer et al. [88] analyzed the effect of a $2%$ 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 [108] to improve the effort in Kogut and Etsion [93] to account for the strain hardening effect. Interestingly, Kogut and Etsion [93] 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 [108] proposed a new empirical formulation similar to the one developed by Kogut and Komvopoulos [80] for three different values of the tangent modulus: $Et=0.025E$⁠, $Et=0.05$E, and $Et=0.1E$⁠. They showed that the interference ratio at which the contact changes to fully plastic also changes based on the ratios $Sy/E$ and $Et/E$⁠. Later, Sahoo et al. [109] 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 (⁠$Et/E=0.02$⁠) 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) [48]. 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 [111] performed such analysis with a technique similar to that in the work of Chang et al. [92] on spherical contacts, which relied on a volume conservation approach for the elastic–plastic contact regime. Jeng and Wang [112] followed a similar methodology but with a template inspired by Zhao et al. [73]. Jamari and Schipper [113] also followed a similar approach, and their model was verified with the experimental data of Chaudhri et al. [81]. Lin and Lin [114] 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.

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. [117] 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 $Ac=Lλ2Sy/(5E′Δ)$⁠. 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. [117] 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 [116] 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. [118] 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 [119]. The experimentally observed phenomenon of asperity persistence [120] can also be explained by this mechanism.

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 < S_y < 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, $E/Sy$⁠. Figure 13 shows the effect of $E/Sy$ (not $E′/Sy$⁠) on the dimensionless contact pressure–area relation for $E=100$ GPa and $Δ/λ=0.0125$⁠. It is clear that as $E/Sy$ 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 $pep*$ 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. [127] 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 $3Sy$⁠, which is in contradiction with the results in Ref. [128]. To illustrate the possible increase in pressure with relatively little hardening, Eq. (33) is plotted in Fig. 15 as a function of $E′/Sy·Δ/λ$⁠, while normalizing the pressure by the yield strength. Here, typical properties for a steel alloy were assumed (⁠$E′=110$ GPa and $Sy=200$ MPa), and the amplitude over the wavelength was varied in the range of $0.00027<Δ/λ<0.1$⁠. 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 [129], 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.

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 [116].

Note that this p* for the geometry defined by Eq. (39) differs from that for other geometries (Eqs. (28) and (33c)).

where $C2$ 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 [116]. In their analysis, the range of material properties and geometries was somewhat limited (⁠$0.00261<Sy/E′<0.0102$ and $0.01<Δ/λ<0.04$⁠). 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 [72]. However, this trend does not occur in sinusoidal contacts. For the fully plastic case, Manners [4] 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.

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 [133] in 1938, Mindlin in 1949 [134], and Mindlin and Deresiewicz in 1953 [135]. 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. [136] followed Mindlin's approach, obtaining the region for sliding of elastic bodies in contact. Hamilton [137] found the yield inception by using Hertz contact pressure and the Mindlin shear stress distribution. Sackfield and Hills [138] modified the stress distribution by considering the effect of the shear stress on surface displacements for two dissimilar cylinders.

Bowden and Tabor [139] 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. [134]. 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 [140] 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 [141] 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. [142] 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. [146] used a similar methodology as Ref. [143] to account for larger deformations on the asperity level by incorporating the effects of the Jackson and Green model [52]. Their model considered the statistical rough surface contact in the static friction model.

where $μlocal$ 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. [150]. 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. [100]. 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 [93]. The rig was also used to confirm the theory of unloading of elastic–plastic spherical contact [106]. Using this test rig, measurements of friction force and contact area at the instant of sliding inception were reported in Ref. [151]. 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 [152]. 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. [146] 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 [153]. Junction growth was observed to increase the contact area by up to 45%. Based on the measurements in Ref. [151], Etsion [154] 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.

This section described original results from finite element analysis of the sinusoidal geometry described by Eq. (38) rather than the typically assumed spherical geometry [156]. 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 $τc$⁠. 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, $Ft/Fn$⁠, as a function of the dimensionless tangential displacement of the flat, $ux/ω0$⁠, for several different combinations of material properties, geometrical parameters, normal pressure, and critical shear strength. In each case, as the tangential loading progresses, $Ft/Fn$ increases gradually until a constant value of $Ftmax/Fn$ 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 $Sy$⁠. 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, $φ=(E′/Sy)(Δ/λ)$⁠, proposed by Gao et al. [157] 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. [149], the value of the critical shear stress $τc$ was set to $Sy/3$⁠, representing the shear strength of two identical metal surfaces bonded together. However, $τc$ often differs from this value as a function of temperature conditions, the presence of contaminants, and lubrication. The effect of $τc$ 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.

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 [24] 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 [158]. Their study neglected the recovery phase of the interaction and therefore could not be related to friction. Jackson et al. [159] 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. [160] studied the sliding and rolling interaction of two-dimensional asperities. Vijaywargiya and Green [161] concurrently studied the sliding interaction and energy loss occurring between cylinders in contact. Mulvihill et al. [162] 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. [163] 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. [164] 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.

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. [168]. 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 [169] 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 $w$ is the work of adhesion, $σ0$ is the maximum adhesive tensile stress, $p0$ is the average elastic–plastic pressure, $a0$ is the elastic–plastic contact radius, and $p0$ and $a0$ are both approximated by the hardness (corresponding to the fully plastic contact case). Generally, as $χ$ increases, more plasticity is expected, while $S$ 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 [172] followed in the footsteps of Chang et al. [168] 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 [150].

Sahoo and Banerjee [173] 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 [174]. 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. [175] 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 $S$ proposed by Mesarovic and Johnson [170] could be used to predict which pull-off mechanism would occur. In a follow-up work, Du et al. [176] 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 $S=1.2$⁠. Note that their calculations assumed that $p0$ was equal to the hardness.

Later, Kadin et al. [177] 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 $p0$ in the expressions for S, and that the well-known Tabor parameter [178] 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.

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 [197]). 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 [220].

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. [179] 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. [26] 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. [221] and Sainsot et al. [222] 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. [223] 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 [224].

Pei et al. [211] 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. [225] also used finite elements to consider the elastic–plastic contact of rough surfaces in a microelectromechanical system. Sahoo and Ghosh [226] 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. [229] attempted to accelerate the elastic–plastic deterministic solution time by using a multilevel approach [230]. Megalingam and Mayuram [231] 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 [139], but its provenance is uncertain.

where $P¯z-d$ is the contact force asperity model. For $A¯z-d$ and $P¯z-d$⁠, 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. [191].

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 [240]. 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 [245]. 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.

where A_r is the real area of contact, $A¯$ is the contact area of a single asperity on a certain scale of roughness, η is the real asperity density, A_n is the nominal contact area, and the subscript i denotes a specific asperity scale level, with i_max denoting the smallest scale level considered. Since each scale bears the same load, $Fi¯$⁠, the single asperity load at the $ith$ scale can be related to the total load by $F=F¯iηiAi-1$$F=F¯iηi Ai−1$⁠. 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 [157] 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.