## 2 Stochastic Inequality Constrained Optimization Problems  
### Mathematical Formulation  
A typical SICOP can be stated formally as shown below [1].  
\[\min f(\mathbf{x}) \tag{1}\]  
subject to  
\[E[g_{j}(\mathbf{x})]\geq b_{j},\ j=1,\ldots,J \tag{2}\]  
\[\mathbf{V}\leq\mathbf{x}\leq\mathbf{U} \tag{3}\]where \(\mathbf{x}=[x_{1},\ldots,x_{K}]\) is a \(K\)-dimensional solution vector, \(f(\mathbf{x})\) denotes the objective function, \(b_{j}\) is a service level, \(E[g_{j}(\mathbf{x})]\) is the expected value of the \(j\)th constrained function, \(J\) is the number of constraints, \(\mathbf{V}=[V_{1},\ldots,V_{K}]\) is the lower bound of a solution vector, and \(\mathbf{U}=[U_{1},\ldots,U_{K}]\) denotes the upper bound of a solution vector. Sufficient simulation replications should be carried out to assure a precise estimate of \(E[g_{j}(\mathbf{x})]\). Let \(L\) denote the number of simulation replications and \(g_{j}^{\ell}(\mathbf{x})\) represent the measure of the \(\ell\)th simulation replication. A straightforward approach to estimate \(E[g_{j}(\mathbf{x})]\) is the sample mean that is computed as follows.  
\[\overline{g}_{j}(\mathbf{x})=\frac{1}{L}\sum_{\ell=1}^{L}g_{j}^{\ell}( \mathbf{x}) \tag{4}\]  
When the value of \(L\) increases, \(\overline{g}_{j}(\mathbf{x})\) obtains a better estimate of \(E[g_{j}(\mathbf{x})]\). There are two important issues that have been concerned in efficiency of stochastic simulation optimization. (i) \(L\) must be large if a sound estimate of \(E[g_{j}(\mathbf{x})]\) is required, which has a decisive influence on the optimum that an optimization technique can find. (ii) \(\overline{g}_{j}(\mathbf{x})\) must be evaluated for various solutions to determine the best solution.  
Since the stochastic inequality constraints are soft, the constrained problem (1)-(2) can be converted to an unconstrained one by addition of a quadratic penalty function [13].  
\[\min F(\mathbf{x})=\eta\times f(\mathbf{x})+(1-\eta)\times\sum_{j=1}^{J}PF_{j} (\mathbf{x}) \tag{5}\]  
where \(\eta\in(0,1)\) denotes an integrated factor, \(F(\mathbf{x})\) represents an integrated objective function, and \(PF_{j}(\mathbf{x})\) is a steep penalty function. The integrated objective function of (5) comprises two terms, the objective function (1) and the violated constraint functions in (2).  
\[PF_{j}(\mathbf{x})=\left\{\begin{array}{c}0,if\ \overline{g}_{j}(\mathbf{x}) \geq b_{j}\\ \left(\overline{g}_{j}(\mathbf{x})\ -\ b_{j}\right)^{2}\times 10^{4},else \end{array}\right. \tag{6}\]  
The penalty function is developed with a steep jump to ensure that the sample mean satisfies the service level. Since the objective function is penalized by addition of a penalty function, a large value of amplification \(10^{4}\) for violated service level constraints is to guarantee an obvious discontinuity of the quadratic penalty function. Let \(L_{a}\) indicate the large value of \(L\), and the accurate evaluation of (4) is referred to as \(L=L_{a}\). An accurate evaluation is defined as an evaluation that is tolerant of a small noise and modeling error. For simplicity, let \(F_{a}(\mathbf{x})\) represent the integrated objective value of a \(\mathbf{x}\) using accurate evaluation.  
As pointed out by the OO method [11], order of a solution is likely kept even it is evaluated using a metamodel. To choose \(N\) significant solutions from solution space within an acceptable computation time, we should construct a metamodel to approximate the integrated objective function of a solution more rapidly and adopt an efficient search method assisted by this metamodel. The metamodel is based on the RMTSs [12], and the search method is the IEHO.