## 3 Novel Approach to Computing Infill Criteria for MOEGO
### A Brief Description of the Monte Carlo Approach  
For the sake of easy implementation, the Monte Carlo approach is considered to calculate the HVI or EHVI in the present work. According to the principle of the theory [35], the HVI (see Figure 1 and Equation (10)) could be calculated via the hit-or-miss method (also called the rejection method) [36]. The basic structure of the Monte Carlo-based HVI (MCHVI) is presented for obtaining the HVI, as shown in Algorithm 1. ```
Input: Current set of Pareto front points \(\boldsymbol{\mathcal{P}}_{n}\), Predicted updating point \(\mathbf{y}\), reference point \(\boldsymbol{\mathcal{R}}\). Output: Hypervolume improvement \(I_{H}\left(\mathbf{y}^{\prime},\boldsymbol{\mathcal{P}}_{n}\right)\)
1 Initialize the Monte Carlo sampling hypercube: \(\left[\min(y_{i}^{\prime},y_{z}^{\prime},\cdots y_{s}^{\prime},r),\max(y_{i }^{\prime},y_{z}^{\prime},\cdots y_{s}^{\prime},r)\right]\), \(i=1,2,\cdots,m\);
3 Generate a set of \(n_{p}\) sample points: \(\boldsymbol{P}\)=\(



work, the value of the MCHVI indicator is approximated by counting the number \(n_{q}\) of the sampling points located in the domain of dependence (red area in Figure 3) by

\[I_{\mathcal{H}}\big{(}\mathbf{y}^{\prime},\ \mathcal{P}_{n}\big{)}\approx\frac{n_{q} }{n_{p}}\times V \tag{14}\]

where \(n_{p}\) is the number of Monte Carlo sampling points (see line 3 of Algorithm 1), and \(V\) is the volume of the hypercube. Such counting-like simple operations are light in programed implementation, and it can be imagined that the value of the MCHVI indicator is proportional to the size of the Monte Carlo sampling points. However, in order to avoid misguidance caused by the approximate MCHVI indicator, the number of Monte Carlo sapling points is preferred to be as large as possible, which in turn requires massive arithmetic operations associated with predicting the value of the MCHVI indicator will be required. This situation could be alleviated through GPU accelerating, which will be mainly considered and addressed in Section 3.2. Once the value of HVI indicator \(I_{\mathcal{H}}\) is obtained, the value of EHVI indicator \(EI_{\mathcal{H}}\) can be approximated by the generalized form of the Monte Carlo estimator [37], which can be expressed as

\[EI_{\mathcal{H}}(\mathbf{x})\approx\frac{1}{N}\sum_{i=1}^{N}I_{\mathcal{H}}( \mathbf{y}^{mc}_{i},\ \mathcal{P}_{n}) \tag{15}\]

where \(\{\mathbf{y}^{mc}_{1},\ \mathbf{y}^{mc}_{2},\cdots,\mathbf{y}^{mc}_{N}\}\) are the \(N\) random samples of the Gaussian distribution with mean \(\hat{y}\) and standard deviation \(s\). The number \(N\) is a user-defined parameter, and in this paper, following the work of [38], a fixed number of 1000 is used for all test cases presented in Section 4. The \(I_{\mathcal{H}}(\mathbf{y}^{mc}_{i},\ \mathcal{P}_{n})\) is predicted by the MCHVI algorithm (see Algorithm 1), which is the main time-consuming part of the EHVI indicator. and hence will be implemented in GPU architecture, as discussed in the next section. ### GPU-Accelerated Infill Criterion for the MOEGO Algorithm

In this section, GPU implementation is presented for computing the value of the MCHVI indicator. As described in Algorithm 1, a mass of independent arithmetic operations associated with the dominated or non-dominated identifications (see Algorithm 1, line 6) are proved to be time-consuming. Fortunately, such tasks are mostly weak-dependent compute-intensive and are very suitable for GPU parallel architecture [39; 40; 41]. Therefore, such a kind of computation is implemented on the GPU to achieve acceleration. In order to obtain a global view of our implementation, the corresponding general framework of the program is given in Figure 4. It should be noted that necessary computing data, including the coordinates of all Monte Carlo sample points, Pareto front points, and updating points, are designed to be prepared on the CPU side and then sent to the GPU side before starting

Figure 3: Example of HVI calculation by the Monte Carlo-based method with two objectives.