*3.4. Optimization Results*

Figure 13 shows the contours of the baseline and optimized shape displayed on the duct's symmetry cross-section. The optimization affects the geometry mostly in the transition from the channel to the AIP. The sharp shape change from the inner wall to the compressor entrance is moderated in the final solution. Moreover, the optimizer reduced the curvature of the inner and outer walls; however, the effect on the concave side is more apparent.

**Figure 13.** Contours of the duct's symmetry plane for reference and optimized shape.

Table 5 groups the optimized objective function values for all considered conditions and compares them with the reference levels. The results are presented on an absolute scale and as a percentage of improvement with respect to the baseline shape. Undoubtedly, the optimizer was successful in finding a solution that improves both objectives for all considered design points.


**Table 5.** Quantitative optimization results.

The greatest improvement level in both objectives is attained for the nominal cruise, which corresponds with the predefined preferences in targets (see Table 3). Furthermore, comparable gains in the pressure loss and the distortion coefficient coincide with the equality in targets set for both objectives. However, a slight consistent preference towards refinement in the distortion coefficient is noticeable. Such an outcome confirms the ability of the ASF formulation to guide optimization toward the presumed target, although a margin has to be considered for potential moderate departures.

Figure 14 shows a cumulated best value of the ASF and the convergence criterion value along the EGO steps. A rapid convergence towards the optimized solutions is evident. The optimizer finds a case significantly improving the initial ASF value (marked by the dashed red line) already in the first iteration. Such superior performance results from a fine design space sampling and high-quality surrogate fit.

**Figure 14.** Cumulated minimum value of the ASF (**a**) and convergence criterion evolution (**b**) in the course of optimization.

The objective space with DoE and EGO evaluations is displayed in Figure 15 for the three considered design points. For all conditions, the instances from the optimization stage converge towards an imaginary line connecting the reference and target points. The vast majority of the Pareto front is not resolved, which is intentional and prioritizes the use of computational resources in the target's neighborhood. Such an observation shows the algorithm's ability to follow the direction defined by a predetermined target.

**Figure 15.** Optimization solutions in the objective function space for the three design points. Markers indicate optimization target (target), evaluation of initial geometry (reference), and best obtained solution (best).

Figure 16 gives insights into the improvements in the distortion coefficient. Figure 16a contrasts radially averaged total pressure values for the reference and optimized solutions. The pressure value is normalized against the average dynamic pressure over the AIP and plotted against the azimuthal position at the compressor face. For the optimized solution, an apparent reduction in the total pressure's lowest peak value manifests for all considered design points. The source of refinement is located in the 150–210◦ sector, corresponding with the transition area from the duct's inner wall to the AIP. Such an improvement translates directly to an enhancement in the distortion coefficient. A slight shift in the total pressure curves for the optimized solution toward higher values is an

effect of reduction in the duct's pressure loss; however, it does not directly impact the distortion metric.

**Figure 16.** Details of the distortion coefficient improvements: (**a**) radially averaged total pressure distribution at AIP. Solid and dashed lines indicate total pressure levels averaged over the AIP and worst 60◦ sector, respectively; (**b**) total pressure at AIP normalized by the average dynamic head.

The distortion coefficient value can be visualized intuitively as an area between a plane-averaged total pressure (solid lines in Figure 16a) and total pressure averaged over the worst 60◦ sector (dashed lines in Figure 16a). The areas are marked with shaded fields with a color corresponding to the reference (gray) and optimized (blue) solutions. The area of the fields representing the improved solution is visibly smaller than the corresponding regions for the baseline case. Such an observation confirms the quantitative results given above in Table 5.

Figure 16b groups the maps of total pressure normalized by the average dynamic head and plotted at the AIP. The pressure field homogenization is visible for all conditions, although the most significant effect is noticeable for Design Point 1. This corresponds with the highest *DC*<sup>60</sup> improvement observed in the quantitative data. The reduction in the low-pressure zone visible in the six o'clock position matches the improvement in the 150–210◦ sector discussed above. A slight improvement is observable in the twelve o'clock region corresponding to the convex wall transition to the AIP. This local enhancement

can also be seen in Figure 16a for sector ~360 ± 30◦. Even though such improvement does not influence the distortion metric, an increase in the pressure field uniformity at the compressor face is a positive effect.

Figure 17 shows details of the flow field with regard to the pressure loss coefficient improvements. Figure 17a displays the flow loss evolution along the duct through a local pressure loss coefficient contour. This metric is conceptually similar to the pressure loss objective; however, it quantifies the drop from the duct's inlet to each spatial location (Equation (33)).

$$dP^\* = \frac{P\_t^{IN} - P\_t}{P\_t^{IN}} \tag{33}$$

**Figure 17.** Details of the pressure loss coefficient improvements: (**a**) maps of local pressure loss coefficient; (**b**) flow streamlines colored by the vertical velocity component at cross−section plane downstream of the AIP.

The most apparent enhancement is the reduction in pressure loss in the transition region from the duct's convex wall to the AIP. Adjusting the shape of the wall's top sector and the above-discussed reduction in the low-pressure zone on the concave side transition removes a secondary flow motion downstream of the AIP. This improvement is visualized in Figure 17b, in which the vertical velocity component values color flow streamlines. The two counter-rotating vortices, visible at the twelve o'clock position in the baseline case, are evidently removed in the optimized solution.

The reduction in the concave wall curvature and expansion of the duct's cross-section results in a diffusing shape, bringing additional pressure recovery from the flow kinetic energy. This effect is subtle but still contributes to the duct's overall performance.

#### **4. Conclusions**

The proposed aerodynamic optimization framework, constructed from state-of-the-art components in advanced surrogates, mesh morphing, and distance-based scalarization, was applied to the multi-point optimization problem of an I-31T airplane air-intake duct. The study aimed to simultaneously improve the duct's pressure loss and flow distortion under three flight conditions: nominal cruise, low-altitude climbing, and high-altitude cruise.

The optimizer obtained both objective values superior to the reference configuration for all considered design points. The consistent level of improvements in the pressure loss and flow distortion confirms the capacity of the ASF-based formulation to guide optimization toward the presumed target. The study results prove the methodology's potential for optimizing complex multi-objective air-intake duct problems in multiple flight conditions while saving substantial computational resources.

Moreover, FANOVA-based sensitivity analysis was recognized as a valuable tool for assessing the importance of particular design variables. Application of this technique is particularly beneficial in Kriging-based frameworks where the use of surrogate predictions balances high computational costs related to the variance-based techniques.

This research, however, may be subject to potential limitations regarding the generalization of its results. Primarily, the number of design variables was relatively moderate, which favors the use of the Kriging surrogate. Although literature sources report a considerable margin for increased design space dimensions, studies on larger problems should be executed to prove the framework's efficiency. Moreover, the results revealed that all considered flight conditions were coherent in the direction of shape improvements. Although justified for this particular application, this scenario might only represent part of the class of problems. Future studies should cover scenarios with design points of contradicting performance requirements to generalize the results. Finally, the optimization was initiated from an already well-designed solution, which aimed to set an ambitious task for the optimizer but also resulted in a relatively low level of required deformations. Optimization studies demanding significant geometry adjustments may be subject to numerical errors resulting from an excessively distorted mesh. Such problems may require the implementation of a solution for handling erroneous data samples.

Optimization studies using the proposed framework could be expanded to more holistic multi-disciplinary problems. The most natural development would include solid body mechanics addressing the structural targets, although involvement of the intuitively more distant performance and economic objectives seems feasible. No fundamental reasons were identified that might inhibit the ASF-based formulation in handling such heterogeneous goals, although this would need to be proved in further studies.

**Author Contributions:** Conceptualization, P.S.D., S.K. and J.Z.; methodology, P.S.D. and S.K.; soft- ˙ ware, P.S.D.; validation, P.S.D.; formal analysis, P.S.D.; investigation, P.S.D.; resources, P.S.D. and J.Z.; ˙ data curation, P.S.D.; writing—original draft preparation, P.S.D.; writing—review and editing, P.S.D. and S.K.; visualization, P.S.D.; supervision, P.S.D., S.K. and J.Z. All authors have read and agreed to ˙ the published version of the manuscript.

**Funding:** This work was co-funded by the Ministry of Education and Science as part of the "Implementation Doctorate first edition" program under the agreement no 3/DW/2017/01/1.

**Data Availability Statement:** Not applicable.

**Acknowledgments:** We would like to express our great appreciation to Karolina Dr˛ezek for her ˙ support in the preparation of figures.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

#### **Appendix A. Leave-One-Out Cross-Validation**

The leave-one-out cross-validation (LOO-CV) technique [91] sequentially omits one sample from the set of *N* observations, and the metamodel is constructed from the remaining *N* − 1 points. Afterward, the objective is predicted at the left-out sample's location, and the estimated value is compared with the known evaluation. The concept of LOO-CV is demonstrated schematically in Figure A1. The operation repeats *N* times, resulting in an assessment of prediction accuracy for the whole observation set.

**Figure A1.** Schematic representation of LOO-CV concept for one-dimensional design space.

The measure of difference between the true objective value *y*(*xi*) and the corresponding LOO-CV prediction *<sup>y</sup>*ˆD−*i*(*xi*) is termed a residual (Equation (A1)).

$$\epsilon\_i = \mathcal{y}(\mathbf{x}\_i) - \mathcal{Y}\_{\mathfrak{D}-i}(\mathbf{x}\_i) \tag{A1}$$

The residuals are standardized using the variance estimator *s* -2 D−*i* (*xi*) available in the Kriging surrogate (Equation (A2)).

$$
\epsilon\_i^z = \frac{\epsilon\_i}{\sqrt{\overline{s\_{\mathfrak{D}}^2}}} \tag{A2}
$$

If the residuals approximately follow the Gaussian distribution, the objective evaluation locates within three standard deviations from the predicted Kriging mean with a confidence level of 99.7%. Thus, a surrogate with good prediction quality should be characterized by standardized residuals bounded in a ±3 interval [23]. The probability of a residual being further away is less than 0.3%. Detecting such an outlier may suggest poor prediction quality in the vicinity of the corresponding observation.

#### **Appendix B. Functional Analysis of Variance (FANOVA)**

Consider that the design space X forms a *<sup>p</sup>*-dimensional hypercube and *<sup>x</sup>* <sup>∈</sup> X is a vector of independent random variables normalized to a range [0, 1]. The surrogate model is described by a square-integrable function *<sup>Y</sup>* <sup>=</sup> *<sup>f</sup>*(*x*) defined in <sup>X</sup>, which can be decomposed using the FANOVA representation [89]:

$$Y(\mathbf{x}) = f\_0 + \sum\_{i=1}^p f\_i(\mathbf{x}\_i) + \sum\_{1 \le i < j \le p} f\_{\bar{\mathbf{y}}}(\mathbf{x}\_i, \mathbf{x}\_j) + \dots + f\_{1, \dots, p}(\mathbf{x}\_1, \dots, \mathbf{x}\_p) \tag{A3}$$

In the above equation, the centered and orthogonal terms denote:

$$\text{mean value } f\_0 = \mathbb{E}[Y(\mathbf{x})] \tag{A4}$$

$$\text{main effects}\, f\_i(\mathbf{x}\_i) = \mathbb{E}[Y(\mathbf{x})|\mathbf{x}\_i] - f\_0 \tag{A5}$$

second-order interactions *fij xi*, *xj* = E *Y*(*x*) % %*xi*, *xj* − *f*<sup>0</sup> − *fi*(*xi*) − *fj xj* (A6)

The interactions of higher orders can be constructed accordingly as conditional expected values.

A similar technique serves for a decomposition of the model output's variance:

$$Var(Y(\mathbf{x})) = \sum\_{i=1}^{p} Var(f\_i(\mathbf{x}\_i)) + \sum\_{1 \le i < j \le p} Var\left(f\_{ij}(\mathbf{x}\_i, \mathbf{x}\_j)\right) + \dots + Var\left(f\_{1, \dots, p}(\mathbf{x}\_1, \dots, \mathbf{x}\_p)\right) \tag{A7}$$

The individual contribution of variable *xi* to the variance in the model's output is quantified by the *main effect Sobol index*:

$$S\_i = \frac{\operatorname{Var}(f\_i(\mathbf{x}\_i))}{\operatorname{Var}(\mathbf{Y}(\mathbf{x}))} \tag{A8}$$

The influence of interaction between any two variables *xi* and *xj* is described by the second-order Sobol' index:

$$S\_{ij} = \frac{\operatorname{Var}\left(f\_{ij}\left(\mathbf{x}\_{i\prime}\mathbf{x}\_{j}\right)\right)}{\operatorname{Var}\left(Y(\mathbf{x})\right)}\tag{A9}$$

All higher-order interactions can be assessed by Sobol' indices, constructed following a similar concept of ratios between decomposed higher-order terms and the overall output variance. The values of main and higher-order indices for all variables sum to unity (Equation (A10)).

$$\sum\_{i=1}^{p} S\_i + \sum\_{1 \le i < j \le p} S\_{ij} + \dots + S\_{1,\dots,p} = 1 \tag{A10}$$

The aggregated contribution to the model's output variance of the *i*-th variable is measured by the total effect Sobol' index:

$$S\_i^T = \frac{\mathbb{E}[Var(Y(x)|x\rangle\langle x\_i\rangle)]}{Var(Y(x))}\tag{A11}$$

In practice, only main and total indices are assessed for economic reasons. Knowing them allows for an estimation of the combined influence of all-order interactions.

The values of the Sobol' indices are to be computed numerically. For this purpose, this study employs the FAST algorithm [92] available in the R package sensitivity [93].
