Frictional Energy Dissipation in Partial Slip Contacts of Axisymmetric Power-Law Graded Elastic Solids under Oscillating Tangential Loads: Effect of the Geometry and the In-Depth Grading

: Due to the rapid development of additive manufacturing, a growing number of components in mechanical engineering are made of functionally graded materials. Compared to conventional materials, they exhibit improved properties in terms of strength, thermal, wear or corrosion resistance. However, because of the varying material properties, especially the type of in-depth grading of Young’s modulus, the solution of contact problems including the frequently encountered tangential fretting becomes signiﬁcantly more difﬁcult. The present work is intended to contribute to this context. The partial-slip contact of axisymmetric, power-law graded elastic solids under classical loading by a constant normal force and an oscillating tangential force is investigated both numerically and analytically. For this purpose, a ﬁctitious equivalent contact model in the mathematical space of the Abel transform is used since it simpliﬁes the solution procedure considerably without being an approximation. For different axisymmetric shaped solids and various elastic inhomogeneities (types of in-depth grading), the hysteresis loops are numerically generated and the corresponding dissipated frictional energies per cycle are determined. Moreover, a closed-form analytical solution for the dissipated energy is derived, which is applicable for a breadth class of axisymmetric shapes and elastic inhomogeneities. The famous solution of Mindlin et al. emerges as a special case. on the normalized tangential force as a function of the normalized tangential displacement in the case of a conical initial gap function.


Introduction
In numerous mechanical engineering applications, machine components are subjected to periodic oscillating loads, causing microslip within the contact interfaces. On the one hand, this cyclic microslip is detrimental since it is associated with the wear formation and may generate fretting fatigue cracks [1]. A classic example is turbine blades, which are mounted on rotating disks by dovetail or fir-tree type roots. Due to the aerodynamic and centrifugal forces as well as blade vibration, the blade-disk joints must withstand considerable loads and relative motion on parts of their contact interfaces is inevitable [2]. Press-fitted shaft hub joints and bolted joints are affected by fretting fatigue in a similar way [3].
On the other hand, the frictional energy dissipation due to partial or gross slip contact can be beneficial since it contributes to effective structural damping. Remaining in the field of turbomachinery applications, underplatform friction dampers are worth mentioning. They are mounted underneath the platforms of adjacent turbine blades with the purpose of limiting the blade vibration amplitude by frictional energy dissipation. In this way they significantly reduce the risk of a high-cycle fatigue failure [4].
Numerous theoretical studies on such elastic contacts subjected to oscillating loads have been carried out in the past. Usually, they are based on the half-space assumption. The solution of the classical contact problem of two parabolically curved bodies subjected to an oscillating tangential force at constant normal force originates from Mindlin et al. [5] and predicts cyclic microslip in a ring surrounding a central stick zone. Later, Johnson [6] verified the resulting fretting damage experimentally and showed that the extension of the damaged region in the contact area agreed well with the theoretically predicted one. The influence of an oscillating force with a fixed oblique line of action (relative to the surface) was treated in a seminal work by Mindlin and Deresiewicz [7]. However, of more practical importance are the loading histories involving both periodically varying tangential and normal forces. They have been studied in more recent works considering fairly general profiles (not just spheres), including rough surfaces [8,9]. The results comprise the evolution of stick-slip patterns and the frictional energy dissipation, which depends significantly on the relative phase between the oscillations in normal and tangential direction.
All of the above-mentioned works have in common that they presuppose elastically homogeneous isotropic solids. However, owing to the increased demands on the performance and service life of mechanical power transmission systems, more and more components with functionally graded material properties are being used. Such materials are characterized by a gradually changing composition or microstructure, according to a predefined law to meet a specific functionality. For example, gears must be tough enough in the interior to withstand fracture, whereas its surface must be hard and wear resistant. Likewise with turbine blades, except that a high heat resistance on the surface is additionally required. Due to the rapid development in additive manufacturing, the production of functionally gradient materials (FGMs) is no longer a challenge and will become less expensive over time [10]. Here, we focus on FGMs whose Young's modulus varies continuously perpendicular to the solid surface. Contact problems of such solids characterized by a predefined in-depth grading of Young's modulus can usually only be solved numerically with great computational effort [11]. The only exceptions are FGMs, whose Young's modulus changes with the depth according to a power law where c 0 denotes a characteristic depth. For these so-called power-law graded elastic solids, several closed-form analytical solutions to various types of contact problems have been derived in the past. Even though this law of in-depth grading is, strictly speaking, of a purely academic nature due to the infinitely large or vanishing Young's modulus at the surface, the basic contact behavior deviates slightly from that obtained for a more realistic piecewise defined law [12]. Significant contributions to the solution of normal contact problems of power-law graded elastic solids have been made by Booker et al. [13] as well as by Giannakopoulos and Suresh [14]. The solutions of normal contact problems with adhesion, according to Johnson, Kendall and Roberts and according to Maugis, respectively, were developed in later years (see especially Refs. [15][16][17]). A rigorous theory for the treatment of general three-dimensional tangential contacts of power-law graded elastic solids has been derived just very recently [18]. An elegant method for solving tangential contacts of axisymmetric power-law graded elastic solids was presented somewhat earlier by Heß and Popov [19,20]. In the present work, this method is applied to study the energy dissipation in partial slip tangential contacts of two axisymmetric power-law graded elastic solids. Figure 1 schematically shows this type of contact, which is subjected to a constant normal force and an oscillating tangential force specified by F z (t) = F N = const. , F x (t) = F A sin(ωt).
Mathematics 2022, 10, 3641 3 of 19 graded elastic solids. Figure 1 schematically shows this type of contact, which is subjected to a constant normal force and an oscillating tangential force specified by ( ) ( ) ( ) (2) Figure 1. Partial slip contact between two axisymmetric power-law graded elastic solids (tangential fretting mode).
The amplitude of the tangential force A F is assumed to be smaller than the maximum frictional force, in order to achieve the partial slip regime and to avoid gross slip. Due to the oscillating tangential loading, cyclic microslip is induced in an annulus surrounding a central stick zone, resulting in frictional energy dissipation which contributes to the effective structural damping of the dynamic systems. Here, we focus on investigating the influence of the elastic inhomogeneity (type of in-depth grading) and the axisymmetric shape of the contacting bodies on the hysteresis loops and thus the energy dissipation. We explicitly emphasize that our theoretical investigations are based on the same simplifying assumptions that were presupposed in the original theory of Mindlin et al. [5,7]. These include in particular: The bodies are assumed to be elastically similar half-spaces with ideally smooth surfaces, incomplete contact occurs, and the local contact obeys Amontons-Coulomb's friction law with a constant friction coefficient. Although the Cattaneo-Mindlin theory is widely accepted and its results are in good qualitative agreement with the experimental ones, its inherent limitations due to the simplified assumptions are evident. Already, the early experimental work of Johnson [21] as well as Goodman and Brown [22] revealed that in the range of small tangential forces, deviations can be observed. Johnson attributed the deviations to the inadequacy of the theoretical assumption of a constant coefficient of friction as a result of the breakdown of contaminant and oxide films by oscillating slip. The practical applicability of the further simplified assumptions, in particular that the energy dissipation is exclusively the result of interfacial slip, was questioned by Etsion [23]. Neither the experimentally observed junction growth [24], nor the plastic yield at the contact interface (or sub-surface) including its consequences in terms of damage can be accounted for within the framework of the applied theory [25,26]. A detailed study on the influence of the elastic mismatch and plasticity, as well as a phase-shifted alternating normal force, was conducted by Patil and Eriten [26]. Despite the aforementioned works that legitimately question the simplifying assumptions, the Cattaneo-Mindlin theory is one of the most cited, widely accepted, and still frequently applied contact theories to date. So, we make no apology applying the corresponding assumptions as well for investigating the influence of the material gradient on the dissipated energy in a partial slip contact of power-law graded elastic solids. For some aspects concerning the suppression of plastic yielding and fretting wear of such materials, the reader is referred to [27]. The amplitude of the tangential force F A is assumed to be smaller than the maximum frictional force, in order to achieve the partial slip regime and to avoid gross slip. Due to the oscillating tangential loading, cyclic microslip is induced in an annulus surrounding a central stick zone, resulting in frictional energy dissipation which contributes to the effective structural damping of the dynamic systems. Here, we focus on investigating the influence of the elastic inhomogeneity (type of in-depth grading) and the axisymmetric shape of the contacting bodies on the hysteresis loops and thus the energy dissipation.
We explicitly emphasize that our theoretical investigations are based on the same simplifying assumptions that were presupposed in the original theory of Mindlin et al. [5,7]. These include in particular: The bodies are assumed to be elastically similar half-spaces with ideally smooth surfaces, incomplete contact occurs, and the local contact obeys Amontons-Coulomb's friction law with a constant friction coefficient. Although the Cattaneo-Mindlin theory is widely accepted and its results are in good qualitative agreement with the experimental ones, its inherent limitations due to the simplified assumptions are evident. Already, the early experimental work of Johnson [21] as well as Goodman and Brown [22] revealed that in the range of small tangential forces, deviations can be observed. Johnson attributed the deviations to the inadequacy of the theoretical assumption of a constant coefficient of friction as a result of the breakdown of contaminant and oxide films by oscillating slip. The practical applicability of the further simplified assumptions, in particular that the energy dissipation is exclusively the result of interfacial slip, was questioned by Etsion [23]. Neither the experimentally observed junction growth [24], nor the plastic yield at the contact interface (or sub-surface) including its consequences in terms of damage can be accounted for within the framework of the applied theory [25,26]. A detailed study on the influence of the elastic mismatch and plasticity, as well as a phase-shifted alternating normal force, was conducted by Patil and Eriten [26]. Despite the aforementioned works that legitimately question the simplifying assumptions, the Cattaneo-Mindlin theory is one of the most cited, widely accepted, and still frequently applied contact theories to date. So, we make no apology applying the corresponding assumptions as well for investigating the influence of the material gradient on the dissipated energy in a partial slip contact of power-law graded elastic solids. For some aspects concerning the suppression of plastic yielding and fretting wear of such materials, the reader is referred to [27].
The paper is structured as follows. Section 2 introduces the applied methods, including the fundamentals of the method of dimensionality reduction for contact problems of powerlaw graded elastic materials. The normal contact and the tangential contact for those contact problems are introduced. Furthermore, the implementation with the method of dimensionality, especially for power-law graded materials, is thematized. Section 3 analyzes the results of the tangential contact of the contact problem. Hysteresis curves and especially the dissipated energy during one cycle of oscillation are investigated for different initial gap functions. Additionally, various materials with corresponding exponents of the elastic inhomogeneity k and characteristic depth c 0 are examined. Section 4 introduces the closed-form analytical solution for the dissipated energy during one cycle of oscillation. The solution is given for any choice of the exponent of the power function n and the exponent of the elastic inhomogeneity k. The last chapter contains the final discussion of the results and the conclusion.

Methods
A powerful mathematical tool for solving boundary value problems in the classical theory of elasticity (and other physical fields) is the integral transform [28]. Examples from contact mechanics are the Fourier transform for solving plane contact problems and the Hankel transform for solving axisymmetric problems. Here, the partial differential equations are first transferred to ordinary equations in the frequency domain by means of the integral transform. The ordinary differential equations are then to be solved with conventional methods including an adaptation to the also transferred boundary conditions. A subsequent back transform then provides a solution of the original boundary value problem.
The Abel transform is another transform for solving contact problems of axisymmetric solids. Whereas all integral transforms have in common that the calculation in the transformed domain requires typically simpler mathematical operations in comparison to those in the original domain, the Abel transform has a decisive advantage: The corresponding mathematics in the transformed domain can be linked to an extremely simple, equivalent one-dimensional contact problem [29]. The latter is directly used in the so-called method of dimensionality reduction (MDR) to solve axisymmetric contact problems of elastically homogeneous solids [30,31]. Beneficial is that the relationships between global quantities, such as normal force, indentation depth and contact radius can be calculated from the equivalent one-dimensional problem, hence no back transform is required. If one is also interested in the local quantities, such as the pressure distribution within the contact area or the surface displacement outside, the inverse transforms must be applied. No information is lost. Taking into account a slightly modified Abel transform, the method is also applicable to the exact mapping of contact problems between axisymmetric power-law graded elastic solids [16,19].

Fundamentals of the MDR for Mapping Contact Problems of Power-Law Graded Elastic Solids
In the following, the application of the MDR will be briefly explained. Figure 2 serves this purpose, contrasting the contact configurations of the original axisymmetric contact problem and the equivalent one-dimensional contact problem. On the left (Figure 2a), the initial contact of two axisymmetric, power-law graded elastic solids is illustrated, where the bodies just contact each other in a single point and no forces are acting. Young's modulus of both bodies changes with the depth according to the in-depth law given by Equation (1), where equal exponents k of the elastic inhomogeneity are assumed. The initial gap function is specified by f (r). On the right-hand side (Figure 2b), however, the extremely simple equivalent contact problem is shown, which physically interprets the mathematics in the Abel transformed domain. It consists of a rigid plane indenter with profile g(x) and a so-called one-dimensional Winkler foundation which can be regarded as independent, linear elastic springs arranged along a line having both normal and tangential stiffness. If both systems are now loaded with a normal force, for example, both lead to exactly the same relationships between the global quantities normal force, contact radius and indentation depth. The same holds for the tangential loading. Due to its simplicity, solving the equivalent one-dimensional contact problem can be accomplished by anyone. It doesn't need a higher knowledge in contact mechanics. However, two preliminary steps are required [19] 2. Transform of the material behavior. The independent normal and tangential stiffnesses of the springs must be defined. Respectively, after the transition from discrete springs to a continuous one-dimensional Winkler foundation, the normal and tangential foundation modulus must be specified. The latter are given by [19].
The coefficients N h and T h , which depend on both Poisson's ratio  and the exponent k of the elastic inhomogeneity are listed in Appendix A. From Equations (4) and (5), it is evident that the elastic foundation moduli depend on the coordinate x . They vary with the lateral distance from the initial contact point to exactly the same power law as specified by the law of in-depth grading. However, we strongly emphasize that this correlation may not be generalized to other laws of in-depth grading. Nevertheless, there is no question that contact problems between solids of arbitrary elastically layered or gradient material can always be mapped to equivalent onedimensional ones [32]. The only difficulty lies in finding the associated mapping rules for If both systems are now loaded with a normal force, for example, both lead to exactly the same relationships between the global quantities normal force, contact radius and indentation depth. The same holds for the tangential loading. Due to its simplicity, solving the equivalent one-dimensional contact problem can be accomplished by anyone. It doesn't need a higher knowledge in contact mechanics. However, two preliminary steps are required [19]: 1.
Profile transform. The profile g(x) of the plane rigid indenter must be calculated by the (modified) Abel transform of the initial gap function f (r) according to

2.
Transform of the material behavior. The independent normal and tangential stiffnesses of the springs must be defined. Respectively, after the transition from discrete springs to a continuous one-dimensional Winkler foundation, the normal and tangential foundation modulus must be specified. The latter are given by [19].
The coefficients h N and h T , which depend on both Poisson's ratio ν and the exponent k of the elastic inhomogeneity are listed in Appendix A.
From Equations (4) and (5), it is evident that the elastic foundation moduli depend on the coordinate x. They vary with the lateral distance from the initial contact point to exactly the same power law as specified by the law of in-depth grading. However, we strongly emphasize that this correlation may not be generalized to other laws of indepth grading. Nevertheless, there is no question that contact problems between solids of arbitrary elastically layered or gradient material can always be mapped to equivalent one-dimensional ones [32]. The only difficulty lies in finding the associated mapping rules for the initial gap function and the material properties. Except for a few cases, their calculation must be carried out in a numerical way.

Normal Contact of Equal Power-Law Graded Elastic Solids with a Power-Law Initial Gap Function
The stated objective of this work comprises the study of hysteresis and energy dissipation in tangential (mode I) fretting of power-law graded axisymmetric solids. For this purpose, the numerical computations and the analytical enhancement will be conducted later based on the equivalent one-dimensional contact problem. Equal power-law graded elastic materials are assumed, hence the tangential and normal contact are uncoupled. Since the normal force is kept constant during the oscillating tangential loading, let us solve the normal contact analytically, separately in advance. In order to represent a wide class of axisymmetric shaped solids, an initial gap function in the form of a power law is assumed Application of the modified Abel transform according to Equation (3) leads to the rigid plane indenter profile of the equivalent one-dimensional model It results from a vertical stretch of the initial gap function by the factor where B(x, y) denotes the complete beta function. By pressing the indenter with the profile according to Equation (7) into the elastic foundation, only those springs undergo a displacement that come into direct contact with the indenter. Consequently, the displacement of the elastic foundation is simply given by where the indentation depth δ z represents the maximum displacement located at the center of contact (x = 0), hence δ z (a) = κ(n, k)A n a n .
In addition to the vertical displacement of the Winkler foundation, a linear vertical force density can be defined representing the distribution of the vertical spring forces per unit length and due to equilibrium, the sum of all spring forces must be equal to the externally applied normal force Taking into account Equations (9)-(11) as well as the normal foundation modulus according to Equation (4), Equation (12) yields It should be emphasized once again that the indentation depth and the normal force as a function of the contact radius, according to Equations (10) and (13), agree exactly with those from the original problem, i.e., the axisymmetric contact of two power-law graded elastic solids with a power-law gap function in the initial undeformed state.

Contact of Equal Power-Law Graded Elastic Solids under a Constant Normal Force and a Monotonically Increasing Tangential Force
The equivalent 1D contact problem for the case where a constant normal force is followed by a monotonically increasing tangential force is illustrated in Figure 3a. Amontons-Coulomb's friction law is applied to each individual spring which means that as long as the tangential spring force is smaller than the normal spring force, multiplied by the coefficient of friction, the spring is in a state of stick, otherwise it slips. Figure 3b exemplarily shows a sticking spring, which therefore undergoes the whole tangential rigid body translation.
It should be emphasized once again that the indentation depth and the normal force as a function of the contact radius, according to Equations (10) and (13), agree exactly with those from the original problem, i.e., the axisymmetric contact of two power-law graded elastic solids with a power-law gap function in the initial undeformed state.

Contact of Equal Power-Law Graded Elastic Solids under a Constant Normal Force and a Monotonically Increasing Tangential Force
The equivalent 1D contact problem for the case where a constant normal force is followed by a monotonically increasing tangential force is illustrated in Figure 3a. Amontons-Coulomb's friction law is applied to each individual spring which means that as long as the tangential spring force is smaller than the normal spring force, multiplied by the coefficient of friction, the spring is in a state of stick, otherwise it slips. Figure 3b exemplarily shows a sticking spring, which therefore undergoes the whole tangential rigid body translation. If we introduce a tangential line load analogous to the linear vertical force density, the following applies [19].
Note that the tangential displacement of a spring at position x is directly proportional to the tangential line load The relationship between the tangential rigid body displacement and the stick radius can be calculated from the requirement that the tangential line load must be continuous at the border of the stick domain.
The externally applied tangential force must in turn be equal to the sum of all the tangential spring forces It should be emphasized once again, that by simply applying the Amontons-Coulomb's friction law to each individual spring of the equivalent 1D model, it is ensured that the relationships between the tangential force x F , the tangential rigid body displacement If we introduce a tangential line load analogous to the linear vertical force density, the following applies [19].
Note that the tangential displacement of a spring at position x is directly proportional to the tangential line load q x (x) = c T (x)u 1D (x). The relationship between the tangential rigid body displacement and the stick radius can be calculated from the requirement that the tangential line load must be continuous at the border of the stick domain.
The externally applied tangential force must in turn be equal to the sum of all the tangential spring forces It should be emphasized once again, that by simply applying the Amontons-Coulomb's friction law to each individual spring of the equivalent 1D model, it is ensured that the relationships between the tangential force F x , the tangential rigid body displacement δ x , as well as the stick radius c, are exactly those being present in the original tangential contact between the two axisymmetric, power-law graded elastic solids. For the tangential contact of two power-law graded elastic solids with the power-law initial gap function defined in the last section, the application of Equations (14)- (16) leads to the following results Let us discuss these results considering a simple conical initial gap (n = 1). Figure 4 illustrates, in this case, the dependence of the normalized tangential force on the normalized tangential displacement for three different materials: An elastically homogeneous material (k = 0) and two power-law graded solids, one with a Young's modulus that decreases with depth (k = −0.5) and another one with an increasing Young's modulus (k = 0.5). The same normal force was assumed in all cases, resulting in different contact radii. The tangential displacement was normalized to the tangential displacement for homogeneous material at the onset of gross slip, δ max x,0 . The corresponding contact radius is denoted by a 0 . The characteristic depth c 0 as a parameter in the material law according to Equation (1) was set equal to this contact radius. Note that the Young's modulus at a depth of z = c 0 is equal to E 0 irrespective of which value the exponent of elastic inhomogeneity k has.
Let us discuss these results considering a simple conical initial gap = ( 1) n . Figure 4 illustrates, in this case, the dependence of the normalized tangential force on the normal ized tangential displacement for three different materials: An elastically homogeneous material = ( 0) k and two power-law graded solids, one with a Young's modulus that de creases with depth =− ( 0.5) k and another one with an increasing Young's modulus = ( 0.5) k . The same normal force was assumed in all cases, resulting in different contac radii. The tangential displacement was normalized to the tangential displacement for ho mogeneous material at the onset of gross slip, max  Two main results emerge from Figure 4. First, the tangential contact stiffness de creases with anincreasing exponent of elastic inhomogeneity. This means that the materia characterized by a soft surface and a stiff core = ( 0.5) k reacts more compliant than tha characterized by a stiff surface and a soft core =− ( 0.5) k . Furthermore, the onset of gros slip is reached for positive k at larger tangential displacements and for negative k a smaller ones, compared to the homogeneous material. These trends hold for other initia gap functions (e.g., parabolic) as well and are similar in behavior to a substrate with sof or stiff coatings, respectively [33]. Two main results emerge from Figure 4. First, the tangential contact stiffness decreases with anincreasing exponent of elastic inhomogeneity. This means that the material characterized by a soft surface and a stiff core (k = 0.5) reacts more compliant than that characterized by a stiff surface and a soft core (k = −0.5). Furthermore, the onset of gross slip is reached for positive k at larger tangential displacements and for negative k at smaller ones, compared to the homogeneous material. These trends hold for other initial gap functions (e.g., parabolic) as well and are similar in behavior to a substrate with soft or stiff coatings, respectively [33]. Figure 5 demonstrates the dependency of the characteristic depth on the normalized tangential force as a function of the normalized tangential displacement. It indicates that a power-law graded material whose elastic modulus increases with depth (k = 0.5) behaves softer when the characteristic depth increases (Figure 5a). For a power-law graded material with a stiff surface but whose Young's modulus decreases with depth (k = −0.5), the opposite effect is observed (Figure 5b). Figure 5 demonstrates the dependency of the characteristic depth on the normalized tangential force as a function of the normalized tangential displacement. It indicates that a power-law graded material whose elastic modulus increases with depth = ( 0.5) k behaves softer when the characteristic depth increases (Figure 5a). For a power-law graded material with a stiff surface but whose Young's modulus decreases with depth =− ( 0.5) k , the opposite effect is observed (Figure 5b).

Numerical Implementation of the MDR for Studying the Tangential Contact of Two Power-Law Graded Elastic Solids under General Loading Scenarios
Due to the enormous reduction of degrees of freedom associated with the mapping of the original problem into an equivalent problem in the Abel space, the method is particularly suitable for fast numerical calculations of contact interfaces. For a direct calculation with the equivalent model, only simple linear spring laws have to be implemented and Amontons-Coulomb's law has to be applied locally for each spring. However, in this case, a displacement-controlled loading procedure is preferred to a force-controlled one and in order to map arbitrary loading scenarios, an incremental version of the stick-slip rule is defined. For a small displacement Δ x δ of the plane rigid indenter, then holds By tracking the incremental difference of the plane rigid indenter position, the tangential displacements of all springs within the contact area can be uniquely determined, thus yielding the values of all tangential spring forces, whose sum corresponds to the total transmitted tangential force. Even for the calculation of the dissipated frictional energy, the equivalent one-dimensional model is convenient. No back transformation into the original space is necessary for this! The tangential force density and the local slip in the equivalent one-dimensional model are used to calculate the dissipated energy. The latter is defined as follows

Numerical Implementation of the MDR for Studying the Tangential Contact of Two Power-Law Graded Elastic Solids under General Loading Scenarios
Due to the enormous reduction of degrees of freedom associated with the mapping of the original problem into an equivalent problem in the Abel space, the method is particularly suitable for fast numerical calculations of contact interfaces. For a direct calculation with the equivalent model, only simple linear spring laws have to be implemented and Amontons-Coulomb's law has to be applied locally for each spring. However, in this case, a displacement-controlled loading procedure is preferred to a force-controlled one and in order to map arbitrary loading scenarios, an incremental version of the stick-slip rule is defined. For a small displacement ∆δ x of the plane rigid indenter, then holds , in a state of slip (19) By tracking the incremental difference of the plane rigid indenter position, the tangential displacements of all springs within the contact area can be uniquely determined, thus yielding the values of all tangential spring forces, whose sum corresponds to the total transmitted tangential force. Even for the calculation of the dissipated frictional energy, the equivalent one-dimensional model is convenient. No back transformation into the original space is necessary for this! The tangential force density and the local slip in the equivalent one-dimensional model are used to calculate the dissipated energy. The latter is defined as follows and especially its incremental change ∆s 1D does not depend on the coordinate. Using the MDR for the calculation of the dissipated energy becomes extremely simple. In this way, Hanisch et al. [34] numerically investigated the energy dissipation in the tangential contact of parabolic homogeneous elastic solids subjected to oscillations in normal and tangential directions with the same frequency, but a phase shift. There is just one method that is similarly simple to the MDR and is excellent for analyzing tangential contacts under arbitrary loading histories. This is the method of memory diagrams (MMD) developed by Aleshin et al. [35]. It is based on replacing the complex traction distribution inside the contact area by a simple internal function containing the same memory information. In a more recent work [36], the MMD has applied to calculate the frictional energy dissipation in the contact of two convex axisymmetric solids subject to arbitrarily varying oblique loading. However, to date the MMD has only been designed for tangential contacts of elastically homogeneous solids. According to the information provided in this section, the MDR has been implemented numerically to study the tangential contact of axisymmetric, power-law graded elastic solids under classical conditions of a constant normal force and an oscillating tangential force. In the following the numerical results are presented.

Numerical Results for Tangential (Mode I) Fretting of Two Axisymmetric Power-Law Graded Elastic Solids
In the following chapter, the numerical results for a tangential (mode I) fretting of two axisymmetric power-law graded elastic solids are investigated. The normal force was kept constant and additionally an oscillating tangential force was applied, as illustrated in Figure 1. The previous chapter describes the numerical implementation of the MDR for such a case.
The object of the study are the hysteresis curves and especially the dissipated energy during one cycle of oscillation. They are examined for different initial gap functions and corresponding indenter forms, various exponents of the elastic inhomogeneity k and quotients of the characteristic depth c 0 /a 0 .
The energy dissipation during a complete cycle |∆W| shall be investigated numerically. Two approaches for the calculation can be considered. The energy dissipation in one period of oscillation T is computed with the numerical integration computing the enclosed area of the hysteresis curve during one cycle. The second approach for the calculation of the energy dissipation is explained in more detail in Section 4. Energy dissipation due to friction appears only in the slip region Figure 6a,b show two sample graphs of the energy dissipation W(t), which reaches its maximum W max at the maximal tangential displacement δ x,max . Independent of the material, all curves show the qualitatively same behavior for various exponents of the elastic inhomogeneity k. The same applies to the indenter form, both curves for a conical and parabolic indenter show the qualitatively same behavior. It is noticeable that on the initial branch during the first tangential displacement up to δ x,max the curve is monotonically non-increasing, but not in the following cycles, because of the deflection of the springs. This is additionally documented in the Supplementary Materials. In the present study, both approaches to calculate the dissipated energy were implemented, which receive the same results.
Discussing the influence of specific parameter on the dissipation energy during one cycle of oscillation, the normal force F N was kept constant in all simulations to ensure comparability. The contact radius in the case of homogeneous material (k = 0) and a power-law initial gap function results from Equation (13) It is used for the normalization and to specify the characteristic depth c 0 by the characteristic quotient c 0 /a 0 . At first, it is chosen as c 0 /a 0 = 1 and varied later. With fixed values E 0 and c 0 /a 0 , the power law graded materials can be compared, varying the exponent of the elastic inhomogeneity k. The dissipated energy during one cycle is normalized with the quotient , which is composed of constants and depends only on the exponent of the power function .

Numerical Results for the Conical Indenter
Starting with the conical indenter, Figure 7a shows the hysteresis curve of the tangential force over the tangential displacement. The tangential displacement is normalized by the maximal tangential displacement of the homogenous case max ,0 x δ . The exponent of the elastic inhomogeneity is varied. The gradient of the hysteresis is smaller for positive and higher for negative exponents than the homogenous case. The enclosed area is larger for positive and smaller for negative exponents than the homogenous case. Since the enclosed area is equal to the dissipated energy, increasing the exponent k increases the dissipated energy. This is also depicted in Figure 7b, the higher the exponent , the higher the dissipated energy. It shows the behavior of the dimensionless dissipated energy in one cycle |∆ | over the traction ratio . The force is the maximal tangential force during one cycle. Increasing the traction ratio causes a power-law increase of the dissipated energy.

Numerical Results for the Conical Indenter
Starting with the conical indenter, Figure 7a shows the hysteresis curve of the tangential force over the tangential displacement. The tangential displacement δ x is normalized by the maximal tangential displacement of the homogenous case δ max x,0 . The exponent of the elastic inhomogeneity k is varied. The gradient of the hysteresis is smaller for positive and higher for negative exponents k than the homogenous case. The enclosed area is larger for positive and smaller for negative exponents k than the homogenous case. Since the enclosed area is equal to the dissipated energy, increasing the exponent k increases the dissipated energy. This is also depicted in Figure 7b, the higher the exponent k, the higher the dissipated energy. It shows the behavior of the dimensionless dissipated energy in one cycle |∆W| over the traction ratio F A µF N . The force F A is the maximal tangential force during one cycle. Increasing the traction ratio causes a power-law increase of the dissipated energy.
The Figure 8a-d demonstrate the influence of the characteristic depth c 0 a 0 . For a negative exponent of the elastic inhomogeneity k = −0.5, the gradient of the hysteresis increases with an increase of the characteristic depth, as shown in Figure 8a. Figure 8b illustrates that the dissipated energy |∆W| rises exponentially with the traction ratio F A µF N . The higher the characteristic depth, the smaller the dissipated energy. The Figure 8c,d show the same dependencies for a positive exponent of the elastic inhomogeneity k = 0.5. The influence of the characteristic depth is reversed. Increasing the characteristic depth, the gradient of the hysteresis decreases and the dissipated energy increases. (c) Tangential force over the tangential displacement for = 0.5; (d) Energy dissipation over the traction ratio for = 0.5.

Numerical Results for the Parabolic Indenter
Considering the parabolic indenter, Figure 9a shows the hysteresis curve of the tangential force over the tangential displacement, varying the exponent of the elastic inhomogeneity . The gradient of the hysteresis is smaller for positive and higher for negative exponents than the homogenous case. Figure 9b shows the behavior of the dimensionless dissipated energy in one cycle |∆ | over the traction ratio . Increasing the trac- Figure 8. Hysteresis curve and energy dissipation of a conical indenter for a fixed elastic inhomogeneity k and various quotients of the characteristic depth c 0 /a 0 : (a) Tangential force over the tangential displacement for k = −0.5; (b) Energy dissipation over the traction ratio F A µF N for k = −0.5; (c) Tangential force over the tangential displacement for k = 0.5; (d) Energy dissipation over the traction ratio F A µF N for k = 0.5.

Numerical Results for the Parabolic Indenter
Considering the parabolic indenter, Figure 9a shows the hysteresis curve of the tangential force over the tangential displacement, varying the exponent of the elastic inhomogeneity k. The gradient of the hysteresis is smaller for positive and higher for negative exponents k than the homogenous case. Figure 9b shows the behavior of the dimensionless dissipated energy in one cycle |∆W| over the traction ratio F A µF N . Increasing the traction ratio increases the dissipated energy. The higher the exponent k, the higher the dissipated energy.
(c) (d) Figure 8. Hysteresis curve and energy dissipation of a conical indenter for a fixed elastic inhomogeneity and various quotients of the characteristic depth 0 / 0 : (a) Tangential force over the tangential displacement for = −0.5; (b) Energy dissipation over the traction ratio for = −0.5; (c) Tangential force over the tangential displacement for = 0.5; (d) Energy dissipation over the traction ratio for = 0.5.

Numerical Results for the Parabolic Indenter
Considering the parabolic indenter, Figure 9a shows the hysteresis curve of the tangential force over the tangential displacement, varying the exponent of the elastic inhomogeneity . The gradient of the hysteresis is smaller for positive and higher for negative exponents than the homogenous case. Figure 9b shows the behavior of the dimensionless dissipated energy in one cycle |∆ | over the traction ratio . Increasing the traction ratio increases the dissipated energy. The higher the exponent , the higher the dissipated energy.
(a) (b) Figure 9. Hysteresis curve and energy dissipation of a parabolic indenter for various elastic inhomogeneities : (a)Tangential force over the tangential displacement; (b) Energy dissipation over the traction ratio .
The Figure 10a-d demonstrate the influence of the characteristic depth 0 / 0 , which has the same qualitative influence on the conical and the parabolic indenter. For a negative exponent of the elastic inhomogeneity = −0.5 the gradient of the hysteresis rises for The Figure 10a-d demonstrate the influence of the characteristic depth c 0 /a 0 , which has the same qualitative influence on the conical and the parabolic indenter. For a negative exponent of the elastic inhomogeneity k = −0.5 the gradient of the hysteresis rises for larger characteristic depths, as shown in Figure 10a. Figure 10b depicts an exponential growth of the dissipated energy |∆W| in dependence of the traction ratio F A µF N . Increasing the characteristic depth decreases the dissipated energy. The Figure 10c-d show the same dependencies for a positive exponent of the elastic inhomogeneity k = 0.5 with a reversed influence of the characteristic depth. Increasing the characteristic depth, the gradient of the hysteresis decreases, and the dissipated energy grows.  (c) Tangential force over tangential displacement for = 0.5; (d) Energy dissipation over the traction ratio for = 0.5.

Analytical Study on Frictional Energy Dissipation
Since the relationships between the externally applied forces, the relative rigid body displacements, and the stick and contact radius of the original problem, exactly emerge from the (fictional) equivalent contact problem in the Abel domain, it is possible to determine the energy dissipated per cycle from it, as well. One way is to evaluate the work carried out by the tangential force on the tangential rigid body displacement of the plane rigid indenter during a complete cycle. Another, perhaps less elaborate way, is to consider the energy dissipation over a complete cycle for each single spring contact and subsequently sum up the energy loss of all spring contacts. Figure 11 illustrates the variation of

Analytical Study on Frictional Energy Dissipation
Since the relationships between the externally applied forces, the relative rigid body displacements, and the stick and contact radius of the original problem, exactly emerge from the (fictional) equivalent contact problem in the Abel domain, it is possible to determine the energy dissipated per cycle from it, as well. One way is to evaluate the work carried out by the tangential force on the tangential rigid body displacement of the plane rigid indenter during a complete cycle. Another, perhaps less elaborate way, is to consider the energy dissipation over a complete cycle for each single spring contact and subsequently sum up the energy loss of all spring contacts. Figure 11 illustrates the variation of the spring force per unit length and the slip for a point of contact at position x during one cycle.   If we start from the initial, tangentially unloaded state, the spring will initially stick until the tangential line load has reached the limiting value µq z (x). Then, the spring begins to slip while the tangential line load remains constant. At the point when the loading phase of the cycle stops and the unloading starts, the spring instantaneously sticks. During the unloading, the spring remains in stick until the tangential line load reaches the limiting value −µq z (x). Then, the reverse slip begins at constant tangential line load, and so on.
Since the hysteresis loop represents a rectangle, the energy dissipated at each spring contact is simply given by ∆W s (x) = 4µq z (x)s max 1D (x) (24) and the total dissipated energy is obtained by summing up the energy losses of all spring contacts Using the definitions in Equations (9), (11), (15) and (20), Equation (25) can be rewritten as follows where we have already assumed power-law graded solids of equal material. According to Equation (26), apart from the normal stiffness of the springs, only the indenter profile of the equivalent one-dimensional contact problem is required to calculate the full energy loss.
In the following, we will calculate the dissipated energy per cycle for an axisymmetric contact with an initial gap function in the form of a power-law under tangential fretting conditions. The gap function is stated in Equation (6) and the rigid plane indenter profile of the equivalent one-dimensional model is given in Equation (7). After the insertion of the profile and considering the normal stiffness, according to Equation (4), the integral can be easily performed. The result is Herein, we have taken into account the normal force determined in Equation (13) and substituted the ratio of the stick radius to the contact radius according to Equation (18). F A denotes the amplitude of the oscillating tangential force. The partial slip regime is present as long as F A < µF N . Using the analytically determined formula for the loss of energy, both the influence of the geometry of the solids and the influence of the elastic inhomogeneity can be discussed. In the case of a parabolic initial gap function (n = 2) and elastically homogeneous material (k = 0) Equation (27) yields where h T (0, ν) = 2/[(1 + ν)(2 − ν)] and E 0 = 2(1 + ν)G 0 have been adopted. The frictional energy dissipation, according to Equation (28) is the famous result from Mindlin et al. [5].
For an arbitrary exponent of the elastic inhomogeneity but still a parabolic initial gap function (n = 2), Equation (27) agrees with a result, derived in [37] based on the Ciavarella-Jäger theorem applied to the original contact problem. This case will be exemplarily discussed below. The dissipated energy in the normalized form is illustrated in Figure 12.
The same normalization was adopted as in the chapter on numerical results, i.e., the same normal force was always applied irrespective of k. Furthermore, the contact radius a 0 that arises in the classical Hertzian contact was taken for normalization. It is evident that the curves for k = −0.5, 1, 0.5 agree well with those from the numerical calculation depicted in Figure 9b (computed, e.g., by the evaluation of the area of the hysteresis curves). Please note that the numerical results were plotted only up to F A = 0.99µF N , while the analytical function was plotted up to F A = µF N . The influence of the characteristic depth can be studied in an analogous way and other shapes than only the conical or parabolic one can be investigated with the help of the analytically derived energy dissipation, as well. νG have been adopted. The frictional energy dissipation, according to Equation (28) is the famous result from Mindlin et al. [5].
For an arbitrary exponent of the elastic inhomogeneity but still a parabolic initial gap function = ( 2) n , Equation (27) agrees with a result, derived in [37] based on the Ciavarella-Jäger theorem applied to the original contact problem. This case will be exemplarily discussed below. The dissipated energy in the normalized form is illustrated in Figure 12. The same normalization was adopted as in the chapter on numerical results, i.e., the same normal force was always applied irrespective of k. Furthermore, the contact radius 0 a that arises in the classical Hertzian contact was taken for normalization. It is evident that the curves for =−0.5, 1, 0.5 k agree well with those from the numerical calculation depicted in Figure 9b (computed, e.g., by the evaluation of the area of the hysteresis curves). Please note that the numerical results were plotted only up to

Discussion and Conclusions
Tangential fretting between two axisymmetric power-law graded elastic solids has been studied by using the simple equivalent one-dimensional contact model in the mathematical space of the Abel transform, both numerically and analytically. One critical question could be why the method was implemented numerically at all, although a closed-form analytical solution exists, which was subsequently derived. Well, on the one hand it should be demonstrated that the equivalent contact problem is ideally qualified as a module for the numerical computation of contact interfaces in dynamic systems, such as multibody systems. The transform of the contact problem is associated with an enormous reduction of degrees of freedom, which means a considerable saving of computation time. Since only linear spring laws need to be implemented, the contact forces can be determined in real time and transferred to the macrodynamic of the system at each time step [38,39].
On the other hand, we intend to show that the equivalent one-dimensional model is well suited to derive the closed-form analytical solutions of contact problems in a simple way. For this purpose, no deeper knowledge of contact mechanics is required. A central role of the method represents the mapping rules for normal contact problems [16,29]. This can be explained by the fact that most of the other contact problems, such as the tangential contact, normal contact with adhesion or rolling contact (even viscoelastic contacts) can be represented as superposition of normal contact problems.
Another critical point in this work certainly forms the assumed power-law of in-depth grading, which serves more as an academic example than a real material law due to the infinite or vanishing Young's modulus, at the surface. To go into this in more detail, let us consider again the elastic inhomogeneities mainly used in the numerical calculations depicted in Figure 13.
For a positive exponent, Young's modulus increases with depth, starting from an (infinitely) compliant surface to an (infinitely) stiff core. For negative exponents the inverse characteristic holds. At the characteristic depth z = c 0 , Young's modulus takes the value E 0 regardless of the exponent k. Our numerical and analytical studies concluded that the tangential contact stiffness decreases with increasing exponent of the elastic inhomogeneity. This means that the material characterized by a soft surface and a stiff core (k = 0.5, k = 0.7) reacts more compliant than a material characterized by a stiff surface and a soft core (k = −0.5, k = −0.7). In addition, the dissipated energy per cycle increases by increasing the exponent of the elastic inhomogeneity k. Wang et al. [11] studied tangential fretting of the contact between an elastic sphere and a substrate of the same elastic material but coated with a functionally graded layer. Depending on the material gradient of the coating, the same qualitative results were obtained with respect to the tangential stiffness and the dissipated energy. Moreover, the influence of the coating thickness was investigated and in fact it correlates completely with the influence of the characteristic depth of the power-law graded elastic material. Increasing the characteristic depth c 0 , the dissipated energy decreases and the tangential stiffness increases for negative exponents of the elastic inhomogeneity k. For positive exponents of the elastic inhomogeneity k, both trends are reversed. mined in real time and transferred to the macrodynamic of the system at each time step [38,39].
On the other hand, we intend to show that the equivalent one-dimensional model is well suited to derive the closed-form analytical solutions of contact problems in a simple way. For this purpose, no deeper knowledge of contact mechanics is required. A central role of the method represents the mapping rules for normal contact problems [16,29]. This can be explained by the fact that most of the other contact problems, such as the tangential contact, normal contact with adhesion or rolling contact (even viscoelastic contacts) can be represented as superposition of normal contact problems.
Another critical point in this work certainly forms the assumed power-law of indepth grading, which serves more as an academic example than a real material law due to the infinite or vanishing Young's modulus, at the surface. To go into this in more detail, let us consider again the elastic inhomogeneities mainly used in the numerical calculations depicted in Figure 13. For a positive exponent, Young's modulus increases with depth, starting from an (infinitely) compliant surface to an (infinitely) stiff core. For negative exponents the inverse characteristic holds. At the characteristic depth = 0 zc , Young's modulus takes the value 0 E regardless of the exponent k . Our numerical and analytical studies concluded that the tangential contact stiffness decreases with increasing exponent of the elastic inhomogeneity. This means that the material characterized by a soft surface and a stiff core ( = 0.5, = 0.7) reacts more compliant than a material characterized by a stiff surface and a soft core ( = −0.5, = −0.7). In addition, the dissipated energy per cycle increases by increasing the exponent of the elastic inhomogeneity k. Wang et al. [11] studied tangential fretting of the contact between an elastic sphere and a substrate of the same elastic material It is obvious that the spring forces do not coincide with the normal and tangential stresses of the original contact problem. The same applies to the displacements and the relative slip. However, these can be determined by the back transform from the image into the original space. At this point, it should be emphasized that the relative displacements are only accurate in the sense of the Cattaneo-Mindlin approximation, hence the model is not suitable for capturing the surface density of the energy dissipation.
Supplementary Materials: The Supplementary Materials contain videos for illustration purpose. Those videos show the values of the tangential displacement of all springs in contact and the stick and slip regions during oscillation. The simulation videos show the contact of a parabolic indenter for various exponents of the elastic inhomogeneity k (k = −0.5, k = 0, k = 0.5) for the traction ratio F A µF N = 0.9 and the quotient of the characteristic depth c 0 a 0 = 1. The hysteresis curve and the energy dissipation are illustrated. The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/math10193641/s1, Video S1: Tangen-tial_fretting_negative_exponent, Video S2: Tangential_fretting_homogenous_material, Video S3: Tangential_fretting_positive_exponent.

Conflicts of Interest:
The authors declare no conflict of interest.