## 3 Multibody Dynamics Analysis  
To estimate the reaction forces influencing wear at a joint within a multibody system, it is necessary to conduct a multibody dynamics analysis. This process involves formulating and solving a system of Differential Algebraic Equations (DAEs) of motion that characterizes the system's dynamics. A concise overview of this analysis is provided in this section.  
A multibody system comprises interconnected bodies capable of undergoing significant displacements and rotations. In this study, the bodies are assumed to be rigid and confined to planar motion. Each body, denoted as \(i\), is associated with a fixed local reference frame (\(x_{i}\)-\(y_{i}\)) and is tracked within a global coordinate system (\(x\)-\(y\)). The position of each body within the system is represented in the global coordinate system as follows:  
\[\mathbf{q}_{i}=[x,y,\phi]_{i}^{T}, \tag{9}\]  
where \(\phi\) represents the orientation of the body-fixed reference frame relative to the global frame. For a multibody system, a set of generalized coordinates that uniquely define both the position and orientation of the system's bodies can be expressed as follows:  
\[\mathbf{q}=\Big{[}[x,y,\phi]_{1}^{T},[x,y,\phi]_{2}^{T}\ldots,[x,y,\phi]_{i}^ {T}\Big{]}^{T}. \tag{10}\]  
The bodies in the system are connected by joints that constrain their motion relative to one another. These constraints, expressed as algebraic equations dependent on time and the generalized coordinates, are known as holonomic kinematic constraints and are formulated as follows:  
\[\mathbf{\Phi}(\mathbf{q},t)=\mathbf{0}. \tag{11}\]  
For a properly constrained system (where the number of constraints equals the number of generalized coordinates), Equation (11) can be solved to determine the position \(\mathbf{q}\) of all system components. Additionally, the velocity \(\dot{\mathbf{q}}\) and acceleration \(\ddot{\mathbf{q}}\) of the components can be calculated by solving Equations (12) and (13), which represent the first and second derivatives of Equation (11), respectively.  
\[\mathbf{\Phi}_{\mathbf{q}}\dot{\mathbf{q}}=-\mathbf{\Phi}_{t} \tag{12}\]  
\[\mathbf{\Phi}_{\mathbf{q}}\ddot{\mathbf{q}}=-\left(\mathbf{\Phi}_{\mathbf{q}} \ddot{\mathbf{q}}\right)_{\mathbf{q}}\dot{\mathbf{q}}-2\mathbf{\Phi}_{\mathbf{ q}t}\dot{\mathbf{q}}-\mathbf{\Phi}_{tt}\equiv\mathbf{\gamma} \tag{13}\]  
The process of determining the position, velocity, and acceleration by solving Equations (11)-(13) is known as kinematic analysis. To analyze the system's response to external loads, a Differential Algebraic Equation (DAE) of motion must be formulated and solved. The derivation of the DAE for a constrained multibody system is beyond the scope of this work and is presented directly in Equation (14) without derivation. A comprehensive explanation of such a derivation, using the Lagrangian approach, can be found in [43].  
\[\mathbf{M}\ddot{\mathbf{q}}\mathbf{+}\mathbf{\Phi}_{\mathbf{q}}^{T}\mathbf{ \lambda}=\mathbf{Q}^{A}. \tag{14}\]  
In Equation (14), \(\mathbf{M}\) represents the mass matrix and the moment of inertia for the system components, \(\mathbf{\lambda}\) is a vector of Lagrange multipliers, and \(\mathbf{Q}^{A}\) is a vector of the applied loads. Equations (13) and (14) are combined to formulate the DAE of motion, which is expressed as follows:  
\[\begin{bmatrix}\mathbf{M}&\mathbf{\Phi}_{\mathbf{q}}^{T}\\ \mathbf{\Phi}_{\mathbf{q}}&\mathbf{0}\end{bmatrix}\begin{bmatrix}\ddot{ \mathbf{q}}\\ \mathbf{\lambda}\end{bmatrix}=\begin{bmatrix}\mathbf{Q}^{A}\\ \mathbf{\gamma}\end{bmatrix}, \tag{15}\]  
where \(\mathbf{\gamma}=\mathbf{\Phi}_{\mathbf{q}}\ddot{\mathbf{q}}=-\left(\mathbf{ \Phi}_{\mathbf{q}}\dot{\mathbf{q}}\right)_{\mathbf{q}}\dot{\mathbf{q}}-2 \mathbf{\Phi}_{\mathbf{q}t}\dot{\mathbf{q}}-\mathbf{\Phi}_{tt}\). Equation (15) can then be solved to determine the response of the system to applied loads. Complete details of the multibody dynamics analysis, as well as the derivations of Equation (9) through Equation (15), can be found in [43].