## 2 Materials and Methods
### Mechanical Engineering Analysis
#### 2.2.5 Discretization  
The next stage in the analysis process is the discretization of the model. Numerous publications detail the discretization stage as a fundamental part of the finite element method [24; 25; 26]. Based on these publications, some of its fundamentals are explained below for better reader comprehension. This stage consists of dividing the model into multiple elements bounded by nodes, forming a mesh. This nodal division is crucial for simplifying calculations, as the analysis variables will be computed at these points. However, the results provided by the software are obtained through the interpolation of each element's nodes, yielding a displacement value per element. This significantly affects the quality of the results. If the analysis is performed with a very large element size, the result gradients will not be accurately visible, potentially omitting relevant information. Conversely, a small element size results in a greater number of nodes, allowing for more precise value determination. Therefore, a balance between computational cost and result quality must be achieved through element size optimization. The elements constituting the mesh can be of various types depending on the model's geometry to obtain better results. For three-dimensional elements, these can be tetrahedra, prisms, or hexahedra, differentiated by the number of nodes, which results in distinct shape functions affecting displacement interpolation. These polynomial shape functions ensure displacement continuity within the discretized domain. Furthermore, it is worth mentioning that the ideal element geometry, regardless of type, should have proportional sides. An element with a long side compared to a short side will have poor quality. This is measurable using the Jacobian matrix, which defines the element's shape functions. Additionally, the differential equations governing the system are discretized, obtaining an algebraic system of equations for each element, such that the global stiffness matrix can be determined by assembling each local stiffness matrix when resolved for an element, correlating nodal displacements with nodal forces. The equilibrium equations for a continuous system are formulated as partial differential equations:  
\[\
abla\sigma+f=0 \tag{2}\]  
where external applied forces are given by the vector \(f\) and the internal stress tensor is \(\sigma\). Thus, once the boundary conditions, loads, and global stiffness matrix are applied, the system of linear equations can be solved, obtaining as a result the vector of nodal displacements, the system being as follows:  
\[f=Ku \tag{3}\]  
where the unknown nodal displacements (\(u\)) result from the nodal forces (f) and the system's stiffness, represented by the global stiffness matrix (K). Subsequently, the presentation and interpretation of the results proceed through post-processing tools, allowing the display of displacements, deformations, and stresses in the discretized model, thereby determining the most affected areas and ensuring an integrity criterion regarding the design. In the case of the present research, the element type employed for discretization was the first-order tetrahedron, as it only has nodes at the vertices of the element (Figure 14). Furthermore, the software generates the mesh automatically, with larger element sizes in larger components and smaller element sizes for smaller components or surfaces, to conform to the geometry. Furthermore, the configuration for automatic mesh generation involves various factors (with the default values offered by the software, which have been adopted in the present research, in parentheses): average element size as a fraction of the bounding box length (0.1), minimum element size as a fraction of the average size (0.2), element growth gradient factor (1.5), and maximum turn angle (\(60^{\circ}\)). This configuration provides significant control in mesh generation in contact areas between components.