**2. Methods**

A 3D model of a 32-mm diameter CoCr head with a 12/14 CoCr neck (as per dimensions reported by Rehmer et al. [12]) was first analyzed under normal walking load profiles [24,25] in order to determine the most critical plane of the taper junction in terms of two important fretting wear parameters of normal contact stress (contact pressure) and micro-motion. An elastic-plastic material model was used for CoCr (ISO 5832-12) with a Young's modulus of 210 GPa, Poisson's ratio of 0.30, yield stress of 910 GPa, ultimate tensile stress of 1350 GPa, and tensile elongation of 15%. The contact pressure and displacement were retrieved for all the nodes in the contacting regions of the head and neck, using a Python script, under the maximum force and moment of the loading profile. The relative displacements (micro-motions) of the contacting nodes were then determined using a MATLAB code. The middle plane passing through the superolateral region of the neck was found as the most critical plane to have the largest area of contact pressure (as shown in Figure 1) together with micro-motions. This critical plane was then employed for the main 2D fretting wear model that will be described next.

**Figure 1.** (**a**) A three-dimensional model of a hip joint implant and (**b**) distribution of contact pressure (in Pa) in both the head and neck under a normal walking gait loading, which indicates that the super-lateral region has the largest contacting area with both contact pressures and micro-motions.

#### *2.1. Fretting Wear Model Development*

The Archard wear formulation (Equation (1)) [26] was used in the fretting wear model of this work:

$$\frac{V}{S} = k \frac{F}{H} \tag{1}$$

where *V* (m3) is the lost volume, *S* (m) is the amplitude of sliding, *k* (Pa−1) is the wear coefficient, *F* (N) is the normal force, and *H* is the hardness of the material [26]. The Archard formulation can be localized and applied to the points of a contact region, which makes it suitable for FE simulations. In addition, the Archard equation has been previously used for fretting wear simulations with success [27–30].

Dividing both sides of the Archard wear equation (Equation (1)) by the area, yields the following.

$$h = K \cdot \mathbb{S} \cdot p \cdot \Delta N \tag{2}$$

where *h* (m) is the depth of wear, *K* is the wear coefficient-to-hardness ratio (*k*/*H*), *p* (Pa) is the normal contact stress, and Δ*N* is the load cycle update interval. The main reason to re-write Equation (1) in the form of Equation (2) was to apply the Archard formula to the contacting nodes in the FE model. For the fretting wear model in the present study, a FORTRAN code was developed to trace and determine the positioning of the contacting nodes through the ABAQUS UMESHMOTION subroutine within an adaptive meshing constraint.

## 2.1.1. Verification

McColl et al. [29] and Ding et al. [28] developed an algorithm based on the Archard equation in order to simulate fretting wear for a pin-on-disc testing system. They reported surface profiles of the disc after various cycles of fretting wear. In order to verify the UMESHMOTION code developed for the head-neck junction in this study, a pin-on-disc model was first generated to replicate the Ding's model. This model had a very similar configuration (materials, geometry, element sizes, meshing structure, normal force, and sliding amplitude and frequency). The surface of the disc after the fretting wear process was evaluated (Figure 2) and compared with the results reported by Ding et al. [28]. Table 1 provides a comparison between the results of this study and those presented by Ding et al. [28] in terms of the width and height of the wear profile for the disc, which shows a very good level of agreemen<sup>t</sup> and verifies the UMESHMOTION code and its accuracy used in this study.



#### 2.1.2. Fretting Wear Model for the Head-Neck Taper Junction

The most critical plane of the head-neck junction that was previously identified by the 3D FE analysis was used to develop a 2D fretting wear model for the taper junction. Figure 2 illustrates the mesh structure of the 2D head-neck junction model and the profile of the corresponding force components (from the normal walking gait cycle) applied to this plane.

For the Archard wear equation (Equation (2)), the wear coe fficient-to-hardness ratio ( *K*) for the CoCr/CoCr head-neck combination was determined from a set of experimental results reported by Maruyama et al. [31]. They employed CoCr/CoCr pin-on-disc experiments under various normal

contact stresses and sliding cycles in a phosphate bu ffered saline (PBS) condition. From their results and using the Archard equation, the wear coe fficient-to-hardness ratio was determined for nine cases tested in their study. The nine *K* values were found to be very close with a maximum di fference of 10% from which an average of *K* = 1.7 × 10−<sup>15</sup> Pa−<sup>1</sup> was calculated and used for the CoCr/CoCr taper junction model.

The coe fficient of friction between CoCr and CoCr in the PBS condition was also obtained from the results reported by Maruyama et al. [31]. Their results for di fferent contact stresses and cycles showed that the friction coe fficient becomes constant at 0.60 after approximately 5000 cycles, which was used in the FE simulations of this work.

The authors' previous work [23] showed that head-neck taper junctions with distal angular mismatches have generally a better resistance to fretting wear when compared to junctions with proximal angular mismatches. Hence, in this paper, a small ye<sup>t</sup> realistic distal angular mismatch (0.01◦) was chosen for all the cases in order to investigate the influence of the assembly load.

The adaptive FE simulation was used to simulate the fretting wear process for one million loading cycles. An adaptive time stepping [28] was used in the simulations with an assumption of constant wear rate during a certain number of cycles ( ΔN). After several preliminary simulations, it was found that ΔN should not be assumed the same for all the periods of loading cycles. Due to the existence of a very small mismatch angle in the geometry of the interface (distal contact type with an angular mismatch of 0.01◦), large contact pressures were induced over the small contacting area at the first loading cycles, which showed that care should be taken for the selection of ΔN. During the fretting wear process, the contacting area expanded gradually, which then reduced the contact pressure. Therefore, ΔN was carefully changed from 50 to 800 loading cycles over the entire fretting wear process. The size of the elements was refined several times and 0.10 mm was found as the most suitable length of the element edge in the contact area, which could provide mesh-independent results. Figure 3a shows that the first layers of the head and neck materials at the interface were meshed with very small structured quadratic elements. These elements need to be small enough to correctly model the contact pressure and relative displacement over the contact area. The sublayers were then meshed by free-quad elements, which allow increased element sizes away from the contact area. The third part of the head and neck models was again meshed by relatively large structured elements. This meshing structure considerably reduced the solution time while providing accurate results. To simulate the interaction between the head and neck, both normal and tangential contact behaviors were defined. Normal contact was simulated using a surface-to-surface contact algorithm within ABAQUS via the "hard" contact option. The tangential interaction was modelled with a classical isotropic Coulomb friction model that was implemented with a sti ffness (penalty) method.

The 2D fretting wear model of the CoCr/CoCr taper junction was assembled with four di fferent assembly forces of 2000 N, 3000 N, 4000 N, and 5000 N. A PYTHON code and a MATLAB code were developed to report the contact pressures and relative micro-motions at the contact interface, and to find the material loss in the form of worn area from the surface at various cycles (up to 1,025,000 cycles) of normal walking gait loading.

**Figure 3.** (**a**) Mesh structure of the 2D model of the most critical plane of the head-neck junction and (**b**) profile of corresponding forces in *Y* and *Z* axes from the normal walking gait cycle.
