## 4 Design Optimization  
### Optimization Methodology  
For contact geometry design optimization, we developed and implemented an Integrated Multi-stage Adaptive Sequential Optimization (IMASO) algorithm. Figure 8 presents the comprehensive workflow of this algorithm, showing the systematic progression from initial sampling through metamodel construction to final optimization and validation. To efficiently explore the design space, we employed a 3\({}^{3}\) LHS design consisting of 27 systematically distributed sampling points. We carefully selected the design variables and their ranges, as detailed in Table 5, based on preliminary analysis and manufacturing constraints, ensuring practical implementability of the optimized design. Figure 8: Flowchart of the optimization algorithm for contact geometry design. \begin{table}
\begin{tabular}{c c c} \hline \hline
**Design Variable** & **Minimum Value** & **Maximum Value** \\ \hline \(R_{p}\) (mm) & 5.25 & 6.25 \\ \(R_{s}\) (mm) & 4.25 & 5.25 \\ \(e\) (mm) & 0.85 & 1.15 \\ \hline \hline \end{tabular}
\end{table} 
Table 5: Design variable ranges for optimization study. We subjected each of the 27 design configurations to detailed FEA to evaluate two critical performance metrics: maximum contact pressure at low pressure (2.0 MPa) and maximum equivalent stress at high pressure (87.5 MPa). Table 6 presents these analysis results, revealing significant performance variations across the design space and highlighting the potential for substantial improvement through optimization. Based on these results, we constructed a kriging metamodel using the established statistical framework commonly employed in the design optimization literature, expressed as follows [14; 15]:  
\[\hat{y}(x)=\mu+Z(x) \tag{13}\]  
where \(\hat{y}(x)\) is the predicted response at design point \(x\), \(\mu\) is the global trend function, and \(Z(x)\) represents local deviations. The spatial correlation necessary for computing \(Z(x)\) is quantified by the following:  
\[R\big{(}x_{i},x_{j}\big{)}=e^{-\sum_{k=1}^{3}\theta_{k}|x_{i,k}-x_{j,k}|^{2}} \tag{14}\]  
where \(\theta_{k}\) are correlation parameters that determine how rapidly correlation decreases with distance in each dimension. To determine optimal values for \(\theta_{k}\), we maximized the log-likelihood function:  
\[L(\theta)=-0.5\Bigg{[}n\text{ln}\Big{(}\sigma^{2}\Big{)}+\text{ln}|R|+\frac{( Y-\mu)^{T}R^{(-1)(Y-\mu)}}{\sigma^{2}}\Bigg{]} \tag{15}\] \begin{table}
\begin{tabular}{c c c c c} \hline \hline
**LHS** & & & & \\
**Sample No.** & \(R_{p}\) **(mm)** & \(R_{s}\) **(mm)** & \(e\) **(mm)** & **Max. Contact Pressure (MPa)** & **Max. Equivalent Stress (MPa)** \\ \hline \(1\) & 5.30 & 4.30 & 0.90 & 123.61 & 167.17 \\ \(2\) & 5.32 & 4.68 & 0.91 & 142.36 & 225.96 \\ \(3\) & 5.35 & 4.85 & 1.13 & 135.34 & 265.00 \\ \(4\) & 5.38 & 4.42 & 1.08 & 139.98 & 221.37 \\ \(5\) & 5.40 & 5.20 & 0.95 & 148.22 & 263.01 \\ \(6\) & 5.45 & 4.55 & 0.92 & 151.79 & 181.43 \\ \(7\) & 5.50 & 4.40 & 1.10 & 134.35 & 203.14 \\ \(8\) & 5.52 & 4.72 & 0.94 & 145.83 & 194.38 \\ \(9\) & 5.55 & 5.15 & 1.07 & 132.55 & 239.38 \\ \(10\) & 5.60 & 5.00 & 0.93 & 131.82 & 219.20 \\ \(11\) & 5.65 & 4.33 & 0.96 & 130.56 & 155.01 \\ \(12\) & 5.7 & 4.60 & 1.05 & 143.94 & 174.29 \\ \(13\) & 5.74 & 4.76 & 1.00 & 137.85 & 175.93 \\ \(14\) & 5.78 & 4.38 & 1.03 & 115.34 & 158.03 \\ \(15\) & 5.8 & 4.90 & 0.88 & 116.12 & 179.10 \\ \(16\) & 5.85 & 5.05 & 1.11 & 116.53 & 224.38 \\ \(17\) & 5.9 & 4.35 & 1.12 & 108.69 & 155.75 \\ \(18\) & 5.92 & 5.17 & 0.98 & 114.49 & 198.04 \\ \(19\) & 5.95 & 4.65 & 0.89 & 104.16 & 157.70 \\ \(20\) & 5.97 & 5.03 & 1.15 & 141.39 & 216.12 \\ \(21\) & 6.00 & 4.80 & 0.97 & 142.80 & 164.65 \\ \(22\) & 6.04 & 4.46 & 1.06 & 99.414 & 163.11 \\ \(23\) & 6.08 & 4.58 & 0.86 & 108.10 & 164.63 \\ \(24\) & 6.10 & 4.50 & 1.02 & 83.649 & 152.95 \\ \(25\) & 6.15 & 4.95 & 1.04 & 124.53 & 175.04 \\ \(26\) & 6.18 & 4.88 & 1.14 & 136.59 & 169.76 \\ \(27\) & 6.20 & 5.10 & 0.87 & 133.79 & 158.29 \\ \hline \hline \end{tabular}
\end{table}