# Efficient initialization for multi-fidelity surrogate-based optimization

Jeroen Wackers, Riccardo Pellegrini, Andrea Serani, Michel Visonneau,

Matteo Diez

► To cite this version:

Jeroen Wackers, Riccardo Pellegrini, Andrea Serani, Michel Visonneau, Matteo Diez. Efficient initialization for multi-fidelity surrogate-based optimization. Journal of Ocean Engineering and Marine Energy, In press, 10.1007/s40722-022-00268-5. [https://doi.org/10.1007/s40722-022-00268-5](https://doi.org/10.1007/s40722-022-00268-5). [hal-03847005](https://HAL.inan逆)
# Efficient Initialization for Multi-Fidelity Surrogate-Based Optimization

Jeroen Wackers<fnref n="1" />*†, Riccardo Pellegrini<fnref n="2" />†, Andrea Serani<fnref n="2" />, Michel Visonneau<fnref n="1" /> and Matteo Diez<fnref n="2" />

*Corresponding author(s). E-mail(s): jeroen.wackers@ec-nantes.fr;
†These authors contributed equally to this work.

## Abstract

The performance of surrogate-based optimization is dependent on the surrogate training set, certainly for realistic optimizations where the high cost of computing the training set data imposes small training set sizes. This is especially true for multi-fidelity surrogate models, where different training sets exist for each fidelity. Adaptive sampling methods have been developed to improve the fitting capabilities of surrogate models, adding training points only where necessary or most useful to the optimization process (i.e., providing the highest knowledge gain) and avoiding the need for an a priori design of experiments. Nevertheless, the efficiency of the adaptive sampling is highly affected by its initialization. The paper presents and discusses a novel initialization strategy with a limited training set for adaptive sampling. The proposed strategy aims to reduce the computational cost of evaluating the initial training set. Furthermore, it allows the surrogate model to adapt more freely to the data. In this work, the proposed approach is applied to single- and multi-fidelity stochastic radial basis functions for an analytical test problem and the shape optimization of a NACA hydrofoil. Numerical results show that the results of the surrogate-based optimization are improved, thanks to a more effective and efficient domain space exploration and a significant reduction of high-fidelity evaluations.

**Keywords:** Multi-fidelity, Surrogate-based optimization, Stochastic Radial Basis Functions, Initial Training Set

<footnote n="1">¹LHEEA Lab, CNRS UMR 6598, Centrale Nantes, 1 rue de la Noë, B.P. 92101, 44321 Nantes cedex 3, France.</footnote>

<footnote n="2">²National Research Council - Institute of Marine Engineering (CNR-INM), Via di Vallerano 139, 00128 Rome, Italy.</footnote>
# 1 Introduction

Automatic shape optimization and uncertainty quantification offer rigorous and effective mathematical approaches to the design and performance assessment of modern ships and ships' subsystems (Karlberg et al., 2013; Harries and Abt, 2019). These methods generally require a large number of evaluations of one or more merit factors: in (Furcas et al., 2020) thousands of geometries of a wake equalizing duct have been simulated before converging to an optimal solution. Similarly, in (Demo et al., 2021; Serani et al., 2022) although computational cost reduction methods are used, hundreds of numerical simulations are still required to converge towards a global optimum. If these performance metrics are evaluated via high-fidelity computations, the computational cost can become prohibitively high and unaffordable for most users (Serani et al., 2022). Surrogate models give a solution to this problem: computations are only performed in a few design points or conditions and an interpolatory or regressive model is built based on these computations (Jin et al., 2001). Shape optimization and uncertainty quantification are then performed using the surrogate model, which is inexpensive to evaluate. Among other surrogate models, radial basis functions have demonstrated their efficacy in analytical tests (Jin et al., 2001) and have been successfully used for design optimization (Volpi et al., 2015), also in high-dimensional spaces (Regis, 2020). An even greater gain in efficiency is obtained through multi-fidelity surrogate models, which combine computationally inexpensive low-fidelity simulations with high-fidelity computations (Peherstorfer et al., 2018). Two or more fidelities can be defined by combining computations performed with different physical models (e.g. Reynolds averaged Navier-Stokes equations and potential flow, Pellegrini et al. 2018) or with different accuracy (e.g. computations performed varying the grid size, Bonfiglio et al. 2018). Radial basis functions have also been successfully used for multi-fidelity optimizations (Nuñez et al., 2018).

The performance of multi-fidelity surrogate models depends on several factors (Fernández-Godino et al., 2019), such as the presence of nonlinearities, the problem dimensionality, the noisy or smooth behavior of the function, and the approach used for the definition of the training set. Numerical experiments show that there is no unique optimal multi-fidelity approach: the best choice of the surrogate model depends on the data being modeled. For example, a significant local variation of the merit factors may require a high density of training data, the presence of noise may require heavy filtering, and finally low-fidelity corrections are only useful when the low- and high-fidelity data are sufficiently correlated (similar). The first two issues led the authors to adopt solution-adaptive sampling for the training points (Serani et al., 2019) and automatically tuned noise filters (Wackers et al., 2020). Such procedures, which automatically adapt themselves to the behavior of the data, are a necessity to achieve efficient multi-fidelity surrogate modeling.

However, adaptive multi-fidelity sampling starts from initial datasets for each fidelity, defined in a non-adaptive manner. Computing these initial datasets may be expensive, while the data may not be optimal for the surrogate
model accuracy. For example, in the authors' previous work, the multi-fidelity model is based on a surrogate model of the low-fidelity solutions and a discrepancy surrogate model based on the difference between high- and low-fidelity solutions. For each fidelity level, the initial training set for starting the surrogate-based optimization is defined using $2D + 1$ samples, where $D$ is the design space dimension. The training points are placed in the center of the domain and the center of the boundaries. This initialization approach requires a considerable computational cost, especially for the highest fidelity. Furthermore, unless the optimum is located on a domain boundary, most of these points are of low importance for the optimization.

The objective of the present work is to introduce a new approach for the definition of the initial training set. The new approach uses only one initial point for all the fidelities except for the lowest. When a single training point is available, the surrogate model prediction is an extrapolation based on that single training point. Thus, the proposed approach reduces the computational cost of the initialization making available a larger budget for the adaptive sampling of the training points. Furthermore, less information is initially provided to the surrogate model in order to have it freely adapt to the data.

The surrogate model is based on stochastic radial-basis functions (SRBF) with a power-law kernel. The power kernel used in the SRBF lacks a compact support, as a consequence the extrapolated prediction may not be well correlated with the desired function behavior, negatively affecting the adaptive sampling. Although other kernels exist with compact support, the SRBF with power kernels is robust, showing good results for several applications (Wackers et al., 2020) and is preferred here. Therefore, a constraint is imposed on both the surrogate model prediction and the associated uncertainty when an extrapolation is performed, to improve the adaptive sampling.

The proposed initialization approach is assessed for an analytical test problem and a simulation-based optimization problem: the drag-coefficient minimization of a NACA hydrofoil. For both problems, one and three fidelities are used and the results are compared with the previous initialization approach of Wackers et al. (2020). The simulations are performed with the unsteady RANS solver ISIS-CFD (Queutey and Visonneau, 2007), developed at Ecole Centrale de Nantes/CNRS and integrated in the FINE/Marine simulation suite from Cadence Design Systems. Mesh deformation and adaptive grid refinement are adopted, different fidelity levels are defined by increasing the grid refinement.

## 2 Multi-fidelity approach

The multi-fidelity surrogate model is based on a summation of surrogate models. Consider $\mathbf{x} \in \mathbb{R}^D$ as a design variables vector of dimension $D$. Let the true function $f(\mathbf{x})$ be assessed by $N$ fidelity levels: the highest-fidelity level is $f_1(\mathbf{x})$, the lowest-fidelity is $f_N(\mathbf{x})$, and the intermediate fidelity levels are $\{f_l\}_{l=2}^{N-1}(\mathbf{x})$. Using $\sim$ to denote surrogate model prediction, the multi-fidelity
4 Efficient Initialization for Multi-Fidelity Surrogate-Based Optimization

(MF) approximation $\hat{f}_l(\mathbf{x})$ of $f_l(\mathbf{x})$ ($l = 1, \dots, N-1$) is the sum of the lowest-fidelity surrogate and surrogates of the difference between subsequent levels (inter-level errors or bridge-functions, $\tilde{\varepsilon}$):

$$ \hat{f}_l(\mathbf{x}) = \tilde{f}_N(\mathbf{x}) + \sum_{k=l}^{N-1} \tilde{\varepsilon}_k(\mathbf{x}) \quad (1) $$

The final multi-fidelity approximation of $f$ is $\hat{f} = \hat{f}_1$. For each $l$-th fidelity level the training set is $\mathcal{T}_l = \{\mathbf{x}_i, f_l(\mathbf{x}_i)\}_{i=1}^{J_l}$, with $J_l$ the training set size. The resulting inter-level error training sets are defined as $\mathcal{E}_l = \{\mathbf{x}_i, \varepsilon_l(\mathbf{x}_i)\}_{i=1}^{J_l}$, where

$$ \varepsilon_l(\mathbf{x}_i) = f_l(\mathbf{x}_i) - \hat{f}_{l+1}(\mathbf{x}_i) \quad (2) $$

The SRBF surrogate models provide both a prediction and an associated uncertainty. The uncertainty $U_{\tilde{f}_N}$ of the lowest-fidelity prediction is considered as uncorrelated with the uncertainty $U_{\tilde{\varepsilon}_l}$ of the inter-level error predictions. Therefore, the uncertainty $U_{\hat{f}_l}$ of the MF prediction can be evaluated as ($l = 1, \dots, N-1$)

$$ U_{\hat{f}_l}(\mathbf{x}) = \sqrt{U_{\tilde{f}_N}^2(\mathbf{x}) + \sum_{k=l}^{N-1} U_{\tilde{\varepsilon}_k}^2(\mathbf{x})} \quad (3) $$

## 2.1 Adaptive sampling method

The multi-fidelity surrogate model is dynamically updated by adding new training points following a deterministic approach. First, a new training point $\mathbf{x}^*$ is identified based on an aggregate-criteria adaptive sampling method, specifically a modified version of the lower-confidence bounding from Cox and John (1992) with equally weighted contributions of $\hat{f}$ and $U_{\hat{f}}$, presented in Serani et al. (2019). It aims to find points with large prediction uncertainty and small objective function value, see Figure 1. Accordingly, the sampling method identifies a new training point by solving the single-objective minimization $\mathbf{x}^* = \operatorname*{argmin}_{\mathbf{x}} [\hat{f}(\mathbf{x}) - U_{\hat{f}}(\mathbf{x})]$ using a deterministic version of the particle swarm optimization algorithm introduced in (Serani et al., 2016). Once $\mathbf{x}^*$ is identified, the fidelity to be evaluated needs to be selected. To achieve this, the surrogate model prediction uncertainty vector is defined as

$$ \mathbf{U} = \{U_{\hat{\varepsilon}_1}/\beta_1, \dots, U_{\hat{\varepsilon}_{N-1}}/\beta_{N-1}, U_{\tilde{f}_N}/\beta_N\} \quad (4) $$

where $\beta_l = c_l/c_1$ with $c_l$ the computational cost associated to the $l$-th level and $c_1$ the computational cost of the highest fidelity. Then, the fidelity level to sample is $l^* = \max\{\mathbf{U}\}$, and the new training point is added to the $l^*$-th training set $\mathcal{T}_{l^*}$ and to the lower-fidelity sets from $l^* + 1$ up to $N$.
<figure><img src="image_1.png" /><figcaption>Fig. 1: Example of the aggregate-criterion adaptive sampling method using one fidelity: (a) shows the initial surrogate model with the associated prediction uncertainty and training set; (b) shows the position of the new training point and the new surrogate model prediction and its uncertainty</figcaption></figure>

## 2.2 Stochastic Radial Basis Functions

Given a (single-fidelity) training set $\mathcal{T} = \{\mathbf{x}_i, f(\mathbf{x}_i)\}_{i=1}^J$, the SRBF surrogate model prediction $\tilde{f}(\mathbf{x})$ is computed as the expected value (EV) over a stochastic tuning parameter of the surrogate model (Volpi et al., 2015), $\tau \sim \text{unif}[1, 3]$

$$ \begin{gathered} \tilde{f}(\mathbf{x}) = \text{EV}[g(\mathbf{x}, \tau)]_{\tau} \quad \text{with} \\ g(\mathbf{x}, \tau) = \text{EV}[\mathbf{f}] + \sum_{j=1}^{M} w_j \|\mathbf{x} - \mathbf{c}_j\|^{\tau} \end{gathered} \qquad (5) $$

where $w_j$ are unknown coefficients, $\|\cdot\|$ is the Euclidean norm, and $\mathbf{c}_j$ are the RBF centers, with $j = 1, \dots, M$ and $M \le J$. If the training set is not affected by numerical noise then exact interpolation of the training set is imposed and the coefficients $w_j$ are computed by solving $\mathbf{A}\mathbf{w} = (\mathbf{f} - \text{EV}[\mathbf{f}])$, with $\mathbf{c}_j = \mathbf{x}_j$ (yielding $M = J$), $a_{ij} = \|\mathbf{x}_i - \mathbf{c}_j\|^{\tau}$, and $\mathbf{f} = \{f(\mathbf{x}_i)\}_{i=1}^J$. If numerical noise affects the training set then noise reduction is achieved by choosing a number of RBF centers $M$ smaller than the number of training points $J$, and $\mathbf{c}_j$ coordinates are defined via $k$-means clustering (Lloyd, 1982) of the training points. Hence, $w_j$ are determined with least squares regression by solving $\mathbf{w} = (\mathbf{A}^{\top}\mathbf{A})^{-1}\mathbf{A}^{\top}(\mathbf{f} - \text{EV}[\mathbf{f}])$. The optimal number of stochastic RBF centers ($M^*$) is defined by minimizing a leave-one-out cross-validation (LOOCV) metrics (Wackers et al., 2020).
The uncertainty $U_{\tilde{f}}(\mathbf{x})$ associated with the SRBF prediction is quantified by the 95%-confidence band using the cumulative density function of $g(\mathbf{x}, \tau)$ (Volpi et al., 2015).

## 3 Initial training set and bounded surrogate model

In this section a new approach for defining the initial training set is proposed. The new method uses a reduced training set (RS) with only one datapoint for the higher fidelities, as opposed to the authors' previous work where a full training set (FS) was used: the central composite design without factorial points. Both the RS and FS approaches can be used with single- or multi-fidelity methods. Table 1 summarizes the RS and FS approaches for the single- and multi-fidelity cases, respectively.

<table><caption>**Table 1:** Comparison between the new reduced training set (RS) and the full training set (FS)</caption><thead><tr><th>Approach</th><th>N</th><th>Fidelity level</th><th>Number of training points</th><th>Training points placement</th></tr></thead><tbody><tr><td rowspan="2">RS</td><td>1</td><td>1</td><td>1</td><td>Center of domain</td></tr><tr><td>&gt; 1</td><td>from 1 to N - 1<br>N</td><td>1<br>2D+1</td><td>Center of domain<br>Center of domain and of the boundaries</td></tr><tr><td rowspan="2">FS</td><td>1</td><td>1</td><td>2D+1</td><td>Center of domain and of the boundaries</td></tr><tr><td>&gt; 1</td><td>from 1 to N - 1<br>N</td><td>2D+1<br>2D+1</td><td>Center of domain and of the boundaries<br>Center of domain and of the boundaries</td></tr></tbody></table>

The challenge for the RS approach is to create a SRBF surrogate model which can handle extrapolation. Using the RS approach, during the first iterations of the adaptive sampling, the surrogate model prediction is an extrapolation based on the limited training points available. If the SRBF surrogate model gives reasonable predictions in this situation, the rest of the adaptive sampling procedure from section 2 can be used unchanged. However, the SRBF with power kernel has a low accuracy when extrapolating, since the positive power laws are unbounded away from the kernel center. Therefore, a bounded surrogate model prediction (identified with the *B* subscript) is adopted, as described in Algorithm 1.

When one training point is available, the surrogate model prediction is set equal to the function value in the training point $f(\mathbf{x}_1)$. This approach is produced naturally by Eq. 5, where the expected value of the training set is added to the radial basis functions thus providing a non-zero prediction when only one training point is available. When more training points are available, the surrogate model prediction is bounded only in regions of the domain far
from these training points. In Algorithm 1 a sigmoid-like function $s(r)$ is used to provide a smooth transition between the SRBF prediction and the bounded prediction

$$s(r) = \frac{1}{1 + e^{\alpha(r-\gamma)}} \qquad (6)$$

where, for the present work, $\alpha = -75$ and $\gamma = 0.2$. To define $r$, the smallest hyperrectangle (whose edges are parallel to the Cartesian coordinated axis) containing the training points is defined and $r$ is the Euclidean distance of $\mathbf{x}$ from the hyperrectangle boundaries.

Regarding the uncertainty, the definition of $U_{\tilde{\epsilon}_{Bl}}(\mathbf{x})$ in Algorithm 1 stems from the consideration that the error surrogates represent errors in the multi-fidelity approximation $\hat{f}$. Thus, for want of something better, the average of the data itself can be used as a surrogate model prediction uncertainty when an extrapolation is performed. For the single-fidelity $U_{\hat{f}_B}$, the same estimation is adopted. It may be noted that when multiple training points are available, bounding the uncertainty according to $U_{\tilde{\epsilon}_{Bl}}(\mathbf{x}) = \min(U_{\hat{f}}, \text{EV}[\mathbf{f}])$ instead of $U_{\tilde{\epsilon}_{Bl}}(\mathbf{x}) = \min(U_{\hat{f}}, 2\text{EV}[\mathbf{f}])$ would be more consistent with the bounding applied when one training point is available. Such bounding will be adopted in the future.

<table><caption>**Algorithm 1** Bounding of the SRBF prediction and associated uncertainty.</caption><tr><td>if N = 1 then</td><td>▷ Single-fidelity case</td></tr><tr><td>if J = 1 then</td><td>▷ One training point available</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;f̂<sub>B</sub>(x) = f(x<sub>1</sub>)</td><td></td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;U<sub>f̂<sub>B</sub></sub>(x) = f(x<sub>1</sub>)</td><td></td></tr><tr><td>else if J &gt; 1 then</td><td>▷ J training points available</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;f̂<sub>B</sub>(x) = f̂(x) [1 - s(r)] + EV[f]s(r)</td><td></td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;U<sub>f̂<sub>B</sub></sub>(x) = min(U<sub>f̂</sub>, 2EV[f])</td><td></td></tr><tr><td>end if</td><td></td></tr><tr><td>else if N &gt; 1 then</td><td>▷ Multi-fidelity case</td></tr><tr><td>if J<sub>l</sub> = 1, l = 1, ..., N - 1 then</td><td>▷ One training point available</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;ε̃<sub>Bl</sub>(x) = ε<sub>l</sub>(x<sub>1</sub>)</td><td></td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;U<sub>ε̃<sub>Bl</sub></sub>(x) = ε<sub>l</sub>(x<sub>1</sub>)</td><td></td></tr><tr><td>else if J<sub>l</sub> &gt; 1, l = 1, ..., N - 1 then</td><td>▷ J<sub>l</sub> training points available</td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;ε̃<sub>Bl</sub>(x) = ε̃<sub>l</sub>(x) [1 - s<sub>l</sub>(r)] + EV[ε<sub>l</sub>]s<sub>l</sub>(r)</td><td></td></tr><tr><td>&nbsp;&nbsp;&nbsp;&nbsp;U<sub>ε̃<sub>Bl</sub></sub>(x) = min(U<sub>ε̃<sub>l</sub></sub>, 2EV[ε<sub>l</sub>])</td><td></td></tr><tr><td>end if</td><td></td></tr><tr><td>end if</td><td></td></tr></table>

Figure 2 shows the first iteration of the RS and FS approaches for a monodimensional example, see Eqs. 7 where $f_1$ and $f_2$ are introduced in Clark et al. (2016) and $f_3$ in Wackers et al. (2020). In the single-fidelity case the RS approach yields a constant surrogate model prediction and associated uncertainty, see Figure 2a. In the multi-fidelity case, the RS approach produces a multi-fidelity prediction that only relies on the low-fidelity surrogate model to
8 Efficient Initialization for Multi-Fidelity Surrogate-Based Optimization

approximate the trend of the desired function. Differently, the FS approach takes advantage also of the medium- and high-fidelity evaluations. The RS approach produces an obviously less accurate prediction in the first iteration but allows to preserve a reasonable trend while using only one evaluation of the medium- and high-fidelity.

<figure><img src="image_2.png" /><figcaption>Fig. 2: Example of the two initialization strategies for the mono-dimensional function. The dots indicate initial data points, the lines indicate the surrogate models, and the gray zone is the uncertainty band</figcaption></figure>

$$
\begin{align}
f_1(x) &= \sin(30(x - 0.9)^4) \cos(2(x - 0.9)) + (x - 0.9)/2 \nonumber \\
f_2(x) &= (f_1(x) - 1 + x)/(1 + 0.25x) \nonumber \\
f_3(x) &= \sin(20(x - 0.87)^4) \cos(2(x - 0.87)) + \nonumber \\
&\quad (x - 0.87)/2 - (2.5 - (0.7x - 0.14)^2) + 2x \tag{7}
\end{align}
$$
# 4 Test problems

## 4.1 Analytical test problem

The analytical test problem is based on the Rosenbrock function with two variables, three fidelities are considered ($N = 3$). The Rosenbrock function is a function defined using two parameters $a$ and $b$. The high-fidelity function $f_1$ is a Rosenbrock function with $a = 1$ and $b = 100$. The medium-fidelity function $f_2$ is the additive combination of a Rosenbrock function, with $a = -2$ and $b = 50$, with a second order function, see Eq. 8. Finally, the low-fidelity function $f_3$ is a transformation of $f_1$ by addition and multiplication of first order functions (Rumpfkeil and Beran, 2020):

$$ 
\begin{aligned}
f_1(\mathbf{x}) &= 100(x_2 - x_1^2)^2 + (1 - x_1)^2 \\
f_2(\mathbf{x}) &= 50(x_2 - x_1^2)^2 + (-2 - x_1)^2 - 80 - 0.5x_10.5x_2 \\
f_3(\mathbf{x}) &= (f_1(x_1, x_2) - 4 - 0.5x_1 - 0.5x_2)/(10 + 0.25x_1 + 0.25x_2)
\end{aligned} 
\quad (8) $$

where $\mathbf{x} \in [-2.048, 2.048] \otimes [-2.048, 2.048]$. The functions are represented in Figure 3, they are challenging from the optimization viewpoint since the minimum is located in a narrow and flat valley. Specifically, the minimum of $f_1$ is in $\mathbf{x} = \{1, 1\}$ and equal to 0.

A synthetic computational cost is associated with the evaluation of each fidelity, based on typical cost ratios for CFD simulations in naval applications. The selected computational cost ratios are equal to $\beta = \{1, 0.1, 0.05\}$ (Mainini et al., 2022).

<figure><img src="image_3.png" /><figcaption>Fig. 3: Test problem, the Rosenbrock function with three fidelity levels</figcaption></figure>

## 4.2 NACA hydrofoil optimization

This problem addresses the drag coefficient minimization of a NACA four-digit airfoil. The following minimization problem is solved

$$ \text{minimize } f(\mathbf{x}) = C_D(\mathbf{x}) \text{ subject to } C_L(\mathbf{x}) = 0.6 \text{ and to } \mathbf{l} \le \mathbf{x} \le \mathbf{u} \quad (9) $$
10 *Efficient Initialization for Multi-Fidelity Surrogate-Based Optimization*

where **x** is the design variables vector; **l** and **u** are the lower and upper bounds of the design space, respectively; $C_D$ and $C_L$ are the drag and lift coefficient, respectively. The equality constraint on the lift coefficient is necessary in order to compare different geometries at the same lift force (equal to the weight of the object), since the drag depends strongly on the lift. The simulation conditions are: velocity $U = 10$ m/s, chord $c = 1$ m, and fluid density $\rho = 1,026$ kg/m$^3$, with a chord based Reynolds number $Re = 8.41 \cdot 10^6$. The hydrofoil shape is defined by the general equation for four-digit NACA foils (Moran, 2003). In this work the design variables vector is defined as $\mathbf{x} = \{t, m\}$, where $t \in [0.030, 0.120]$ is the maximum thickness and $m \in [0.025, 0.065]$ is the maximum camber value with the maximum camber position fixed at $p = 0.4$. Tests are run with one and three fidelity levels ($N = 1, 3$), see Figure 4.

<figure><img src="image_4.png" /><figcaption>**Fig. 4:** NACA hydrofoil computational grids for ISIS-CFD: (a) Fine grid, 12.8k cells, (b) Medium grid, 5.7k cells, and (c) Coarse grid, 3.6k cells</figcaption></figure>

Numerical simulations are performed with the ISIS-CFD incompressible unstructured finite-volume RANS solver for multi fluid flow. In ISIS-CFD, the velocity field is obtained from the momentum conservation equations and the pressure field is extracted from the mass conservation constraint transformed into a pressure equation. These equations are similar to the SIMPLE method (Rhie and Chow, 1983), but have been adapted for flows with discontinuous density fields. The unstructured discretization is face-based. While all unknown state variables are cell-centered, the systems of equations used in the implicit time stepping procedure are constructed face by face. Therefore, cells with an arbitrary number of arbitrarily-shaped constitutive faces are accepted. The code is fully parallel using the message passing interface protocol. A detailed description of the solver is given by Queutey and Visonneau (2007).

Computational grids are created through adaptive grid refinement (AGR) (Wackers et al., 2014, 2017), to optimize the efficiency of the solver and to simplify the automatic creation of suitable grids. The AGR method adjusts the computational grid locally, during the computation, by dividing the cells of an original coarse grid. The decision where to refine comes from a refinement criterion, a tensor field $C(x, y, z)$ computed from the flow. The tensor is based on second derivatives of pressure and velocity, which gives an approximate indication of the local truncation errors. The grid is refined until the dimensions
$d_{p,j}$ ($j = 1, 2, 3$) of each hexahedral cell $p$ satisfy $||C_p d_{p,j}|| = T_r$. The refinement criterion based on the second derivatives of the flow is not very sensitive to grid refinement (Wackers et al., 2017), so the cell sizes everywhere are proportional to the constant threshold $T_r$.

For the MF optimization, grid adaptation is used to take into account the need for several fidelities. The interest of this procedure is that different fidelity results can be obtained by running the same simulations and simply changing the threshold $T_r$. Thus, it is straightforward to automate the MF simulations. Highest- to lowest-fidelity simulations require about 17, 9, and 5 minutes of wall-clock time to converge, respectively. The resulting computational cost ratios are equal to $\beta = \{1, 0.5, 0.3\}$.

# 5 Numerical results

The performance of the RS approach is compared with the performance of the FS approach at a fixed computational budget equivalent to 45 high-fidelity simulations, for both problems. Even though the analytical problem is cheap and could be iterated to convergence, this limited budget is chosen to resem- ble realistic simulation-driven optimization problems where the training set evaluation is expensive.

## 5.1 Analytical test problem

The performance of the method applied for the analytical test is assessed by the convergence of the predicted minimum value ($\hat{f}_1(\mathbf{x}_{min})$), its validation by a high-fidelity evaluation ($f_1(\mathbf{x}_{min})$), the prediction uncertainty in the minimum point ($U_{\hat{f}_1}(\mathbf{x}_{min})$), and the quantification of the effectiveness of the surrogate model in identifying the location of the minimum in the variable space by the computation of the location error ($E_x$):

$$E_x = \sqrt{\frac{1}{D} \sum_{j=1}^{D} \left( \frac{x_{min,j} - \tilde{x}_j}{u_j - l_j} \right)^2} \quad (10)$$

where $\mathbf{x}_{min}$ is the position found by the surrogate model-based optimization, $\tilde{\mathbf{x}}$ if the position of the true reference minimum, and $l_j$ and $u_j$ (for $j = 1, \dots, D$) are the lower and the upper bounds of the variables domain, respectively.

The evaluation of the analytical test is not affected by numerical noise, therefore an interpolative formulation of the SRBF is used. The computational costs associated with the initial RS training sets are 1 for $N = 1$ and 1.35 for $N = 3$, while for FS the costs are 5 and 5.75, respectively.

Figure 5 shows the response surfaces with the training sets for the two initialization approaches using one fidelity. The response surface of the RS method is a constant flat surface whereas for the FS approach a concavity is already evident, see Figures 5a and 5b. Figures 5c and 5d show the response surface for a normalized computational cost equal to 10: the RS method has
<figure><img src="image_5.png" /><figcaption>**Fig. 5:** Analytical problem, results for the single fidelity case. Response surfaces and training sets for the two initialization approaches at the first iteration, for a computational cost (CC) = 10, and at the last iteration of the adaptive sampling</figcaption></figure>

already performed some exploration of the neighborhood of the region of the minimum, whereas the FS method distributed new samples in a much smaller area. At the last iteration, the sampling method with RS has explored the entire valley where the minimum is located whereas with FS the sampling
<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>**Fig. 6:** Analytical problem, results for the single fidelity case. From top to bottom: convergence of $\hat{f}_1(\mathbf{x}_{\min})$, $U_{\hat{f}_1}(\mathbf{x}_{\min})$, and $f_1(\mathbf{x}_{\min})$ for FS (a) and RS (b) approaches (error bars: $U_{\hat{f}_1}(\mathbf{x}_{\min})$); convergence of $E_x$ (c)</figcaption></figure>

method has partially explored the valley and has then focused on the opposite side, which creates a less reliable surrogate model (see Figures 5c and 5d).

Figure 6 shows the convergence of the surrogate model prediction and uncertainty $\hat{f}_1(\mathbf{x}_{\min})$, $U_{\hat{f}_1}(\mathbf{x}_{\min})$, and the true value in the predicted optimum $f_1(\mathbf{x}_{\min})$, and the position error $E_x$, for the two initialization approaches using one fidelity. It is worth noting that for both approaches, the uncertainty estimation is reliable: the validated minimum almost always lies within the
<figure><img src="image_9.png" /><figcaption>Fig. 7: Analytical problem, results for the single fidelity case. Training set size</figcaption></figure>

multi-fidelity model prediction uncertainty. Furthermore, although the RS is initially less accurate than the FS approach (since it is not really intended for single-fidelity use), it eventually converges better: the predicted and validated optima are more accurate, with lower uncertainty, and more stable behavior when new samples are added. Figure 7 shows that the sampling method with RS performs more sampling iterations than with FS: the reduced size of the startset makes available more simulations for the adaptive procedure.

Figure 8 shows the response surfaces with the training sets for the two initialization approaches using three fidelities. At the first iteration the surrogate model prediction with the RS approach shows the same trend as with the FS approach, see Figs. 8a-b. Figures 8c and 8d show the response surface for a normalized computational cost equal to 10, where the RS method has already started to explore the region of the minimum, whereas the FS method with its more expensive startset has only distributed a few new samples in a much smaller area. At the last iteration, the sampling method with RS (Figure 8c) used low-fidelity data to explore the entire valley where the minimum is located. Thanks to this exploration, the high-fidelity points are placed intelligently, in the center of the valley. Differently, the sampling method with FS (Figure 8d) has only partially explored the valley and did not perform much exploration in the rest of the domain; the HF points are of little use since most points are placed far away from the optimum or scattered throughout the sampled region. Finally, sampling occurred on only one side of the optimum position, which reduces the accuracy of the surrogate model. Therefore, RS has produced a more effective sampling.

Figure 9 shows the convergence of the surrogate model prediction and uncertainty $\hat{f}_1(\mathbf{x}_{\min})$, $U_{\hat{f}_1}(\mathbf{x}_{\min})$, the true value in the predicted optimum $f_1(\mathbf{x}_{\min})$, and the position error $E_x$, for the two initialization approaches using three fidelities. While several oscillations are present, the RS approach achieves better results than FS: it predicts an optimum earlier and it converges deeper, because of the more efficient exploration. Furthermore, the difference between the predicted and the validated minimum is always within the prediction uncertainty. Also, comparing with Figure 6 shows that the three-fidelity results are
<figure><img src="image_10.png" /><figcaption>**Fig. 8:** Analytical problem, results for the three fidelity case. Response surfaces and training sets for the two initialization approaches at the first iteration, for a computational cost (CC) = 10, and at the last iteration of the adaptive sampling</figcaption></figure>

better than the single-fidelity optimization. Figures 10a-c show that the sampling method with RS performs a significantly higher number of iterations than with the FS approach, thanks to the smaller computational cost spent for the evaluation of the initial training set. Furthermore, the sampling method
<figure><img src="image_11.png" /><figcaption>Fig. 10: Analytical problem, results for the three fidelity case. Training sets size</figcaption></figure>

is used to quantify the accuracy of the minimum predicted by the surrogate model $\hat{f}_1(\mathbf{x}_{\text{min}})$ and the validation error ($E_v$) is used to quantify the error in the identification of the true minimum $f_1(\tilde{\mathbf{x}})$ in the function space:

$$E_p = \left| \frac{\hat{f}_1(\mathbf{x}_{\text{min}}) - f_1(\mathbf{x}_{\text{min}})}{R_1} \right| \quad (11)$$
$$E_v = \left| \frac{f_1(\mathbf{x}_{\min}) - f_1(\tilde{\mathbf{x}})}{R_1} \right| \qquad (12)$$

where $\mathbf{x}_{\min}$ is the position found by the surrogate model-based optimization, and $R_1$ is the reference high-fidelity function range.

The challenge of this problem for the surrogate model optimization is that the data are strongly affected by numerical noise (Wackers et al., 2020), especially the low-fidelity evaluations. Figure 11 provides an illustration of the typical noise levels. Therefore a least-squares regressive approach, as described in section 2.2, is used to fit the surrogate models. The data pertaining to the NACA hydrofoil optimization with the sampling method and the FS approach are taken from Wackers et al. (2020). The computational costs associated with the initial RS training sets are 1.0 for $N = 1$ and 3.0 for $N = 3$, while for FS the costs are 5.0 and 9.2, respectively. The $R_1$ value is the high-fidelity function range of the FS initial training set and is equal to $1.523E-3$.

Figure 12 shows the response surfaces with the training sets for the two initialization approaches, using one fidelity. Again, the response surface of the RS method is a constant flat surface whereas the FS approach produces a concavity, see Figs. 12a and 12b. Figures 12c and 12d show the response surface for a computational cost equal to 10, the FS approach has identified the region of the minimum whereas the RS approach does not. Figures 12e and 12f show the last iteration of the adaptive sampling: the RS approach has explored around the region of the minimum with only one training point in the opposite side of the domain (in $\{1, 1\}$), whereas the FS approach has also sampled the upper corners, which are less relevant for the optimization.

Figure 13 presents the response surfaces for three fidelities. The two response surfaces at the first iteration are significantly different, see Figures 13a and 13b: the RS surrogate model is reasonable thanks to the LF contribution, while the FS model misinterprets the data. Figures 13c and 13d show the response surface for a computational cost equal to 10, the FS approach has added three low-fidelity samples without actually changing the prediction due to the influence of the high-fidelity initial training set. Differently, the RS

<figure><img src="image_12.png" /><figcaption>Fig. 11: NACA hydrofoil optimization problem. Training sets in the neighborhood of $x_2 = 0$</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>(a) RS, first iteration</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>(b) FS, first iteration</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>(c) RS, CC = 10</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>(d) FS, CC = 10</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>(e) RS, last iteration</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>(f) FS, last iteration</figcaption></figure>

**Fig. 12:** NACA hydrofoil optimization problem, results for the single fidelity case. Response surfaces and training sets for the two initialization approaches at the first iteration, for a computational cost (CC) = 10, and at the last iteration of the adaptive sampling

approach identified the region of the minimum. Figures 13e and 13f show the last iteration of the adaptive sampling method. The sampling method with the RS approach did not sample the upper region of the domain, where the drag coefficient is higher. However, like for the analytical test, the region around
<figure><img src="image_19.png" /><figcaption>(a) RS, first iteration</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>(b) FS, first iteration</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>(c) RS, CC = 10</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>(d) FS, CC = 10</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>(e) RS, last iteration</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>(f) FS, last iteration</figcaption></figure>

**Fig. 13:** NACA hydrofoil optimization problem, results for the three fidelity case. Response surfaces and training sets for the two initialization approaches at the first iteration, for a computational cost (CC) = 10, and at the last iteration of the adaptive sampling

the optimum is thoroughly explored using low fidelity data. Differently, the sampling method with the FS approach sampled also the upper corners of the domain, while most of the effort is spent on medium-fidelity points very close to the perceived optimum. Although the response surface is mostly different
<table><caption>**Table 2:** NACA hydrofoil optimization problem, summary of the results</caption><thead><tr><th>Initial set</th><th>N</th><th>x<sub>1</sub></th><th>x<sub>2</sub></th><th>f&#x0302;(x<sub>min</sub>)</th><th>U<sub>f</sub>(x<sub>min</sub>)/R<sub>1</sub>%</th><th>f(x<sub>min</sub>)</th><th>E<sub>p</sub>%</th><th>E<sub>v</sub>%</th><th>E<sub>x</sub>%</th></tr></thead><tbody><tr><td>Reduced</td><td>1</td><td>0.4173</td><td>0.0000</td><td>7.213E-3</td><td>0.06</td><td>7.219E-3</td><td>0.41</td><td>0.52</td><td>2.81</td></tr><tr><td>Full</td><td>1</td><td>0.3799</td><td>0.0000</td><td>7.213E-3</td><td>0.03</td><td>7.214E-3</td><td>0.01</td><td>0.03</td><td>0.16</td></tr><tr><td>Reduced</td><td>3</td><td>0.4354</td><td>0.0000</td><td>7.050E-3</td><td>2.91</td><td>7.224E-3</td><td>11.4</td><td>0.85</td><td>4.09</td></tr><tr><td>Full</td><td>3</td><td>0.3615</td><td>0.0000</td><td>7.162E-3</td><td>1.90</td><td>7.218E-3</td><td>0.78</td><td>0.09</td><td>1.14</td></tr><tr><td>Reference</td><td>1</td><td>0.3776</td><td>0.0000</td><td>-</td><td>-</td><td>7.212E-3</td><td>-</td><td>-</td><td>-</td></tr></tbody></table>

for the two approaches, the region of the minimum is correctly identified by both, which is essential for optimization where the only goal is to get closest to the optimum. The RS approach only adds one high-fidelity simulation to the initial training set, evaluating it in the center of the region of the minimum. Differently, the FS approach adds two high-fidelity simulations to the initial training set, evaluating them in the two corners of the lower region of the domain.

Table 2 summarizes the minimum position, value, and associated uncertainty and the $E_p$ and $E_v$ metrics. While the single-fidelity RS approach performs well, the table shows that the multi-fidelity RS result is problematic. Even though all the verified minima have errors below 1% with respect to the reference, the prediction error for multi-fidelity RS is high: the surrogate model values in the optimum are not reliable. This is also confirmed by the value of the normalized prediction uncertainty ($U_f/R_1$) which is higher for the multi-fidelity RS.

The reason for this is shown in Figure 13c. Compared to the single-fidelity result 12c, the multi-fidelity surrogate model is more irregular. This implies that, despite the noise filtering, the adaptive surrogate model has overfitted the low-fidelity training data, which are the most affected from the numerical noise. Table 3 summarizes the training sets size and the number of RBF centers at the last iteration. The RS approach uses a larger number of RBF centers for the low-fidelity training set in comparison to the FS approach. This results in the RS approach being more affected from the numerical noise by overfitting. Therefore, the surrogate model includes noise so it is not accurate. For FS, the low-fidelity surrogate model has ended up with fewer RBF centers (see section 2.2): 20, compared with 89 for RS. Thus, the noise filtering is more effective for the FS approach, which explains the better results.

<table><caption>**Table 3:** NACA hydrofoil optimization problem, summary of the training sets and RBF centers at the final iteration</caption><thead><tr><th>Initial set</th><th>N</th><th>J<sub>1</sub>(M<sub>1</sub><sup>*</sup>)</th><th>J<sub>2</sub>(M<sub>2</sub><sup>*</sup>)</th><th>J<sub>3</sub>(M<sub>3</sub><sup>*</sup>)</th></tr></thead><tbody><tr><td>Reduced</td><td>1</td><td>45 (32)</td><td>-</td><td>-</td></tr><tr><td>Full</td><td>1</td><td>45 (36)</td><td>-</td><td>-</td></tr><tr><td>Reduced</td><td>3</td><td>2 (2)</td><td>22 (19)</td><td>105 (89)</td></tr><tr><td>Full</td><td>3</td><td>7 (9)</td><td>19 (18)</td><td>96 (20)</td></tr></tbody></table>
However, this problem with the noise filtering is not directly related with the startsets. Qualitative assessment of the sample points placement shows that the RS initialization has the desired effect of good exploration with sparing but intelligent use of high-fidelity samples.

# 6 Conclusions

Multi-fidelity adaptive sampling is generally started from a non-adaptive initial set of data, which can be expensive to compute for the higher fidelities. Furthermore, in a surrogate model-based optimization problem high-fidelity points far from the optimum do not provide useful information. Therefore, a new initialization approach is defined here with a reduced number of high-fidelity data points. The authors’ previous initialization approach uses a full training set (FS) of $2D + 1$ samples for all the fidelities, where $D$ is the design space dimension. In the multi-fidelity case, the new reduced training set (RS) approach uses $2D + 1$ points as training set only for the lowest (and cheapest to evaluate) fidelity and one single point as initial training set for the other fidelities. A single-fidelity case is also defined which uses only one point as the initial training set. Reducing the computational cost of the initialization means that a larger computational budget is available for an adaptive sampling of the design space.

The proposed startset can be used with the authors’ existing adaptive sampling method, as long as the surrogate model is able to extrapolate from a limited number of data points. This capability is developed here for stochastic radial basis functions with a power kernel (SRBF). Since the power kernel used in standard SRBF lacks a compact support, a bound is imposed on both the surrogate model prediction and the associated uncertainty when an extrapolation is performed, to improve the adaptive sampling. For the bounding, the prediction is gradually blended towards the mean of the data, which is also used as an upper bound of the uncertainty (this is deemed acceptable since the multi-fidelity discrepancy surrogate models represent errors in the prediction of the function).

While single-fidelity performance is not the main objective of the RS startset, the different tests show that the single-fidelity RS approach works well. The distribution of the data points is sensible, with a good balance between exploration and exploitation, while the optimization results are comparable with those obtained from FS startsets. This proves that the adaptive sampling is capable of creating a good surrogate model from a single datapoint, without needing a full startset.

For multi-fidelity data, the first effect of RS is to reduce the costs of the startset and to limit the number of HF points with little practical value. In the tests, the RS approach reduced the computational cost of evaluating the initial training set by 76% for the analytical problem and by 67% for the NACA hydrofoil optimization. Furthermore, it is shown that RS provides a first useful estimate of the optimum at a far lower cost than FS and that, in
all cases, more computations were available for the adaptive sampling. Finally, in the data, the number of HF points far away from the optimum is reduced. Thus, this first objective is met.

A second observation from the tests is that RS leads to more effective exploration of the domain. Moreover, this exploration is mainly performed with low- and medium-fidelity data. This is due to the bounding of the uncertainty, which makes the method less interested in exploring the extrema of the highest-fidelity error surrogate model. On the contrary, HF points are added where they are most useful: near the perceived optimum. For the Rosenbrock function, the HF points fill the minimum valley, while for the NACA hydrofoil optimization problem, the adaptive sampling method adds only one high-fidelity sample in the region of the minimum. This property is of great importance for optimization in higher dimensions, where the exploration of the large design spaces is very costly. Therefore, exploration with cheap low-fidelity data is crucial.

However, the NACA case also shows that the startset is not the only relevant aspect of a multi-fidelity surrogate model. Indeed, the high errors for this case are likely due to the procedure used for the selection of the RBF centers, which is not effective in filtering out the numerical noise, negatively affecting the adaptive sampling procedure and leading to overfitting.

Future work aims to investigate the effects of the RS approach on the adaptive sampling using several analytical test problems. Specifically, the analytical test set defined within the NATO STO Research Task Group AVT-331 on “Goal-driven, multi-fidelity approaches for military vehicle system-level design” will be considered. Furthermore, an improved strategy to perform regression with the SRBF in presence of numerical noise will be investigated in order to avoid overfitting. Still, the current tests indicate that the reduced startset and bounded surrogate model form a useful component of an optimal surrogate modeling strategy.

**Acknowledgments.** CNR-INM is partially supported by the Office of Naval Research through NICOP grant N62909-18-1-2033, administered by Dr. Woei-Min Lin, Dr. Elena McCarthy, and Dr. Salahuddin Ahmed of the Office of Naval Research and Office of Naval Research Global. The work is conducted in collaboration with NATO STO Research Task Group AVT-331 on “Goal-driven, multi-fidelity approaches for military vehicle system-level design”.

## References

Bonfiglio, L., P. Perdikaris, S. Brizzolara, and G. Karniadakis. 2018. Multi-fidelity optimization of super-cavitating hydrofoils. *Computer Methods in Applied Mechanics and Engineering* 332: 63–85.

Clark, D.L., H.R. Bae, K. Gobal, and R. Penmetsa 2016. Engineering design exploration utilizing locally-optimized covariance kriging. In *18th AIAA Non-Deterministic Approaches Conference*, pp. 0428.
Cox, D.D. and S. John 1992. A statistical method for global optimization. In *IEEE International Conference on Systems, Man, and Cybernetics*, Chicago, IL, pp. 1241–1246.

Demo, N., M. Tezzele, A. Mola, and G. Rozza. 2021. Hull shape design optimization with parameter space and model reductions, and self-learning mesh morphing. *Journal of Marine Science and Engineering* 9(2).

Fernández-Godino, G., M. Chanyoung Park, N.H. Kim, and R.T. Haftka. 2019. Issues in deciding whether to use multifidelity surrogates. *AIAA Journal* 57(5): 2039–2054.

Furcas, F., G. Vernengo, D. Villa, and S. Gaggero. 2020. Design of wake equalizing ducts using RANSE-based SBDO. *Applied Ocean Research 97*: 102087.

Harries, S. and C. Abt. 2019. Faster turn-around times for the design and optimization of functional surfaces. *Ocean Engineering* 193: 106470.

Jin, R., W. Chen, and T.W. Simpson. 2001. Comparative studies of meta-modelling techniques under multiple modelling criteria. *Structural and multidisciplinary optimization* 23(1): 1–13.

Karlberg, M., M. Löfstrand, S. Sandberg, and M. Lundin. 2013. State of the art in simulation-driven design. *International Journal of Product Development* 18(1): 68–87.

Lloyd, S. 1982. Least squares quantization in PCM. *IEEE transactions on information theory* 28(2): 129–137.

Mainini, L., A. Serani, M.P. Rumpfkeil, E. Minisci, D. Quagliarella, H. Pehli-van, S. Yildiz, S. Ficini, R. Pellegrini, F. Di Fiore, et al. 2022. Analytical benchmark problems for multifidelity optimization methods. *arXiv preprint arXiv:2204.07867*.

Moran, J. 2003. *An Introduction to Theoretical and Computational Aerodynamics*. Courier Corporation.

Nuñez, L., R.G. Regis, and K. Varela. 2018. Accelerated random search for constrained global optimization assisted by radial basis function surrogates. *Journal of Computational and Applied Mathematics* 340: 276–295.

Peherstorfer, B., K. Willcox, and M. Gunzburger. 2018. Survey of multifidelity methods in uncertainty propagation, inference, and optimization. *SIAM Review* 60(3): 550–591.
Pellegrini, R., A. Serani, R. Broglia, M. Diez, and S. Harries 2018. Resistance and payload optimization of a sea vehicle by adaptive multi-fidelity metamodeling. In *2018 AIAA/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference*, pp. 1904.

Queutey, P. and M. Visonneau. 2007. An interface capturing method for free-surface hydrodynamic flows. *Computers & Fluids* 36(9): 1481-1510.

Regis, R.G. 2020. Large-scale discrete constrained black-box optimization using radial basis functions. In *2020 IEEE Symposium Series on Computational Intelligence (SSCI)*, pp. 2924-2931.

Rhie, C.M. and W.L. Chow. 1983. A numerical study of the turbulent flow past an isolated airfoil with trailing edge separation. *AIAA Journal* 17: 1525-1532.

Rumpfkeil, M.P. and P.S. Beran 2020. Multi-fidelity, gradient-enhanced, and locally optimized sparse polynomial chaos and kriging surrogate models applied to test problems. In *AIAA Scitech Forum*, online.

Serani, A., C. Leotardi, U. Iemma, E.F. Campana, G. Fasano, and M. Diez. 2016. Parameter selection in synchronous and asynchronous deterministic particle swarm optimization for ship hydrodynamics problems. *Applied Soft Computing* 49: 313-334.

Serani, A., R. Pellegrini, M. Diez, J. Wackers, C.E. Jeanson, P. Queutey, and M. Visonneau. 2019. Adaptive multi-fidelity sampling for CFD-based optimization via radial basis function metamodels. *International Journal of Computational Fluid Dynamics* 33(6-7): 237-255.

Serani, A., F. Stern, E.F. Campana, and M. Diez. 2022. Hull-form stochastic optimization via computational-cost reduction methods. *Engineering with Computers* 38: 2245-2269.

Volpi, S., M. Diez, N. Gaul, H. Song, U. Iemma, K.K. Choi, E.F. Campana, and F. Stern. 2015. Development and validation of a dynamic metamodel based on stochastic radial basis functions and uncertainty quantification. *Structural and Multidisciplinary Optimization* 51(2): 347-368.

Wackers, J., G.B. Deng, E. Guilmineau, A. Leroyer, P. Queutey, and M. Visonneau. 2014. Combined refinement criteria for anisotropic grid refinement in free-surface flow simulation. *Computers & Fluids* 92: 209-222.

Wackers, J., G.B. Deng, E. Guilmineau, A. Leroyer, P. Queutey, M. Visonneau, A. Palmieri, and A. Liverani. 2017. Can adaptive grid refinement produce grid-independent solutions for incompressible flows? *Journal of Computational Physics* 344: 364-380.