Time- and Space-Dependent Reliability-Based Design With Envelope Method

Wu, Hao; Du, Xiaoping

doi:10.1115/1.4056599

Abstract

Deterministic optimization may lead to unreliable design results if significant uncertainty exists. Including reliability constraints in reliability-based design (RBD) can solve such a problem. It is difficult to use current RBD methods to deal with time- and space-dependent reliability when responses vary randomly with respect to time and space. This study employs an envelope method for time- and space-dependent reliability for the optimal design. To achieve high accuracy, we propose an inverse envelope method that converts a time- and space-dependent limit-state function into a time- and space-independent counterpart and then use the second-order saddlepoint approximation to compute the probability of failure. The strategy is to find an equivalent most probable point for a given permitted probability of failure for each reliability constraint. To achieve high efficiency, we use a sequential optimization process to decouple the double-loop structure of RBD. The overall optimization is performed with a sequence of cycles consisting of deterministic optimization and reliability analysis. The constraints of the deterministic optimization are formulated using the equivalent most probable points. The accuracy and efficiency are demonstrated with four examples, including one mathematical problem and three engineering problems.

1 Introduction

Engineers always encounter uncertainty in material properties, part dimensions, manufacturing processes, and operational environments [1–3] in all stages of product design and development. Reliability-based design (RBD) is a typical methodology to manage uncertainty by identifying optimal design variables and ensuring satisfied reliability in the design stage [4–6].

RBD minimizes a cost-type objective function while satisfying reliability constraints. If responses are static, meaning that they are time and space invariant, we have static RBD problems, for which many mature methodologies are available, such as double-loop methods [7,8], single-loop single variable approaches [9], sequential optimization method [5,6], and safety factor approaches [10,11]. Du and Chen performed a RBD method by sequential optimization and reliability analysis so that the search of design variables and reliability analysis are executed with a serial of cycles of deterministic optimization and reliability analysis. This helps reduce the computational time [5]. Liang et al. proposed a computationally efficient RBD approach using a single-loop process where the search of the optimum design variables and the reliability analysis are performed simultaneously [9]. Yin and Du developed a modified RBD approach to mechanical component design so that the traditional safety factor design can be used without optimization and complex reliability analysis [10]. Kharmanda demonstrated that the safety factor-based RBD approach is efficient and robust with a new concept of the sequential loop procedure [11]. Tu et al. use a performance measure approach to main robustness and efficiency for evaluating reliability constraints [12].

Many responses are also time-dependent due to time varying stochastic operation conditions and product aging [13]. For instance, the function generator mechanism [14] involves time-dependent motion output. Static RBD methods are not able to handle time-dependent problems. They are extended to time-dependent RBD, and new time-dependent RBD methods have been investigated. In Refs. [15,16], a nested extreme response surface approach accurately carries out time-dependent reliability analysis and determines the optimal designs with good efficiency. The sequential optimization and reliability analysis is extended to time-dependent problems with both stationary stochastic process loads and random loads, and it effectively solves design optimization with dependent reliability constraints [17]. The equivalent most probable point method is proposed to transform the original time-variant RBD problem into an equivalent time-invariant RBD problem formulated by the performance measure approach [18]. The time-dependent concurrent reliability-based design optimization method is developed to improve the confidence of design results with reduced experimental cost and increased efficiency [19]. In Ref. [20], a sequential Kriging modeling approach is introduced to deal with the reliability constraint for time-variant RBD involving stochastic processes.

Although the time-dependent RBD methods have been developed, the most general problems should be addressed where the limit-state functions may take input of stochastic processes, random fields, and tempo-spatial variables. It is still a challenge to address the time- and space-dependent RBD as the research about this topic is limited. Several time- and space-dependent reliability analysis methods have been proposed. Hu and Mahadevan developed a surrogate modeling approach for reliability analysis of a multidisciplinary system [21]. Shi et al. presented a method for the moment estimation of the extreme response using two strategies. One strategy is combining the sparse grid technique and the fourth-moment method while the other is combining the dimensional reduction with the maximum entropy method [22]. Shi and Lu proposed an active learning Kriging method for dealing with dynamic reliability analysis for structure with temporal and spatial parameters [23]. Wei and Du combined first-order reliability method (FORM) and second-order reliability method (SORM) for the time- and space-dependent reliability analysis [24,25]. Yu and Wang developed a general decoupling approach with simulation-based method addressing reliability assessment for time- and space-variant system reliability-based design optimization [26]. Wu and Du extended the envelope method into time- and space-dependent reliability analysis [27,28]. Using the same strategy of the envelope method, we integrate the time- and space-dependent reliability analysis method with optimization in this study [29].

This study extends time- and space-dependent reliability analysis to reliability-based design. The main contributions cover both aspects of analysis and design. (1) Analysis: The new method requires an inverse reliability analysis, whose task is to find an equivalent most probable point (MPP) for a given design and reliability requirement. This is much more challenging than the current direct reliability analysis since the inverse analysis needs to call the direct analysis repeatedly. An efficient analysis method is developed that reduces the number of calls for the direct analysis. (2) Design: We integrate the inverse analysis and optimization by a more efficient means. The new design method performs deterministic optimization and analysis sequentially without embedding the later into the former. During the interactive process, the analysis algorithm performs only partial inverse reliability analysis and converges to the solution of the full inverse analysis at the end of the process. The efficiency is therefore much higher than the direct integration of optimization and inverse analysis.

The remainder of the paper is organized as follows. Section 2 gives a brief introduction to the sequential RBD. Section 3 introduces a time- and space-dependent RBD model with the envelope method and associated procedure. Section 4 presents four examples, followed by conclusions in Sec. 5.

2 Problem Statement

In this section, we define the problem that this study addresses. We also briefly review the sequential RBD, based on which the new method is developed.

2.1 Problem Statement.

This study addresses the most general RBD which includes time- and space-dependent reliability. The RBD model is defined as

{\begin{matrix} min_{d, μ_{X}} f (d) \\ \begin{matrix} s . t . Pr {y_{i} = g_{i} (d, X, z) \leq 0, \exists z \in Ω} \leq {\tilde{p}}_{f i} \\ i = 1, 2, \dots, n_{g} \\ d^{L} \leq d \leq d^{U} \end{matrix} \end{matrix}

(1)

where

f (d)

is the objective function, and

d

is the vector of design variables with their lower and upper bounds

d^{L}

and

d^{U}

⁠, respectively.

Pr {y_{i} = g_{i} (d, X, z) \leq 0, \exists z \in Ω}

is the probability of failure for the ith response y_i. The associated reliability constraint is that the probability of failure should be smaller than or equal to its allowable value

{\tilde{p}}_{f i}

or

1 - {\tilde{R}}_{i}

⁠, where [R_i] is the desired reliability. n_g is the number of constraints. In the constraint, the limit-state function is defined by

\begin{matrix} y_{i} = g_{i} (d, X, z) \end{matrix}

(2)

in which

X = [X_{1}, \dots, X_{n}]^{T}

are input random variables. The time is given by

z_{1} \in [{\underline{z}}_{1}, {\bar{z}}_{1}]

⁠, and the spatial variables are given by

z_{k} \in [{\underline{z}}_{k}, {\bar{z}}_{k}] (k = 2, \dots, m)

⁠, where m = 4. (z₂, z₃, and z₄ represent the x, y, and z coordinates, respectively.) Then,

z = [z_{1}, z_{2}, \dots, z_{m}]^{T}

is a vector of the temporal and spatial variables bounded on

Ω = [{\underline{z}}_{k}, {\bar{z}}_{k}]

⁠. Note that the time and spatial variables are quite different physically because they represent different physical aspects. However, they are mathematically the same variables, and it is the reason we present both in the vector

z

⁠.

Ω

is a rectangular domain.

Accordingly, the probability of failure for a general response y is defined by

\begin{matrix} p_{f} = \Pr {y = g (d, X, z) \leq 0, \exists z \in Ω} \end{matrix}

(3)

where

\exists

means “there exists at least one.”

Solving the time- and space-dependent RBD is time consuming since evaluating p_f in Eq. (3) is computationally expensive and p_f must be evaluated repeatedly during the optimization.

2.2 Review of Sequential Reliability-Based Design.

Sequential RBD methods in general are more efficient than double-loop methods. It decouples the optimization loop and reliability analysis loop and performs the two loops sequentially. Usually the FORM is used for reliability analysis. FORM searches for the MPP, from which the probability of failure is easily calculated.

The optimal design point is first found from the deterministic optimization loop and then FORM is performed to search for the MPP at this optimal point in the reliability analysis loop. The MPP is then used to modify a reliability constraint for the next deterministic optimization, which is followed by the next reliability analysis. This process repeats until convergence.

The deterministic optimization is formulated by

\begin{matrix} {\begin{matrix} min_{d} & f (d) \\ s . t . & g_{i} (d, T (u_{X}), z) \leq 0, i = 1, 2, \dots, n_{g} \end{matrix} \end{matrix}

(4)

where

u_{X}

is the MPP in the standard normal space for the ith reliability constraint from the reliability analysis. T(·) is the transformation operator for transforming random variables

X

to independent and standard normal variables

U

⁠. The result of the optimization is the optimal design point

\tilde{d}

⁠.

After the deterministic optimization, the reliability analysis or the MPP search is performed at

\tilde{d}

for each constraint. The MPP is obtained through an optimization problem given by

{\begin{matrix} min_{d, μ_{X}} & g_{i} (d, T (u_{X}), z) \\ s . t . & ‖ u_{X} ‖ = β_{i} \end{matrix}

(5)

where β_i is the desired reliability index. It is calculated by

β_{i} = - Φ^{- 1} ({\tilde{p}}_{f i})

⁠, and Φ⁻¹(·) is inverse cumulative density function.

The final solution can be found after a few cycles of deterministic optimization and reliability analysis. As a result, the efficiency in general is higher than solving the original RBD model directly. Since FORM may not be accurate for highly nonlinear limit-state functions, several studies employ the SORM with higher accuracy and lower efficiency [30,31].

3 Sequential Reliability-Based Design With the Envelope Method

3.1 Overview.

The aim of this work is to include time- and space-dependent reliability constraints in optimization. To achieve high accuracy, we use the second-order saddlepoint approximation (SOSPA) to calculate the reliability. To achieve high efficiency, we use sequential RBD. In the original sequential RBD, the MPP is directly related to the permitted probability of failure by

{\tilde{p}}_{f i} = Φ (- ‖ u_{XMPP} ‖)

⁠. In the present study, the MPP is not directly related to the permitted time- and space-dependent probability of failure, and the relationship is unknown and nonlinear. Note that in order to make our method applicable, the MPP is not coupled with design variables. The challenge is to find an equivalent MPP

{\tilde{u}}_{XMPP}

⁠, which satisfies

\begin{matrix} \Pr {g_{i} (d, T ({\tilde{u}}_{XMPP}), z) \leq 0, \exists z \in Ω} = {\tilde{p}}_{f i} \end{matrix}

(6)

The model of searching for the equivalent MPP becomes

{\begin{matrix} \begin{matrix} min & g_{i} (d, T (u_{X}), z) \end{matrix} \\ \begin{matrix} s . t . & ‖ u_{X} ‖ = β_{i} \end{matrix} \\ p_{f i} = \Pr {g_{i} (d, T (u_{X}), z) \leq 0, \exists z \in Ω} = {\tilde{p}}_{f i} \end{matrix}

(7)

This gives the solution ${\tilde{u}}_{XMPP} = u_{X}$ ⁠.

The envelope method can be used for reliability analysis when the equivalent MPP ${\tilde{u}}_{XMPP}$ is available from Eq. (7). To solve Eq. (7), we at first search for the MPP using FORM at the optimal point of $z$ parameter that minimizes g_i(·), and then the probability of failure p_fi is calculated using the envelop method. We update β_i iteratively until $p_{f i} = {\tilde{p}}_{f i}$ ⁠. The sequential RBD with the envelope method involves cycles of deterministic optimization and equivalent MPP search (reliability analysis).

A design variable could be a distribution parameter of a random variable. For instance, the design task is to determine the length l of a shaft, which follows a normal distribution with the mean μ_l and standard deviation σ_l. μ_l is to be determined, and σ_l is known from the tolerance of l. Then we treat μ_l as a design variable; namely, d = μ_l and l = d + X, where X is a random variable that follows a normal distribution with mean of o and standard deviation σ_l. Design variables are therefore always deterministic in this study.

3.2 Time- and Space-Dependent Reliability Analysis.

This subsection discusses the time- and space-dependent reliability using the envelope method. We have already developed a direct envelope method, which calculates the MPP and probability of failure for a given design. The new method in this study is the inverse envelope method, whose task is to find the equivalent MPP for a given probability of failure. The new method is more difficult and time consuming.

There are two cases encountered in the optimization that this study addresses.

Case 1: Calling the limit-state function once returns only one response at a specific time instant and spatial point, or a point of $z \in Ω$ ⁠. Solving for the equivalent MPP in Eq. (5) needs a double-loop procedure. This double procedure will be discussed in Sec. 3.2.1.
Case 2: Calling the limit-state function once returns all responses at all specific time instants and spatial points, or all points of $z \in Ω$ ⁠. For instance, if we call a computational fluid dynamics simulation, we obtain all 4D pressure and velocity fields with respect to time and space [3,32]. For this case, solving for the equivalent MPP in Eq. (5) needs only a single loop.

Since case 1 is much more difficult than case 2, we focus our discussions primarily on case 1. The probability of failure can be evaluated by the extreme value of the limit-state function.

\begin{aligned} p_{f} & = Pr {g (d, T (u), z) \leq 0, \exists z \in Ω} \\ = Pr {min_{z \in Ω} g (d, T (u), z) \leq 0} \end{aligned}

(8)

Equation indicates that a failure occurs if the minimum response is negative. The function of the extreme response is equivalent to the envelope function

\begin{matrix} G (X) = min_{z \in Ω} g (d, T (u), z) = g (d, T (u), \tilde{z}) \end{matrix}

(9)

where

G (X)

is the envelope function and it is the global minimum of

g (d, T (u), \tilde{z})

⁠.

G (X)

is time- and space-independent and only depends on

X

⁠.

3.2.1 Search for the Equivalent Most Probable Point Using FORM.

There are two constraints in Eq. (7), and directly computing them is too expensive due to a double-loop procedure. We propose a sequential procedure to find the equivalent MPP.

At first, the MPP search is performed by giving the initial reliability index β at

{\tilde{z}}^{(0)} = [z_{1}^{0}, z_{2}^{0}, \dots, z_{m}^{0}]

with the following model:

{\begin{matrix} min_{u} & g (d, T (u), z) \\ s . t . & ‖ u ‖ = β \end{matrix}

(10)

Equation (10) produces the initial MPP $u_{XMPP}^{(1)}$ ⁠. The next analysis is to find the optimal time and space parameter ${\tilde{z}}^{(1)}$ by fixing the random variables at its realization $u_{XMPP}^{(1)}$ ⁠. The next optimal ${\tilde{z}}^{(1)}$ can be obtained by another optimization ${\tilde{z}}^{(1)} = \underset{z \in Ω}{argmin} g (d, T (u_{XMPP}^{(1)}), z)$ ⁠. The details of solving ${\tilde{z}}^{(1)}$ are illustrated in Sec. 3.2.3. The process repeats until convergence. This sequential procedure produces the worst-case MPP $u_{XMPP}$ and associated optimal $\tilde{z}$ ⁠. After the worst-case MPP is found, the probability of failure with FORM is estimated by $p_{f} = Pr {g (d, T (u), z) \leq 0, \exists z \in Ω} = Φ (- ‖ u_{XMPP} ‖)$ ⁠. However, FORM may not be accurate enough. As a result, we use second-order envelope method to achieve high accuracy.

3.2.2 The Envelope Method.

The envelope function is in general nonlinear, and it is tangent to all the instantaneous limit-state functions. It can be approximated by the second-order approximation method as a quadratic function at the MPP. As indicated in Ref. [14], the MPP of the envelope function is the worst-case MPP of the limit-state function. We also need the gradient

\nabla G

and Hessian matrix

H

at the MPP of the envelope function. The quadratic function is formed as follows:

\begin{matrix} G (U) = a + b^{T} U + U^{T} C U \end{matrix}

(11)

where

{\begin{matrix} a = \frac{1}{2} {(u_{XMPP})}^{T} H u_{XMPP} - \nabla G {(u_{XMPP})}^{T} u_{XMPP} \\ b = \nabla G (u_{XMPP}) - H u_{XMPP} = ({\tilde{b}}_{1}, {\tilde{b}}_{2}, \dots, {\tilde{b}}_{n}) \\ C = \frac{1}{2} H = diag ({\tilde{c}}_{1}, {\tilde{c}}_{2}, \dots, {\tilde{c}}_{n}) \end{matrix}

(12)

\nabla G (u_{XMPP}) = {[{\frac{\partial G}{\partial U_{1}} |}_{u_{XMPP}}, \dots ., {\frac{\partial G}{\partial U_{n}} |}_{u_{XMPP}}]}^{T}

is the gradient of the envelope function, and

H

is the Hessian matrix given by

\begin{matrix} H = {[\begin{matrix} \frac{\partial^{2} G}{\partial U_{1}^{2}} & \dots & \frac{\partial^{2} G}{\partial U_{1} \partial U_{n}} \\ ⋮ & ⋱ & ⋮ \\ \frac{\partial^{2} G}{\partial U_{n} \partial U_{1}} & \dots & \frac{\partial^{2} G}{\partial U_{n}^{2}} \end{matrix}]}_{u_{XMPP}} \end{matrix}

(13)

The envelope function

G (X)

at

u_{XMPP}

is given by

\begin{matrix} G (U) = min_{z \in Ω} g (d, U, z) = g (U, \tilde{z}) |_{u_{XMPP}} \end{matrix}

(14)

where

\tilde{z} = [{\tilde{z}}_{1}, \dots, {\tilde{z}}_{m}]

is the optimal point where the global minimum of function

g (d, U, z)

occurs. For an easier demonstration without loss of generality, we use a two-parameter case

\tilde{z} = [{\tilde{z}}_{1}, {\tilde{z}}_{2}]

as an example to derive the Hessian matrix and gradient. The envelope function satisfies the following equations:

\begin{matrix} {\begin{matrix} \dot{g} (d, U, {\tilde{z}}_{1}, z_{2}) = 0 \\ \dot{g} (d, U, z_{1}, {\tilde{z}}_{2}) = 0 \end{matrix} \end{matrix}

(15)

where

\dot{g}

is the derivative of g with respect to z_i.

The first derivative of

G (U)

at

u_{XMPP}

is

\begin{matrix} \frac{\partial G}{\partial U_{i}} = \frac{\partial g}{\partial U_{i}} + \frac{\partial g}{\partial {\tilde{z}}_{1}} \frac{\partial {\tilde{z}}_{1}}{\partial U_{i}} + \frac{\partial g}{\partial {\tilde{z}}_{2}} \frac{\partial {\tilde{z}}_{2}}{\partial U_{i}} \end{matrix}

(16)

Plugging Eq. into Eq. yields

\begin{matrix} \frac{\partial G}{\partial U_{i}} = \frac{\partial g}{\partial U_{i}} \end{matrix}

(17)

The gradient of the envelope function

\nabla G

is equal to the gradient of the limit-state function

\nabla g

at the MPP. Subsequently, the second derivative of

G (U)

with respect to input random variables U_j at

u_{XMPP}

is

\begin{matrix} \frac{\partial^{2} G}{\partial U_{i} \partial U_{j}} = \frac{\partial}{\partial U_{j}} (\frac{\partial G}{\partial U_{i}}) = \frac{\partial}{\partial U_{j}} (\frac{\partial g}{\partial U_{i}}) \\ = \frac{\partial^{2} g}{\partial U_{i} \partial U_{j}} + \frac{\partial^{2} g}{\partial U_{i} \partial {\tilde{z}}_{1}} \frac{\partial {\tilde{z}}_{1}}{\partial U_{j}} + \frac{\partial^{2} g}{\partial U_{i} \partial {\tilde{z}}_{2}} \frac{\partial {\tilde{z}}_{2}}{\partial U_{j}} \end{matrix}

(18)

Taking the derivative of Eq. with respect to U_j yields

\begin{matrix} \frac{\partial \dot{g}}{\partial U_{j}} + \frac{\partial \dot{g}}{\partial {\tilde{z}}_{1}} \frac{\partial {\tilde{z}}_{1}}{\partial U_{j}} = 0 \end{matrix}

(19)

Rearranging the equation, we have

\begin{matrix} \frac{\partial {\tilde{z}}_{1}}{\partial U_{j}} = - \frac{\partial \dot{g}}{\partial U_{j}} / \frac{\partial \dot{g}}{\partial {\tilde{z}}_{1}} = - \frac{\partial^{2} g}{\partial {\tilde{z}}_{1} \partial U_{j}} / \frac{\partial^{2} g}{\partial {\tilde{z}}_{1}^{2}} \end{matrix}

(20)

Similarly,

\begin{matrix} \frac{\partial {\tilde{z}}_{2}}{\partial U_{j}} = - \frac{\partial \dot{g}}{\partial U_{j}} / \frac{\partial \dot{g}}{\partial {\tilde{z}}_{2}} = - \frac{\partial^{2} g}{\partial {\tilde{z}}_{2} \partial U_{j}} / \frac{\partial^{2} g}{\partial {\tilde{z}}_{2}^{2}} \end{matrix}

(21)

By plugging Eqs. and into Eq. , we obtain the Hessian matrix at

u_{X}

and

{\tilde{z}}_{i}

⁠.

\begin{matrix} {\frac{\partial^{2} G}{\partial U_{i} \partial U_{j}} |}_{u^{*}, {\tilde{z}}_{k}} = \frac{\partial^{2} g}{\partial U_{i} \partial U_{j}} - \sum_{k = 1}^{2} \frac{\partial^{2} g}{\partial U_{i} \partial {\tilde{z}}_{k}} \frac{\partial^{2} g}{\partial U_{j} \partial {\tilde{z}}_{k}} / \frac{\partial^{2} g}{\partial {\tilde{z}}_{k}^{2}} \end{matrix}

(22)

We use the forward finite difference method with step size $δ = max (| u | / 1000, ε)$ ⁠, where $ε = 10^{- 4}$ ⁠, to calculate the derivations in Eq. (22).

The worst-case MPP $u_{MPP}$ ⁠, gradient $\nabla G$ ⁠, and the Hessian matrix are now available. The SOSPA is then used to estimate the probability of failure p_fi after the envelope function is approximated. SOSPA is in general more accurate than FORM because it yields an accurate probability estimation especially in the tail area of a distribution. The details of the implementation of SOSPA are given in Refs. [31,33]. If p_fi is not equal to ${\tilde{p}}_{f i}$ ⁠, we should update the reliability index β_i [34], and the MPP search is executed again using Eq. (10). The process is repeated until $p_{f i} = {\tilde{p}}_{f i}$ ⁠. When the probability of failure is equal to the required one, it will produce the equivalent MPP ${\tilde{u}}_{XMPP}$ ⁠. The detailed procedure is given as follows.

Step 1: Set the initial reliability index β^(h) = β⁽¹⁾, the initial optimal point $z^{(h - 1)} = z^{(0)}$ ⁠, and the initial MPP $u_{XMPP}^{(h - 1)} = u_{0}$ ⁠. The iteration counter h = 1.
Step 2: Search for the inverse MPP using Eq. (10), and obtain the MPP $u_{MPP}^{(h)}$ by solving the following optimization model:
${\begin{matrix} min_{u} & g (d, T (u), z) \\ s . t . & ‖ u ‖ = β^{(h)} \end{matrix}$
(23)
Step 3: Determine the optimal point ${\tilde{z}}^{(h)}$ ⁠. The optimal point ${\tilde{z}}^{(h)}$ minimizes the limit-state function.
$\begin{matrix} {\tilde{z}}^{(h)} = \underset{z \in Ω}{argmin} g (d, T (u_{XMPP}^{(h)}), z) \end{matrix}$
(24)
The optimization method we use in this study is global efficient optimization (EGO) [35].
Step 4: Check the convergence criterion.
$\begin{matrix} ‖ u_{XMPP}^{(h)} - u_{XMPP}^{(h - 1)} ‖ \leq ε_{1} \end{matrix}$
(25)
The tolerance $ε_{1}$ takes a small positive value, for example, 10⁻⁴. If $‖ u_{XMPP}^{(h)} - u_{XMPP}^{(h - 1)} ‖ \leq ε_{1}$ ⁠, terminate the iteration. Otherwise, set h = h + 1, and return to step 2.
Step 5: Calculate the gradient $\nabla G$ and Hessian matrix $H$ of the envelope function at $u_{XMPP}^{(h)}$ ⁠. Calculate the probability of failure p_f using SOSPA from the above information $u_{XMPP}^{(h)}$ ⁠, gradient $\nabla G$ ⁠, and Hessian matrix $H$ ⁠.

3.2.3 Find the Optimal $\tilde{z}$ ⁠.

The global minimum value of

g (d, T (u_{XMPP}), z)

occurs at

\tilde{z} = [{\tilde{z}}_{1}, {\tilde{z}}_{2}]

⁠, which is given by

\begin{matrix} \tilde{z} = \underset{z \in Ω}{argmin} g (d, T (u_{XMPP}), z) \end{matrix}

(26)

Finding the optimal point still needs sequential loops. The first partial derivatives of the limit-state function with respect to z_i at MPP are obtained as follows:

\begin{matrix} {\begin{matrix} \frac{\partial g (d, T (u_{MPP}), z_{1}, z_{2})}{\partial z_{1}} = 0 \\ \frac{\partial g (d, T (u_{MPP}), z_{1}, z_{2})}{\partial z_{2}} = 0 \end{matrix} \end{matrix}

(27)

The optimal point $\tilde{z} = [{\tilde{z}}_{1}, {\tilde{z}}_{2}]$ can be obtained by solving Eq. (27). For an explicit and simple limit-state function, the solution of the derivative equations can be obtained analytically. For an implicit and complicated limit-state function, EGO can be used. EGO has been widely used in various areas because it can efficiently search for the global optimum [25].

EGO generates the training points of the input

z^{in} = [z_{1}^{(i)}, z_{2}^{(i)}]_{i = 1, 2, \dots, h} = [z^{(i)}]_{i = 1, 2, \dots, h}

⁠, where h is the number of initial training points, and the training points of the output dataset are

y^{in} = [g (T (u_{XMPP}), z^{(i)})]_{i = 1, 2, \dots, h}

⁠. Once the training dataset

(z^{in}, y^{in})

is ready, the next step is to train the initial model using the Gaussian process regression. The initial surrogate model is

\hat{y} = g (z) = g (T (u_{XMPP}), z) = F (z)^{T} γ + e (z)

⁠, where

F (z)^{T} γ

is a deterministic term,

e (z)

is a vector of regression functions,

γ

is a vector of regression coefficients, and

e (z)

is a stationary Gaussian process with zero mean and a covariance is

Cov (e (z_{i}), e (z_{j})) = σ_{e}^{2} C (z_{i}, z_{j})

⁠, where

σ_{e}^{2}

is process variance and C(·, ·) is the correlation function [36]. The initial model may not be accurate; hence new training points are then added one by one so that the model is continuously refined. EGO selects a new training point

z_{new}

using the expected improvement (EI) metric defined by

\begin{matrix} z_{new} = \underset{z}{argmin EI} (z) \end{matrix}

(28)

where EI is computed by

\begin{aligned} EI (z) & = E [max (y^{*} - y, 0)] \\ = (y^{*} - μ (z)) Φ (\frac{y^{*} - μ (z)}{σ (z)}) + σ (z) ϕ (\frac{y^{*} - μ (z)}{σ (z)}) \end{aligned}

(29)

where

y^{*} = min_{i = 1, 2, \dots, k} g (z^{i})

⁠,

μ (z)

and

σ (z)

are the mean and standard deviation of

\hat{y}

⁠, respectively, and ϕ(·) is the probability density function.

3.3 Extension to Case 2.

After finishing case 1, we now discuss case 2. We have the complete time- and space-dependent responses with respect to

z

in

Ω

from a single call of the limit-sate function. This case has the most general limit-state function

y (z) = g (X, F (z), z)

⁠. The probability of failure is calculated by

\begin{matrix} p_{f} = Pr {g (X, F (z), z) \leq 0, \exists z \in Ω} = Pr {min_{z \in Ω} g (X, F (z), z) \leq 0 \leq 0} \end{matrix}

(30)

Equation indicates that failure happens when the minimum value of the limit-state function

g (X, F (z), z)

is negative. Since calling the limit-function returns a complete hypersurface of the response

y (z)

with respect to time and space, the minimum value

min_{z \in Ω} y (z)

is known. Thus, the MPP can be obtained from the following model:

{\begin{matrix} min & \sqrt{u^{T} u} \\ s . t . & min_{z \in Ω} y (z) = 0 \end{matrix}

(31)

where

min_{z \in Ω} y (z)

is a function of

u

and is obtained by calling the limit-state once at

u

⁠. Therefore, a single single-loop MPP is needed. This is more efficient than the sequential loop approach. The expansion optimal linear estimation (EOLE) method [37] is used to expand the random field response with respect to independent standard Gaussian random variables

ξ = (ξ_{1}, ξ_{2}, \dots, ξ_{r})

⁠, where r is the dimension of

ξ

⁠. Then, the limit-state function becomes

y = g (\tilde{X}, z)

⁠, where

\tilde{X} = (X, ξ)

⁠. Thus, the proposed method can still work. Take EOLE as an example for a two-dimensional random field

F (z)

with

z \in (z_{1}, z_{2})

⁠. z₁ and z₂ are discretized into

n_{z_{1}} n_{z_{2}}

points, and the autocorrelation coefficient matrix is given by

\begin{matrix} Σ = {[ρ (z_{i}, z_{j})]}_{n_{z_{1}} n_{z_{2}} \times n_{z_{1}} n_{z_{2}}} \end{matrix}

(32)

where

ρ (z_{i}, z_{j})

is the correlation between two points

z_{i} (i = 1, 2, \dots, n_{z_{1}} n_{z_{2}})

and

z_{j} (j = 1, 2, \dots, n_{z_{1}} n_{z_{2}})

in the domain of

F (z)

⁠.

Then the EOLE expansion is given by

\begin{matrix} F (ξ, z) \approx μ (z) + σ (z) \sum_{k = 1}^{r} \frac{ξ_{k}}{\sqrt{λ_{k}}} ϕ_{k}^{T} Σ (:, z), k = 1, 2, \dots, r \end{matrix}

(33)

where

μ (z)

is the mean of

F (z)

⁠, and

σ (z)

is the standard deviation of

F (z)

⁠. ξ_k (k = 1, 2, …, r) are independent standard normal variables,

λ = (λ_{1}, λ_{2}, \dots, λ_{r})

is the eigenvalue vector, and ϕ₁, ϕ₂, …, ϕ_r are the corresponding eigenvectors obtained from autocorrelation coefficient matrix

Σ

⁠. Note that r is determined as the smallest integer that meets the following criterion:

\begin{matrix} \frac{\sum_{j = 1}^{r} λ_{k}}{\sum_{j = 1}^{n_{z_{1}} n_{z_{2}}} λ_{k}} \geq η \end{matrix}

(34)

where η is a hyperparameter determining the accuracy of the expansion. It takes a value close to, but not larger than 1. The smaller is η, the less accurate is the expansion. If η = 1, the expansion is exact at the points of discretization. η could be set to 0.9–0.99 [37].

3.4 Procedure of Sequential Reliability-Based Design With Envelope Method.

After discussing all the details, we now summarize the procedure of the proposed method.

Step 1: Initialize parameters.
Step 2: Perform deterministic optimization. For k = 1, solve deterministic optimization at means of random variable. For k > 1, perform the following deterministic optimization using the equivalent MPP ${\tilde{u}}_{i, XMPP}^{(k - 1)}$ obtained from the (k − 1)th cycle. The solution is $d^{(k)}$ ⁠.
${\begin{matrix} min_{d, μ_{X}} & f (d) \\ s . t . & g_{i} (d, T (u_{i, XMPP}^{(k - 1)}), z) \leq 0, i = 1, 2, \dots, n \\ d^{L} \leq d \leq d^{U} \end{matrix}$
(35)
Step 3: Perform time- and space-dependent reliability analysis at $d^{(k)}$ for each constraint. At first, search for the equivalent MPP given β^(k). Obtain ${\tilde{u}}_{XMPP}^{(k)}$ ⁠, gradient $\nabla G ({\tilde{u}}_{XMPP}^{(k)})$ ⁠, and Hessian matrix $H_{{\tilde{u}}_{XMPP}^{(k)}}$ ⁠. Note that if the inputs of limit-state function are random variables and random fields, the method in Sec. 3.3 is used to find the ${\tilde{u}}_{XMPP}^{(k)}$ ⁠. Next, calculate the probability of failure $p_{f i}^{(k)}$ using SOSPA.
Step 4: Check the convergent criterion by
$\begin{matrix} ε = | \frac{p_{_{f i}}^{(k)} - [P_{f i}]}{[P_{f i}]} | \leq ε_{2} \end{matrix}$
(36)
If convergence is reached, the optimal solution is found at $d^{(k)}$ ⁠, and the process stops. Otherwise update the reliability index β^(k+1) by
$\begin{matrix} {\tilde{u}}_{XMPP}^{(k + 1)} = β^{(k + 1)} \frac{\nabla g ({\tilde{u}}_{XMPP}^{(k)})}{‖ {\tilde{u}}_{XMPP}^{(k)} ‖} \end{matrix}$
(37)
The flowcharts of the proposed method are given in Figs. 1 and 2.
The new method has the same limitations as other MPP-based reliability methods. Its accuracy will deteriorate if there are multiple MPPs, which occur when the envelope function is piecewise and indifferentiable. We could use the method of composite limit-state function for time-dependent reliability to solve this problem [38]. Its extension to time- and space-dependent reliability needs further investigation.
The overall process can typically converge within five cycles. Since the process starts from deterministic optimization in the first cycle with mean values of random variables, the reliability of an active constraint is approximately 50% at the optimal point, much lower than the required reliability. In the following reliability analysis, the MPP is found for the given design point and the required reliability. Then in the deterministic optimization of the second cycle, the MPP is plugged into the corresponding constraint function, which guarantees that the reliability will be improved and much closer to the target reliability, but still smaller than the target reliability. Then the reliability analysis is conducted again. More cycles are performed until convergence. Hence the design process always starts from the infeasible region, and the design point will gradually approach constraint boundaries and the final optimal point. On other words, the reliability will gradually be improved until reaching the target. The entire design process will likely converge if both deterministic optimization and reliability analysis converge. The method, however, may not work if the reliability analysis does not identify a global MPP; for instance, a local MPP is found in one and a global MPP is found in another cycle.

Fig. 1

View large Download slide

Sequential RBD with the envelope method

Fig. 2

View large Download slide

Inversed time- and space-dependent reliability analysis

4 Examples

A mathematical example is provided to show the feasibility of the proposed method. Three engineering examples are then used to demonstrate the computational efficiency and accuracy of the proposed method compared with double-loop method using the direct second-order reliability method (SORM/DL) and double-loop method with FORM (FORM/DL). The accuracy is assessed by the probability of failure obtained from Monte Carlo simulation (MCS) at the optimal points while the efficiency is measured by the number of function calls. The percentage error is computed by

\begin{matrix} ε = \frac{| p_{f} - p_{f}^{MCS} |}{p_{f}^{MCS}} \times 100 % \end{matrix}

(38)

where p_f is the result from a non-MCS method while

p_{f}^{MCS}

is from MCS.

Since case 1 is more complicated by case 2 and case 1 is the focus of this study, all the examples discussed here belong to case 1.

4.1 Example 1: A Mathematical Problem.

This example is a mathematical problem. Two independent random variables X₁ and X₂ are normally distributed with

X_{1} \sim N (μ_{X_{1}}, 0.6)

and

X_{2} \sim N (μ_{X_{2}}, 0.6)

⁠. The time t varies over the interval [0,1], and the spatial variable s changes over the interval [7]. The design variables are

μ_{X_{1}}

and

μ_{X_{2}}

⁠. The limit-state function is defined by

\begin{matrix} g (X, s, t) = \frac{80}{(s^{2} - s + X_{1}^{2} + 8 X_{2} + t - \sin (t) + 5)} - 1 \end{matrix}

(39)

The RBD model is defined as follows:

{\begin{matrix} min_{μ_{X_{1}} μ_{X_{2}}} & f = - (μ_{X_{1}} + μ_{X_{2}}) \\ s . t . & \Pr {g (X, s, t) > 0} \geq Φ (β) \\ - 5 \leq μ_{X_{1}} \leq 10 \\ - 5 \leq μ_{X_{2}} \leq 10 \end{matrix}

(40)

The allowable reliability index β is 3. The optimal results are achieved by following the detailed steps of the proposed method. The proposed method produces the following results: design variables are w = 2.7177 and h = 4.5691; the objective function is f = −7.2867. The results by SORM/DL and FORM/DL are provided in Table 1. The proposed method and SORM/DL have the same accuracy as their errors are the same. This shows that the former method converges to the same solution from the latter method; the former method is much more efficient than the latter method as indicated by the number of function calls N_call in Table 1. As expected, the proposed method is more accurate than FORM/DL that uses a first-order approximation. The proposed method is also more efficient than FORM/DL although the former and latter methods use second and first approximations, respectively. The convergence history is shown in Fig. 3.

Fig. 3

View large Download slide

Iterative process of the objective function and design variables

Table 1

Results of Example 1

Method	New method	SORM/DL	FORM/DL
Obj.	−7.2867	−7.2867	−7.3096
μ	(2.7177, 4.5691)	(2.7177, 4.5691)	(2.7275, 4.5821)
$p_{f}^{MCS}$ (10⁻³)	1.3679	1.3666	1.4942
Error (%)	1.33	1.23	10.68
N_calls	341	4317	747

4.2 Example 2: Three Bar System.

A truss structure is shown in Fig. 4. Each bar of the system has its cross-sectional area A_i and modulus of elasticity, E_i, i = 1, 2, 3. The coefficient of thermal expansion of all bars is $α = 12 \times 10^{- 6}^{\circ} C^{- 1}$ ⁠. The temperature change is related to the installation height of the truss structure and is given by $Δ T = T e^{- 0.01 {(Δ h^{2} + 2 Δ h - 2)}^{2}}$ ⁠, where Δh ∈ [1, 6] m is the difference between two different installation heights. A downward force P = P₀(0.9 + 0.1cos(0.2t)) is applied at joint A, where t ∈ [0, 10] years. All the random variables are given in Table 2. The design variables are the cross-sectional areas of the bars $μ_{A_{1}}$ ⁠, $μ_{A_{2}}$ ⁠, and $μ_{A_{3}}$ ⁠.

Fig. 4

View large Download slide

A truss structure

Table 2

Random variables of Example 2

Variable	Mean	Standard deviation	Distribution
A₁(mm²)	$μ_{A_{1}}$	0.6	Normal
A₂(mm²)		0.6	Normal
A₃(mm²)	$μ_{A_{3}}$	0.6	Normal
E₁(GPa)	200	20	Normal
E₂(GPa)	200	20	Normal
E₃(GPa)	200	20	Normal
P₀(KN)	40	6	Normal
L_AD(mm)	200	2	Normal
L_AB(mm)	231	2.31	Normal
L_AC(mm)	283	2.83	Normal
$T (^{\circ} C)$	35	7	Normal
σ_yield(GPa)	7.5 × 10⁸	4 × 10⁷	Normal

Variable	Mean	Standard deviation	Distribution
A₁(mm²)	$μ_{A_{1}}$	0.6	Normal
A₂(mm²)		0.6	Normal
A₃(mm²)	$μ_{A_{3}}$	0.6	Normal
E₁(GPa)	200	20	Normal
E₂(GPa)	200	20	Normal
E₃(GPa)	200	20	Normal
P₀(KN)	40	6	Normal
L_AD(mm)	200	2	Normal
L_AB(mm)	231	2.31	Normal
L_AC(mm)	283	2.83	Normal
$T (^{\circ} C)$	35	7	Normal
σ_yield(GPa)	7.5 × 10⁸	4 × 10⁷	Normal

The objective function is to minimize the weight given by

f = μ_{A_{1}} μ_{L_{AD}} + μ_{A_{2}} μ_{L_{AB}} + μ_{A_{3}} μ_{L_{AC}}

There are two failure modes for this truss structure. The first failure mode is that the perpendicular displacement of joint A, denoted as Δδ, is greater than the allowable displacement δ₀, and the limit-state function is defined by

\begin{matrix} g_{1} (X, s, t) = δ_{0} - Δ δ \end{matrix}

(41)

where δ₀ = 0.64, and the perpendicular displacement of joint A is calculated by

Δ δ = \frac{A}{B}

where

\begin{aligned} A & = L_{AD} (P A_{1} E_{1} L_{AC} \cos θ_{1}^{2} + P A_{2} E_{2} L_{AB} \cos θ_{2}^{2} + A_{1} A_{3} E_{1} E_{3} L_{AC} T α \cos θ_{1}^{2} + A_{2} A_{3} E_{2} E_{3} L_{AB} T α \cos θ_{2}^{2} \\ + A_{1} A_{2} E_{1} E_{2} T α (L_{AB} \sin α_{1} \cos θ_{2}^{2} + L_{AC} \sin θ_{2} \cos θ_{1}^{2} + L_{AC} \sin θ_{1} \cos θ_{2} \cos θ_{1} + L_{AB} \sin θ_{2} \cos θ_{2} \cos θ_{1})) \end{aligned}

\begin{aligned} B = A_{1} A_{3} E_{1} E_{3} L_{AC} \cos θ_{1}^{2} + A_{2} A_{3} E_{2} E_{3} L_{AB} \cos θ_{2}^{2} + A_{1} A_{2} E_{1} E_{2} L_{AD} (\sin θ_{2}^{2} \cos θ_{1}^{2} + \sin θ_{1}^{2} \cos θ_{2}^{2} + 2 \sin θ_{1} \sin θ_{2} \cos θ_{1} \cos θ_{2}) \end{aligned}

θ_{1} = \arctan (\frac{L_{AD}}{\sqrt{L_{AB}^{2} - L_{AD}^{2}}})

θ_{2} = \arctan (\frac{L_{AD}}{\sqrt{L_{AC}^{2} - L_{AD}^{2}}})

The second failure mode occurs when the stress at the joint is greater than the yield strength, and the limit-state function is defined by

\begin{matrix} g_{2} (X, s, t) = σ_{y} - σ \end{matrix}

(42)

where

σ = \frac{C}{D}

\begin{aligned} C & = A_{1} A_{2} A_{3}^{2} E_{1} E_{2} E_{3} T α (L_{AC} \cos α_{1}^{2} \sin θ_{2} + L_{AB} \cos α_{2}^{2} \sin θ_{1} - 2 L_{AD} \sin θ_{1} \sin θ_{2} \cos θ_{1} \cos θ_{2} \\ + L_{AC} \cos θ_{1} \cos θ_{2} \sin θ_{1} + L_{AB} \cos θ_{1} \cos θ_{2} \sin θ_{2} - L_{AD} \cos θ_{2}^{2} \sin θ_{1}^{2} \\ - L_{AD} \cos θ_{1}^{2} \sin θ_{2}^{2}) + A_{2} A_{3} E_{2} E_{3} F L_{AB} L_{AB} \cos θ_{2}^{2} + A_{1} A_{3} E_{1} E_{3} F L_{AC} \cos θ_{1}^{2} \end{aligned}

\begin{aligned} D & = A_{1} A_{2} E_{1} E_{2} L_{AD} (\cos θ_{1}^{2} \sin θ_{2}^{2} + \cos θ_{2}^{2} \sin θ_{1}^{2}) + A_{1} A_{3} E_{1} E_{3} L_{AC} \cos θ_{1}^{2} \\ + A_{2} A_{3} E_{2} E_{3} L_{AB} \cos θ_{2}^{2} + 2 A_{1} A_{2} E_{1} E_{2} L_{AD} \sin θ_{1} \sin θ_{2} \cos θ_{1} \cos θ_{2} \end{aligned}

The RBD model is given by

{\begin{matrix} min_{μ_{A_{1}} μ_{A_{2}} μ_{A_{3}}} f = μ_{A_{1}} μ_{L_{AC}} + μ_{A_{2}} μ_{L_{AD}} + μ_{A_{3}} μ_{L_{AC}} \\ s . t . Pr {g_{1} (X, s, t) = δ_{0} - Δ δ > 0} \geq Φ (β_{1}) \\ Pr {g_{2} (X, s, t) = σ_{y} - σ > 0} \geq Φ (β_{2}) \\ 10 \leq μ_{A_{1}} \leq 60 \\ 10 \leq μ_{A_{1}} \leq 80 \\ 10 \leq μ_{A_{1}} \leq 60 \end{matrix}

(43)

The allowable reliability indexes are β_i = 2.3, i = 1, 2. The results from three are given in Table 3. The solution of the proposed method is μ = [53.41, 38.83, 60] mm², which satisfy all the design constraints. The probabilities of failure at the optimal design points are $p_{f 1}^{MCS} = 0.0107$ and $p_{f 2}^{MCS} = 0.0107$ ⁠, which is the same as the probabilities of failure with SORM/DL method. Besides, the proposed method has higher efficiency than SORM/DL in terms of the N_calls. Consequently, the proposed method achieves the most accurate and efficient results.

Table 3

Results of Example 2

Method	New method	SORM/DL	FORM/DL
Obj. ( × 10⁴ mm³)	3.6616	3.6624	3.7189
μ (mm)	(53.41, 38.83,60)	(53.31, 38.95,60)	(53.05, 38.57,60)
$p_{f 1}^{MCS}$	0.0107	0.0107	0.0148
$ε_{p_{f 1}}$ (%)	0.25	0.46	4.89
$p_{f 2}^{MCS}$	0.0107	0.0107	0.0112
$ε_{p_{f 2}}$ (%)	0.24	0.46	4.89
N_calls	2194	23,500	7264

Method	New method	SORM/DL	FORM/DL
Obj. ( × 10⁴ mm³)	3.6616	3.6624	3.7189
μ (mm)	(53.41, 38.83,60)	(53.31, 38.95,60)	(53.05, 38.57,60)
$p_{f 1}^{MCS}$	0.0107	0.0107	0.0148
$ε_{p_{f 1}}$ (%)	0.25	0.46	4.89
$p_{f 2}^{MCS}$	0.0107	0.0107	0.0112
$ε_{p_{f 2}}$ (%)	0.24	0.46	4.89
N_calls	2194	23,500	7264

4.3 Example 3: A Cantilever Beam.

Figure 5 shows that the end of cantilever beam is subjected to two forces F₁ and F₂. The length of the cantilever beam L is 100 in. The objective is to minimize the volume f = μ_wμ_hL, where w and h represent the width and height of the beam cross section, respectively. There are two failure modes. The first failure mode is that the stress at the fixed end is greater than the allowable yield stress S_y, and the second failure mode is that the tip displacement of the beam is greater than the allowable displacement D₀ = 2.5 in. The two limit-state functions are given by

\begin{matrix} g_{1} (X, s, t) = S_{y} - \frac{6 L}{w h} (\frac{F_{1}}{h} - \frac{F_{2}}{w}) \end{matrix}

(44)

\begin{matrix} g_{2} (X, s, t) = D_{0} - \frac{4 L^{3}}{E w h} \sqrt{{(\frac{F_{1}}{h^{2}})}^{2} + {(\frac{F_{2}}{w^{2}})}^{2}} \end{matrix}

(45)

Fig. 5

View large Download slide

A cantilever beam

The distributions of the random design variables and random parameters are given in Table 4. The force F₁ is a non-stationary Gaussian random field, whose mean is

μ_{F_{1}} = 500 e^{0.01 ({(s - (1 / 2))}^{2} + {(t - 6)}^{2})} lbs

and standard deviation is

σ_{F_{1}} = 50 lbs

⁠. The spatial variable is s ∈ [0, 1] in. and temporal variable is t ∈ [0, 10] years. The autocorrelation coefficient function is given by

\begin{matrix} ρ_{F_{1}} (s_{1}, t_{1}; s_{2} t_{2}) = \exp [- {(\frac{s_{1} - s_{2}}{10})}^{2} - {(\frac{t_{1} - t_{2}}{10})}^{2}] \end{matrix}

(46)

Table 4

Distributions of variables in Example 3

Variable (unit)	Mean	Standard deviation	Distribution
w(in.)	μ_w	0.01	Normal
h(in.)	μ_h	0.01	Normal
F₂(lb)	1 × 10³	1 × 10²	Normal
E(psi)	2.9 × 10⁷	1 × 10⁵	Normal
S_y(psi)	3.9 × 10⁴	500	Weibull

The RBD model is formulated as

\begin{matrix} {\begin{matrix} min_{μ_{w} μ_{h}} f = μ_{w} μ_{h} L \\ s . t . \Pr {g_{1} (X, s, t) > 0} \geq Φ (β_{1}) \\ \Pr {g_{2} (X, s, t) > 0} \geq Φ (β_{2}) \\ 1 \leq μ_{w} \leq 4 \\ 1 \leq μ_{h} \leq 4 \end{matrix} \end{matrix}

(47)

The allowable reliability indexes are β_i = 3, i = 1, 2. Table 5 shows that the optimal design variables are w = 3.9541 in. and h = 2.2531 in. in, and the objective function is f = 890.9152 in.³ by the proposed method. The probabilities of failure obtained at the optimal design variables by MCS are $p_{f 1}^{MCS} = 1.3633 \times 10^{- 4}$ and $p_{f 2}^{MCS} = 1.3645 \times 10^{- 4}$ ⁠. The results are more accurate than those of SORM/DL and FORM/DL methods. The proposed method is the most efficient method as the number of function calls is 1358 compared with those of SORM/DL and FORM/DL, which are 12,904 and 2098, respectively. In general, the proposed method is the best method for this example due to its high accuracy and efficiency.

Table 5

Results of Example 3

Method	New method	SORM/DL	FORM/DL
Obj. (in.³)	889.3705	889.3766	890.0996
μ (in.)	(3.9375, 2.2587)	(3.9373, 2.2589)	(3.9388, 2.2573)
$p_{f 1}^{MCS} (\times 10^{- 4})$	1.3435	1.3429	1.7410
$ε_{p_{f 1}}$ (%)	0.4740	0.5184	0.918
$p_{f 2}^{MCS}$ ( × 10⁻⁴)	1.3399	1.3348	1.4195
$ε_{p_{f 2}}$ (%)	0.7407	1.1185	5.1561
N_calls	1358	12,904	2098

Method	New method	SORM/DL	FORM/DL
Obj. (in.³)	889.3705	889.3766	890.0996
μ (in.)	(3.9375, 2.2587)	(3.9373, 2.2589)	(3.9388, 2.2573)
$p_{f 1}^{MCS} (\times 10^{- 4})$	1.3435	1.3429	1.7410
$ε_{p_{f 1}}$ (%)	0.4740	0.5184	0.918
$p_{f 2}^{MCS}$ ( × 10⁻⁴)	1.3399	1.3348	1.4195
$ε_{p_{f 2}}$ (%)	0.7407	1.1185	5.1561
N_calls	1358	12,904	2098

4.4 Example 4: Thermal Deflection of a Bimetallic Beam.

Thermal expansion or contraction of a bimetallic beam occurs due to temperature change. The temperature in the operating room varies during the day and night and is given by ΔT = T(0.01sin(0.1t) + 1), where t ∈ [0, 24] h. The typical bimetallic beam consists of two materials bonded together: copper and invar. E_C is Young’s modulus of copper, and E_I is Young’s modulus of invar. The length of the beam depends on the location of the installation, which is given by L = L₀(− s² + s + 1), where s ∈ [0, 1] m. When the temperature change as a thermal load applies on the beam, the beam will deflect in the perpendicular direction at the right end side shown in Fig. 6. The design variables are $d = (h, w)$ ⁠, where h and w represent the height and width of the cross-sectional area of the beam, and their means are μ_h and μ_w, which are to be determined. All the random variables are listed in Table 6.

Fig. 6

View large Download slide

Deflection of the bimetallic beam

Table 6

Random variables

Variable	Mean	Standard deviation	Distribution
w(m)	μ_w	5 × 10⁻⁴	Normal
h(m)	μ_h	5 × 10⁻⁴	Normal
L(m)	1 × 10⁻¹	1 × 10⁻³	Normal
E_C(Pa)	1.37 × 10¹¹	1.37 × 10⁷	Lognormal
E_I(Pa)	1.30 × 10¹¹	1.3 × 10⁷	Lognormal
$T (^{\circ} C)$	130	13	Lognormal

Variable	Mean	Standard deviation	Distribution
w(m)	μ_w	5 × 10⁻⁴	Normal
h(m)	μ_h	5 × 10⁻⁴	Normal
L(m)	1 × 10⁻¹	1 × 10⁻³	Normal
E_C(Pa)	1.37 × 10¹¹	1.37 × 10⁷	Lognormal
E_I(Pa)	1.30 × 10¹¹	1.3 × 10⁷	Lognormal
$T (^{\circ} C)$	130	13	Lognormal

The failure mode is that the deflection exceeds δ = 8 × 10⁻³. The limit-state function is given by

\begin{matrix} g (X, s, t) = δ - Δ (d, E_{C}, E_{I}, Δ T) \end{matrix}

(48)

where Δ(d, E_C, E_I, ΔT) is solved by the finite element method, and matlab partial differential equation toolbox is used to solve this example.

The objective is to minimize the weight of this beam. The RBD model is defined by

\begin{matrix} {\begin{matrix} min_{μ_{w} μ_{h}} f = μ_{w} μ_{h} \\ s . t . \Pr {g (X, s, t) = δ - Δ (d, E_{C}, E_{I})} \geq Φ (β) \\ 8 \times 10^{- 4} \leq μ_{w} \leq 2 \times 10^{- 3} \\ 1 \times 10^{- 4} \leq μ_{h} \leq 1 \times 10^{- 2} \end{matrix} \end{matrix}

(49)

All the results are listed in Table 7. The proposed method is the most accurate method among the three methods as the error is only 3.04% as it uses second-order saddlepoint to find the optimal design results. It is the most efficient method as well with the evidence that the number of function calls is only 332 since it uses the sequential loop to improve the computation efficiency.

Table 7

Results of Example 4

Method	New method	SORM/DL	FORM/DL
Obj. ( × 10⁻⁶m²)	1.9197	1.9181	1.9312
μ (10⁻⁴)	(8, 2.3996)	(8, 2.3975)	(8, 2.4319)
$p_{f}^{MCS}$	0.0013	0.0014	0.0012
Error (%)	3.04	5.71	10.43
N_calls	332	3645	573

5 Conclusions

This study develops a new sequential RBD with the envelope method for time- and space-dependent reliability. The challenge in this work is to search for the equivalent most probable point (MPP), which can be found by iterating the MPP search and updating the equivalent reliability index. When limit-state functions have input with random variables and the temporal/spatial domain Ω, a single-loop MPP search is used which is much more efficient than the sequential loop approach. Once the equivalent MPP is available, the time- and space-dependent problem is transformed into a static counterpart and the second-order saddlepoint approximation is used to estimate the reliability with higher accuracy. The equivalent MPP assures that the overall optimization is performed sequentially in cycles of deterministic optimization and reliability analysis. The proposed strategy has been proven to be effective in four examples.

The proposed method is more accurate than the first-order methods since it uses the second-order saddlepoint approximation to estimate the reliability. The new method, however, shares the same limitations as the MPP-based methods. Its accuracy will deteriorate if multiple MPPs exist. The other limitation is that the method cannot handle the case where the MPP occurs on boundaries of the time and space domain. In this case, the accuracy of reliability prediction will also deteriorate. The future work is to use a system reliability method to accommodate multiple MPPs and derive a new algorithm for an MPP on the boundaries of the time and space domain.

Acknowledgment

The support from National Science Foundation under grant under Grant No. 1923799 is acknowledged.

Conflict of Interest

There are no conflicts of interest.

Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request.

References

1.

Frangopol

,

D. M.

, and

Maute

,

K.

,

2004

, “Reliability-Based Optimization of Civil and Aerospace Structural Systems,”

Engineering Design Reliability Handbook

,

CRC Press

,

Boca Raton, FL

, pp.

559

–

590

.

2.

Aoues

,

Y.

, and

Chateauneuf

,

A.

,

2010

, “

Benchmark Study of Numerical Methods for Reliability-Based Design Optimization

,”

Struct. Multidiscipl. Optim.

,

41

(

2

), pp.

277

–

294

.

Google Scholar

Crossref

3.

Yu

,

H.

,

Khan

,

M.

,

Wu

,

H.

,

Zhang

,

C.

,

Du

,

X.

,

Chen

,

R.

,

Fang

,

X.

,

Long

,

J.

, and

Sawchuk

,

A. P.

,

2022

, “

Inlet and Outlet Boundary Conditions and Uncertainty Quantification in Volumetric Lattice Boltzmann Method for Image-Based Computational Hemodynamics

,”

Fluids

,

7

(

1

), p.

30

.

Google Scholar

Crossref

4.

Wang

,

Z.

,

Huang

,

H.-Z.

, and

Du

,

X.

,

2010

, “

Optimal Design Accounting for Reliability, Maintenance, and Warranty

,”

ASME J. Mech. Des.

,

132

(

1

), p.

011007

.

Google Scholar

Crossref

5.

Du

,

X.

, and

Chen

,

W.

,

2002

, “

Sequential Optimization and Reliability Assessment Method for Efficient Probabilistic Design

,”

ASME 2002 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference

,

Montreal, Quebec, Canada

,

Sept. 29–Oct. 2

.

Google Scholar

Crossref

6.

Du

,

X.

,

Guo

,

J.

, and

Beeram

,

H.

,

2008

, “

Sequential Optimization and Reliability Assessment for Multidisciplinary Systems Design

,”

Struct. Multidiscipl. Optim.

,

35

(

2

), pp.

117

–

130

.

Google Scholar

Crossref

7.

Dering

,

M. L.

, and

Tucker

,

C. S.

,

2017

, “

A Convolutional Neural Network Model for Predicting a Product’s Function, Given Its Form

,”

ASME J. Mech. Des.

,

139

(

11

), p.

111408

.

Google Scholar

Crossref

8.

Yang

,

R.

, and

Gu

,

L.

,

2004

, “

Experience With Approximate Reliability-Based Optimization Methods

,”

Struct. Multidiscipl. Optim.

,

26

(

1

), pp.

152

–

159

.

Google Scholar

Crossref

9.

Liang

,

J.

,

Mourelatos

,

Z. P.

, and

Nikolaidis

,

E.

,

2007

, “

A Single-Loop Approach for System Reliability-Based Design Optimization

,”

ASME. J. Mech. Des.

,

129

(

12

), pp.

1215

–

1224

.

Google Scholar

Crossref

10.

Yin

,

J.

, and

Du

,

X.

,

2021

, “

A Safety Factor Method for Reliability-Based Component Design

,”

ASME J. Mech. Des.

,

143

(

9

), p.

091705

.

Google Scholar

Crossref

11.

Kharmanda

,

G.

, and

Olhoff

,

N.

,

2007

, “

Extension of Optimum Safety Factor Method to Nonlinear Reliability-Based Design Optimization

,”

Struct. Multidiscipl. Optim.

,

34

(

5

), pp.

367

–

380

.

Google Scholar

Crossref

12.

Tu

,

J.

,

Choi

,

K. K.

, and

Park

,

Y. H.

,

1999

, “

A New Study on Reliability-Based Design Optimization

,”

ASME. J. Mech. Des.

,

121

(

4

), pp.

557

–

564

.

Google Scholar

Crossref

13.

Wu

,

Y.-T.

,

Shin

,

Y.

,

Sues

,

R.

, and

Cesare

,

M.

,

2001

, “

Safety-Factor Based Approach for Probability-Based Design Optimization

,”

19th AIAA Applied Aerodynamics Conference

,

Anaheim, CA

,

June 11–14

, p.

1522

.

Google Scholar

Crossref

14.

Hu

,

Z.

, and

Du

,

X.

,

2016

, “

Reliability-Based Design Optimization Under Stationary Stochastic Process Loads

,”

Eng. Optim.

,

48

(

8

), pp.

1296

–

1312

.

Google Scholar

Crossref

15.

Wu

,

H.

,

Hu

,

Z.

, and

Du

,

X.

,

2021

, “

Time-Dependent System Reliability Analysis With Second-Order Reliability Method

,”

ASME J. Mech. Des.

,

143

(

3

), p.

031101

.

Google Scholar

Crossref

16.

Wu

,

H.

, and

Du

,

X.

,

2020

, “

Time-Dependent System Reliability Analysis With Second Order Reliability Method

,”

ASME 2020 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference

,

Virtual, Online

,

Aug. 17–19

, p.

V11BT11A045

.

Google Scholar

Crossref

17.

Wang

,

Z.

, and

Wang

,

P.

,

2012

, “

A Nested Extreme Response Surface Approach for Time-Dependent Reliability-Based Design Optimization

,”

ASME J. Mech. Des.

,

134

(

12

), p.

121007

.

Google Scholar

Crossref

18.

Fang

,

T.

,

Jiang

,

C.

,

Huang

,

Z.

,

Wei

,

X.

, and

Han

,

X.

,

2018

, “

Time-Variant Reliability-Based Design Optimization Using an Equivalent Most Probable Point

,”

IEEE Trans. Reliab.

,

99

(

1

), pp.

1

–

12

.

Google Scholar

Crossref

19.

Wang

,

Z.

,

Cheng

,

X.

, and

Liu

,

J.

,

2018

, “

Time-Dependent Concurrent Reliability-Based Design Optimization Integrating Experiment-Based Model Validation

,”

Struct. Multidiscipl. Optim.

,

57

(

4

), pp.

1523

–

1531

.

Google Scholar

Crossref

20.

Li

,

M.

,

Bai

,

G.

, and

Wang

,

Z.

,

2018

, “

Time-Variant Reliability-Based Design Optimization Using Sequential Kriging Modeling

,”

Struct. Multidiscipl. Optim.

,

58

(

3

), pp.

1051

–

1065

.

Google Scholar

Crossref

21.

Hu

,

Z.

, and

Mahadevan

,

S.

,

2016

, “

A Single-Loop Kriging Surrogate Modeling for Time-Dependent Reliability Analysis

,”

ASME J. Mech. Des.

,

138

(

6

), p.

061406

.

Google Scholar

Crossref

22.

Shi

,

Y.

,

Lu

,

Z.

,

Zhang

,

K.

, and

Wei

,

Y.

,

2017

, “

Reliability Analysis for Structures With Multiple Temporal and Spatial Parameters Based on the Effective First-Crossing Point

,”

ASME J. Mech. Des.

,

139

(

12

), p.

121403

.

Google Scholar

Crossref

23.

Shi

,

Y.

, and

Lu

,

Z.

,

2019

, “

Dynamic Reliability Analysis for Structure With Temporal and Spatial Multi-Parameter

,”

Proc. Inst. Mech. Eng. Part O J. Risk Reliab.

,

233

(

6

), pp.

1002

–

1013

.

24.

Wei

,

X.

, and

Du

,

X.

,

2019

, “

Robustness Metric for Robust Design Optimization Under Time- and Space-Dependent Uncertainty Through Metamodeling

,”

ASME J. Mech. Des.

,

142

(

3

), p.

031110

.

Google Scholar

Crossref

25.

Wei

,

X.

, and

Du

,

X.

,

2019

, “

Uncertainty Analysis for Time- and Space-Dependent Responses With Random Variables

,”

ASME J. Mech. Des.

,

141

(

2

), p.

021402

.

Google Scholar

Crossref

26.

Yu

,

S.

, and

Wang

,

Z.

,

2019

, “

A General Decoupling Approach for Time- and Space-Variant System Reliability-Based Design Optimization

,”

Comput. Methods Appl. Mech. Eng.

,

357

, p.

112608

.

Google Scholar

Crossref

27.

Wu

,

H.

, and

Du

,

X.

,

2022

, “

Envelope Method for Time- and Space-Dependent Reliability Prediction

,”

ASCE-ASME J. Risk Uncertain. Eng. Syst. B: Mech. Eng.

,

8

(

4

), p.

041201

.

Google Scholar

Crossref

28.

Wu

,

H.

,

2022

, “

Probabilistic Design and Reliability Analysis with Kriging and Envelope Method

,” Ph.D. dissertation, Purdue University, West Lafayette, IN.

29.

Wu

,

H.

, and

Xiaoping

,

D.

,

2022

, “

Envelope Method for Time- and Space-Dependent Reliability-Based Design

,”

Proceedings of the ASME 2022 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference

,

St. Louis, MO

,

Aug. 14–17

.

Google Scholar

Crossref

30.

Hu

,

Z.

, and

Du

,

X.

,

2018

, “

Saddlepoint Approximation Reliability Method for Quadratic Functions in Normal Variables

,”

Struct. Saf.

,

71

, pp.

24

–

32

.

Google Scholar

Crossref

31.

Wu

,

H.

, and

Du

,

X.

,

2020

, “

System Reliability Analysis With Second-Order Saddlepoint Approximation

,”

ASCE-ASME J. Risk Uncertain. Eng. Syst. B: Mech. Eng.

,

6

(

4

), p.

041001

.

Google Scholar

Crossref

32.

Yu

,

H.

,

Khan

,

M.

,

Wu

,

H.

,

Du

,

X.

,

Chen

,

R.

,

Rollins

,

D. M.

,

Fang

,

X.

,

Long

,

J.

,

Xu

,

C.

, and

Sawchuk

,

A. P.

,

2022

, “

A New Noninvasive and Patient-Specific Hemodynamic Index for the Severity of Renal Stenosis and Outcome of Interventional Treatment

,”

Int. J. Numer. Methods Biomed. Eng.

,

7

, p.

e3611

.

Google Scholar

Crossref

33.

Wu

,

H.

, and

Du

,

X.

,

2019

, “

System Reliability Analysis With Second Order Saddlepoint Approximation

,”

ASME 2019 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference

,

Anaheim, CA

,

Aug. 18–21

, p.

V02BT03A039

.

Google Scholar

Crossref

34.

Hu

,

Z.

, and

Du

,

X.

,

2019

, “

Efficient Reliability-Based Design With Second Order Approximations

,”

Eng. Optim.

,

51

(

1

), pp.

101

–

119

.

Google Scholar

Crossref

35.

Jones

,

D. R.

,

Schonlau

,

M.

, and

Welch

,

W. J.

,

1998

, “

Efficient Global Optimization of Expensive Black-Box Functions

,”

J. Global Optim.

,

13

(

4

), pp.

455

–

492

.

Google Scholar

Crossref

36.

Wu

,

H.

,

Zhu

,

Z.

, and

Du

,

X.

,

2020

, “

System Reliability Analysis With Autocorrelated Kriging Predictions

,”

ASME J. Mech. Des.

,

142

(

10

), p.

101702

.

Google Scholar

Crossref

37.

Li

,

C.-C.

, and

Der Kiureghian

,

A.

,

1993

, “

Optimal Discretization of Random Fields

,”

J. Eng. Mech.

,

119

(

6

), pp.

1136

–

1154

.

38.

Mourelatos

,

Z. P.

,

Majcher

,

M.

,

Pandey

,

V.

, and

Baseski

,

I.

,

2015

, “

Time-Dependent Reliability Analysis Using the Total Probability Theorem

,”

ASME J. Mech. Des.

,

137

(

3

), p.

031405

.

Google Scholar

Crossref

Time- and Space-Dependent Reliability-Based Design With Envelope Method

Abstract

1 Introduction

2 Problem Statement

2.1 Problem Statement.

2.2 Review of Sequential Reliability-Based Design.

3 Sequential Reliability-Based Design With the Envelope Method

3.1 Overview.

3.2 Time- and Space-Dependent Reliability Analysis.

3.2.1 Search for the Equivalent Most Probable Point Using FORM.

3.2.2 The Envelope Method.

3.2.3 Find the Optimal $\tilde{z}$ ⁠.

3.3 Extension to Case 2.

3.4 Procedure of Sequential Reliability-Based Design With Envelope Method.

4 Examples

4.1 Example 1: A Mathematical Problem.

4.2 Example 2: Three Bar System.

4.3 Example 3: A Cantilever Beam.

4.4 Example 4: Thermal Deflection of a Bimetallic Beam.

5 Conclusions

Acknowledgment

Conflict of Interest

Data Availability Statement

References

Get Email Alerts

Cited By

Journals

Conference Proceedings

eBooks

Resources

Opportunities

Time- and Space-Dependent Reliability-Based Design With Envelope Method

Abstract

1 Introduction

2 Problem Statement

2.1 Problem Statement.

2.2 Review of Sequential Reliability-Based Design.

3 Sequential Reliability-Based Design With the Envelope Method

3.1 Overview.

3.2 Time- and Space-Dependent Reliability Analysis.

3.2.1 Search for the Equivalent Most Probable Point Using FORM.

3.2.2 The Envelope Method.

3.2.3 Find the Optimal z~⁠.

3.3 Extension to Case 2.

3.4 Procedure of Sequential Reliability-Based Design With Envelope Method.

4 Examples

4.1 Example 1: A Mathematical Problem.

4.2 Example 2: Three Bar System.

4.3 Example 3: A Cantilever Beam.

4.4 Example 4: Thermal Deflection of a Bimetallic Beam.

5 Conclusions

Acknowledgment

Conflict of Interest

Data Availability Statement

References

Get Email Alerts

Cited By

Related Articles

Related Proceedings Papers

Related Chapters

Journals

Conference Proceedings

eBooks

Resources

Opportunities

This Feature Is Available To Subscribers Only

3.2.3 Find the Optimal $\tilde{z}$ ⁠.