Elastoplastic Indentation Response of Sigmoid/Power Functionally Graded Ceramics Structures

Due to the applicability of new advanced functionally graded materials (FGMs) in numerous tribological systems, this manuscript aims to present computational and empirical indentation models to investigate the elastoplastic response of FG substrate under an indention process with spherical rigid punch. The spatial variation of the ceramic volume fraction through the specimen thickness is portrayed using the power law and sigmoid functions. The effective properties of two-constituent FGM are evaluated by employing a modified Tamura–Tomota–Ozawa (TTO) model. Bilinear hardening behavior is considered in the analysis. The finite element procedure is developed to predict the contact pressure, horizontal displacement, vertical deformation, and permanent deformation of FG structure under the rigid cylindrical indentation. The empirical forms for permanent deformation were evaluated and assigned. Model validation with experimental works was considered. The convergence of the mesh and solution procedure was checked. Numerical studies were performed to illustrate the influence of gradation function, gradation index, and indentation parameters on the contact pressure, von Mises stresses, horizontal/vertical displacements, and permanent plastic deformation. The present model can help engineers and designers in the selection of an optimum gradation function and gradation index based on their applications.


Introduction
At the end of the 20th century, the problems of singular stress distribution, stress peeling effects, and delamination observed in laminates' structure were solved by Japanese's scientists [1], who developed functional graded materials (FGMs). FGMs have continuous variation, from ceramic/metal phase to metal/ceramic content. Therefore, their material properties can be tailorable through applications [2]. FGMs can be used in various applications such as aerospace, marine, automobile, biomedical, defense, nuclear, turbine and MEMS/NEMS [3].
Different elastoplastic models were proposed to describe the elastoplastic response of FGMs by many researchers. In 1993, Williamson et al. and Drake et al. studied the elastoplastic response of FG cylindrical specimen geometries during a cooling process by using the finite element method (FEM). Akis [4] investigated the yielding and elastoplastic deformation of FG spherical vessels under internal pressure by using Tresca's yield criteria and small deformation theory. Vaghefi et al. [5] implemented the von Mises yield criterion and isotropic strain hardening rule to consider the nonlinear elastoplastic behavior of FGM thick plate under the combination of mechanical and thermal loads using the meshless local Petrov-Galerkin method. Burzyński et al. [6] studied the elastoplastic nonlinear behavior of FGM shells under compression via Cosserat-type kinematics and FEM. Hasrati et al. [7] explored the elastoplastic postbuckling of thick plate by using Prandtl-Reuss and shear deformation Hencky constitutive relations. The stress-strain relation of the material was controlled by the Ramberg-Osgood elastoplastic model. Vakil and Zajkani [8] developed the lower order strain gradient model for analyzing the plastic behavior of FG crystalline microbeam structures. Vaghefi [9] studied the 3D thermo-elastoplastic bending response of FG skew plates subjected to combined thermal and mechanical loads by using von Mises isotropic hardening theory and the Prandtl-Reuss flow rule. Zhang et al. [10] developed the micromechanical model with experimental validation to predict the elastoplastic response of FGM manufactured by the vibration sedimentation process. Saeedi et al. [11] investigated the thermo-elastoplastic response of a thick FG cylindrical shell under internal pressure and temperature gradient by using the successive approximation method and differential quadrature method. Liu et al. [12] used a mesh-free reproducing kernel particle and penalty methods to study the elastoplastic behavior of the 2D structure with FGM composition.
Based on the Tamura-Tomota-Ozawa (TTO) elastoplastic model, Huang et al. [13,14] evaluated the buckling response of elastoplastic FGM cylindrical shells under torsional loads. Zhang et al. [15] studied the buckling behavior of elastoplastic FG cylindrical shells under combined compression and pressure. Nguyen et al. [16] examined the nonlinear bending of elastoplastic FGM beams under the influence of nonuniform distributed loads. Moita et al. [17] presented a nonlinear elastoplastic FG axisymmetric shell under thermal load based on a conical frustum FE model. Zhang et al. [18] investigated the thermal buckling of elastoplastic FGM Euler-Bernoulli beams under transversely nonuniform temperature rise. Zhu and Cai [19] developed an enhanced strain rate-dependent model to predict the elastoplastic response of FG composites under impact loads by using the Taylor dislocation model and Johnson-Cook model. All these studies demonstrated the efficiency of the TTO model to simulate the elastoplastic behavior of FGM.
Contact mechanics is presented in a wide variety of tribological systems whenever two bodies contact and roll over each other. Ma et al. [20] presented a normal contact stress, electric displacement and magnetic induction of 2D half-plane under a rigid flat/cylindrical punch by using the Fourier transform technique. Liu et al. [21] examined the contact pressure and deformation of FG substrate indented by rigid cylindrical/spherical punch by using the Hankel integral transform technique and the transfer matrix. Vasiliev [22] investigated the indentation of an electroelastic piezoelectric FG semi-infinite substrate by a rigid indenter under applied force and electric charge. Hou et al. [23] analytically investigated the stress gradient in the coated structure under tilted circular flat punch by using 3D elasticity theory. Asiri et al. [24] developed a predictive model to study the elastoplastic response during loading and unloading of nanocomposite structures indented by spherical punch. Chen and Yue [25] assessed the normal contact pressure between two elastic FGM spheres with different values of shear modulus and Poisson's ratio. Wagih et al. [26] developed a numerical FEM to study the elastoplastic indentation of FGM substrate under a frictionless contact condition by using the TTO model. Melaibari et al. [27] presented a comprehensive numerical model to predict the deformation and contact force generated by indentation testing on the orthotropic microplates. All these studies were performed over either an orthotropic or isotropic FG plate following the power law gradation. Despite the importance of power law function in the gradation of FGM plates, the sigmoid function showed better stress transfer and a reduction in stress concentration within the thickness of the plate.
Considering the discussed state of the art, it was found that the gradation sigmoid function and elastoplastic contact have not been considered before. Hence, this article aimed to present a comprehensive analysis on the indentation behavior of sigmoid FG elastoplastic material, considering the gradation of all material parameters. Moreover, an empirical relation was developed to predict the permanent deformation in FG elastoplastic materials after indentation by knowing only the gradient index of the material, which will assist the engineer and scientist in their design. This article is organized as follows: Section 2 illustrates the gradation functions, rule of mixture, modified TTO model, bilinear hardening model, and homogenization procedure. Description of the problem, a finite element model and discretization techniques are presented in Section 3. Problem validation and mesh convergence are presented and discussed in Section 4. The numerical parametric studies and empirical permanent indentation form are illustrated and studied in Section 5. Section 6 summarizes the main points and conclusions of the present study

Elastoplastic Functionally Graded Material
Consider an FGM specimen in the (x, y) coordinate system with width b and thickness h, as shown in Figure 1a. The material at the top and bottom surfaces of the specimen at y = ∓h/2 are made of pure metal (ductile phase) and pure ceramic (brittle phase), respectively. In this study, the spatial variation of the ceramic volume fraction, through the specimen thickness, is described using the power law (P-FGM) and sigmoid functions (S-FGM), respectively. The sigmoid-law distribution of FGM results in smooth stress distributions [28,29]. The volume fraction for S-FGM can be described by two power law functions [30]: For P-FGM, the volume fraction is given as [30]: where the subscripts "m" and "c" refer to the metal and ceramic constituents, respectively, and V i is the volume fraction of the constituting phase i. k is the volume fraction index describing the material distribution. The effective elastic properties of two-constituent FGM are evaluated employing the TTO model, i.e., modified rule of mixtures. The effective values of the modulus of elasticity E(y) and Poisson's ratio ν(y) can be expressed as [31]: where E i and ν i are, respectively, the Young's modulus and Poisson's ratio of the constituting material i. The stress transfer parameter q depends on the microstructure and material properties of FGM. To approximate the actual effects of microstructure interaction in FGMs, a nonzero finite value of q is used [32]. Following the TTO model, failure due to plastic deformation of a nonhomogeneous material having ductile and brittle phases is governed by the ductile constituent. During the plastic deformation of FGMs, the stress concentrations around the cracks and flaws of the ceramic phase are considerably reduced by the ductility and good shear strength of the metal phase and, thus, FGMs exhibit no brittle failure [32,33]. On the basis of this assumption, the effective yield strength S Y (y) and tangent modulus H(y) of an elastoplastic FGM are evaluated using the stress transfer parameter as follows: The overall elastoplastic behavior of FGMs is predicted by assuming a bilinear hardening behavior of the metal phase [33][34][35] (Figure 1b). The material properties of titaniumtitanium boride (Ti/TiB) FGM, in which Ti and TiB represent the ductile metal and hard ceramic, respectively, are as follows: E m = 107 GPa, E c = 3 75 GPa, v m = 0.34, v c = 0.14, S Ym = 450 MPa, H m = 14 GPa, and q = 4.5 GPa [13].
Variations in the effective modulus of elasticity, yield strength, and tangent modulus of the Ti/TiB FGM specimen based on the sigmoid and power law functions are illustrated in Figure 2. It is depicted that for P-FGM, isotropic homogeneous specimens of pure metal and pure ceramic materials are recovered by setting k = 0 and k → ∞, respectively. On the other hand, the S-FGM specimen behaves as an isotropic homogeneous specimen at k = 0, with properties between those of Ti and TiB constituents. When k → ∞, the specimen exhibits a bi-material phase with metal and ceramic phases at the lower and upper half portions of the specimen, respectively. Thus, the transformation of the S-FGM from an isotropic homogeneous material phase to a bi-material phase is controlled by varying the gradient index k, as depicted in Figure 2a,c,e.

Finite Element Modeling
A 2D axisymmetric FE model on ANSYS software was implemented to simulate the indentation test experiment on functional graded structures. The FG substrate was defined as a rectangle of 12 mm width and 5 mm thickness, as shown in Figure 3a. Due to the axisymmetry of the material and indentation tests, half of the substrate was simulated with 6 mm width and 5 mm thickness and the axisymmetric boundary conditions were applied to the model at the left side of the substrate and the indenter geometries, as shown in Figure 3b. The spherical indenter with 1 mm radius was defined as a rigid body to apply the indentation load. A frictionless contact pair was deployed to define the contact between the indenter and the specimen. The 2D structural plane elements (PLANE182) were used to mesh the substrate with graded mesh, where fine mesh was applied below the indenter and coarse mesh was applied to the rest of the substrate. The contact between the indenter and the substrate was defined as a surface-to-surface algorithm using contact (CONTA172) and target elements (TARGET169) for the contact pairs.
The augmented Lagrangian method was used as the contact algorithm. The elastic and plastic functional graded properties were applied using macro files where the constitutive Equations (1)-(6) were implemented including a loop to update the material properties of each element based on its position. The two gradation distributions with different gradation indices were implemented in the program using material code.
Fixed boundary conditions were defined for the nodes attached to the lower edge of the specimen. Displacement was applied to the pilot node of the indenter in the negative y-direction. The output of the simulation of the force-displacement response was recorded at the pilot node of the indenter. Four different meshes with different densities were deployed to determine the best mesh density to achieve accurate results with reasonable simulation time.

Validation and Convergence
The validation of the simulation procedure of loading and unload indentation test with experimental results obtained by Artian et al. [36] is presented as shown in Figure 4. As illustrated, the numerical results obtained by FE are very close to the experimental test and a good agreement has been proven. During loading, the contact force was increased in the linear form to the maximum value of 300 N and indentation of 90 µm. After that, unloading was observed, and the contact force decreased linearly to 80 N and 82 µm indentation depth. The slope of the unloading was higher than the slope of the loading test. During unloading, after contact force of 80 N and indentation depth of 82 µm, the nonlinear phenomenon was observed. The permanent indention was recorded at 65 µm. To demonstrate the accuracy and effectiveness of the finite element model, four different FE grids were constructed and compared. The proposed methodology was applied to investigate the FG elastoplastic indentation behaviors (i.e., indentation pressure, horizontal, vertical, permanent horizontal, and permanent vertical displacements) for an elastoplastic P-FGM Ti/TiB at gradation index k = 0.5, considering FE grids with different element size, as illustrated in Figures 5 and 6.
As shown in Figure 5, four meshes are considered in the analysis: mesh 1 (60/5462) for 60 elements in the vertical direction and total elements in the domain of 5462; mesh 2 (70/6372); mesh 3 (80/8882); and mesh 4 (85/10,892). The maximum contact pressures for meshes 1, 2, 3 and 4 are 9.9664 × 10 −3 MPa, 9.9836 × 10 −3 MPa, 9.9617 × 10 −3 MPa, and 9.5596 × 10 −3 MPa, respectively, and corresponding contact areas are 0.4375, 0.4464, 0.44 and 0.4412. Therefore, the results of coarse meshes (meshes 1 and 2) are very close to each other, but deviate by 4.5% from the minimum indentation pressure observed in the case of mesh 4. The effects of mesh size on the horizontal and vertical displacements are portrayed in Figure 6. It is concluded that by refining the mesh, the horizontal and vertical displacements are changed. The horizontal displacement changed from 6.0071 µm at mesh 1 to 6.1623 µm at mesh 2; the permanent horizontal displacement changed from 3.139 µm to  Figures 5 and 6, it can be proved that the response predicted from the present numerical simulation is very close to the experimental data, and the numerical results are convergent.

Numerical Results and Discussions
In this section, the applicability of the proposed numerical methodology to investigate the nonlinear indentation response of elastoplastic FG materials is demonstrated. The effects of material characteristics distribution throughout the thickness direction are investigated, considering both S-FGM and P-FGM. The mechanical indentation behavior is obtained for the two different distributions at different values of the material gradation index. The normalized penetration depth is defined as ∆ = U y U , with U y being the resulting vertical displacement while U is the applied prescribed vertical displacement.
Variations of the resultant indentation force versus the normalized penetration depth for both S-FGM and P-FGM at different material gradation index are illustrated in Figure 7. It is demonstrated that the indentation force required to attain a certain indentation penetration depth is governed by the material gradation index for both S-FGM and P-FGM distributions. In general, increasing the material gradation index decreases the material flexibility, which increases the required indentation force for a certain indentation depth. Additionally, increasing the material gradation index increases the ceramic content of FGM, which decreases the enclosed loop between loading and unloading paths, leading to elastic response at higher values of k (k = 8). It is observed that larger values of the resultant indentation force are required to attain a certain indentation penetration depth for P-FGM compared with the corresponding S-FGM distributions, especially at k = 2 and k = 4. Moreover, enclosed loops with a larger size between loading and unloading paths are observed for S-FGM compared with the corresponding cases for P-FGM due to its higher flexibility.
The indentation pressure profiles throughout the indentation interface for both S-FGM and P-FGM distributions for different material gradation indices are shown in Figure 8. It is revealed that the indentation pressure increases when increasing the material gradation index for both material distributions due to the increase in material stiffness. Additionally, for the same indentation depth, larger values of the resulting indentation pressure are produced with P-FGM compared with S-FGM. This effect is most significant at a material gradation index for k > 1.   The map of the residual von Mises equivalent stress after unloading within the elastoplastic Ti/TiB specimen at different gradient indices for both S-FGM, and P-FGM distributions is shown in Figure 11. It is depicted that the residual von Mises equivalent stress after unloading decreases when increasing the material gradation index due to the increase in the overall system stiffness for both material distributions. Additionally, larger values of the residual von Mises equivalent stress are detected for P-FGM at k = 0.5 compared to those obtained by S-FGM. For k = 2 and 4, the response is reversed, and larger values of the residual von Mises equivalent stress are induced in S-FGM compared with the corresponding cases of P-FGM.  Both normalized horizontal and vertical displacement profiles along the upper surface of the elastoplastic FG indented specimen at different gradient indices for both S-FGM and P-FGM distributions at full load are shown in Figures 12 and 13, respectively. It is demonstrated that the absolute value of the induced displacement increases with increase in the material gradient parameter for both material distributions. Due to the nature of the applied load and prescribed vertical displacement, the effect of the material gradient parameter on the absolute value of the vertical displacement was insignificant. Comparing the two material distributions, larger values of the induced displacements are produced by the P-FGM due to its lower compliance.  After unloading, the induced horizontal and vertical permanent displacements were significantly affected by the material distributions. Profiles for the normalized horizontal and vertical permanent displacements throughout the upper surface of the functionally graded indented specimen at different gradient indices for both S-FGM, and P-FGM distributions after unloading are illustrated in Figures 14 and 15, respectively. It is observed that for 0 < k ≤ 1, the absolute value of the normalized permanent horizontal displacement increases when increasing the material gradient index for both material distributions while, for k > 1, the absolute value of these permanent horizontal displacements decreases with increasing k. On the other hand, the absolute value of the normalized permanent vertical displacement increases with the increase in material gradient index.  Comparison of the two material distributions revealed that for the same vertical indentation depth, larger values of the absolute permanent displacements are detected by the S-FGM due to its higher compliance while, at higher values of the material gradient index, k = 8 complete recovery is attained by the P-FGM and reached zero permanent displacements.
The investigated mechanical behaviors of the elastoplastic FG indentation behavior for both S-FGM and P-FGM distributions at indentation depths 0.1 and 0.2 mm are summarized in Table 1. It is demonstrated that both the total and elastic strain energies increase with increases in material gradient index and indentation depth for both material distributions. Additionally, larger values of these strain energies are produced for P-FGM compared with the corresponding values obtained for S-FGM. On the other hand, the dissipation energy due to the plastic deformation is significantly affected by the material distributions. For k = 0.25, 0.5, 1.0, and 2.0, the dissipated energy increases with the increase in material gradient index and decreases for k = 4 and 8. Additionally, for all material gradient indices, the S-FGM distribution produces larger values for dissipative energy than those induced in the P-FGM distribution.   Based on the obtained numerical results, the following demonstrate empirical relations for the permanent depth with the material gradient index. where m and n are material-dependent parameters, which depend on the properties of the considered elastoplastic functionally graded material [37], incorporating the gradation of all the material properties. The dependency of the normalized maximum permanent penetration depth on the material gradient index is investigated using the proposed numerical procedure and the empirical formula (Equations (7) and (8)) for two different indentation depths: U = 0.1 and 0.2 mm. As shown in Figure 16, it can be observed that a good correlation is obtained between the predictions using the empirical formulas illustrated in Equations (7) and (8) with material parameters m = 0.92 and n = 0.6 and the FE results. Consequently, using the obtained empirical formulas and evaluating the material parameters m and n, the indentation depth of an elastoplastic functionally graded material with any gradation index could be predicted with no need for computational experimental efforts.
To illustrate the effectiveness of the obtained empirical formulas, the maximum and minimum errors between the investigated normalized maximum permanent penetration depth by the empirical equations and the corresponding results obtained by the FE are illustrated in Figure 17. It is seen that the maximum error, where m and n are material-dependent parameters, which depend on the properties of the considered elastoplastic functionally graded material [37], incorporating the gradation of all the material properties. The dependency of the normalized maximum permanent penetration depth on the material gradient index is investigated using the proposed numerical procedure and the empirical formula (Equations (7) and (8)) for two different indentation depths: U = 0.1 and 0.2 mm. As shown in Figure 16, it can be observed that a good correlation is obtained between the predictions using the empirical formulas illustrated in Equations (7) and (8) with material parameters m = 0.92 and n = 0.6 and the FE results. Consequently, using the obtained empirical formulas and evaluating the material parameters m and n, the indentation depth of an elastoplastic functionally graded material with any gradation index could be predicted with no need for computational experimental efforts. To illustrate the effectiveness of the obtained empirical formulas, the maximum and minimum errors between the investigated normalized maximum permanent penetration depth by the empirical equations and the corresponding results obtained by the FE are illustrated in Figure 17. It is seen that the maximum error, ℇ, is bound between −8 < ℇ < 4, validating the effectiveness of the obtained empirical formulas for both S-FGM and P-FGM distributions.
, is bound between −8 < where m and n are material-dependent parameters, which depend on the properties of the considered elastoplastic functionally graded material [37], incorporating the gradation of all the material properties. The dependency of the normalized maximum permanent penetration depth on the material gradient index is investigated using the proposed numerical procedure and the empirical formula (Equations (7) and (8)) for two different indentation depths: U = 0.1 and 0.2 mm. As shown in Figure 16, it can be observed that a good correlation is obtained between the predictions using the empirical formulas illustrated in Equations (7) and (8) with material parameters m = 0.92 and n = 0.6 and the FE results. Consequently, using the obtained empirical formulas and evaluating the material parameters m and n, the indentation depth of an elastoplastic functionally graded material with any gradation index could be predicted with no need for computational experimental efforts. To illustrate the effectiveness of the obtained empirical formulas, the maximum and minimum errors between the investigated normalized maximum permanent penetration depth by the empirical equations and the corresponding results obtained by the FE are illustrated in Figure 17. It is seen that the maximum error, ℇ, is bound between −8 < ℇ < 4, validating the effectiveness of the obtained empirical formulas for both S-FGM and P-FGM distributions.

Conclusions
The nonlinear mechanical behavior of elastoplastic FGM under spherical indentation is studied and analyzed considering two different material distributions. Both the power and sigmoid functions are exploited to express the spatial variation of the ceramic volume fraction throughout the indented FGM thickness direction. The modified Tamura-Tomota-Ozawa (TTO) model is adopted to describe the elastoplastic FGM distribution. A computational finite element procedure is exploited to investigate the nonlinear mechanical indentation behavior. The computational methodology is checked and verified. Numerical results are obtained and discussed. Based on the obtained numerical results, the dependency of the normalized maximum permanent indentation depth on the material gradient index is empirically obtained for the considered FGM. The results obtained revealed that: • S-FGM distribution results in more compliant material than P-FGM distributions. Thus, a smaller indentation force and a larger enclosed loop size between the loading and unloading path are predicted, especially for gradient index k > 2.

•
The von Mises equivalent stress at full was significantly affected by the material gradient index; it is increased by increasing the material gradient index for both S-FGM and P-FGM distributions. In general, the FGM characteristics' distribution throughout the thickness direction significantly affects the mechanical behavior of the elastoplastic material under spherical indentation. This mechanical behavior could be controlled by selecting the suitable material distribution law with a controlled material gradient index.

•
The study can be extended for different elastoplastic materials with a wide range of properties to generalize the developed empirical equation, so that the response of the general elastoplastic sigmoid FGM became well understood and predictable using these empirical relations.