A decoupling algorithm with first-order asymptotic integration for reliability-based design optimization

Abstract

Conventional decoupling approaches usually employ first-order reliability method to deal with probabilistic constraints in a reliability-based design optimization problem. In first-order reliability method, constraint functions are transformed into a standard normal space. Extra non-linearity introduced by the non-normal-to-normal transformation may increase the error in reliability analysis and then result in the reliability-based design optimization analysis with insufficient accuracy. In this article, a decoupling approach is proposed to provide an alternative tool for the reliability-based design optimization problems. To improve accuracy, the reliability analysis is performed by first-order asymptotic integration method without any extra non-linearity transformation. To achieve high efficiency, an approximate technique of reliability analysis is given to avoid calculating time-consuming performance function. Two numerical examples and an application of practical laptop structural design are presented to validate the effectiveness of the proposed approach.

Introduction

Most engineering structures inevitably involve many uncertainties, such as geometric sizes, material parameters, loads, and work conditions.^1,2 Although these uncertainties are usually small, the combined effects can lead to large fluctuations and even failure in the structural performance.^1,3 Reliability-based design optimization (RBDO)^4,5 can assess the reliability of design option under uncertain factors in the optimization process; therefore, it can achieve a reliable optimal result. A two-layer nested optimization is involved in RBDO, which contains the design optimization in the outer layer and the reliability analysis in the inner layer. For most practical problems, the nested optimization is computationally expensive. Two aspects to improve accuracy and efficiency for solving RBDO are: (1) to select an appropriate reliability analysis method and (2) to combine it effectively with the design optimization in the outer layer. ⁶ At present, the first-order reliability method (FORM),^7–11 such as the reliability index approach (RIA)^7,8 or the performance measurement approach (PMA) ⁹ has been adopted by many important RBDO approaches.^12–14 The FORM can perform reliability analysis efficiently by establishing a linear approximation at the most probable point (MPP) for the constraint function transformed into the standard normal space. Yu et al. ¹⁵ introduced the RIA into the overall framework which contained a set of interactive design steps with trade-off analysis and what-if study. Cheng et al. ¹⁶ performed RBDO by solving a sequence of sub-programming problems that consist of an approximate objective function subjected to a set of approximate constraint functions. In the approach, the reliability index and its sensitivity can be obtained approximately. In the sequential optimization and reliability assessment (SORA) proposed by Du et al., ¹⁷ the PMA was used to conduct reliability analysis and formulate the approximate equivalent feasible domain of the probabilistic constraints, and whereby the reliability analysis is decoupled from the design optimization. Due to the high efficiency of SORA framework, a series of important variants are further developed to enhance the performance in terms of efficiency,^18–22 accuracy,^23,24 or convergence.^25–27 In the above RBDO methods, the efficiency is improved by employing FROM, but the accuracy may be not enough because of the limitation of FORM.^28,29 In the FORM, each non-normal variable is required to be transformed into the standard normal variable, and the performance function is linearized at MPP in the standard normal space. The transformation is non-linear and may increase the non-linearity of a performance function. As a result, the error of FORM will also increase, and the overall accuracy of RBDO will be reduced. Therefore, it is significant for enhancing the accuracy of RBDO to conduct the constraint reliability analyses in the original space.

At present, there are few researches^30,31 on such RBDO algorithms. In the work by Chan et al., ³² each probability constraint is transformed into a deterministic constraint in the original space, and then an adaptive sequential linear programming algorithm is proposed for the solution of the updated optimization problem. This algorithm exhibits high efficiency, but it is not suitable for problems with non-normal random variables. The SORA-SPA algorithm ³⁰ adopted the existing decoupling framework of SORA; ¹⁷ the difference is that the constraint reliability analyses are carried out by first-order saddle point approximation (FOSPA) ³¹ in the original space to improve the accuracy of RBDO. Therefore, it is named as SORA-SPA. Because the original FOSPA method is only for reliability calculation, but performing SORA requires the percentile performance of each constraint. For calculating the percentile performance, the original FOSPA process is modified into a sequence steps with solving three optimization problems. Objectively speaking, when the non-normal-to-normal transformation increases the non-linearity of constraint functions, SORA-SPA can exhibits a higher efficiency than the FORM-based RBDO method. However, there are also some disadvantages. First, the flowchart of modified FOSPA is complex; comparing with the FORM, the FOSPA involves several optimization processes and is more prone to convergence problems. Second, the entire iterative process of SORA-SPA still needs to conduct reliability analysis many times, and each reliability analysis requires repeatedly calculating the constraint functions. The constraint functional evaluation is generally very time-consuming for practical engineering problems. Therefore, there is still a gap between the efficiency of SORA-SPA and actual demand in many cases.

To improve accuracy and maintain efficiency, the proposed approach in this article performs the constraint reliability analysis in the original space, and develops an efficient decoupling framework for RBDO problems. The rest of this article is organized as follows. The general concepts and decoupling strategies of RBDO are reviewed in section “RBDO problem and decoupling strategy.” The formulation of proposed approach is given in section “Formulation of the proposed approac.” The effectiveness of the proposed approach is demonstrated by two numerical examples and an application of practical laptop structural design. A conclusion is made in section “Conclusion.”

RBDO problem and decoupling strategy

Statement of RBDO problem

A typical RBDO problem can be expressed as follows ¹⁷

$$\begin{array}{c}\underset{\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}}}{\min }f(\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}},{\mathbf{\mu }}_{\mathbf{P}})\hfill \\ \mathrm{s}.\mathrm{t}.R_j=\mathrm{Pr}(g_j(\mathbf{d},\mathbf{X},\mathbf{P})\ge 0)\ge R_j^t,j=1,2,\dots ,n_g\hfill \\ {\mathbf{d}}^l\le \mathbf{d}\le {\mathbf{d}}^u,{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{l}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{u}}\hfill \end{array}$$

(1)

where d represents the n_d-dimensional deterministic design vector, X represents the n_X-dimensional random design vector, P represents the n_P-dimensional random parameter vector; μ denotes the mean vector of a random vector. d and μ_X can be determined by the designer; superscript l and u indicate the upper and lower bounds of the value rang, respectively. f is the objective function; for the design point $(\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}})$ , the value of f can be calculated by substituting μ_X and μ_P. There are n_g numbers of probabilistic constraints; g_j denotes the jth constraint performance function; $g_j\ge 0$ means the structure can meet performance requirements. Pr denotes probability calculation; R_j is the probability of $g_j\ge 0$ , namely the reliability of the jth constraint. R_j should be greater than or equal to the target reliability $R_j^t$ . For a certain constraint, its reliability requirements can also be given by the target reliability index ${\beta }^t$ , the relationship between the them is as follows

$$R^t=\Phi \left({\beta }^t\right),{\beta }^t={\Phi }^{-1}\left(R^t\right)$$

(2)

where $\Phi$ and ${\Phi }^{-1}$ denote standard normal cumulative distribution function (CDF) and its inverse function, respectively.

In the optimization process, reliability analysis needs to be conducted for any design point to determine whether it satisfies a probabilistic constraint. Currently, FORM is the most commonly used method of handling constraint reliability analysis. Among them, PMA ¹⁰ usually exhibits a higher efficiency in RBDO analysis because of its fixed search scope. ³³

In PMA, the constraint function is first mapped from the original random space into the normal standard space. Second, an optimization problem is formulated to search the MPP. If the functional value at MPP is greater than or equal to zero, the probabilistic constraint is satisfied.

Since there is no need to distinguish the random design variables from random parameters, the n_Z-dimensional random vector Z is used to represent all the random variables, namely $\mathbf{Z}=[\mathbf{X},\mathbf{P}]$ . The mapping relationship between the original random space (Z space) and the normal standard space (U space) can be written as

$$U_i={\Phi }^{-1}\left(F_Z\left(Z_i\right)\right),Z_i=F_Z^{-1}\left(\Phi \left(U_i\right)\right),i=1,2,\dots ,n_Z$$

(3)

where $Z_i$ denotes the ith component of $\mathbf{Z}$ . $F_Z$ and $F_Z^{-1}$ represent the CDF and inverse function of the random variable, respectively. For a constraint, the functional value at the MPP is the minimum, and the distance from the MPP to the original point is equal to ${\beta }^t$ . Thus, the optimization problem for searching MPP can be written as

$$\begin{array}{c}\underset{\mathbf{U}}{\min }G\left(\mathbf{U}\right)\hfill \\ \mathrm{s}.\mathrm{t}.\Vert \mathbf{U}\Vert ={\beta }^t\hfill \end{array}$$

(4)

where G is the constraint function g in the U space. The advanced mean value method ¹⁰ can be used to solve the above equation to obtain the MPP ${\mathbf{U}}^{*}$ and the percentile performance $G^R$ . $G^R\ge 0$ is equivalent to $R\ge R^t$ , which means the design point satisfies the probabilistic constraint.

Through the above analysis, performing RBDO involves a two-level nested optimization with the outer level being the optimization for design variables and the inner level being reliability analysis. In practical engineering problems, performance functions are usually based on time-consuming numerical simulation technology, such as finite element models (FEMs) ³⁴ and multi-body dynamics models. ³⁵ Nested optimization with simulation models usually result in extremely low computational efficiency. ³⁶

Decoupling approach

Decoupling approach has now become the most effective method to solve the nested optimization in the RBDO. And the SORA ¹⁷ is a widely used decoupling strategy. Its basic idea is to convert the nested optimization into a sequence cycles of reliability analysis and deterministic optimization. In each iteration step, the reliability analysis by using the PMA is first conducted for the precious design point to search the MPP at each probabilistic constraint. Second, according to the MPP, each probabilistic constraint is equivalent to the conventional deterministic constraint, and whereby the updated deterministic optimization is formulated. After solving it, the current solution is obtained.

The detailed steps of SORA are described below. In the first iteration step $(k=1)$ , the deterministic optimization problem is formulated as follow

$$\begin{array}{c}\underset{\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}}}{\min }f(\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}},{\mathbf{\mu }}_{\mathbf{P}})\hfill \\ \mathrm{s}.\mathrm{t}.g_j(\mathbf{d},{\mathbf{\mu }}_{\mathbf{Z}})\ge 0,j=1,2,\dots ,n_g\hfill \\ {\mathbf{d}}^l\le \mathbf{d}\le {\mathbf{d}}^u,{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{l}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{u}}\hfill \end{array}$$

(5)

Here, the random variables in equation (1) are replaced by the corresponding mean values. And the solution $({\mathbf{d}}^{\left(1\right)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{\left(}1\mathbf{\right)}}\mathbf{)}$ of the first iteration step is obtained by solving it. In the second iteration step, the reliability analysis of PMA is first carried out for each probabilistic constraint at the precious design $({\mathbf{d}}^{(k-1)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}\mathbf{k}\mathbf{-}1\mathbf{)}}\mathbf{)}$ . For a certain constraint, the MPP ${\mathbf{U}}^{*\left(k\right)}$ can be searched by solving equation (4), and it can be transformed to the MPP in the original space, namely ${\mathbf{Z}}^{*\left(k\right)}$ . Second, the shifting vector is calculated

$${\mathbf{S}}^{\left(k\right)}={\mathbf{\mu }}_{\mathbf{Z}}^{\mathbf{(}\mathbf{k}\mathbf{-}1\mathbf{)}}\mathbf{-}{\mathbf{Z}}^{\mathbf{*}\mathbf{\left(}\mathbf{k}\mathbf{\right)}}$$

(6)

By moving the boundary $g(\mathbf{d},{\mathbf{\mu }}_{\mathbf{Z}})=0$ toward the feasible domain by ${\mathbf{S}}^{\left(k\right)}$ , the approximate boundary of the actual probabilistic constraint boundary is obtained, which is written as $g_j(\mathbf{d},{\mathbf{\mu }}_{\mathbf{Z}}-{\mathbf{S}}_j^{\left(k\right)})=0$ . Correspondingly, the updated deterministic optimization problem can be formulated

$$\begin{array}{c}\underset{\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}}}{\min }f(\mathbf{d},{\mathbf{\mu }}_{\mathbf{X}},{\mathbf{\mu }}_{\mathbf{P}})\hfill \\ \mathrm{s}.\mathrm{t}.g_j(\mathbf{d},{\mathbf{\mu }}_{\mathbf{Z}}-{\mathbf{S}}_j^{\left(k\right)})\ge 0,j=1,2,\dots ,n_g\hfill \\ {\mathbf{d}}^l\le \mathbf{d}\le {\mathbf{d}}^u,{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{l}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}\mathbf{\le }{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{u}}\hfill \end{array}$$

(7)

The above equation is solved to get the solution $({\mathbf{d}}^{\left(k\right)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{\left(}\mathbf{k}\mathbf{\right)}}\mathbf{)}$ in the current iteration step. The reliability analysis and the deterministic optimization are performed alternatively until convergence.

In SORA, the non-linear transformation for the constraint function in PMA introduces the extra error into the reliability analysis. To enhance accuracy, the SORA-SPA ³⁰ algorithm was developed as a variant of SORA. It employs FOSPA ³¹ to perform the constraint reliability analyses in the original space without converting the constraint performance functions to the standard normal space. In the kth iteration step, the reliability analysis for a constraint at the current design point $(\mathbf{d}{\mathbf{d}}^{(k-1)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{\left(}\mathbf{k}\mathbf{\right)}}\mathbf{,}{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}\mathbf{k}\mathbf{-}1\mathbf{)}}\mathbf{)}$ mainly includes the three steps as follows:

1. To establish linear approximation $L^{(k)}(\mathbf{Z})$ for the performance function of $g(\mathbf{Z})$ at the most possible point (MLP), which is the point with the maximum joint probability density at the boundary of $g(\mathbf{Z})=0$ ;

2. To calculate the performance measure value for the target reliability index by solving the nested optimization problem $\mathrm{Pr}(L^{(k)}(\mathbf{Z})\ge g^{*(k)})=R^t$ , and the detailed solution process refers to the literature; ³¹

3. To search for another MLP (namely ${\mathbf{Z}}^{*(k)}$ ) along the boundary $L^{(k)}(\mathbf{Z})=g^{R(k)}$ for the formulation of the shifting vector ${\mathbf{S}}^{(k)}$ as equation (6).

Each of these steps needs to solve the optimization problem. In the first step, the performance functional evaluation is involved; and the second step contains the nested optimization with the outer layer being the root-finding problem and the inner layer being the reliability analysis.

Formulation of the proposed approach

The SORA-SPA is more accurate than the FORM-based RBDO methods when the non-normal-to-normal transformation increases the non-linearity of probabilistic constraints. However, it can be observed that the constraint reliability analysis of SORA-SPA is complicated, and there is no guarantee of convergence due to the three optimization problems involved in it. Furthermore, each constraint reliability analysis contains the performance functional evaluation, which is time-consuming for most engineering problems. The whole optimization process needs a number of reliability analyses, especially for problems with several constraints. Therefore, the issues of efficiency seem still challenging for SORA-SPA when dealing with more complex engineering problems.

In this section, a decoupling approach is formulated to provide an alternative tool to enhance accuracy and efficiency for performing RBDO. Similar to the SORA, the framework of proposed approach converts the nested optimization into the sequential reliability analysis and design optimization steps. The calculation of the shifting vector is the key point for the proposed approach and also the fundamental difference from the SORA ¹⁷ and its varieties.^18–26,32 On one hand, the shifting vector determines the deviation between the actual boundary of the probabilistic constraint and the equivalent boundary of the deterministic constraint in the iterative process. Due to the equivalent boundaries driving the design point to approach the optimal solution, the small deviation between the two boundaries can speed up the convergence. On the other hand, the formulation of the shifting vector also directly affects the overall efficiency of RBDO. The proposed approach calculates the shifting vector based on the result of reliability analysis. The features of the process can be summarized as: (1) the first-order asymptotic integration method is employed for the reliability analysis in the original space to enhance the accuracy; (2) an approximate technology is given for avoiding calculate the time-consuming constraint function to improve efficiency; and (3) an increment is used for fine-tuning the value of the shifting vector between successive iterations to ensure the convergence. The procedure of proposed approach is detailed below.

Constraint reliability analysis

In the first iteration, the deterministic optimization as equation (5) is solved to obtain the design point $({\mathbf{d}}^{(1)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}1\mathbf{)}}\mathbf{)}$ . In each subsequent iteration step, the reliability analysis and deterministic optimization are sequentially performed. The constraint reliability analysis is conducted by the first-order asymptotic integration method,^37,38 which concludes the three steps as follows. First, the most likelihood point (MLP) is searched along the boundary of constraint function in the original space. Second, a linear approximation for constraint function is established at the MLP. Finally, the constraint reliability is calculated by the asymptotic integration algorithm.

The MLP is the point with the maximum joint probability density at the limit-state boundary of $g(\mathbf{Z})=0$ . For a certain constraint, the searching process can be formulated as the optimization problem as follow

$$\begin{array}{c}\underset{\mathbf{Z}}{\max }{f}_{\mathbf{Z}}\left(\mathbf{Z}\right)\hfill \\ \mathrm{s}.\mathrm{t}.g\left(\mathbf{Z}\right)=0\hfill \end{array}$$

(8)

where $f_{\mathbf{Z}}(\mathbf{Z})$ represents the joint probability density function (PDF) of Z, and can be written as $f_{\mathbf{Z}}(\mathbf{Z})=\underset{1}{\overset{n_Z}{\Pi }}{f}_{Z_i}(Z_i)$ , $Z_i$ is the ith component of the vector Z, $f_{Z_i}$ denotes the PDF of $Z_i$ .

The linear approximation of the constraint function is established at ${\mathbf{Z}}^{*(k)}$

$$L^{(k)}\left(\mathbf{Z}\right)=g\left({\mathbf{Z}}^{*(k)}\right)+{(\mathbf{Z}-{\mathbf{Z}}^{*(k)})}^{\mathrm{T}}·\nabla g\left({\mathbf{Z}}^{*(k)}\right)$$

(9)

Here, $\nabla g({\mathbf{Z}}^{*(k)})$ represents the gradient vector of $g(\mathbf{Z})$ at ${\mathbf{Z}}^{*(k)}$ .

For $L^{\left(k\right)}\left(\mathbf{Z}\right)$ , the constraint reliability can be calculated as follow

$$\begin{array}{c}\multicolumn{1}{c}{R^{\left(k\right)}=\underset{g\left(\mathbf{Z}\right)\ge 0}{\int }{f}_{\mathbf{Z}}\left(\mathbf{Z}\right)d\mathbf{Z}}\\ \approx \underset{L^{\left(k\right)}\left(\mathbf{Z}\right)\ge 0}{\int }{f}_{\mathbf{Z}}\left(\mathbf{Z}\right)d\mathbf{Z}\hfill \\ =\underset{L^{\left(k\right)}\left(\mathbf{Z}\right)\ge 0}{\int }\exp \left(\ln f_{\mathbf{Z}}\left(\mathbf{Z}\right)\right)d\mathbf{Z}\hfill \end{array}$$

(10)

Let $h(\mathbf{Z})=\ln f_{\mathbf{Z}}(\mathbf{Z})$ , and the second-order Taylor expansion of $h(\mathbf{Z})$ at ${\mathbf{Z}}^{*(k)}$ is established

$$\begin{array}{c}h\left(\mathbf{Z}\right)\approx h\left({\mathbf{Z}}^{*\left(k\right)}\right)+{(\mathbf{Z}-{\mathbf{Z}}^{*\left(k\right)})}^{\mathrm{T}}·\nabla h\left({\mathbf{Z}}^{*\left(k\right)}\right)\hfill \\ +\frac{1}{2}{(\mathbf{Z}-{\mathbf{Z}}^{*\left(k\right)})}^{\mathrm{T}}·{\nabla }^{2}h\left({\mathbf{Z}}^{*\left(k\right)}\right)·(\mathbf{Z}-{\mathbf{Z}}^{*\left(k\right)})\hfill \\ =h\left({\mathbf{Z}}^{*\left(k\right)}\right)+\frac{1}{2}{\left({\mathbf{v}}^{\left(k\right)}\right)}^{\mathrm{T}}·{\mathbf{B}}^{\left(k\right)}·{\mathbf{v}}^{\left(k\right)}\hfill \\ -\frac{1}{2}{(\mathbf{Z}-{\mathbf{Z}}^{*\left(k\right)}-{\mathbf{B}}^{\left(k\right)}·{\mathbf{v}}^{\left(k\right)})}^{\mathrm{T}}·{\left({\mathbf{B}}^{\left(k\right)}\right)}^{-1}·\hfill \\ (\mathbf{Z}-{\mathbf{Z}}^{*\left(k\right)}-{\mathbf{B}}^{\left(k\right)}·{\mathbf{v}}^{\left(k\right)})\hfill \end{array}$$

(11)

where

$${\mathbf{v}}^{\left(k\right)}=\nabla h\left({\mathbf{Z}}^{*\left(k\right)}\right),{\mathbf{B}}^{\left(k\right)}=-{\left({\nabla }^{2}h\left({\mathbf{Z}}^{*\left(k\right)}\right)\right)}^{-1}$$

(12)

The mean and variance of the constraint functional values are calculated by

$$\begin{array}{c}{\mu }_g^{\left(k\right)}={(\nabla g\left({\mathbf{Z}}^{*\left(k\right)}\right))}^{\mathrm{T}}·{\mathbf{B}}^{\left(k\right)}·{\mathbf{v}}^{\left(k\right)}\hfill \\ {\left({\sigma }_g^{\left(k\right)}\right)}^2={(\nabla g\left({\mathbf{Z}}^{*\left(k\right)}\right))}^{\mathrm{T}}·{\mathbf{B}}^{\left(k\right)}·\nabla g\left({\mathbf{Z}}^{*\left(k\right)}\right)\hfill \end{array}$$

(13)

The constraint reliability index obtained by the first-order second-moment method ⁹ is

$${\beta }^{\left(k\right)}=\frac{{\mu }_g^{\left(k\right)}}{{\sigma }_g^{\left(k\right)}}$$

(14)

So far, the MLP ${\mathbf{Z}}^{*(k)}$ and the constraint reliability index ${\beta }_j^{(k)}$ have been obtained; theoretically, the constraint reliability analysis in the kth iteration has been completed. As mentioned above, however, the first step of the asymptotic integration method is to solve the optimization problem of equation (8). The solving process requires repeated calculation of the constraint function. In many cases, the constraint functions of engineering problem are based on the time-consuming simulation models. Repeatedly calling simulation models will likely lead to an expensive computational costs. ³⁹ Therefore, an approximate algorithm is given to avoid calculating the constraint function when searching the MLP.

Compared with performing the reliability analysis by equation (10) based on the linear approximation at other points, the error of that at the MLP is smallest. ³⁸ Moreover, the solutions’ difference between the two adjacent iteration steps is generally small, and thus that of corresponding MLPs of a constraint is small. With the iteration conducting, the difference becomes smaller. As shown in Figure 1, in the proposed algorithm, $g(\mathbf{Z})$ in equation (8) is replaced with $L^{(k-1)}(\mathbf{Z})$ which is the linear approximation at the precious MLP of ${\mathbf{Z}}^{*(k-1)}$

$$\begin{array}{c}\underset{\mathbf{Z}}{\min }{f}_{\mathbf{Z}}\left(\mathbf{Z}\right)\hfill \\ \mathrm{s}.\mathrm{t}.L^{(k-1)}\left(\mathbf{Z}\right)=g\left({\mathbf{Z}}^{*(k-1)}\right)\hfill \\ +{(\mathbf{Z}-{\mathbf{Z}}^{*(k-1)})}^T·\nabla g\left({\mathbf{Z}}^{*(k-1)}\right)=0\hfill \end{array}$$

(15)

Figure 1.

The approximate algorithm for searching the MLP.

The calculation of the constraint function $g(\mathbf{Z})$ is no longer necessary when searching the MLP in the current iteration step by equation (15), and thus the calculation cost is almost negligible. Noted that, if and only if $k=2$ , $g(\mathbf{Z})$ cannot be replaced due to the lack of $L^{(1)}(\mathbf{Z})$ .

The calculation of the shifting vector

The SORA and its variants determine the equivalent constraint boundaries by recalculating the shifting vectors in each iteration step. For some complex cases, such as the constraint with multiple MLPs, the differences between the shifting vectors of adjacent iteration steps can be large, which may lead to numerical oscillations and whereby impact the convergence of the iteration process. ³⁰ To improve convergence, the concept of increment ⁴⁰ recently proposed by the authors is introduced, which is calculated in the kth iteration step for fine-tuning the shifting vector in the (k-1)th iteration step. This fine-tuning process can be expressed as ${\mathbf{S}}^{(k)}={\mathbf{S}}^{(k-1)}+\Delta {\mathbf{S}}^{(k)}$ . The updated shifting vector drives the equivalent constraint boundary of $g({\mathbf{\mu }}_{\mathbf{Z}}-({\mathbf{S}}^{(k-1)}+\Delta {\mathbf{S}}^{(k)}))=0$ to gradually approach the actual probabilistic constraint boundary of $\mathrm{Pr}(g(\mathbf{Z})\ge 0)\ge R_j^t$ .

Figure 2 is used to illustrate the calculation process of the increment. For ease of description, the constraint function contains only two random variables. $D=\{g:g({\mathbf{\mu }}_{\mathbf{Z}})\ge 0\}$ represents the feasible region obtained by the original boundary of $g({\mathbf{\mu }}_{\mathbf{Z}})=0$ , and $S=\{R:\mathrm{Pr}(g(\mathbf{Z})\ge 0)\ge R^t\}$ represents that of the probabilistic constraint boundary. Because of the uncertainty of the random vector Z, $S\subset D$ . The deterministic boundary of $g({\mathbf{\mu }}_{\mathbf{Z}}-{\mathbf{S}}^{(k-1)})=0$ is established to be equal to the probabilistic constraint boundary in the $(k-1)\mathrm{th}$ iteration step. Under driven by the equivalent boundary, the previous design point ${\mathbf{\mu }}_{\mathbf{Z}}^{\mathbf{(}\mathbf{k}\mathbf{-}1\mathbf{)}}$ is located in the region of D. As described in section “Constraint reliability analysis,” the constraint reliability for ${\mathbf{\mu }}_{\mathbf{Z}}^{\mathbf{(}\mathbf{k}\mathbf{-}1\mathbf{)}}$ is performed to get the MLP ${\mathbf{Z}}^{*(k)}$ and the reliability index ${\beta }^{(k)}$ . Assuming that ${\beta }^{(k)}$ is less than the target reliability index ${\beta }^t$ , it means that the probabilistic constraint is still unsatisfied; the gap between them is $\Delta {\beta }^{(k)}={\beta }^t-{\beta }^{(k)}$ . To enhance the reliability of the design point, the equivalent boundary needs to be further shifted toward the feasible region in the kth iteration step. Specifically, the boundary $g({\mathbf{\mu }}_{\mathbf{Z}}-{\mathbf{S}}^{(k-1)})=0$ is shifted by $\Delta {\beta }^{(k)}$ along the MLP gradient direction. In addition, the standard deviations of components in the random vector Z are different. The larger is the standard deviation of the component $Z_i$ , the greater is the increment in this dimension. To sum up, the components of $\Delta {\mathbf{S}}^{(k)}$ can be expressed as

$$\Delta S_i^{\left(k\right)}={\sigma }_{\mathrm{Zi}}·({\beta }^t-{\beta }^{(k)})·{\alpha }_i^{\left(k\right)},i=1,2,\dots ,n_Z$$

(16)

Figure 2.

The calculation for the increment of shifting vector.

where the subscript i denotes the ith component of the vector; ${\mathbf{\alpha }}^{\mathbf{(}\mathbf{k}\mathbf{)}}$ is the normalized gradient vector which can be written as ${\mathbf{\alpha }}^{\mathbf{(}\mathbf{k}\mathbf{)}}=\nabla g({\mathbf{Z}}^{*(k)})/||\nabla g({\mathbf{Z}}^{*(k)})||$ ; ${\mathbf{\sigma }}_{\mathbf{Z}}$ is the standard deviation vector of Z. The increments $\Delta {\mathbf{S}}_j^{(k)},j=1,2,\dots ,n_g$ are calculated for probabilistic constraints, and then the corresponding equivalent constraints can be obtained.

After the above equivalent process, the updated deterministic optimization problem as equation (7) is also formulated. The current solution $({\mathbf{d}}^{(k)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}\mathbf{k}\mathbf{)}}\mathbf{)}$ is obtained by solving it. The reliability analysis and the deterministic optimization are carried out alternately until the convergence is achieved.

The computational procedure

As shown in Figure 3, the proposed approach is summarized as the following steps:

• Step1: After analyzing the practical engineering problem, the RBDO model is formulated as equation (1).

• Step2: Set the iteration step $k=1$ and solve equation (5) to obtained the solution $({\mathbf{d}}^{(1)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}1\mathbf{)}}\mathbf{)}$ .

• Step3: Set $k:=k+1$ and search the MLP ${\mathbf{Z}}_j^{*(k)}$ for each probabilistic constraint; if $k=1$ , equation (8) is used, else equation (15) is used.

• Step4: Calculate the reliability index ${\beta }_j^{(k)}$ by equations (10)–(14).

• Step5: Calculate the increment $\Delta {\mathbf{S}}_j^{(k)}$ by equation (16) for probabilistic constraints and convert them into the equivalent deterministic constraints.

• Step6: Formulate the optimization problem as equation (7) and solve it to get updated solution $({\mathbf{d}}^{(k)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}\mathbf{k}\mathbf{)}}\mathbf{)}$ .

• Step7: Repeat Step 3 to Step 6 until the following conditions are achieved

$$\{\begin{array}{c}\hfill R_j^{\left(k\right)}\ge R_j^t,j=1,2,\dots ,n_g\hfill \\ \multicolumn{1}{c}{\left|\frac{f^{\left(k\right)}-f^{(k-1)}}{f^{\left(k\right)}}\right|\le {\epsilon }_r}\end{array}$$

(17)

Figure 3.

The flowchart of the proposed approach.

where ${\epsilon }_r$ is the minimal error limit; $f^{(k)}=f({\mathbf{d}}^{(k)},{\mathbf{\mu }}_{\mathbf{X}}^{\mathbf{(}\mathbf{k}\mathbf{)}}\mathbf{,}{\mathbf{\mu }}_{\mathbf{P}}\mathbf{)}$ denotes the objective functional value in the $k\mathrm{th}$ iteration step.

It can be found that all of the above steps are carried out in the original space. Because there is no non-linear transformation for constraint functions, extra error will not be introduced into the reliability analysis compared to other FORM-based approach. It can be also found that the calculation of constraint functions is not necessary for the reliability analysis in all iterations except the second iteration. This can further improve efficiency. These characteristics can be exhibited by example analysis in the nest section.

Numerical examples and discussions

Three examples with different complexity are used to test the proposed approach, which include a standard test case, a structural design problem, and an application of laptop structural design. By comparing the results of SORA, the performance of proposed approach is investigated in terms of accuracy, efficiency, and convergence. To verify the accuracy of all results, the double-loop strategy is used to get reference solutions, in which the sequential quadratic programming ⁴¹ is employed to solve the design optimization in outer level, and the Monte-Carlo method ⁴² is adopted to conduct the reliability analysis in inner level. The number of constraint functional computation is counted to measure the efficiency of the methods. In addition, the same solution parameters and convergence tolerances are set for all methods to keep objective comparisons as possible.

A standard test case

The example is a widely used standard test case for RBDO methods^18–25

$$\begin{array}{c}\underset{{\mathbf{\mu }}_{\mathbf{X}}}{\min }f\left({\mathbf{\mu }}_{\mathbf{X}}\right)={\mu }_{X_1}+{\mu }_{X_2}\hfill \\ \mathrm{s}.\mathrm{t}.\mathrm{Pr}(g_j\left(\mathbf{X}\right)\ge 0)\ge \Phi \left({\beta }^t\right),j=1,2,3\hfill \\ g_1\left(\mathbf{X}\right)=1-\frac{X_1^2{X}_2}{20}\hfill \\ g_2\left(\mathbf{X}\right)=1-\frac{{(X_1+X_2-5)}^2}{30}-\frac{{(X_1-X_2-12)}^2}{120}\hfill \\ g_3\left(\mathbf{X}\right)=1-\frac{80}{(X_1^2+8X_2+5)}\hfill \\ 0\le {\mu }_{X_i}\le 10,\frac{{\sigma }_{X_i}}{{\mu }_{X_i}}=0.05,i=1,2\hfill \end{array}$$

(18)

As mentioned above, the difference between the FORM-based methods^12,13 and the proposed approach is that the former performs the reliability analysis after mapping the constraint functions into the standard normal space and the latter does it in the original space. To show the impact of this difference on the validity of the results, a variety of situations are considered in solving this problem. First, four distribution situations are set for the two random variables. In the Situations 1–3, their distribution is the same type which respectively is normal, log-normal, ands Gumbel. In the Situation 4, $X_1$ obeys a log-normal distribution, and $X_2$ obeys a Gumbel distribution. Second, the constraint target reliability indexes are set to different values: ${\beta }_{(I)}^t=2.5$ , ${\beta }_{(\mathrm{II})}^t=3.0$ , ${\beta }_{(\mathrm{III})}^t=3.5$ . As a result, there are $4\times 3=12$ situations in total.

All calculation results are listed in Table 1. The results show that the optimal design of proposed approach is very close to the reference result and meets the constraint reliability requirements under each situation; while the optimal design of SORA violates Constraint 1 under Situations 7, 8, 10, and 11. In order to analyze the reason, Figure 4 shows the curves of Constraint 1 transformed into the standard normal space under Situations 1, 7, and 10, namely Curve_0, Curve_1, and Curve_2 respectively. Because $X_1,X_2$ are normal random variables under Situation 1, Curve_0 can represent the non-linearity of the constraint function in the original space. By comparing the three curves, it can be found that the non-linearity of Curve_1 and Curve_2 is higher than that of Curve_0. This shows that the non-linearity of the constraint function increases obviously after transformed into the standard normal space under Situations 7 and 10. In SORA, the non-linearity transformation introduces the extra error into each constraint reliability analysis and thus leads to the unreliable results of the overall RBDO analysis. This also shows to some extent that performing constraint reliability analysis in the original space is important for the validity of RBDO solution.

Table 1.

The computational results under the 12 situations.

Situation	Distribution type	TRI	Reference result	SORA			Proposed approach
			Design variable	Design variable	Reliability index	Times of constraint calculation	Design variable	Reliability index	Times of constraint calculation
		${\beta }^t$	${\mu }_{X_1}^{*},{\mu }_{X_2}^{*}$	${\mu }_{X_1}^{*},{\mu }_{X_2}^{*}$	${\beta }_1,{\beta }_2,{\beta }_3$	nf	${\mu }_{X_1}^{*},{\mu }_{X_2}^{*}$	${\beta }_1,{\beta }_2,{\beta }_3$	nf
1	$X_1,X_2:$ Normal	2.5	3.234,2.569	3.232,2.572	2.5, 2.6,5.0	192	3.230,2.572	2.5, 2.5,5.0	94
2		3.0	3.270,2.676	3.265,2.683	3.0, 3.1,5.0	211	3.264,2.682	3.0, 3.0,5.0	98
3		3.5	3.302,2.790	3.302,2.798	3.5, 3.5,5.0	218	3.300,2.800	3.6, 3.5,5.0	92
4	$X_1,X_2:$ Log-normal	2.5	3.218,2.554	3.222,2.558	2.5, 2.5,5.0	192	3.246,2.651	2.5, 2.5,5.0	101
5		3.0	3.255,2.656	3.250,2.658	3.0, 3.0,5.0	211	3.243,2.671	3.0, 3.0,5.0	115
6		3.5	3.286,2.756	3.280,2.759	3.6, 3.6,5.0	214	3.265,2.781	3.5, 3.7,5.0	117
7	$X_1,X_2:$ Gumbel	2.5	3.307,2.726	3.299,2.736	2.4, 2.5,5.0	210	3.308,2.724	2.5, 2.5,5.0	108
8		3.0	3.410,3.006	3.387,3.016	2.9, 3.0,5.0	228	3.401,3.006	3.0, 3.0,5.0	118
9		3.5	3.542,3.372	3.503,3.400	3.5, 3.6,5.0	255	3.504,3.396	3.5, 3.5,5.0	120
10	X1: Log-normal; X2: Gumbel	2.5	3.189,2.655	3.181,2.646	2.4, 2.5,5.0	202	3.190,2.656	2.5, 2.5,5.0	104
11		3.0	3.030,3.137	3.180,2.850	2.9, 3.0,5.0	207	3.030,3.139	3.0, 3.9,5.0	105
12		3.5	3.191,3.119	3.176,3.123	3.5, 3.5,5.0	216	3.111,3.273	3.5, 3.7,5.0	113

Figure 4.

Constraint function curve in the standard normal space.

A structural design problem for a roof truss

Figure 5 shows a structural design problem for a roof truss. ⁴³ The top boom and compression members are concrete, and the bottom boom is steel. Their cross-sectional areas are $A_1,A_2$ , and the elastic modules are $E_1,E_2$ , respectively. The roof is subjected to the vertical load of $q$ . The design variable are $A_1,A_2$ and the objective is to minimize the material cost. Only one constraint is involved that the deformation of the roof in the vertical direction should not exceed 0.03 m. The RBDO problem is expressed as follows

$$\begin{array}{c}\underset{{\mu }_{A_1},{\mu }_{A_2}}{\min }f=20224{\mu }_{A_1}+364{\mu }_{A_2}\hfill \\ \mathrm{s}.\mathrm{t}.\mathrm{Pr}(g=\frac{{\mathrm{qL}}_0^2}{2}(\frac{3.81}{A_1{E}_1}+\frac{1.13}{A_2{E}_2})-0.03\ge 0)\ge \Phi \left({\beta }^t\right)\hfill \\ {\beta }^t=2,L_0=12\mathrm{m}\hfill \\ 0.02{\mathrm{m}}^2\le {\mu }_{A_1}\le 0.10{\mathrm{m}}^2,0.0006{\mathrm{m}}^2\le {\mu }_{A_2}\le 0.002{\mathrm{m}}^2\hfill \end{array}$$

(19)

As listed in Table 2, two situations are considered for the distribution type of random variables. The SORA-SPA method and proposed approach are used to solve equation (19) and the same initial design $(0.04{\mathrm{m}}^2,0.001{\mathrm{m}}^2)$ is selected. Under both situations, the two methods converge to the stable solution with satisfying the probabilistic constraint, and all computational results are listed in Table 3. First, the results of the two methods are very close to the reference result of the double-loop under each situation. The maximum difference in the objective functional value between proposed approach with double-loop is only 0.3% and that between the SORA-SPA with double-loop is 0.15%. Second, the proposed approach exhibits obvious efficiency advantage than the SORA-SPA. In both cases, the number of function calls of this method is about 23% and 57% of that of SORA-SPA, respectively. Under the two situations, the constraint calculation times of proposed approach are 23% and 57% of the SORA-SPA, respectively. The reason is that the proposed approach only needs to solve one optimization problem with constraint calculations, while the SORA-SPA needs to solve two.

Table 2.

The random variables in the design problem of the roof truss.

Variable	Symbol	Mean	Standard deviation	Distribution type
				Situation 1	Situation 2
Cross-sectional area of concrete beam (m2)	$A_1$	${\mu }_{A_1}$	0.06	Normal	Normal
Cross-sectional area of steel beam (m2)	$A_2$	${\mu }_{A_2}$	0.06	Normal	Normal
Elastic modulus of concrete (MPa)	$E_1$	$1\times {10}^{11}$	0.06	Normal	Log-normal
Elastic modulus of steel (MPa)	$E_2$	$2\times {10}^{10}$	0.06	Normal	Log-normal
Load (MPa)	$q$	$2\times {10}^4$	0.06	Normal	Gumbel

Table 3.

The computational results in the design problem of the roof truss.

Results	Symbol	Situation 1			Situation 2
		Double loop	SORA-SPA	Proposed method	Double loop	SORA-SPA	Proposed method
Design variables (m)	${\mu }_{A_1}^{*}$	3.06E–1	3.05E–1	3.06E–1	3.02E–1	3.00E–1	3.01E–1
	${\mu }_{A_2}^{*}$	1.02E–3	1.03E–3	1.02E–3	1.01E–3	1.01E–3	1.01E–3
Objective functional value $(\$)$	$f^{*}$	31.73	31.86	31.84	31.35	31.40	31.34
Reliability index	$\beta$	2.00	2.00	2.00	2.00	2.02	2.00
Iteration steps	$N_I$	–	3	3	–	3	3
Times of the constraint calculation	$N_F$	–	356	74	–	401	89

SORA-SPA: sequential optimization and reliability assessment—saddle point approximation.

Figure 5.

A structural design problem for a roof truss. ³⁶

An application of laptop structural design

A laptop has become an essential tool for many people due to its portability and rich functionality. In order to enhance the market competitiveness, a laptop is designed to be thinner structurally and more powerful functionally. The issue of heat dissipation is always a challenge in the structural design of laptop. Moreover, the user’s settings directly determine the function of the laptop. Therefore, it is significant for a designer to fully consider the uncertainties involved in the working conditions of the laptop. ⁴⁴

As shown in Figure 6, the structural design problem of a 13-inch laptop is considered based on heat dissipation. The laptop adopts the forced air cooling unit which combining two centrifugal fans with a heat pipe. The heat of key components on the mainboard, such as the central processing unit (CPU) and the graphics processing unit (GPU), is conducted to the heat sink through the heat pipe. Centrifugal fans drive air to enter the laptop from the inlets on sides and pass through the heat sink to take away the heat. In addition, the heat of the mainboard can also be conducted into the ambient through the bracket and cover. In the practical problem, the thickness of the fans, bracket, and cover can be considered as the three design variables, namely $d$ , $X_1$ , and $X_2$ . The design goal is to minimize the thickness of the cooling unit, which is expressed as $\delta =d+X_1+X_2$ . The CPU temperatures under two different cases are considered as design constraints. Case 1: two centrifugal fans work properly, and the heat dissipation mainly relies on the airflow driven by the fans. Under Case 1, the typical power dissipation of CPU is $8.0\mathrm{W}$ . Case 2: both centrifugal fans have failed, and the heat dissipation is achieved by conducting through the bracket and cover. To avoid further damage to the laptop due to overheating, CPU is designed to automatically reduce the performance, and its power dissipation is controlled within $6.0\mathrm{W}$ under Case 2. The design variables $X_1,X_2$ and the heat dissipation of the components are random variables, and the detailed information of distribution is listed in Table 4. The RBDO model for the structural design of the laptop is summarized as follows

Table 4.

The random variables in the design problem of the laptop.

Random Variable	Symbol	Mean	Coefficient of variance	Distribution type
Bracket thickness (mm)	$X_1$	${\mu }_{X_1}$	0.03	Uniform
Cover thickness (mm)	$X_2$	${\mu }_{X_2}$	0.03	Uniform
CPU heat dissipation under Case 1 (W)	$P_1$	8.0	0.05	Normal
GPU heat dissipation (W)	$P_2$	2.0	0.1	Gumbel
Convertor heat dissipation (W)	$P_3$	0.5	0.1	Normal
CPU heat dissipation under Case 2 (W)	$P_4$	6.0	0.05	Normal

CPU: central processing unit.

$$\begin{array}{c}\underset{d,{\mu }_{X_1},{\mu }_{X_2}}{\min }\delta =d+{\mu }_{X_1}+{\mu }_{X_2}\hfill \\ \mathrm{s}.\mathrm{t}.\mathrm{Pr}(g_j=T_0^C-T_j^C\ge 0)\ge \Phi \left({\beta }^t\right),j=1,2\hfill \\ T_0^C={45}^{\circ }C,1.0\mathrm{mm}\le d\le 5.0\mathrm{mm}\hfill \\ 0.5\mathrm{mm}\le {\mu }_{X_i}\le 5.0\mathrm{mm},i=1,2\hfill \end{array}$$

(20)

Here, the constraint functions is based on the FEM, which consists of 35,511 eight-node hexahedral elements and 108,612 four-node tetrahedron elements. The FEM under Case 1 or Case 2 is obtained by setting corresponding boundary condition as shown in Figure 7. To achieve parameterization, the second-order polynomial response surfaces of $T_j^C$ in the constraint function are established by the sampling of corresponding FEM for 100 times. And then they can be written as

$$\begin{array}{c}\multicolumn{1}{c}{T_1^C=0.04P_1^2+2.93P_1+0.26P_2^2+1.09P_3+1.05X_1^2}\\ +0.2X_1{X}_2+0.15X_{1}d-6.5X_1\hfill \\ +0.04X_2^2+0.06X_{2}d+0.18d-1.78X_2\hfill \\ +0.07d^2-2.07d+37.82\hfill \\ \multicolumn{1}{c}{T_2^C=0.45P_2^2+0.02P_2{P}_3-0.05P_2{P}_4+1.62P_3^2}\\ +0.18P_4^2+1.34P_4\hfill \\ -1.42X_1^2+0.14X_1{X}_2-0.27X_2^2+4.16X_2+24.91\hfill \end{array}$$

(21)

To test the accuracy of the response surfaces, the comparison between the results of each response surface and corresponding FEM is carried out at six points which are randomly selected in the design space. The results are listed in Table 5. The maximal error of the response surfaces does not exceed 5%. Usually, it is acceptable for such engineering problems.

Table 5.

Accuracy test of the response surfaces in the design problem of the laptop.

Test points $(d,X_1,X_2,P_1^C,P_2^C,P^G,P^A)$	Relative error from the FEMs
	$T_1^C$ (%)	$T_2^C$ (%)
2.20, 1.78, 2.95, 8.5, 2.2, 0.7	3.6	4.5
1.54, 3.62, 1.42, 7.6, 1.6, 0.2	1.3	0.8
8.13, 3.08, 2.53, 9.5, 1.2, 0.2	2.1	3.6
7.78, 1.81, 3.12, 4.3, 24.9, 0.1	3.9	1.6
7.51, 3.13, 2.69, 10.6, 3.0, 0.3	1.5	0.1
3.49, 2.81, 0.79, 6.0, 1.3, 0.5	2.4	3.3

FEM: finite element model.

Figure 6.

The structural drawing of the laptop.

Figure 7.

The FEMs of the laptop.

To better investigate the practicability of proposed approach, the target reliability index ${\beta }^t$ in equation (20) is given as two situations: ${\beta }_{(I)}^t=2.5$ and ${\beta }_{(\mathrm{II})}^t=3.0$ . The initial design $(d^{(0)},{\mu }_{X_1}^{(0)},{\mu }_{X_2}^{(0)})=(8.0\mathrm{mm},2.0\mathrm{mm},2.0\mathrm{mm})$ is determined according to the designer’s experience, and corresponding objective functional value is ${\delta }^{(0)}=12\mathrm{mm}$ . The constraint reliability analyses are conducted at the initial design to distinguish its validity. The reliability index of two constraints are respectively ${\beta }_1^{(0)}=3.7$ and ${\beta }_2^{(0)}=1.6$ . It can be observed that ${\beta }_2^{(0)}$ is less than ${\beta }_{(I)}^t$ or ${\beta }_{(\mathrm{II})}^t$ and Constraint 2 is violated. In other words, the initial design satisfies the requirement of heat dissipation when the fans working properly, but the laptop may be damaged by overheating when the fan fails. To address the issue, the RBDO is performed by the proposed approach. The computational results are listed in Table 6, and the iterative process is shown as Figure 8. First, after five iteration steps, the stable solutions are converged which differ from the reference results by only 0.2%, 0.2% under the two situations; it shows that the proposed approach has good convergence and accuracy for this problem. Second in efficiency, the calculations of the constraints are 399 times and 402 times, respectively. Although the response surfaces of constraint functions have been formulated, this is not a necessary step of the proposed approach. If the FEMs are called directly when performing RBDO, it will be important for efficiency to reduce the times of constraint functional calculation. The time of each FEM calculation is 0.25 h when using the personal computer with the CPU of i7-4710HQ and the RAM 8G and thus calling the FEM by 402 times means the computational cost of about 100 h. It is feasible for this practical problem in terms of efficiency. Finally, it can be found that, if a higher reliability requirement is set for the cooling unit, the laptop will become thicker and exhibit relatively poor portability. For example, the thickness ${\delta }^{*}$ of the cooling unit when setting ${\beta }_{(\mathrm{II})}^t=3.0$ is 8% more than that of ${\beta }_{(\mathrm{II})}^t=2.5$ . Based on the above analysis, the proposed approach is an efficient tool to help the designer to balance the performance and reliability in the laptop design process, and finally achieve the optimal design. In addition, SORA-SPA is also used to solve equation (21). The iteration progresses of this algorithm and proposed approach are shown as Figure 8. It can be seen that, under the two conditions, proposed approach achieves convergence after less number of iteration steps. Because the increment is adopted for fine-tuning the value of the shifting vector between successive iterations, the proposed approach shows better convergence for this engineering application.

Table 6.

The computational results in the design problem of the laptop.

Result	Symbol	Case 1 $({\beta }_{(I)}^t=2.5)$			Case 2 $({\beta }_{(\mathrm{II})}^t=3.0)$
		Double loop	SORA-SPA	Proposed approach	Double loop	SORA-SPA	Proposed approach
Design variables (mm)	$d^{*}$	5.22	5.25	5.21	6.02	6.07	6.04
	${\mu }_{X_1}^{*}$	2.32	2.32	2.33	2.39	2.38	2.39
	${\mu }_{X_2}^{*}$	2.22	2.21	2.22	2.15	2.12	2.14
Objective functional value (mm)	${\delta }^{*}$	9.76	9.78	9.76	10.55	10.57	10.57
Reliability index	$\beta$	2.5, 2.5	2.5, 2.5	2.5, 2.5	3.0, 3.0	3.0, 3.0	3.0, 3.0
Iterations	$N_I$	9	6	5	11	7	5
Times of constraint calculation	$N_F$	–	619	399	–	958	402

SORA-SPA: sequential optimization and reliability assessment—saddle point approximation.

Figure 8.

The iterative process for the design problem of the laptop.

Conclusion

In this article, a decoupling approach of RBDO is developed which exhibits a good performance in terms of accuracy and efficiency. This approach is particularly applicable when the non-normal-to-normal transformation increases the non-linearity of probabilistic constraints. The good accuracy is achieved by conducting the constraint reliability analysis in the original space to avoid the error from the non-linear transformation. The high efficiency is obtained by proposing the approximate technology to reduce the calculation times of the constraint function. The convergence is ensured by introducing the concept of increment to drive the equivalent constraint boundary to gradually approximate the probabilistic constraint boundary. The validity and practicability of this approach are demonstrated by two numerical examples and the design application for a laptop. Also, it seems promising to extend the approach to deal with some other important problems in the future, such as system reliability design and multidisciplinary reliability design.