Numerical Modeling of the Effect of Randomly Distributed Inclusions on Fretting Fatigue-Induced Stress in Metals

Abstract

The analysis of fretting fatigue plays an important role in many engineering fields. The presence of heterogeneity may affect the performance of a machine or a structure, including its lifetime and stability. In this paper, the effect of randomly distributed micro inclusions on the fretting fatigue behaviour of heterogeneous materials is analysed using the finite element method (FEM) for different sizes, shape and properties of inclusions. The effect of micro inclusions on macroscopic material properties is also considered by representative volume element (RVE). It is shown that the influence of micro inclusions on macroscopic material properties cannot be ignored, and the shape and size of the inclusions have less effect on the macroscopic material properties as compared to the material properties of inclusion and volume ratio. In addition, various parameters of inclusions have little effect on the peak tensile stress, which remains almost the same as homogeneous material. Peak shear stress occurs at many places inside the specimen, which can result in multiple cracking points inside the specimen, as well as at the contact surface. Moreover, the stress band formed by the stress coupling between adjacent inclusions may have an important influence on the direction of crack growth.

1. Introduction

Fretting fatigue can significantly affect structural performance in many engineering applications [1]. The fatigue life of the structure may be reduced by the fretting phenomenon, which occurs due to small oscillatory motion between two contact surfaces [2]. The reduction in fatigue life can reach up to 50% [3]. The failure process is generally characterized by two phases, crack nucleation [4,5] and crack propagation [6,7]. Researchers are concerned with contact stresses because it directly affects the initiation of cracks [8]. The contact stresses are affected by various factors, such as loading magnitude, contact geometry, surface imperfections and slip amplitude, which consequently affect the failure process. On the other hand, material properties also have significant influence on contact stresses. By common consensus, heterogeneity of material has a great influence on the life and stability of the structure [9]. For many engineering components, the required design lifetime can be affected significantly by heterogeneity. Therefore, material heterogeneity has been widely studied in recent years [10,11,12,13,14]. Generally, materials contain fibers [15,16,17,18,19], particles [20,21], precipitates, porosities, or voids/cracks [14,22] at the micro level, which are common causes of the heterogeneity. Generally, scanning electron microscopy (SEM) is used to observe the microstructure of material [23]. There are several kinds of typical inclusions in alloys; Al₂O₃, MgO, Al₂MgO₄, CaSO_4, TiB₂, Al₃Ti and refractory brick (Al, Si, O). Typical geometric form of inclusions are particles, films or group of films and rods [24]. As a result, stress concentration will appear at the interface between different materials or voids, which may give rise to shorter fatigue life. Thus, in addition to contact area, the subsurface area can also be affected by the heterogeneity, and in some cases it is strongly affected [25]. There are many ways to analyze the behavior of heterogeneous materials. Previous researchers have proposed several models to predict the effect of inclusions and defects on metal fatigue strength. Murakami and Endo [25] proposed an engineering guide to predict the fatigue strength of components with heterogeneity. Some experimental results about high-strength steels showed that fatigue failures were mostly caused by the inclusions inside the matrix [26,27,28,29,30,31].

In short, both fretting fatigue and heterogeneity may significantly influence the lifetime and stability of mechanical components. Numerous factors can affect micro stress field, thus affecting the fretting fatigue behavior [8]. The macroscopic fatigue failure occurs due to the distribution of micro-stress. In the literature, homogeneous materials are widely assumed to study fretting fatigue problems [14]. Generally, the micro cracks are observed inside the slip zone and at the contact edge. For experimental methods, it is difficult to get the details of stress field and initiation of crack. Therefore, numerical modeling is an efficient way to solve the problem of fretting fatigue [8,32,33,34,35,36,37,38] and fracture [39,40,41,42,43]. However, only few studies have taken into account the heterogeneity of material under fretting fatigue conditions [14,44,45]. Kumar et al. used numerical analysis to study the influence of heterogeneity on stress distribution in fretting fatigue problems [14]. They studied heterogeneity by considering micro-voids in the material, which is found in metal alloys [46]. They found that the effect of heterogeneity on shear stress was greater than on normal stress. Normally, for homogeneous material, the peak shear stress appears at the contact interface. However, in the case of heterogeneous materials, the peak may shift to the micro-voids. In the practical situation, the materials always have some heterogeneity due to heat treatment process. The inclusions can be regarded as foreign particles, which exist in metals [47]. According to the relative position of the inclusions and the material surface, there can be surface inclusions and internal inclusions [48].

In this study, we focus on two kinds of common inclusions in aluminum alloy 2024-T3. Considering the effect of micro inclusions on macroscopic material properties by means of representative volume element (RVE), a numerical analysis is presented to study the influence of local randomly distributed inclusions on the stress distribution of fretting fatigue specimens. The influence of different variables, namely material properties, volume ratio and shape of the inclusions is investigated. This study presents the variation of different stresses at the contact interface, as well as under the surface, and thus allows predicting probable crack initiation sites. From the result of stress analysis, conclusions are drawn.

2. Line Contact under Partial Slip

Before we study the influence of inclusions on fretting fatigue behavior in heterogeneous material, a brief review of Hertzian contact is presented here, which also can be used to verify the finite element (FE) model for homogenous material. The contact is said to be line contact if it is between a cylinder and a flat body (or two parallel cylinders). We have an analytical solution [49] to get the normal stress of the contact surface, if there is only a normal force F to compress the two bodies:

$$p(x)=-p_{\max }\sqrt{1-{(\frac{x}{a})}^2}$$

(1)

where $p\left(x\right)$ represents the normal stress distribution, $a$ is the semi contact width and $p_{max}$ the peak value of normal contact stress at the centre of the contact zone, and is given by [49]:

$$p_{\max }=\sqrt{\frac{F{E}^{\ast }}{t\pi R^{\ast }}}$$

(2)

where t is the thickness of both cylinder and the flat specimen and $R^{\ast }$ is the equivalent radius and $E^{\ast }$ is the equivalent modulus of elasticity, i.e.,:

$$\frac{1}{R^{\ast }}=\frac{1}{R_1}+\frac{1}{R_2}$$

(3)

$$\frac{1}{E^{\ast }}=\frac{1-v_1^2}{E_1}+\frac{1-v_2^2}{E_2}$$

(4)

$R_1$ is the radius of the pad and $R_2$ is the radius of the specimen, which is flat for our configuration. $E_1$ and $E_2$ are Young’s moduli, and $v_1$ and $v_2$ are Poisson’s ratios for the two bodies, respectively. The parameter $a$ is given by:

$$a=2\sqrt{\frac{F{R}^{\ast }}{t\pi E^{\ast }}}$$

(5)

Figure 1 shows the normal stress distribution of the contact interface of the specimen. On the whole contact surface (between ‒1 ≤ x/a ≤ 1), it is a parabolic distribution.

To illustrate the effects of tangential loading $Q$ a modified formulation is required [50,51]. Near the contact edge, $c\le \left|x\right|<a$, this is the slip zone and there will be relative sliding between the contact interfaces. Therefore, we can express the shear traction as $\mu p\left(x\right)$. Whereas, near the contact centre, $\left|x\right|<c$, there is a centrally symmetrical stick zone (indicated as stick zone 1 in Figure 2), which means that the contact surfaces will move together. So that in this zone the shear traction $q\left(x\right)$ does not reach the critical value $\mu p\left(x\right)$.

The contact shear traction can be modeled by the Coulomb friction law, which is given by:

$$q(x)=q_1(x)+q_2(x)$$

(6)

where $q_1(x)=\mu p(x)$. $q_1$ represents shear traction due to slip. $q_2$ is added as a perturbation in the solution of slip condition, to represent shear traction in the stick zone.

$$q_1(x)=-\mu p_{\max }\sqrt{1-{\left(\frac{x}{a}\right)}^2}$$

(7)

Moreover, in the slip zones, the perturbation $q_2(x)$ is zero.

$$q_2(x)=0,c\le \left|x\right|\le a$$

(8)

However, for the stick zones, it is given by:

$$q_2(x)=\mu p_{\max }\frac{c}{a}\sqrt{1-{\left(\frac{x}{c}\right)}^2},\left|x\right|\le c$$

(9)

where:

$$\frac{c}{a}=\sqrt{1-\frac{Q}{\mu F}}$$

(10)

Finally, when there are normal loads $F$ and tangential loads $Q$ simultaneously, the shear stress distribution on the contact surface is expressed as:

$$q(x)=\{\begin{array}{c}-\mu p_{\max }\sqrt{1-{\left(\frac{x}{a}\right)}^2},c<\left|x\right|<a\\ -\mu p_{\max }\left[\sqrt{1-{\left(\frac{x}{a}\right)}^2}-\frac{c}{a}\sqrt{1-{\left(\frac{x}{c}\right)}^2}\right],\left|x\right|\le c\end{array}$$

(11)

Furthermore, there also may be bulk stresses ${\sigma }_{axial}$ acting on the body, thus the perturbation shear traction $q_2\left(x\right)$ is given by [52]:

$$q_2(x)=\mu p_{\max }\frac{c}{a}\sqrt{1-{\left(\frac{x+\mathrm{e}}{c}\right)}^2},\left|x+\mathrm{e}\right|\le c$$

(12)

In the above formula, there is an eccentric displacement $e=\frac{a{\sigma }_{axial}}{4\mu p_{\max }}$ in the stick zones (the stick zone changes into stick zone 2) as shown in Figure 2. For the slip zone, the shear traction still follows the Coulomb friction model $q(x)=\mu p(x)$. Therefore, the shear stress is given by:

$$q(x)=\{\begin{array}{c}-\mu p_{\max }\sqrt{1-{\left(\frac{x}{a}\right)}^2},-a\le x<-c-e\&c-e<x\le a\\ -\mu p_{\max }\left[\sqrt{1-{\left(\frac{x}{a}\right)}^2}-\frac{c}{a}\sqrt{1-{\left(\frac{x+e}{c}\right)}^2}\right],\left|x+e\right|\le c\end{array}$$

(13)

An example of the above analytical solutions is shown in Figure 3.

3. Finite Element Model and Validation

From Figure 4, we can see the experimental setup by a schematic view [53] for contact of two cylindrical pads and a flat specimen. Under the action of normal load F, the two fretting pads maintain the contact. In addition, the coefficient of friction between the contact surfaces is 0.65 [53]. The cyclic axial load ${\sigma }_{\mathrm{axial}}$ is acting on the right side of the specimen. Two springs are attached to the fretting pad, which will generate the tangential load Q, under the combined effect of these loads, so that fretting fatigue will occur around the contact area.

As geometries and loads are symmetric, a simplified model of the structure can be constructed, which is one pad and half of the specimen. In the same way as previous research [8,14,35,38], we can model the load and boundary conditions as shown in Figure 5.

The radius of the pad and width of the specimen are taken from [53]. As shown in Figure 5, the length of the specimen is 40 mm and the width is equal to 5 mm. In addition, the thickness of two bodies and the radius of the cylindrical pad are 4 mm and 50 mm, respectively. At the top surface of the cylindrical pad, a normal load $F$ is applied. The pad is restrained from both sides in the x-direction. On the right hand-side and left hand-side of the specimen, an axial stress ${\sigma }_{axial}$ and a reaction stress ${\sigma }_{reaction}$ are applied, respectively. The value of ${\sigma }_{reaction}$ is given by [54]:

$${\sigma }_{reaction}={\sigma }_{axial}-\frac{Q}{A_s}$$

(14)

where Q is the tangential force between the contact surfaces. In addition, the cross sectional area of the specimen is expressed as $A_s$. The bottom side of the specimen is fixed in the y-direction. In order to verify our FE model and study the effect of inclusions on fretting fatigue, the experimental data used in this paper is taken from the work of Talemi and Wahab [53]. In this study, both in validation models and in parametric studies, the FF2 [53] load case has been used, with F = 543N, ${\sigma }_{axial}$ = 115 MPa and $Q_{\max }$ = 186.25 N.

Both pad and specimen are made of aluminum alloy 2024-T3, which are widely used in the aviation field. Here, we choose two kinds of inclusions Al₂CuMg and Al₂O₃ that are very commonly embedded in aluminum alloy 2024-T3 [55]. Their SEM observations showed that inclusions in material are highly discrete and randomly distributed. So we model the heterogeneity of materials by representative volume element method using DIGIMAT-FE which is a tool that considers the effects of microstructure on macroscopic material properties as shown in Figure 6. It should be noted that this study does not consider the heterogeneity of cylindrical pad material.

In this way, we can get the macro material properties (equivalent elastic modulus $E^{\ast }$ and equivalent Poisson’s ratio ${\mu }^{\ast }$) that consider the effect of microscopic inclusions [56,57,58,59]. Here, we just consider the elastic material response, because for all specimens with inclusions, von-Mises stress is always below the yield limit under such loading conditions. This is common in fretting fatigue problems. According to the previous literature, the original material properties of aluminum alloy 2024-T3 [53], Al₂CuMg [60] and Al₂O₃ [61] in this paper are given in

Table 1. The original material properties involved in this paper.

Material	Modulus (GPa)	Poisson’s Ratio
Aluminum alloy 2024-T3	72.1	0.33
${\mathsf{A}\mathrm{l}}_2\mathrm{CuMg}$	120.5	0.2
${\mathsf{A}\mathrm{l}}_2{\mathsf{O}}_3$	380	0.2

. This article assumes that the cylindrical pad is a homogeneous aluminum alloy 2024-T3.

According to the SEM study by Merati [55], here 2%, 4%, and 6% volume ratio υ between inclusions and matrix material was chosen. From the experimental observation by Hashimoto et al. [62] and other FEM research about the inclusion [62,63], we consider the inclusions as idealized spherical and ellipsoid with 23 μm to 65 μm diameter, perfectly bonded with matrix material. Due to the randomness and uniformity of the inclusions’ distribution, the RVE is constructed as a cube, which is subjected to periodic boundary conditions. In order to study the effect of inclusion, various sizes and distribution are considered. As an example, an RVE with spherical Al₂O₃ inclusions, 65 μm diameter, and 6% volume ratio has been studied first.

This case has the strongest inclusions, and the convergence of RVE size is studied for it. Four different kinds and sizes of RVEs and corresponding mesh models are shown in Figure 7. For the convergence study nine different sizes have been calculated and the corresponding macro material properties are shown in Figure 8. When the RVE size reaches 325 μm, the macro elastic modulus remains around 78.84 GPa, Poisson’s ratio will be around 0.32 and the relative difference from adjacent results is less than 1%. Therefore, the converged RVE size of this case is 325 μm and it can be applied to all other cases. In order to be safe, each case has been calculated three times, and then the average of the results is taken. In this way, we can get the macroscopic material properties of a specific heterogeneous aluminum alloy.

The fretting fatigue FE model of heterogeneous material can be built in two ways: (a) use the equivalent homogenized material in the whole specimen, or (b) model a small area near the contact region using the heterogeneous material with inclusions and use equivalent homogenized material in the rest of the specimen. Since inclusion will cause stress concentration [25,30,48,55], and fretting fatigue has maximum stress near the contact area [8], the second way is chosen in order to study the effect of inclusion on the stress distribution near the contact area. The partition diagram of the specimen is shown in Figure 9.

After that, a parametric 2D finite element model is created in ABAQUS [64] using Python script. A higher-order element always causes an instability in the stress value on the contact surfaces [65]. Therefore, we chose the CPE4R element (plane strain element, 2D, four nodes) instead of eight-node elements to mesh both parts. As shown in Figure 3, the stress distribution in the contact region is very complicated and stress amplitude is large. In particular, near the border between slip and stick zone, it changes very rapidly. In order to get precisely the stress distribution and contact stresses, the model is refined near the region of contact and inclusions. The boundary and loading conditions are as described at the beginning of this section. The contact behavior is described by the master-slave algorithm. A Lagrange multiplier is used to establish the contact between the pad and the specimen. The other contact methods can also be applied as shown by other researchers [66]. The slave surface is defined on the top surface of the specimen and the master surface is defined on the bottom surface of the pad. For homogeneous material, the stress distribution at the contact interface can be obtained analytically, if the assumptions of the Hertz solution are met. The most important two assumptions are: (a) pure elasticity and (b) the size of the contact area is small enough relative to both contact bodies. The first assumption is met as only linear elasticity is considered in this study. The second one is also called the half-space assumption [52]. The contact width for all load cases in the experiment [53] is 0.47 mm, which is less than one tenth of the height (5 mm) of the sample. A comparison between the analytical solution of Equation (13) and simulation results with different mesh sizes, for the case of homogeneous material, is shown in Figure 10a. Element sizes of 5 μm, 3 μm and 2 μm around the contact zone have been chosen for the convergence study. The mesh refinement showed convergence towards the peak analytical solution for the shear traction. Finally, according to the results, a 2 µm × 2 µm element size is used around the contact zone, which is smaller than in most of the previous numerical studies [8,32,33,34,35,36,37,38]. The simulation results and the theoretical results may not be exactly the same, due to numerical errors and geometric constraints [67]. However, the difference between simulation and theoretical solution is less than 2% (green line and red line in the Figure 10a), which can be considered as good enough to validate our FE contact model.

For comparative analysis, all the cases involved in this study are given in

Table 2. All the cases involved in this study.

Number	Volume Ratio	Type	Size (µm)	Aspect Ratio	${\mathit{E}}^{\ast }\left(\mathbf{GPa}\right)$	${\mathit{\mu }}^{\ast }$
Case 1	4%	Al2CuMg	44	1	73.7005	0.32584
Case 2	2%	Al2O3	44	1	74.269	0.32791
Case 3	4%				76.615	0.32579
Case 4	6%				78.2946	0.3243
Case 5	4%	Al2O3	23	1	76.5524	0.32503
Case 6			65		76.5046	0.32537
Case 7	4%	Al2O3	53.889	1.5	76.318	0.3256
Case 8			62.225	2	76.163	0.32255
Case 9	4%	Al2O3	23 to 65	1	76.2028	0.32532
Case 10	4%	Void	44	1	69.1	0.327

. The aspect ratio describes the evolution of inclusion from sphere to ellipsoid. The size refers to the diameter of the ball or the length of the long axis of the ellipsoid. In case 7 and case 8, we consider the ellipsoid inclusions and in order to control them, they have the same cross-sectional area as the spherical inclusion with a diameter of 44 µm and the long axis length of the ellipsoid converted to 53.889 µm and 62.225 µm, respectively.

From the results of case 3 and case 5 to 9, it can be seen that the macroscopic material properties have little to do with the size and shape of inclusions. However, by comparing case 1 with case 3, it is shown that the material properties of inclusions have a significant impact on macro elastic modulus. Similarly, by comparing case 3, 4 and 5, it is shown that the volume fraction ratio of inclusion also obviously affects the macroscopic material properties. The macroscopic Poisson’s ratio has hardly changed for the different cases.

As the contact width for all load cases in the experiment [53] is 0.47 mm. A 2 mm × 2 mm rectangle is chosen in Area 1 in Figure 9. The numerical study in this paper is divided into two parts. The first part is the simulation of the nine inclusion-structures in Table 2, which have completely randomly distributed inclusions. For the second part, since the inclusion locations of each case are randomly generated, each case cannot have the same inclusion order. This makes it impossible to use the control variable method to explore the effect of the shape, size and position of the inclusions on the contact stress distribution. Therefore, we artificially placed four inclusions below the sample contact area to compare the effects of different inclusions on the surface stress distribution. These inclusions vary in size, shape, and material corresponding to different cases. Figure 11 shows a completed finite element model, corresponding to case 3 in Table 2.

4. Numerical Results and Discussion

4.1. Completely Randomly Distributed Inclusions

4.1.1. Stress Peak and Its Location

There are many experimental studies [48,55,68,69] and numerical studies [63] on the fatigue problem of heterogeneous materials. However, limited numerical research on fretting fatigue of heterogeneous materials is available [14]. Therefore, the first part of our research is to study the cases with random distribution inclusions. In real metal materials, inclusions and defects are common and randomly distributed. In the fretting fatigue problem of homogeneous materials, the peak of shear stress appears between the stick zone and the slip zone, as shown in Figure 10b. The peak tensile stress and peak von-Mises stress in the whole specimen appears near the edge of contact [8]. From experimental materials and loading conditions, the point of occurrence of von-Mises stress and tensile stress peaks was observed at x = 0.47 mm, the shear stress peak occurred at x = 0.21 mm on the contact surface. But for the heterogeneous materials with randomly distributed inclusions (Table 2), significant stress concentration inside the specimen is observed, as shown in Figure 12. It is the stress distribution below the contact surface of the specimen. It can be seen, from comparison between Figure 10b and Figure 12c, that the shear stress peak is transferred from the contact surface to the inside of the specimen in the heterogeneous situation. It also showed that there may be multiple high stress concentration points inside the structure, hereby forming an influencing group of intrusions, eventually causing multiple fatigue cracking points inside the material as observed in the experiment [55].

Due to the complete randomness of the inclusion distribution, it is difficult to investigate the effect of particle size, shape and other factors on the stress distribution by the control variable method; such as maintaining the same inclusion material, shape, size but different volume ratio (case 2, 3, 4) to study the effect of volume ratio on the surface stress distribution of the sample. There is no regularity in the results, that is, the stress peak is almost determined by the critical defects in each case. It is similar as the experimental observation, that the size of the particles is not necessarily related to the fatigue life [55].

Based on the performed simulations, peak stress value and locations are determined, as shown in Figure 13.

Although various parameters (type, volume ratio, size, shape) of inclusions are considered in the numerical model, the data from Figure 13a indicates that tensile stress and shear stress are similar to homogeneous materials except for the von-Mises stress, which is significantly higher than the case of homogeneous materials. As shown in Figure 13b, the stress peak location is the same as the homogeneous specimen on the contact surface as mentioned before, the ordinate value is equal to 1, otherwise (below the contact surface), is equal to 2. These figures show that, for a fretting fatigue problem, heterogeneous materials containing randomly distributed inclusions have almost the same tensile and shear stress peaks and locations as homogeneous materials on the contact surface. Due to the presence of hard inclusion, the equivalent elastic modulus of the material becomes larger, resulting in a situation where the stress peak in heterogeneous material is sometimes even lower than that of the homogeneous material. But at the same time, the shear stress inside the specimen is also relatively large, comparable to peak value on the surface. Therefore, for materials with shear stress as the main fatigue index, it is likely that cracks will occur simultaneously on the surface as well as inside (Figure 12c). The von-Mises stress will be higher and appear below the contact surface with a high probability. Moreover, we can notice that for Case 10, which uses voids instead of inclusions, all of the three kinds of stress peak are significantly increased and appear inside the specimen near the edge of voids. This confirms experimental observations, that porosity is the main cause of fatigue damage, followed by oxide [68] and provides more useful results data.

4.1.2. Stress Peak Location Characteristic

From Figure 13b, we can see the shear stress is randomly appearing inside or on the surface of the material. In Figure 12c, it is also indicated that there will be multiple cracking points appearing on the surface or inside the material because both have roughly the same shear stress. For homogeneous materials, the peak value of shear stress appears on the contact surface. Also, in the left side of the sample below the surface, the shear stress is relatively higher as shown in Figure 10b. So when there are inclusions in the material, the peak value of shear stress is likely to transfer to this area (left side of sample inside). The local inclusion alignment of the stress and peak location of case 3 and case 4 are shown in Figure 14. Therefore, this shows that the inclusion in the high stress region is more likely to cause the transfer of stress peaks from surface to inside the specimen. There is also an effect of mutual coupling between the inclusions. Forming a stress band between the inclusions may directly affect the expansion direction of cracks at a later stage.

In all cases, except for the case of void (case 10), case 5 has the highest von-Mises stress. The local inclusion alignment of the stress peak location of case 5 is shown in Figure 15a. Stress coupling between the inclusions also appears here; smaller and denser inclusions form a more pronounced stress band as shown in Figure 15.

4.2. Randomly Distributed and Manually Placed Inclusions

After observing the preliminary analysis of the effects of random phenomena on inclusion, we artificially placed four inclusions below the contact surface of the specimen for each case shown in Figure 11. The geometric center is zero, and the abscissas corresponding to the four inclusions are 0.43, 0, 0.21 and 0.47. These four points correspond to the two peaks of shear stress, the geometric center and the contact edge point. As for the vertical position, their original position is 100 µm below the specimen surface. In order to analyze the effect of inclusion on the contact surface stress distribution, the size and shape of inclusions are corresponding to each case.

4.2.1. Effect of Inclusions Type

As commonly known, there are many kinds of inclusions inside metallic materials depending on the process of smelting impurities. Here, we consider two different inclusions to compare case 1 and case 3 (Table 2).

As we can see in Figure 16a, in both of homogenous and heterogeneous specimens, the von-Mises stress peak is near the contact edge, and it changes very sharply near the edge of the contact. The presence of inclusions makes the peak of the von-Mises stress increase significantly, and Al₂CuMg has a more obvious effect than Al₂O₃. This means that in the heavy load condition, the presence of inclusions may accelerate the distortion of the contact edge of the fretting fatigue contact member. The effect of inclusion on the surface tensile stress is not obvious. However, the inclusion will affect the normal and shear stresses.

4.2.2. Effect of Distance from Surface

The volume ratio (cases 2, 3 and 4) is not easy to measure for a single or several inclusions, so here we refer to the volume ratio as the particle crowding, which is reflected here as the distance between the inclusion and the contact surface. Thus, here for case 2, case 3 and case 4, the distance between the center of particle to the surface is 120 µm, 100 µm and 80 µm, respectively. The results are shown in Figure 17.

In the experimental study [48], the authors believe that the distance between the inclusions or the distance from the surface of the defect is important for fatigue damage. From our results, it is clear from Figure 17 that the distance from the surface has a more significant effect on the shear stress and normal stress. The smaller the distance from the surface, the greater the peak value of the shear stress generated. However, it can also be seen that the influence of the depth on the surface contact stress distribution is different for the inclusions at different lateral positions.

4.2.3. Effect of Inclusion Size

From Figure 18, we can see that the effect of the inclusion size is more pronounced than the influence of the inclusion type on the surface stress distribution (Figure 16).

4.2.4. Effect of Inclusion Shape

Here we consider the effect of different aspect ratios of inclusions. As shown in Figure 19, there is no significant difference between inclusion shape cases. For all cases, the effect of inclusion parameters on the tensile stress is very modest.

For the current work, we considered different stresses for analysis in order to visualize the contribution of each stress. The local stress and stress concentration effects will have a direct impact on the fretting fatigue life of the specimen. The next article shall contain the application of damage parameters (critical plane approach) to compute crack initiation lives, where we shall present the comparison of numerical lives with experimental lives.

5. Conclusions

In this study, in order to consider the effect of inclusions in aluminum alloy 2024-T3 on the fretting fatigue stress distribution, a numerical research method combining RVE with a finite element method is conducted. Based on the research results, we can draw the following conclusions.

For the cases considered in this study, with the same inclusion material to volume ratio, the shape and size of the inclusions have little effect on the macroscopic material properties. However, the material properties of inclusion and volume ratio have a significant impact on macro elastic modulus. The macroscopic Poisson’s ratio has hardly changed with the change of different parameters.

Due to the randomness of the inclusion distribution, the fatigue cracking nucleation depends on the most dangerous inclusion. The von-Mises stress peak very likely increases and transfers from the contact surface to the interior of the specimen. However, various parameters of inclusions have little effect on the peak value and position of the tensile stress; it remains almost same as in the homogeneous material. Among the materials with inclusions, the peak value and position of the shear stress are consistent with the homogeneous material, but sometimes it is also transferred to the inside of the specimen. Moreover, in all cases, higher shear stress occurs in many places inside the sample, which can result in multiple cracking points inside the sample as well as at the contact surface.

Due to the presence of hard inclusions, the equivalent elastic modulus of the material becomes larger, resulting in a situation where the stress peak in heterogeneous material is sometimes even lower than that of the homogeneous material. A void compared to inclusion will cause larger stress and randomness of stress distribution. Therefore, the influence of microstructure on macroscopic material properties should be considered in engineering.

Inclusion in high stress areas is more likely to cause stress concentration leading to peak transfer. Stress coupling between adjacent inclusions will form a high stress band, which has an important effect on the direction of crack growth.

The size of the inclusions and the distance from the surface have a more significant effect on the surface stress distribution than the type and shape of inclusions. The effect of inclusion parameters on the tensile stress is insignificant.