Structural Safety 82 (2020) 101875
Contents lists available at ScienceDirect
Structural Safety journal homepage: www.elsevier.com/locate/strusafe
A combined Importance Sampling and active learning Kriging reliability method for small failure probability with random and correlated interval variables Xiao-Xiao Liua, Isaac Elishakoffb, a b
T
⁎
State Key Laboratory for Strength and Vibration of Mechanical Structures, School of Aerospace Engineering, Xi’an Jiaotong University, Xi’an 710049, China Department of Ocean and Mechanical Engineering, Florida Atlantic University, Boca Raton, FL 33431-0991, USA
ARTICLE INFO
ABSTRACT
Keywords: More general hybrid reliability analysis Random variable Correlated interval variable Active learning Kriging Importance Sampling
The existing hybrid reliability analysis (HRA) method (Yang et al., 2015; Zhang et al., 2015; Yang et al., 2015) is found not suitable for estimating small failure probabilities. Meanwhile, the previous ALK-HRA algorithm (ALKHRA: an active learning HRA method combining Kriging and Monte Carlo simulation) reduces its numerical efficiency when number of uncertain variables increases. Furthermore, the ALK-HRA approach with both random and interval/ellipsoid variables cannot deal with complex “multi-source uncertainty” problems. In order to overcome these issues, therefore, the following strategies is proposed: 1) First, a more general HRA (MGHRA) method with both random and parallelepiped convex variables is developed. Within the MGHRA method, the parallelepiped convex model is employed to describe independent and correlated interval variables in a unified framework. 2) Sequentially, we propose an original and implementable approach called ALK-MGHRA-IS for active learning MGHRA method and Kriging-based Importance Sampling. The MGHRA method, which is capable of handling the complicated “multi-source uncertainty” problems, associates the Kriging metamodel, and its advantageous stochastic property with Importance Sampling to accurately evaluate bounds of small failure probabilities with respect to interval variables. Actually, the calculated failure probability is still a random variable when the approximations of the proposed method are employed. The proposed method enables the correction of the FORM-UUA approximation with only a few function computations. To further improve the efficiency of the proposed ALK-MGHRA-IS, an optimization method based on Karush–Kuhn–Tucker conditions (KKT) is introduce to relieve the burden of searching the extreme values. Four numerical examples are investigated to demonstrate the efficiency and accuracy of the proposed method.
1. Introduction Various uncertainties with reference to structure dimensions, loading, material properties, incomplete information, manufacturing imprecision, etc., always exist in actual engineering problems. They dramatically affect the structural safety and performance. Reliability analysis is the technique which provides the assessment of the safety degree of a system by taking into such uncertainties account. In numerous engineering applications, a frequently encountered circumstance is that [4–8]: some of the uncertainties can be described using specific probabilistic models and others need to be treated as convex models due to lacking of imperfect knowledge. Since the two types of uncertainties mentioned above often coexist, hybrid reliability analysis (HRA) with both random and convex variables has been paid more and more attention in recent years [1,9–14]. A series of the ⁎
corresponding HRA methods have been proposed. For example, Du and his co-authors [15–17] proposed a unified uncertainty analysis (UUA) method based on the first-order reliability method (FORM), which was named as FORM-UUA. Compared with interval arithmetic [18,19], FORM-UUA is capable of applying to black-box performance functions. To improve the robustness of FORM-UUA, Yao et al. [20] developed an enhanced FORM-UUA with single-level optimization. Jiang et al. [21] constructed an equivalent single-layer optimization algorithm, which can only calculate the maximum failure probability. Xiao et al. [22] presented a new method termed as MVFOSPA-UUA to estimate the lower and upper bounds of the system probabilities of failure. Luo et al. [12] proposed the mathematical programming method (MPM) and the single-loop iterative (SLI) method to calculate the minimum reliability index, which was used as constraint in reliability-based design optimization (RBDO) with both random and convex variables. When the
Corresponding author. E-mail address:
[email protected] (I. Elishakoff).
https://doi.org/10.1016/j.strusafe.2019.101875 Received 5 April 2018; Received in revised form 8 July 2019; Accepted 9 July 2019 0167-4730/ Published by Elsevier Ltd.
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
performance functions are highly nonlinear or have multiple design points, these algorithms are very inaccurate. To solve these problems, Yang et al. [1,3] proposed an active learning Kriging model for HAR with both random and convex variables. It was named ALK-HAR. However, ALK-HAR will reduce its numerical efficiency if the number of uncertain variables increases as revealed in [1–3]. Meanwhile, ALKHAR is found not suitable for estimating small failure probabilities in engineering problems. Rare event estimation requires a large number of samples, and so a great amount of searches for the extreme values need to be implemented based on global optimization algorithm. Even though this process has no need of calling a specific performance function, it becomes more difficult and time-consuming. To improve the computational efficiency and solve small probability events, Importance Sampling (IS) [23,24] is the most widely used method in reliability analysis. IS requires introducing the Importance Sampling density function, which is usually taken in the standard space as an n-Gaussian probability density function (PDF) centered on the most probable failure point (MPFP). When IS simulates a population centered on the MPFP, a large number of the simulated samples are located in the failure region where the probability density is relatively larger. When an isolated MPFP exists in the performance function, the estimated failure probability is unbiased by using IS. In addition, the number of samples simulated by IS is significantly less than that simulated by Monte Carlo Simulation (MCS). Although IS can reduce immensely the number of calls to performance functions, it is still not applicable with time-consuming numerical models. To further solve the efficiency problems, the probability reliability analysis (PRA) methods combining Kriging model with IS have been proposed by some researches [25–27]. Echard et al. [28] proposed an efficient reliability method termed as AK-IS for small failure probabilities with time-demanding models. Dubourg and Sudret [29] proposed a so-called metaAK-IS method for the reliability sensitivity analysis. Cadini et al. [30] developed an enhanced method called as metaAK-IS2 for sampling multiple failure regions of low probability. Zhao et al. [31] presented an improvement of the adaptive IS method based on Kriging metamodel and active learning mechanism. Meanwhile, Zhao [32] developed an efficient interval IS method to calculate the bounds of the limit state probability. In these existing literatures, some of convex variables are described with interval models and others are denoted as ellipsoidal models. Nevertheless, as demonstrated in [33], an interval model or an ellipsoidal model cannot deal with complex “multi-source uncertainty” problems. For this purpose, Jiang and his coworkers proposed a more general convex model called multidimensional parallelepiped (MP) model, which can consider independent and correlative interval variables. Without doubt, the existing AK-IS [25–32] and ALK-HRA [1–3] methods cannot deal with complex “multi-source uncertainty” problems as only random and interval/ellipsoid variables are considered in these literatures. Therefore, the following schemes is proposed to overcome these issues. First, a more general HRA (MGHRA) method with both random and parallelepiped convex variables, which can describe independent and correlated interval variables in a unified framework, is constructed. Sequentially, this paper can propose an original and implementable approach called ALK-MGHRA-IS for active learning MGHRA method and Kriging-based Importance Sampling. The MGHRA method, which can handle the complex “multi-source uncertainty” problems, associates the Kriging model, and the advantageous stochastic property involved in the MGHRA with Importance Sampling to evaluate bounds of small failure probabilities with respect to interval variables. Actually, the failure probability calculated through the proposed method is still a random variable when the approximations are used. Therefore, the lower and upper bounds of the failure probability of an MGHRA problem with both random and dependent interval variables are defined as
P Uf = P {min G (X , Y ) < 0|Y Y
}
P fL = P {max G (X , Y ) < 0|Y Y
}
(2)
where P Uf and P fL represent the upper and lower bounds, respectively; min G (·) and max G (·) are the extreme values of the performance Y
Y
function with respect to Y. The proposed method enables the correction of the FORM-UUA approximation with only a few function computations. To make the subsequent MGHRA more convenient, a matrix transformation [34,35] is introduced to transform the correlated intervals into independent ones. Furthermore, an optimization methodbased Karush–Kuhn–Tucker (KKT) conditions is recommended to relieve the burden of searching the extreme values. The remainders of this paper are organized as follows. Various methods for assessing the small failure probabilities with both random and parallelepiped convex variables is presented in Section 2. FORMUUA method based on KKT conditions is introduced in Section 3 which will be the cornerstone of the proposed method. The methodology of ALK-MGHRA-IS method is elaborated in Section 4. Four examples are used to demonstrate the efficiency and accuracy of the proposed method in Section 5. Conclusions are given in Section 6. 2. Assessment of small failure probability with both random and parallelepiped convex variables Actually, aleatory and epistemic uncertainties often in many engineering applications. Generally, aleatory uncertainty is described through precise probabilistic model. Epistemic uncertainty, which stems from the lack of knowledge or incomplete information, may has an important impact on the reliability analysis. Therefore, quantifying the structural reliability with both aleatory and epistemic variables needs to be attracted much in such circumstance. Epistemic uncertainty can be described using various approaches, such as convex model (or interval model), Bayesian approaches, evidence theory and fuzzy set theory. These types of the uncertainty description have been gradually improved. The interval model is a good tool to describe the epistemic uncertainty, but it is very hard to get credible reference of the uncertainty just through the defined boundaries. Inside characteristics of the uncertainties will be completely neglected so that the threshold is roughly estimated. Subsequently, the prediction ranges of the bound values of the failure probability will be overestimated. In fact, although lack of information in the early design stage, frequently, there are some experience from previous work more than thresholds. In this stand of point, the epistemic parameters, which is described using parallelepiped model, should be a better choice since it is able to take full use of the limited experience or data to find out the prediction ranges. Parallelepiped model as description of epistemic uncertainty, is capable of dealing with the coexistence of independent and relevant uncertain variables. As for the choice of the non-probabilistic model, it is recommended that a reasonable model for different uncertainty issues can be chosen by considering the amount of information in early stage. Meanwhile, the choice of the epistemic uncertainty model depends on the maturity of computation tool, mathematical methods and importance of design variables. Therefore, a MGHRA problem, which considers the structural input parameters as both random and correlated interval variables (parallelepiped variables) for the first time, will be addressed for evaluating the small failure probability. First of all, some approaches such as Monte Carlo Simulation (MCS), Important Sampling (IS), and FORM-UUA method based on KKT conditions (FORM-UUA-KKT) are introduced to deal with the small failure probability, which involves the random and correlated interval variables. In order to further improve the computational efficiency and the accuracy, however, the combination of Active Learning Kriging (ALK) and IS can be proposed to deal with the MGHRA problem. The proposed method has be called ALK-MGHRA-IS in the previous section because it is an ALK method combining with IS to deal with the MGHRA problem. When only the IS method is used to
(1) 2
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
solve the MGHRA, for comparison purpose, this type of approach is called MGHRA-IS.
stressed that the failure probability are bounds with respect to interval variables Y . Actually, these calculated failure probabilities are still random variables when the approximations of the proposed method are employed. The lower and upper bounds of the failure probability of the MGHRA problem have been expressed in Eqs. (1) and (2). When the original interval space are transformed into the standard interval space, the upper bound P Uf and the lower bound P Uf will be rewritten as
2.1. Monte Carlo simulation (MCS) with both random and correlated interval variables 2.1.1. The multidimensional parallelepiped (MP) model To solve the complex “multi-source uncertainty” problems, Jiang’s team [33,35,36] has developed the MP model to deal with the independent and correlated interval variables in a unified framework. Generally, the correlation between interval variables should be obtained based on large amounts of experiment samples or simulation samples. In the MP model, a MP is utilized to encircle the given samples for constructing the uncertain domain. Therefore, the variation ranges and the correlation of the intervals can be well exhibited by the size and shape of the uncertain domain. When a two-dimensional problem is taken as an example, the MP model degenerates into a parallelogram. In this model, the marginal intervals of interval variables Y1 and Y2 are expressed as Y11 and Y21, respectively. The two marginal intervals are represented by
Y1I = [Y1L,Y1U ],Y1C = (Y1L + Y1U )/2, Y1W = (Y1U
Y1L )/2
(3)
Y2I = [Y2L,Y2U ],Y2C = (Y2L + Y2U )/2, Y2W = (Y2U
Y2L )/2
(4)
Yi =
YiW | (i , j )|
n j=1
n ik k
+ YiC , i = 1, 2, …, n
Y 1I
a
O
H (Y C2 , Y C1 )
D C
IFU (X ) f (X )dX
(8)
1 min G (X , Y ) Y
0
0 min G (X , Y ) > 0
(9)
Y
In like manner, Eq. (2) can be written as
P fL = where
IFL (X ) f (X )dX
IFL (·)
IFL (X ) =
(10)
can be expressed as
1 max G (X , Y ) Y
0
0 max G (X , Y ) > 0 Y
(11)
To obtain globally optimal solutions, the failure probabilities P Uf and P fL can be obtained by combining MCS with the modified DIRECT algorithm [37]. 2.2. Important sampling (IS) with both random and correlated interval variables Since MCS sometimes is found not suitable for dealing with timedemanding models, IS as an alternative approach is introduced to carry out the MGHRA problem. IS requires an importance sampling density (ISD) function, which is usually taken in the standard space as an nGaussian PDF centered on the MPFP. Hence, the iso-probabilistic transformation T is used to convert random variables X to the standard normal spaceU = T (X ) . Sequentially, Eq. (8) can be rewritten as
P Uf =
IFU (U ) f (U )dU
(12)
where f (U ) is the joint Gaussian PDF of the standard random variables U and IFU (U ) can be rewritten as
IFU (U ) =
1 min g (U , Y ) Y
0
0 min g (U , Y ) > 0
(13)
Y
A
b
B
(7)
IFU (X ) =
2.1.2. Monte Carlo simulation (MCS) When both random and correlated interval variables coexist, the performance function is denoted as G (X , Y ) , where Y = [Y1, Y2, …YnY ]T represents the vector of correlated interval variables, described in Section 2.1.1, and X = [X1 , X2 , …XnX ]T represents the vector of random variables. Due to the presence of the correlative intervals, the failure probability is an uncertain-but-bounded variable. However, it should be
YW 2
cube}
where X = [X1 , X2 , …XnX represents the vector of random variables. IFU (·) is the failure indictor function associated with the minimum response and it can be expressed as
where the lower and upper bounds of Yi are defined as Yi L and Yi U , respectively.
¦ rd ¸
P fL = P {max G (X , ) < 0|
]T
(5)
k=1
(6)
P Uf =
Y1I
Y1U
cube }
When the MGHRA problem is solved by combining with MCS, Eq. (1) can be written as
where and are the lower bound and the upper bound of respectively; Y2L and Y2U are the lower bound and the upper bound of Y2l , respectively; Y1C and Y2C are the midpoints of the two marginal intervals; Y1W and Y2W are interval radii of the two marginal intervals. If Y1 and Y2 are independent of each other, the uncertain domain of Y1 and Y2 will manifest as a rectangular domain [Y1L, Y1U ] × [Y2L, Y2U ], denoted as rd . If Y1 and Y2 are correlated, the uncertain domain will form a parallelogram domain . Fig. 1 shows that the center of coincides with that of rd . The definition of correlation coefficient between Y1 and Y2 is elaborated in Appendix A. For a general n-dimensional problem, the definition of the correlation coefficient between any two interval variables Yi and Yj is elaborated in Appendix B. Sequentially, the variable Y1 involved in n-dimensional interval variables, is expressed as:
Y1L
P Uf = P {min G (X , ) < 0|
where g (U , Y ) represents the performance function with respect to Y within U space. When the original interval space are converted to the standard interval space, Eq. (13) can be rewritten as
YW 1 ¦ ¸
1 min g (U , ) IFU (U ) =
0 min g (U , ) > 0 cube
Y 2I
0
cube
(14)
where g (U , ) represents the performance function with respect to within U space. By introducing the ISD function h (U ) , Eq. (8) can be expressed as
Fig. 1. The MP model with the two-dimensional problem. 3
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff f (U )
PUf IS=
IFUIS (U ) h (U ) h (U )dU
(
f (U )
= E IFUIS (U ) h (U )
)
Failure domain
(15)
where IFUIS (U ) equals to Eq. (13) or Eq. (14). Then, the integral (15) can be estimated by the following expression.
P Uf IS
Pf
UIS
1 = NUIS
NUIS
IFUIS (U (i ) )
i=1
IS samples
(i ) )
f (U h (U
(i ) )
MPFPU
(16) MPFPL
where NUIS is the total number of samples with IS. It should be noted that the samples {U (i), i = 1, …, NUIS} , are simulated according to the ISD function h (U ) . The variance and the coefficient of variation of the upper bound estimator P f UIS are respectively given by
(
)
Var P f UIS =
(
NUIS
1 1 NUIS NUIS
IFUIS (U
f (U h (U
(i ) )
i=1
(i ) ) 2 (i ) )
(P f UIS )2
Y
min G (U ,Y )=0
O
)
P f UIS
Fig. 2. Importance Sampling-based MGHAR.
index. From Fig. 2, it is observed that L is the distance from the origin to the MPFPL and U is the distance from the origin to the MPFPU. The MPFPL and MPFPU correspond to the most probable failure points (MPFPs) on min G (U , Y ) = 0 and max G (U , Y ) = 0 , respectively, as
(18)
Likewise, Eq. (10) using IS can be estimated as
P fLIS
P f LIS =
1 NLIS
N LIS
IFLIS (U
i=1
(i ) )
f (U h (U
Y
(i ) ) (i ) )
(19)
(
)
(
)
Var P f LIS =
1 1 NLIS NLIS
Cov P f LIS =
i=1
IFLIS (U
(i ) )
f (U h (U
(i ) ) (i ) )
or Eqs. (22) and (23) are double-layer nesting optimization problems. Taking Eq. (22) as an example, the process can be decomposed into: Outer-layer optimization for PRA
(20)
L
Inner-layer optimization for interval analysis
G (U ) = min G (U , Y ) s.t. |C
m (U , Y ) =
where
and
L
1 2
U + c G (U , Y )
(29)
U G (U , Y )
c>
(22)
(30)
The search direction dk is denoted as
U
U
(28)
where the constant c should satisfy
= min U Y
(27)
e
where is the step size and can be obtained by minimizing the following merit function
= min U
s.t. max G (U , Y ) = 0
Y|
Uk + 1 = Uk + dk
U
U
Y 1
In this paper, the iHL-RF [15] is used to calculate MPFPs. When performing interval analysis, an optimization method based on KKT conditions is introduce to relieve the burden of searching the extreme values. Sequentially, the computational efficiency can be improved due to the reduction of searching times. Still using Eq. (22) as an example, the FORM-UUA-KKT algorithm can be summarized as follows: Step 1. Construct the initial MP model and map the standard δspace to the original interval space. Initialize a set of U0 , Y0 and the iteration counter k = 0 . Step 2. In iteration k + 1, the MPFPL can be obtained by
FORM-UUA as a gradient-based optimization method is widely used to solve the HRA problem. The important thing to note, however, is that the efficiency of this method depends on the number of uncertain variables and the dimensionality of the performance function. When FORM-UUA is adopted to calculate failure probabilities, Eqs. (8) and (10) can be respectively rewritten as:
Y
(26)
s.t. G (U , Y ) = 0
2.3. FORM-UUA-KKT with both random and correlated interval variables
s.t. min G (U , Y ) = 0
= min U U
(21)
when n-dimensional interval variables exhibit as dependent quantities, the limit state G (U , Y ) = 0 shows a bunch of limit state surfaces in the standard random space. It should be noted that if the chosen ISD function is inapplicable, the computational efficiency cannot be improved and the failure probability estimator will be biased. To avoid such occurrences, the MPFP-based IS method is applied to deal with the MGHRA problem. Moreover, the FORM-UUA method based on KKT conditions can be utilized to obtain an MPFP. It can be named FORMUUA-KKT. Although IS greatly reduces the number of function calls in comparison with MCS, it is still unacceptable for conducting rare-event estimation. Therefore, the use of ALK method is motivated to improve MGHRA-IS. The so-called FORM-UUA-KKT method, presented first, is the benchmark of the ALK-MGHRA-IS method.
L
(25)
h (U ) = f (U + MPP L )
2
(P f LIS )2
(24)
h (U ) = f (U + MPPU )
Var(P f LIS ) P f LIS
Y
shown in Fig. 2. Subsequently, the ISD function can be obtained by coordinate transformation, expressed as follows:
The variance and the coefficient of variation of the lower bound LIS estimator P f are respectively given by N LIS
Y
(17)
Var (P f UIS)
Cov P f UIS =
max G (U ,Y )=0
dk =
(23)
G y (Uk , Yk ) UkT
G (Uk, Yk )
G y (Uk , Yk )
2
G y (Uk, Yk )
Uk
To reduce the computational cost, one needs to find a value
represent the maximum and the minimum reliability 4
(31) such
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
the above, when (X (j ) , Y1 ) keeps away from the vicinity of G (X , Y ) = 0 , the Kriging model can rightly predict the sign of the minimal value. When the point (X (j ) , Y1 ) is located in the vicinity of G (X , Y ) = 0 , Kriging model cannot assure a right prediction for the sign of extreme values. However, in almost all cases, the probability that the extreme points concentrate on the vicinity of the limit state is small. Therefore, even if the Kriging model is difficult to guarantee the accuracy of predicting the sign in the whole space, it can still keep high accuracy for MGHRA-IS. Consequently, when the Kriging model has a prediction error for the sign, Corollary 1 can be largely true. Based on the above corollaries, one can know that the construction of a predicted Kriging model providing a right prediction for the sign of the performance function can meet the accuracy requirement of MGHRA-IS. For example, if only the Kriging model µG (X , Y ) is able to have a right prediction for the sign of G (X , Y ) , the extreme values of G (X , Y ) and µG (X , Y ) will have the same signs.
that the merit function (29) is adequately reduced. The following rule proposed by Du [15] is used to find = bh such that
h = max{b h |m (Uk + bhdk , Yk )
m (Uk, Yk ) < 0} b > 0
(32)
2 Uk G (Uk, Yk )
+ 10 . where b = 0.5 and c = Step 3. If k > 1, check KKT conditions in the following expressions
H (Y ) =
nY
H (Y j )
j=1
H (Y j ) = Y j = Y j L ,
{Y = Y j
{Y j
U
,
j
L
G (Uk + 1, Yk ) Yj
< Yj < Y j U ,
G (Uk + 1, Yk ) Yj
0
0 G (Uk + 1, Yk ) Yj
}
=0
} (33)
If satisfied, Yk + 1 = Yk and go to Step 4; otherwise, go to Step 3. Step 4. Find Yk + 1 by calculating Yk + 1 = argmin G (Uk + 1, Y ).
3.2. Revisiting Kriging model
Y
Step 5. Check convergence. If |G (Uk + 1, Yk + 1)| and 1 L = U Uk + 1 Uk k+1 ; 2 ( 1 and 2 are the margins of error), then otherwise, k = k + 1, go to (2). It should be noted that H (Y ) in Eq. (33) should be replaced by H ( Y ) when solving Eq. (23). If huge number of variables exist, reliability sensitivity index should be computed by developing efficient algorithms. The equations of sensitivity indices of four types with both random and correlated interval variables can be derived by combining FORM with dimensional reduction. The sensitivity indices are elaborated in Appendix D.
For convenience, this subsection use the vector to replace all of uncertain variables (X , Y ) and then the performance function can be written as G ( ) . Kriging model expresses the unknown function G ( ) as [26]: (34)
G( ) = F ( , ) + z( )
where F ( , ) is the polynomial regression part which represents an approximation of the response in mean and is the vector of regression coefficients. The second part z ( ) is a Gaussian random process whose statistical characteristics are denoted as (35)
E [z ( )] = 0
3. Kriging metamodel
Cov [z (c ), z (d )] =
Y
Y
Y
2 [R (
, c , d )]
(37)
where is the standard deviation of z ( ) , Cov [·] represents the covariance between two arbitrary points c and d , R ( , c , d ) is a correlation function and can be defined using several models. The Gaussian correlation function is chosen here, and it can be formulated as
This paper motivates the use of ALK method to improve the efficiency of MGHRA through combining with Important Sampling (IS). The ALK method only locally approximates in the region of interest rather than in the overall uncertain space. Therefore, one can implement the MGHRA with function evaluations as few times as possible when the IS approach is adopted. The core technology of ALK-MGHRAIS is that the predicted Kriging model only providing a right prediction for the sign of G (X , Y ) can meet the requirement of MGHRA-IS in accuracy. The predicted Kriging model is denoted as µG (X , Y ) and the following corollaries can be obtained based on the existing properties [1]. Corollary 1. If sign(µG (X , Y )) = sign(G (X , Y )) , then the sign of min G (X , Y ) is consistent with that of min µG (X , Y ) . Y
(36)
2
Var [z ( )] =
3.1. The methodology of active learning Kriging (ALK) for MGHRA
n
R ( , c , d ) = exp
i (ci i=1
di ) 2
(38)
where n is the dimension of samples, i is the ith component of , ci and di are the ith component of the points c and d respectively. Meanwhile, the predicted value µG ( ) and predicted variance G2 ( ) of Kriging model are elaborated in Appendix C. Kaymaz [38] indicated that the parameter in Eqs. (C.3)–(C.4) (see Appendix C) had a significant impact on the accuracy of the Kriging model. Therefore, the method of maximum likelihood estimation is utilized to determine . Subsequently, the unconstrained optimization problem must be carried out by the following expression [39]:
Corollary 2. If sign(µG (X , Y )) = sign(G (X , Y )) , then the sign of max G (X , Y ) is consistent with that of max µG (X , Y ) . Taking Corollary 1 as an example, the minimal point of G (X , Y ) at a sample X (j) is denoted as Y1 and the minimal point of µG (X , Y ) is denoted as Y1 . When a point (X (j ) , Y1 ) keeps away from the neighborhood of G (X , Y ) = 0 and G (X (j ) , Y1 ) < 0 , Kriging model µG (X , Y ) can rightly predict the sign of the performance function, namely, µG (X (j ) , Y1 ) < 0 . Then there must exist an infinitesimal neighbourhood = {(X (j) , Y )|lim (X (j ) , Y ) (X (j) , Y1 ) < } , such that (X (j ) , Y1)
opt
(
1
= argmin |R|k
2
)
(39)
It is noteworthy that the optimal parameter opt should be searched through the global optimization scheme and still the modified DIRECT algorithm [37] is used in this paper. As mentioned above, the Kriging model provides a predictor µG ( ) for G ( ) . Because a Gaussian random variate G ( )Ñ(µG ( ), G2 ( )) , there have been some uncertainties in this prediction process and µG ( ) is only an approximate value rather than an actual value of G ( ) . Thereby, there exists a risk that the sign of G ( ) is mistakenly predicted. In order to well perform the prediction for the sign of G ( ) , the point at which the sign of G ( ) has the maximum risk value to be mistakenly predicted needs to be identified. Then the aforesaid point should be added into DoE and the sign prediction of G ( ) can be greatly improved. To recognize this point and construct an ALK model, a learning function such as EIF [40], EFF [41]
0
keeps away from the vicinity of G (X , Y ) = 0 and G (X (j ) , Y1 ) < 0 . Therefore, Kriging model can rightly predict the sign of the performance function, namely, µG (X (j ) , Y1 ) < 0 . Subsequently, it must be deduced that µG (X (j ) , Y1 ) µG (X (j ) , Y1 ) < 0 . Moreover, when the point (X (j ) , Y1 ) keeps away from the neighborhood of G (X , Y ) = 0 and G (X (j ) , Y1 ) > 0 , there must exist a point (X (j ) , Y1 ) such that G (X (j ) , Y1 ) G (X (j) , Y1 ) > 0 which is out further away from the vicinity of G (X , Y ) = 0 . Then Kriging model µG (X , Y ) can rightly predict the sign of the performance function, namely, G (X (j ) , Y1 ) > 0 . From 5
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
and the information entropy-based H [42] is required for offering an accurate sign prediction. However, these learning functions are only suitable for PRA. Here, a criterion called Expected Risk Function (ERF) is introduced to measure the degree of the before-mentioned risk, expressed as [43]:
i = 1, …, NIS . (b) Generation of the MP convex samples is composed of three steps. Step 1 is the construction of the MP model. Step 2 is the generation of interval samples Y. Latin hypercube sampling (LHS) is used to generate sampling points uniformly covering the space defined by Eq. (B.4) (see Appendix B) and mapped into the original interval space defined in Eq. (B.2) through Eq. (B.3) and Eq. (3). Meanwhile, the lower and upper bounds of the intervals Y are defined respectively as Y L and Y U . (c) The set of these samples as candidate points is denoted as IS + LHS .
E [R ( )] =
sign(µG ( )) µG ( )
sign(µG ( ))
µG ( ) G(
)
+
G(
)
µG ( ) G(
) (40)
where sign(·) is the sign function, (·) and (·) are the cumulative distribution function (CDF) and the PDF of the standard normal distribution, respectively. ERF demonstrates the risk degree that the sign of a performance function at a sampling point is capable of being expected to be wrongly predicted through the Kriging model. If a point is found capable of maximizing ERF, the sign of G ( ) at this point will have the maximum risk to be mistakenly predicted. Thus, this point must be added into the initial DoE. The advantage of ERF is that the point enriching the initial samples can be identified and the prediction for the sign of a performance function can be enhanced. Next, the procedure description for closely integrating ALK with MGHRA-IS will be given.
4.4. Step 4. Identify the new training point among the set of candidate points The sampling point among IS + LHS which can maximize ERF is selected as a new training point which should be added into the DoE. The point is expressed as (U ( ), Y ( ) ) . 4.5. Step 5. Check the stopping condition If the maximization of ERF is sufficiently small, the sign of G (U , Y ) will have little risk to be mistakenly predicted. Skip to Stage (7). Otherwise, go to Stage (6). The stopping condition adopted in this manuscript is E [R ((U ( ), Y ( ) ))]/(|G (µX , Y¯ )| + ) 10 4 [1], where µX is the mean of random variables X , Y¯ is the nominal value of the in10 6 . terval variables Y and
4. ALK-MGHRA-IS to assess small failure probability ALK-MGHRA-IS requires the calculation of MPFP. For random variables, a population of the IS samples centred at the MPFP can be simulated using the ISD function mentioned in Eqs. (24) and (25). For the MP convex variables, the samples are chosen uniformly distribution in the uncertain space defined by Eq. (3). Sequentially, these samples as candidate points are used to construct the ALK model. It should be specially explained that FORM-UUA might need a great deal of function calls for searching the MPFP. This step is quite time consuming as it is implemented on the true performance function G ( ) . However, a precise MPFP is not required and a few iterations providing an approximated MPFP is enough to meet the requirement of IS in accuracy. The procedure of ALK-MGHRA-IS is described as follows, and its flowchart is shown in Fig. 3.
4.6. Step 6. Update the DoE and construct a new Kriging model If the stopping condition mentioned in Stage (5) cannot be satisfied, the performance function G (U , Y ) at the point (U ( ), Y ( ) ) needs to be calculated. Subsequently, the point (U ( ), Y ( ) ) is added into DoE and a new Kriging model is constructed. Go back to Stage (4). 4.7. Step 7. Calculate the minimum failure probability and the coefficient of variation Calculate the minimum failure probability P f LIS and the coefficient of variationCov (P f LIS ) , respectively defined by Eqs. (19) and (21). If Cov (P f LIS ) is less than a predefined convergence tolerance (5%), P f LIS can be used as the estimator of the minimum failure probability. Otherwise, the population NIS needs to be enlarged with more samples and ALK-MGHRA-IS goes back to Stage (3) to obtain new candidate points.
4.1. Step 1. FORM-UUA approximation (a) The random variables needs to be transformed into the standard normal space U . For the MP convex (correlated interval) variables, two steps are required. Step 1 is the creation of the MP model elaborated in Section 2.1.1. Step 2 is the transformation of the standard interval space into the original interval space defined in Eq. (B.2) (see Appendix B) via Eq. (B.3) (see Appendix B) and Eq. (11). (b) The MPFPU corresponding to max G (U , Y ) = 0 can be calculated
4.8. Step 8. Calculate the maximum failure probability and the coefficient of variation
Y
using FORM-UUA. Because an approximated MPFPU is sufficient to satisfy the precision requirement of IS, 1 and 2 involved in FORM-UUAKKT can be set to larger values.
Analogous to Stage (7), if the coefficient of variation Cov (P f UIS) defined by Eq. (18) is less than a predefined convergence tolerance, P f UIS can be used as the estimator of the maximum failure probability. Otherwise, the MPFPL corresponding to min G (U , Y ) = 0 should be Y
4.2. Step 2. Define the initial DoE
calculated using FORM-UUA and ALK-MGHRA-IS repeats Stages (1–8) until the stopping condition is satisfied again. It should be noted that the variation range of marginal intervals are sometimes relatively small in practical engineering applications. Therefore, the short distance between MPFPU and MPFPL can be observed. In most cases, a population of the IS samples centred at one of the two MPFPs can be used to estimate the lower and upper of failure probabilities. ALK-MGHRA-IS utilizes this character to calculate one of the two MPFPs and the minimum failure probability. Subsequently, one tentatively employs the same IS samples to estimate the maximum failure probability. Through this treatment, the number of performance function computations can be drastically reduced because FORM-UUAKKT might execute only once.
The initial DoE is formed with the samples at which the function calls have already been performed during FORM-UUA approximation of ALK-MGHRA-IS Step 1. These samples are named FORM-UUA DoE and are enough to start the procedure of constructing an initial Kriging model. 4.3. Step 3. Generate a large amount of the samples as candidate points (a) For random variables, a population of the IS samples centred at the MPFPU can be simulated using the ISD function defined by Eq. (24). The samples of this population are noted U (i) = {U1 (i), …, Un (i)} for 6
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
Transform input variables (random and correlated interval variables) into standard normal space MPFPU obtained by FORM-UUA
Generate candidate points based on IS and LHS: random points obtained by IS; MP points obtained by LHS
Define an initial FORM-UUA DoE
Build an initial Kriging metamodel
Execute active learning process and identify the new training point
MPFPL obtained by FORM-UUA Update the previous DoE and construct the new Kriging model
No
Whether the stopping criterion is satisfied Yes
No
Calculate the minimum failure probability and coefficient of variation
Calculate the maximum failure probability and coefficient of variation
Whether the coefficient of variation satisfies the convergence tolerance
Whether the coefficient of variation satisfies the convergence tolerance
Yes
No
Yes End
End
Fig. 3. The flowchart of the proposed ALK-MGHAR-IS.
5. Numerical examples and discussions
et al. [43] is considered. The performance function is defined as
In this section, four examples are investigated to illustrate the efficiency and accuracy of the proposed method. The first mathematical problem is used to demonstrate the difference of the precise MPFP and approximate MPFP, and also can be employed to illustrate the superiority of ALK-MGHRA-IS compared with other approaches. The second example is employed to show that ALK-MGHRA-IS is able to deal with the bounds of small failure probabilities with respect to correlated interval variables. The third example is used to verify the efficiency and accuracy of the proposed method with respect to a more complex application problem. The fourth example as a practical engineering case is used to demonstrate the performance of the advocated method to the black-box performance function. In these examples, the IS method with a precise MPFP is performed to individually examine the results obtained through other numerical methods. Moreover, the HRA problem with both random and independent interval variables is compared with ALK-MGHRA-IS, in order to evaluate the role of the range of the interval variables on the precision of failure probabilities. Finally, the efficiency and accuracy of ALK-MGHRA-IS are checked through comparing with FORM-UUA-KKT.
G (X , Y ) =
x1 y1
N1
+
x2 y2
N2
1
(41)
where X = (x1, x2 is the vector of standard Gaussian random variables; Y = (y1 , y2 )T denotes the correlated interval variables; the marginal intervals are all [4, 5]. Assume that the correlation coefficient between the two interval variables that defined by Eq. (A.1) (see Appendix A) is 12= 21 = 0.2 . The correlation matrix can be expressed as:
)T
=
1 0.2 0.2 1
(42)
Therefore, the shape matrix C can be obtained:
C=
0.4167 0.0833 0.0833 0.4167
(43)
Then, the uncertain domain of the MP model can be constructed using Eqs. (B.2) and (43). Thereby, the MP model can be converted to the standard δ-space. When N1 2 and/or N2 2 , Eq. (41) is a highly nonlinear function in the random space and it includes an isolated MPFP. The realization of G (X , Y ) = 0 for values of Y is shown in Fig. 4, when N1 = 1 and N2 = 3 of Eq. (41) is considered in this section. ALK-MGHRA-IS is suggested to calculate the lower and upper bounds of failure probabilities. The existing investigations [1–3] indicated that when FORM-UUA performed a large amount of function
5.1. Example 1: A mathematical problem A mathematical problem from Cimellaro and Reinhorn [42] and Liu 7
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
calls to the performance function and six iterations. In order to test the performance of the proposed method, ALK-MGHRA, IS2 and IS1 are implemented with 100 different populations of NIS = 104 samples centred on the MPFP identified by FORM-UUA-KKT. The median and coefficient of variation of different methods are listed in Table 1. It can be seen that both IS2 and IS1 obtain very similar results. This finding illustrates that only an approximate MPFP is required when the IS method is performed. The FORM-UUA DoE samples cluster, which corresponds to 15 calls to performance function required for calculating the approximate MPFP, is used as the initial DoE and the population of the IS samples in Fig. 4 is chosen as candidate points. Then, the ALK model should be constructed to have a right prediction for the sign of the performance function. It can be observed from Fig. 4 that 7 training points at which the sign of G (X , Y ) has the maximum risk to be wrongly predicted are added into the DoE and then the ALK model satisfies the precision requirement as expected. The sign of minimizing G (X , Y ) at all of the IS samples should be predicted by the ALK model and then the estimator of the minimum failure probability can be obtained. From Table 1, it can be observed that ALK-MGHRA-IS with NIS = 104 samples is enough to estimate the precise failure probability. Moreover, the proposed method has a coefficient of variation on the failure probability close to 2.2%. However, ALK-MGHRA-MCS needs at least NMCS = 4×107 samples under the same conditions and then the extreme values of the Kriging model at such a plenty of simulated samples need to be calculated by using global optimization algorithm. Although this process needs not call the true performance function, ALK-MGHRA-MCS is still time-consuming. When calculating the maximum failure probability, one observes that the distance between MPFPU and MPFPL is quite close. Therefore, the previously constructed ALK model can rightly predict the sign of maximizing G (X , Y ) at the foregoing IS samples. Sequentially, the coefficient of variation and the maximum failure probability can be obtained. Because the coefficient of variation is sufficiently small (2.5% 5\% ), the above “tentative operation” is reasonable to use P f UIS as the estimate of the maximum failure probability. In this example, ALK-MGHRA-IS only calculates the MPFP once and then a Kriging model is enough to obtain the accurate bounds of the failure probability. However, FORM-UUA must be respectively implemented to obtain the bounds of the MPFPs. Because the performance function is highly nonlinear in the random space, ALK-MGHRA-IS obtains very precise results in comparison to FORM-UUA. In short, the proposed method exhibits much better than other different methods whatever in efficiency or accuracy. Finally, in order to illustrate the role of the range of Y on the precision of estimating failure probabilities, ALK-HRA-IS (ALK-HRA-IS: an existing active learning HRA method combining Kriging and IS method) is compared with ALK-MGHRA-IS. In the method of ALK-HRA-IS, the uncertain-bound variables are described with interval variables rather than parallelepiped convex variables (correlated interval variables). From Table 1, it is observed that the calculated failure probabilities with ALK-MGHRA-IS is more conservative than those with ALK-HRA-IS.
Fig. 4. Classification of the population of NIS points by ALK-MGHAR-IS.
calls, the calculated MPFP could be enough precision to obtain the required failure probability. Therefore, a large enough iteration number is required to obtain the MPFP. At this point, the MPFP with a large enough iteration number is called as the precise MPFP. Before ALKMGHRA-IS is performed, the MPFP is required to construct the ISD function such as Eqs. (24) and (25). However, it should be stressed that FORM-UUA might need a great deal of function calls for searching the MPFP. This step is quite time consuming as it is implemented on the true performance function. Therefore, only a few iteration number is required to obtain an approximate MPFP. Meanwhile, Table 1 have showed that a precise MPFP is not required and a few iterations providing an approximated MPFP is enough to meet the requirement of IS in accuracy. As mentioned earlier, an approximate MPFPU (AMPFPU) as shown in Fig. 4 can be obtained using FORM-UUA, which requires 15 calls to the performance function and four iterations. After a few iteration number is employed to obtain an approximate MPFPU (AMPFPU), the ISD function is constructed through using the AMPFPU. Sequentially, a population of the IS samples with NIS = 104 centred on the AMPFPU can be obtained by the aforementioned ISD function. Once the simulated IS samples are obtained, the lower bound of failure probability can be directly estimated by calling the performance function. Thus, the above method of calculating the low bound of failure probability is named as IS2. A population of the IS samples with NIS = 104 centred on the AMPFPU is observed in Fig. 4. Similarly, a large enough iteration number is employed to obtain a precise MPFPU (PMPFPU), and then the ISD function can be developed by using the PMPFPU. Sequentially, a population of NIS = 104 samples centred on the PMPFPU can be obtained. After the simulated IS samples are obtained through the above process, the lower bound of failure probability can be directly evaluated via the function calls. Therefore, the above approach is called IS1. A precise MPFPU (PMPFPU) requires 23 Table 1 Results of example 1 using several algorithms. Method
PU f
P fL
Function calls
Cov (upper bound/lower bound)
Error (upper bound/lower bound)
IS1 IS2 ALK-MGHRA-IS2 ALK-MGHRA-IS1
6.39 × 10−6 6.55 × 10−6 6.56 × 10−6
3.63 × 10−6 3.52 × 10−6 3.53 × 10−6
6.58 × 10
3.53 × 10
23 + 200 × 104 15 + 200 × 104 15 + 7 23 + 7
0.0217/0.0232 0.0226/0.0246 0.0226/0.0247 0.0223/0.0247
– 2.51%/3.03% 2.66%/3.03% 2.97%/2.75%
FORM-UUA-KKT
5.61 × 10 6 7.68 × 10−6
–
20.19%/15.43%
ALK-HRA-IS 1 2
6
6
3.02 × 10 6 4.19 × 10−6
0.0297/0.0303
22 + 6 23 + 23
Importance Sampling using a precision MPFP. Importance Sampling using an approximate MPFP. 8
12.20%/16.80%
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
F (t) Z (t) C1
minimum failure probability and the coefficient of variation can be computed by substituting the maximal response into Eqs. (19) and (21). The minimum value of G (U , Y ) is obtained by using the aforementioned IS samples and then the maximum failure probability and the coefficient of variation are calculated by substituting the minimal response into Eqs. (16) and (18). This method is also called IS1. Results of different methods are listed in Table 3. Because both the coefficients of variation are small enough, P f UIS and P f LIS are able to be used as the
F1 F(t)
m
C2
t1
t
Fig. 5. Example 2: nonlinear oscillator.
estimates of actual failure probabilities P Uf IS and P fLIS . Moreover, an approximate MPFPU is available by FORM-UUA approximation, which needs 36 calls to the performance function. The bounds of failure probability are calculated once again through IS2, Analogous to the procedure of IS1. From Table 3, it can be observed that the results of both IS2 and IS1 are basically coincide. That demonstrates that it needs not cost too much energy to be able to obtain a very precise result when IS2 is applied to calculate the MPFPU. The sampling points corresponding to the 36 calls to performance function required for IS2 are used as the initial DoE for ALK-MGHRA-IS and the IS samples (NIS = 104) are chosen as candidate points. The subsequent ALK model can be constructed for estimating the bounds of failure probabilities. It can be seen from Table 3 that ALK-HRA-IS requires 29 training points and it obtains the results with 65 function evaluations. The DoE including initial and added training points is shown in Fig. 6. It can been seen that most of added training points are located in the vicinity of the limit state. This illustrates that ALK-MGHRA-IS is capable of only locally approximating the performance function in the region of interest rather than in the whole uncertainty space. It should be noted that some of initial DoE focus on the vicinity of the limit state. That is because FORM-UUA tends to obtain the actual MPFPU which concentrates on the limit state after a couple of iterations. FORM-UUA requires the calculations of the MPFPL and MPFPU with 122 function computations. Because Eq. (44) is not a highly nonlinear function, both FORM-UUA and the proposed method exhibit very well in accuracy. However, the efficiency of ALK-MGHRA is much better than that of FORM-UUA due to this small probability event. For rare event probability estimation, a large amount of samples (at least NMCS = 1011) as candidate points are required if ALK-MGHRA-MCS is implemented. In the meantime, the extreme values of the Kriging model at so many sampling points entails to be computed by using global optimization algorithm. Therefore, this time-consuming task is unthinkable and unacceptable. However, ALKMGHRA-IS is able to produce highly precise results with only 104 sampling points. This demonstrates the efficiency and accuracy of the advocated method. In the end, the comparisons between ALK-HRA-IS and ALK-MGHRA-IS is shown in Table 3. It is derived that only interval variables will result in an overestimation of structural reliability. The negative evaluations will be hard to guarantee the structural safety.
These results indicate that the correlated interval variables have a positive influence on the estimation of structural reliability. 5.2. Example 2: Nonlinear oscillator The second example is a single degree of freedom nonlinear oscillator (Fig. 5), which was used in [44–46]. As shown in Fig. 5, the oscillator is subject to a rectangular pulse load. The performance function is defined as
G (X , Y ) = 3r
2F1 sin m 02
2 0 t1
2
(44)
where 0 = (c1 + c2)/ m ; X = (m , r , F1 denotes the vector of random variables and Y = (c1, c2, t1 )T denotes the vector of the interval variables. c1 and c2 are treated as dependent quantities and t1 is independent to the other interval variables. If the limit state function G (X , Y ) < 0 , the instability of nonlinear oscillator will happen. Sequentially, the nonlinear oscillator will be out of operation. Details are listed in Table 2. Assume that the correlation coefficients between each pair of interval variables are and c1 c 2 = c 2 c1 = 0.3 c1 t 1 = t 1 c1 = c 2 t1 = t1 c 2 = 0 . The correlation matrix can be constructed as:
)T
1 = 0.2 0 c1
0.2 0 1 0 0 1 c2 t1
c1 c2 t1
(45)
Thus, the shape matrix C can be calculated:
0.0833 0.0167 0 C = 0.0017 0.0083 0 0 0 1
(46)
Similarly, the constructed MP model with the shape matrix (46) needs to be transformed into the standard δ-space. When FORM-UUA and ALK-MGHRA-IS are performed, the standard interval space entails to be mapped into the original interval space. FORM-UUA requires the transformation of random variables into the standard U-space. A precise MPFPU can be calculated using FORM-UUA. Then a population of NIS = 104 samples centred on this point is obtained through the ISD function. The optimization algorithm is used to solve the maximum value of G (U , Y ) at each simulated sample. Subsequently, the
5.3. Example 3: Crank–Slider mechanism A crank-slider mechanism is adopted in a construction machine as depicted in Fig. 7. The length of the crank a, the length of the coupler b, the internal diameter d1 and the external diameter d2 are treated as correlated interval variables. The external force P and the Young’s modulus of the material of the coupler E are considered as random variables. Due to the rugged environment of the construction site, the coefficient of friction µ between the slider C and the ground NN is nonprobabilistic uncertainty. Because of installation error, the installation positions of the slider and the crank have a certain error, denoted as e. The coefficient of friction µ and the error e are regarded as independent interval variables. All of interval variables can be represented by using the MP model. Details including random variables and the MP variables are provided in Table 4. The performance function is defined as [16]:
Table 2 Uncertain variables of example 2. Variables
Variables type
Parameter 1
Parameter 2
m r F1 c1 c2 t1
Normal Normal Normal MP MP MP
1.01 0.51 1.01 0.9* 0.09* 0.5*
0.052 0.12 0.22 1.1** 1.0** 0.6**
* = lower bound of marginal interval. ** = upper bound of marginal interval. 1 = mean of normal distribution. 2 = standard deviation.
G (X , Y ) =
9
3E (d 4 2
64b2
d14 )
P (b
a)
(a + b ) 2
e2
µe
(47)
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
Table 3 Results of example 2 using several algorithms. Method
PU f
P fL
Function calls
Cov (upper bound/lower bound)
Error (upper bound/lower bound)
IS1 IS2 ALK-MGHRA-IS2 ALK-HRA-IS2
1.179 × 10−7 1.158 × 10−7 1.155 × 10−7
1.156 × 10−9 1.139 × 10−9 1.136 × 10−9
61 + 300 × 10−4 36 + 300 × 104 36 + 29 35 + 28
0.0259/0.0367 0.0283/0.0352 0.0281/0.0352 0.0307/0.0319
– 1.78%/1.47% 2.04%/1.73% 15.60%/15.57%
0.995 × 10 7 1.293 × 10−7
FORM-UUA-KKT 1 2
0.976 × 10 9 1.202 × 10−9
61 + 61
–
9.67%/3.98%
Importance Sampling using a precision MPFP. Importance Sampling using an approximate MPFP.
independent to the other parameters. The correlation matrix of the uncertain parameters can be created as follows:
1 0.7 0.2 0.1 0 0 a 0.7 1 0.5 0.6 0 0 b 0.2 0.5 1 0.8 0 0 d1 = 0.1 0.6 0.8 1 0 0 d2 0 0 0 0 1 0 e 0 0 0 0 0 1 µ a b d1 d2 e µ
Then the shape matrix C can be calculated by substituting Eqs. (48) into (A.4). The subsequent MP model is constructed with the shape matrix C , and it needs to be converted into the standard δ-space. When FORM-UUA is carried out and a large amount of candidate points needs to be generated, the standard δ-space can be mapped into the original interval space. Results of different methods are shown in Table 5. It is observed that IS1 obtains the precise MPFPU with 57 function computations and IS2 obtains the approximate MPFPU with 48 calls to performance function. Furthermore, both calculation results are very close. Again, IS does not require a very accurate MPFPU. The FORM-UUA DoE corresponds to the 48 function evaluations required for IS2 and which are used as the initial DoE and the population of IS samples with NIS = 104 are selected as candidate points. In creating Kriging model, 26 training points at which the sign of the performance function has the maximum risk to be wrongly predicted need to be added into DoE. ALKMGHRA-IS needs a total of 74 calls to the performance function. The DoE obtained by ALK-MGHRA-IS is shown in Fig. 8. It is observed that most of added training points and some of initial DoE still concentrates on the neighborhood of the limit state. That is because ALK-MGHRA-IS largely approximates the neighboring region of the limit state G (X , Y ) = 0 and a calculated MPFPU is prone to the actual one after several iterations of FORM-UUA. Based on the Kriging model, the bounds of failure probability and the corresponding coefficients of variation can be effectively solved by the IS method. 104 sampling points are used and the extreme values of the performance function are explored by global optimization algorithm. From Table 5, it is seen that ALK-MGHRA-IS obtains the very precise bounds of failure probability with only 104 samples. However, implementing MCS based on ALK-MGHRA needs at least 108 sampling points in the same circumstances. In order to ensure the accuracy of ALK-MGHRA-MCS, the minimal value of the Kriging model at so many simulated samples should be searched using global optimization algorithm. Therefore, the realization of ALK-MGHRA-MCS is a time-consuming process. This perfectly illustrates that the advocated method exhibits much better than other different methods whatever in accuracy or in efficiency. Meanwhile, Table 5 also provides the comparisons of the proposed method with ALK-HRA-IS. It can be summarized that uncorrelated interval variables have a negative impact on the precision of evaluating the failure probabilities, especially rare events.
Fig. 6. DoE of example 2 obtained with ALK-MGHAR-IS.
B
d2
b
a
d1
M
B
M
A
e N
P N
C
Fig. 7. Example 3: crank–slider mechanism. Table 4 Uncertain variables of example 3. Variables
Variables type
Parameter 1
Parameter 2
a/mm b/mm d1/mm d2/mm e/mm μ/mm P/KN E/Gpa
MP MP MP MP MP MP Normal Normal
98* 398* 50* 58* 100* 0.1* 2801 2001
102** 402** 54** 62** 150** 0.3** 282 202
(48)
* = lower bound of marginal interval. ** = upper bound of marginal interval. 1 = mean of normal distribution. 2 = standard deviation.
where X denotes the vector of random variables and Y denotes the vector of the MP variables. Assume that correlatively exists between the components of crank-slider and the diameters of M-M section, although the components of crank-slider and the diameters of M-M section might be independent. The coefficient of friction µ and the error e are
5.4. Example 4: three-bay five-story frame A three-bay five-story frame [47] shown in Fig. 9 is investigated in this subsection. From Fig. 9, it is seen that five concentrated loads are 10
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
Table 5 Results of example 3 using several algorithms. Method
PU f
P fL
Function calls
Cov (upper bound/lower bound)
Error (upper bound/lower bound)
IS1 IS2 ALK-MGHRA-IS2 ALK-HRA-IS2
9.285 × 10−6 9.535 × 10−6 9.393 × 10−6
5.179 × 10−6 5.126 × 10−6 5.198 × 10−6
57 + 600 × 104 48 + 600 × 104 48 + 26 47 + 25
0.0208/0.0212 0.0216/0.0223 0.0217/0.0229 0.0298/0.0305
– 2.69%/1.02% 1.16%/0.36% 5.25%/10.17%
FORM-UUA-KKT 1 2
8.797 × 10 6 9.662 × 10−6
4.652 × 10 6 5.325 × 10−6
57 + 57
–
4.06%/2.82%
Importance Sampling using a precision MPFP. Importance Sampling using an approximate MPFP.
information is listed in Table 6. These random variables are assumed to be mutually independent. Five concentrated loads are described by using the MP model. Assume that the correlation coefficients between any two loads are Pi Pj = Pj Pi = 0.5, i j, i = 1, 2. ..,5; j = 1, 2. ..,5. Therefore, the correlation matrix can be constructed as follows:
1 0.5 = 0.5 0.5 0.5
0.5 1 0.5 0.5 0.5
0.5 0.5 1 0.5 0.5
0.5 0.5 0.5 1 0.5
0.5 0.5 0.5 0.5 1
P1 P2 P3 P4 P5
(49)
P1 P2 P3 P4 P5
Sequentially, the MP model can be constructed by using Eq. (B.2). To guarantee the serviceability of the frame, the floor displacement should not exceed a threshold. Therefore, the failure criterion is defined by an implicit performance function:
G (X , Y ) = Zlim
where Zlim = 0.2 is the allowable value of the horizontal deformation; (X , Y ) is the horizontal top floor displacement which can be calculated by the FEA platform; the vector X represents all random variables; the vector Y denotes all the correlated interval variables shown in Table 7. The FE model here is developed by using the platform OpenSees [48–50]. Nonlinear beam-column elements are used to model column and beam members. All the beam-column connections are modeled as fully rigid. The floor slabs are modeled with a rigid diaphragm MPCs (multi-point constraints). The horizontal top floor displacement is obtained through OpenSees. Results calculated by different methods are presented in Table 8. For this complicated engineering problem, the results of both IS2 and IS1 access to each other. It can be observed from Table 8 that IS1 obtains the precise MPFPU with calling the platform OpenSees 69 times and IS2 obtains the approximate MPFPU with calling the FEA platform only 53 times. This illustrates that when performing the IS method it does not need to compute a very accurate MPFPU. The sampling points that corresponds to the 53 function calls required for IS2 are taken as the initial DoE and the population of the IS samples are chosen as candidate points. When constructing the Kriging model, 32 training points at which the sign of the implicit performance function has the largest risk to be wrongly predicted should be added into DoE. The presented method requires 85 calls to the FEA platform. From Table 8, it is seen that ALK-MGHRA-IS obtains the very precise results with calling OpenSees only 85 times. Based on the constructed ALK
Fig. 8. DoE of example 3 obtained with ALK-MGHAR-IS.
12'
P5
12'
P4
12'
P3
12'
P2
16'
P1
25'
30'
25'
Fig. 9. Example 4: three-bay five-story frame. Table 6 Random variables of example 4. Variables
Variables type
Mean
Standard deviation
Young’s modulus of columns, E1 Young’s modulus of beams, E2 Moment of inertia of columns, I1 Moment of inertia of beams, I2 Cross section area of columns, A1 Cross section area of beams, A2
Normal Normal Normal Normal Normal Normal
497,000 454,000 3.00 2.69 6.00 4.00
40,000 40,000 0.35 0.65 1.20 1.20
(50)
(X , Y )
Table 7 The MP variables of example 4.
Note: the units of Ei, I, and Ai are kips/ft2, ft4, and ft2, respectively.
exerted on the frame with four columns and five beams. The Young’s modulus, the cross sectional area and the moment of inertia of each member are considered as random variables. Their statistical 11
Variables
Variables type
Lower bound of marginal interval
Upper bound of marginal interval
P1/kips P2/kips P3/kips P4/kips P5/kips
MP MP MP MP MP
12.0 16.0 20.0 24.0 28
15.5 19.5 23.5 27.5 32
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
Table 8 Results of example 4 using several algorithms. Method
PU f
P fL
Function calls
Cov (upper bound/lower bound)
Error (upper bound/lower bound)
IS1 IS2 ALK-MGHRA-IS2 ALK-HRA-IS2
1.736 × 10−5 1.708 × 10−5 1.704 × 10−5
1.26810−5 1.255 × 10−5 1.253 × 10−5
69 + 500 × 104 53 + 500 × 104 53 + 32 52 + 31
0.0228/0.0253 0.0237/0.0251 0.0237/0.0251 0.0299/0.0322
– 1.61%/1.03% 1.84%/1.18% 36.75%/24.21%
FORM-UUA-KKT 1 2
1.098 × 10 5 1.969 × 10−5
0.961 × 10 5 1.375 × 10−5
69 + 69
–
11.83%/8.44%
Importance Sampling using a precision MPFP. Importance Sampling using an approximate MPFP.
based Importance Sampling. ALK-MGHRA-IS can evaluate the bounds of small failure probabilities of complicated structures involving in “multisource uncertainty” problems and performance functions linked to time-consuming FE models. The proposed method enables the correction of the FORM-UUA approximation with only a few function computations. To make ALK-MGHRA-IS more efficient, an optimization method based on Karush–Kuhn–Tucker conditions (KKT) is introduce to relieve the burden of searching the extreme values. Firstly, we create the ISD function centered on the MPFP found using FORM-UUA approximation. Sequentially, a population of IS samples centered on the MPFP is simulated via the constructed ISD function and can be chosen as candidate points. To minimize the function calls, ALK-MGHAR-IS uses a learning function, named as ERF originated from the predicted Kriging model. The sign of minimizing response function at each IS sample is then predicted by using the created ALK model. Through this treatment, the minimum failure probability is obtained by IS method. Similarly, the maximum failure probability can be estimated using the identical Kriging model and the same IS samples. If the coefficient of variation is sufficiently small, P f UIS defined in Eq. (16) can be taken as an estimator of the maximal failure probability. Three academic examples and one engineering example validate the computation benefit of the proposed method. From four cases, it is demonstrated that ALK-MGHRA-IS exhibits much better than FORMUUA whatever in accuracy or in efficiency. The proposed method only requires an approximated MPFP which satisfies the demand of IS in accuracy after several iterations. Besides, ALK-MGHRA-IS makes FORM-UUA-KKT execute only once for searching one MPFP. When FORM-UUA is employed to estimate the bounds of failure probability, FORM-UUA-KKT needs to be respectively implemented for exploring bounds of MPFPs. Meanwhile, ALK-MGHRA-IS can obtain very precise failure probabilities with only 104 samples. Therefore, the MGHRA method closely integrating IS and the Kringing model can commendably meet the required precision in comparison to other methods. However, it should be noted that the proposition of FORM-UUA stage of the presented algorithm is based on the assumption that there exists a single failure region with a unique MPFP. Considering the existence of multiple-disconnected failure regions, ALK-MGHAR-IS may cause biased estimates of the bounds on failure probability. To eliminate this drawback, an adaptive strategy by coupling multi-modal IS with the ALK model will be investigated in our further work. An effective algorithm should be exploited to automatically explore the existing multiple MPFPs corresponding to multiple, non-connected failure domains.
Fig. 10. DoE of example 4 obtained with ALK-MGHAR-IS.
model, the extreme values of the performance function at each simulated IS sample are searched by global optimization algorithm. The sign of minimal function value at each sample is predicted using the Kriging model and then the bounds of failure probability can be calculated. However, the implementation of MCS based on ALK-MGHRA requires at least 108 samples. Apparently, this procedure becomes very timeconsuming. That is why the combination of IS method and ALK model should be proposed in MGHRA. The DoE originated from iterative process of ALK-MGHRA-IS is depicted in Fig. 10. From which, again it is observed that the present method concentrates more energy on the neighboring region of G (X , Y ) = 0 . That is the reason why the advocated method can maintain the efficiency. Furthermore, the results with ALK-HRA-IS are also shown in Table 7. From the comparisons of the ALK-HRA-IS method with the proposed approach, it can be found that the uncorrelated interval variables have an adverse effect on the precision of evaluating the small failure probabilities, involving in the complex engineering applications. Therefore, it is demonstrated from this numerical example that ALK-MGHRA-IS can deal with complex “multisource uncertainty” engineering problems with implicit performance functions, especially small probability events. 6. Conclusion The existing ALK-HRA approach is found incapable of assessing small failure probabilities, and it is also not suitable for dealing with “multisource uncertainty” problems. Therefore, a MGHRA problem with both random and MP convex variables can be proposed in this paper. Then, we can present an original and implementable approach called ALK-MGHRA-IS for active learning MGHRA method and Kriging-
Acknowledgements Our authors acknowledge the support of National Natural Science Foundation of China (grant number 11802224), China Postdoctoral Science Foundation (grant number 2018M633495).
Appendix A The correlation coefficient between Y1 and Y2 is defined by the following expression [34,35]: 12
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
b a b+a
=
Y1 Y2
(A.1)
where a and b represent the semi-axis length in the direction HB and HA , respectively. When the two intervals Y1 and Y2 are linearly correlated, the correlation coefficient Y1 Y2 = 1 or Y1 Y2 = 1. When Y1 and Y2 are positively or negatively correlated, there exists 0 < Y1 Y2 < 1 or 1 < Y1 Y2 < 0 . If Y1 and Y2 are independent of each other, the MP model will degenerates into a traditional interval model with Y1 Y2 = 0 . Then the correlation matrix can be defined as:
=
Y1 Y1
Y1 Y2
Y2 Y1
Y2 Y2
=
1
Y1 Y2
Y2 Y1
1
(A.2)
1. where Y1 Y2 = Y2 Y1 and 1 Y1 Y2 Therefore, the uncertain domain of the interval variables can be explicitly quantified as: 1T 1R 1 (Y
|
Y C )|
1
e = |C
Y|
(A.3)
e
where | | represents the absolute value sign and the matrices T , R , C and the vectors Y , Y C , Y and e are described as:
T = diag(w1,w2), wi =
R=
1 2 j =1
| (i , j )|
i = 1,
2
(A.4)
diag(Y1W ,Y2W ),
C = RT ,
Y = ( Y1 Y2 )T , T
Y C = ( Y1C Y 2C ) ,
Y C,
Y=Y
e = (1 1)T . Appendix B For a general n-dimensional problem, the correlation coefficient between any two interval variables Yi and Yj can be defined as the definition of Eq. (A.2), the correlation matrix is defined as:
=
11
12
1n
1
21
22
2n
21
1
n1
n2
nn
n1
n2
where
|C
1
12
ij .
Analogous to
1n 2n
1
(B.1)
is an n × n symmetric matrix. Subsequently, the uncertain domain can be expressed as follows [34,35]:
Y|
(B.2)
e
where C is defined as an n × n shape matrix and Y = Y Y C . It should be noted that all of the uncertain parameters’ marginal intervals and correlations are required for creating an MP model. Moreover, the process of constructing an MP model has been elaborated in [33,34,36]. To make the proposed reliability more convenient, the correlative interval variables should be converted to independent ones by a matrix transformation [34,35]. Introduce a vector :
= (RT ) 1 (Y
Y C) = C
1
(B.3)
Y
Sequentially, the intervals can be mapped into the standard cube
= { || |
space and the uncertain domain
becomes an n-dimensional cubic domain
cube :
(B.4)
e}
where =( 1, 2, …, n )T , i = [ 1, 1], i = 1, 2, …, n . After the matrix transformation, the uncertain domain becomes a standard multidimensional cube with one-half the length of each side 1. Subsequently, for two interval variables Y1 and Y2 , the standard interval space can be mapped into the original interval space by the following expression.
Y1 = Y2 =
Y1W
w1 1
Y2W
w2 2
+ +
Y1W
w1 Y1 Y2 2
Y2W w2
Y2 Y1 1
+ Y1C + Y2C
where w1 and w2 are defined in Eq. (A.4) and w1 = w2 ,
(B.5) 1
and
2
are independent of each other.
Appendix C Kriging model can provide the predictions for G ( ) at unknown points after a design of experiments (DoE) is selected to determine its statistical parameters. A given DoE is defined as [ (1) , (2), , (k) ]T with (i ) the ith training points, and the vector of corresponding responses is denoted as G = [G ( (1) ), G ( (2) ), , G ( (k ) )]T with G ( (i) ) the ith response. In this paper, ordinary Kriging is adopted to simplify Eq. (42) and then F ( , ) is 13
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
replaced by a constant . Therefore, the predicted value µG ( ) and predicted variance
µG ( ) = 2( G
rT (
+ 2
)=
) R (G
) of Kriging model G ( ) can be respectively expressed as: (C.1)
1 ) 1T R 1r ( )) 2 1T R 11
(1
r T ( ) Rr ( ) +
1
2( G
(C.2)
In Eqs. (C.1)–(C.2), 1 denotes k-dimensional unit vector, r ( ) is the correlation vector between
r ( ) = [R ( ,
(1) ),R (
,
and ith training point, and it can be defined as (C.3)
, x , x (2) ), R ( , x , x (k) )
In Eqs. (C.1)–(C.2), R is the correlation matrix between each pair of training points and it can be expressed as
R( ,
(1) ,
(1))
R( ,
(1) ,
(m ) )
R( ,
(m ) ,
(1) )
R( ,
(m) ,
(m) )
R= The parameters
and
2
(C.4)
in Eqs. (C.1)–(C.2) are respectively given by (C.5)
= (1T R 11) 11T R 1G 2
=
1 (G k
1)T R 1 (G
1)
(C.6)
Appendix D The sensitivity analysis of the correlated interval variables with respect to the failure probability of performance function is derived as follows: The following definition is firstly given p i
= P Uf = (Y i
P fL; P fC = (PUf + PfL )/2 U
L
Yi )/2;
YiC
= (Y i
U
(D.1) (D.2)
L
+ Yi )/2
are the marginal interval width and marginal interval central value of failure probability, respectively; i and are the width and where p and the central value of marginal intervals. According to the existing literature [1], the maximum and minimum failure probability are only reflected in the interval bounds. Therefore, the following four types are described as: Type I sensitivity p/ i p/ i is the sensitivity of the width of the Pf bounds, p , with respect to the interval width of the i th interval variable Yi , i , p is defined by
YiC
P fC
(1) P Uf appears on Yi U , P fL appears on Yi L p
=
(PUf
i
P fL )
=
i
PUf
1 2
Yi
P fL
+
U
Yi L
(D.3)
(2) P Uf appears on Yi L , P fL appears on Yi U p
=
(PUf
i
P fL )
PfL
1 2
=
i
U
Yi
PUf
+
Yi
L
(D.4)
Type II sensitivity i P fC / i is the average probability of failure P fC with respect to YiC
P fC /
(1) P Uf appears on Yi U , P fL appears on Yi L
P fc
=
i
U L 1 (P f + P f ) 1 = 2 4 i
P Uf Yi
P fL
U
Yi L
(D.5)
(2) P Uf appears on Yi L , P fL appears on Yi U
P fc
=
i
U L 1 (P f + P f ) = 2 i
1 4
PUf Yi
U
+
P fL Yi L
Type III sensitivity p/ C p/ Y i is the sensitivity of the width of the probability of failure
(D.6)
YiC
p
with respect to the average of the ith marginal interval variable, YiC .
(1) P Uf appears on Yi U , P fL appears on Yi L p YiC
=
(PUf
PfL ) YiC
=
PUf
PfL
Yi U
Yi
(D.7)
L
14
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
(2) P Uf appears on Yi L , P fL appears on Yi U p YiC
(PUf
=
PfL )
=
YiC
PUf Yi
PfL
L
(D.8)
U
Yi
Type IV sensitivity P fC / YiC is the sensitivity of the average of the probability of failure P fC with respect to the average of the ith marginal interval variable, YiC .
P fC /
YiC
(1) P Uf appears on Yi U , P fL appears on Yi L
P fC YiC
=
U L 1 (P f + P f ) 1 = C 2 2 Yi
PUf Yi
U
P fL
+
Yi
L
(D.9)
(2) P Uf appears on Yi L , P fL appears on Yi U
P fC YiC
=
where
U L 1 (P f + P f ) 1 = 2 2 YiC
PU f Yi U
PU f
and
Yi L
PfL Yi
U
PUf
+
Yi
L
(D.10)
are the marginal interval upper of the probability of failure with respect to the upper and lower of partial derivatives of marginal P fL
interval variables, respectively;
Yi U
P fL
and
are the marginal interval lower of the probability of failure with respect to the upper and lower of
Yi L
partial derivatives of marginal interval variables. The partial derivatives are obtained: PU f
= Yi U
U
U)
(
(µg3 )U
1 3!
( (µ )
3 3 U [(µg2 ) 2 ]U g
1 3!
Yi U
(µg2 )U 5 3 (µg3 )U [(µg2) 2 ]U Y U 2 i
[(µg2 ) 2]U
2(µg4 )U [(µg2 ) 3]U
(µg4 )U
+
1 4!
+
1 4!
((µg4 )U [(µg2) 2]U
+
10 6!
2 (µg3 )U [(µg2) 2 ]U
+
10 6!
2 (µg3 )U [(µg2) 2 ]U
+
10 6!
where
Yi U
(
3
(
3
( (µ )
U
=
(5) (
3)
(µg3 )U
)
Yi U
(µg2 )U Yi U U
U)
Yi U
3
[(µg2) 2 ]U
(6) (
(µg2 )U 5 3 (µg3 )U [(µg2 ) 2 ]U Y U 2 i
)
(7) (
1 (µg1 )U [(µg2) 2 ]U
U)
U
U
U)
+ ...
Yi U
(3) (
U)
(4) (
+ ...
U)
+ ...
+ ...
)
2 3 3 U [(µg2) 2 ]U g
Yi U
(4) (
[(µg2 ) 2 ]U
3
Yi U
)
U)
+ ... (6) (
U)
+ ...
...
Yi U
(D.11)
is the upper of reliability index with respect to the partial derivative of marginal interval upper;
Yi U
upper of performance function with respect to the partial derivative of marginal interval upper; function with respect to the partial derivative of marginal interval upper;
(µg4 )U
the partial derivative of marginal interval upper. PU f
= Yi L
U
U)
(
(µg3 )U
1 3!
Yi L
( (µ )
3 3 U [(µg2) 2 ]U g
1 3!
Yi L
)
(4) (
[(µg2) 2 ]U
(µg2 )U 5 3 (µg3 )U [(µg2 ) 2 ]U Y L 2 i
[(µg2) 2]U
2(µg4 )U [(µg2 ) 3]U
(µg4 )U
3
+
1 4!
+
1 4!
((µg4 )U [(µg2 ) 2]U
+
10 6!
2 (µg3 )U [(µg2 ) 2 ]U
+
10 6!
2 (µg3 )U [(µg2 ) 2 ]U
+
10 6!
Yi L
(
3
(
3
( (µ )
3)
)
(5) (
(µg3 )U Yi U
U) 3
U
Yi L
[(µg2 ) 2 ]U
(µg2 )U Yi L
(6) (
(µg2 )U 5 3 (µg3 )U [(µg2) 2 ]U Y L 2 i
)
(7) (
U)
U
Yi L
U)
(4) (
is the variance
is the third moment upper of performance
is the fourth moment upper of performance function with respect to
+ ...
Yi L
(3) (
Yi U
U)
+ ... + ...
+ ...
)
2 3 3 U [(µg2) 2 ]U g
U
U)
Yi U
(µg3 )U Yi U
(µg2 )U
U)
+ ... (6) (
U)
+ ...
...
(D.12)
15
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff P fL
= Yi U
L
U)
(
(µg3 ) L
1 3!
( (µ ) [(µ ) ] )
1 3!
Yi U
3 L g
(µg4 ) L
(µg2 ) L 5 3 (µg3 ) L [(µg2) 2 ]L Y U 2 i
[(µg2) 2]L
2(µg4 ) L [(µg2 ) 3]L
+
1 4!
+
1 4!
((µg4 ) L [(µg2 ) 2]L
+
10 6!
2 (µg3 ) L [(µg2 ) 2 ]L
+
10 6!
2 (µg3 ) L [(µg2 ) 2 ]L
+
10 6!
P fL
= Yi L
Yi U
(
3
(
3
3 L g
(µg3 ) L
)
Yi U
U) (µg3 ) L Yi L
1 3!
(µg4 ) L
L
L)
L
3 L g
3
3 L 2
2 g
[(µg2) 2]L
2(µg4 ) L [(µg2 ) 3]L
1 4!
+
1 4!
((µg4 ) L [(µg2 ) 2]L
+
10 6!
2 (µg3 ) L [(µg2 ) 2 ]L
+
10 6!
2 (µg3 ) L [(µg2 ) 2 ]L
+
10 6!
(
3
(
3
3)
) )
( (µ ) [(µ ) ] ) 3 L g
2 g
3 L 2 2
(5) (
(µg3 ) L Yi U
(4) (
2 L 5 L (µg ) 2] Yi L
3 (µg3 ) L [(µg2) 2
L) 3
(µg2 ) L Yi L L
Yi L
[(µg2) 2 ]L
L)
L
Yi L
+ ...
L)
+ ... (6) (
L)
+ ... (D.13)
L)
L
Yi L
(3) (
+ ...
L)
(4) (
L)
+ ... + ...
+ ... (6) (
(µg2 ) L 5 3 (µg3 ) L [(µg2) 2 ]L Y L 2 i (7) (
L)
...
Yi U
[(µg2) 2 ]L
+
Yi L
(7) (
+ ...
+ ... (6) (
( (µ ) [(µ ) ] )
1 3!
Yi L
3
[(µg2) 2 ]L
+ ...
L)
(4) (
Yi U
Yi U
L
Yi U
(3) (
(µg2 ) L
L
L)
L)
(µg2 ) L 5 3 (µg3 ) L [(µg2 ) 2 ]L Y U 2 i
)
3 L 2 2
2 g
(5) (
3)
( (µ ) [(µ ) ] ) (
(4) (
[(µg2) 2 ]L
3
Yi U
3 L 2
2 g
L)
+ ... (6) (
L)
+ ...
...
(D.14)
Appendix E. Supplementary data Supplementary data to this article can be found online at https://doi.org/10.1016/j.strusafe.2019.101875.
Methods 2015;12(4):1540006. [14] Du XP, Sudjianto A, Huang B. Reliability-based design with the mixture of random and interval variables. J Mech Des 2005;127(6):1068–76. [15] Du XP. Interval reliability analysis, in: Proceedings of the ASME 2007 design engineering technical conference and computers and information in engineering conference, Las Vegas, Nevada, USA, 2007; 1103–1109. [16] Guo J, Du XP. Sensitivity analysis with mixture of epistemic and aleatory uncertainties. AIAA J 2007;45:2337–49. [17] Du XP. Unified uncertainty analysis by the first order reliability method. J Mech Des 2008:130. 091401-091410. [18] Qiu Z, Wang J. The interval estimation of reliability for probabilistic and nonprobabilistic hybrid structural system. Eng Fail Anal 2010;17(5):1142–54. [19] Wang J, Qiu Z. The reliability analysis of probabilistic and interval hybrid structural system. Appl Math Model 2010;34(11):3648–58. [20] Yao W, Chen X, Huang Y, van Tooren M. An enhanced unified uncertainty analysis approach based on first order reliability method with single-level optimization. Reliab Eng Syst Saf 2013;116:28–37. [21] Jiang C, Lu GY, Han X, Liu LX. A new reliability analysis method for uncertain structures with random and interval variables. Int J Mech Mater Des 2012;8(2):169–82. [22] Xiao NC, Huang HZ, Wang Z, Liu Y, Zhang XL. Unified uncertainty analysis by the mean value first order saddlepoint approximation. Struct Multidiscip Optim 2012;46(6):803–12. [23] Melchers RE. Radial importance sampling for structural reliability. J Eng Mech 1990;116(1):189–203. [24] Cadini F, Gioletta, Zio E. An improvement of a metamodel-based importance sampling algorithm for estimating small failure probabilities. Vulnerability, uncertainty, and risk: quantification, mitigation, and management. ASCE; 2014. p. 2104–14. [25] Balesdent M, Morio J, Marzat J. Kriging-based adaptive importance sampling algorithms for rare event estimation. Struct Saf 2013;44:1–10. [26] Dubourg V, Sudret B, Deheeger F. Metamodel-based importance sampling for structural reliability analysis. Probab Eng Mech 2013;33:47–57. [27] Cadini F, Gioletta A, Zio E. Improved metamodel-based importance sampling for the
References [1] Yang X, Liu Y, Gao Y, Zhang Y, Gao Z. An active learning Kriging model for hybrid reliability analysis with both random and interval variables. Struct Multidiscip Optim 2015;51(5):1003–16. [2] Zhang Y, Liu Y, Yang X, Zhao B. An efficient Kriging method for global sensitivity of structural reliability analysis with non-probabilistic convex model. Proc Instit Mech Eng Part O: J Risk Reliab 2015;229(5):442–55. [3] Yang XF, Liu YS, Zhang YS, Yue ZF. Probability and convex set hybrid reliability analysis based on active learning Kriging model. Appl Math Model 2015;39(14):3954–71. [4] Luo YJ, Li A, Kang Z. Reliability-based design optimization of adhesive bonded steel–concrete composite beams with probabilistic and non-probabilistic uncertainties. Eng Struct 2011;33(7):2110–9. [5] Kang Z, Luo YJ. Reliability-based structural optimization with probability and convex set hybrid models. Struct Multidiscip Optim 2010;42(1):89–102. [6] Qiu Z, Wang J. The interval estimation of reliability for probabilistic and nonprobabilistic hybrid structural system. Eng Fail Anal 2010;17(5):1142–54. [7] Wu J, Luo Z, Zhang N, Zhang Y. A new uncertain analysis method and its application in vehicle dynamics. Mech Syst Sig Process 2015;50:659–75. [8] Guo J, Du XP. Reliability analysis for multidisciplinary systems with random and interval variables. AIAA J 2010;48(1):82–91. [9] Jiang C, Han X, Lu GY. A hybrid reliability model for structures with truncated probability distributions. Acta Mech 2012;223(9):2021–38. [10] Guo J, Du XP. Reliability sensitivity analysis with random and interval variables. Int J Numer Meth Eng 2009;78:1585–617. [11] Jiang C, Long XY, Han X, Tao YR, Liu J. Probability-interval hybrid reliability analysis for cracked structures existing epistemic uncertainty. Eng Fract Mech 2013;112:148–64. [12] Luo YJ, Kang Z, Li A. Structural reliability assessment based on probability and convex set mixed model. Comput Struct 2009;87(21):1408–15. [13] Jiang C, Zheng J, Ni BY, Han X. A probabilistic and interval hybrid reliability analysis method for structures with correlated uncertain parameters. Int J Comput
16
Structural Safety 82 (2020) 101875
X.-X. Liu and I. Elishakoff
[28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38]
performance assessment of radioactive waste repositories. Reliab Eng Syst Saf 2015;134:188–97. Echard B, Gayton N, Lemaire M, Relun N. A combined importance sampling and Kriging reliability method for small failure probabilities with time-demanding numerical models. Reliab Eng Syst Saf 2013;111:232–40. Dubourg V, Sudret B. Meta-model-based importance sampling for reliability sensitivity analysis. Struct Saf 2014;49:27–36. Cadini F, Santos F, Zio E. An improved adaptive kriging-based importance technique for sampling multiple failure regions of low probability. Reliab Eng Syst Saf 2014;131:109–17. Zhao H, Yue Z, Liu Y, Gao Z, Zhang Y. An efficient reliability method combining adaptive importance sampling and Kriging metamodel. Appl Math Model 2015;39(7):1853–66. Zhang H. Interval importance sampling method for finite element-based structural reliability assessment under parameter uncertainties. Struct Saf 2012;38:1–10. Jiang C, Zhang QF, Han X, Qian YH. A non-probabilistic structural reliability analysis method based on a multidimensional parallelepiped convex model. Acta Mech 2014;225(2):383–95. Ni BY, Jiang C, Han X. An improved multidimensional parallelepiped non-probabilistic model for structural uncertainty analysis. Appl Math Model 2016;40(7):4727–45. Jiang C, Fu CM, Ni BY, Han X. Interval arithmetic operations for uncertainty analysis with correlated interval variables. Acta Mech Sin 2015:1–10. Jiang C, Zhang QF, Han X, Liu J, Hu DA. Multidimensional parallelepiped model—a new type of non-probabilistic convex model for structural uncertainty analysis. Int J Numer Meth Eng 2015;103(1):31–59. Tavassoli A, Hajikolaei KH, Sadeqi S, Wang GG, Kjeang E. Modification of DIRECT for high dimensional design problems. Eng Optim 2014;46:810–23. Kaymaz I. Application of kriging method to structural reliability problems. Struct
Saf 2005;27:133–51. [39] Wang H, Zhu X, Du Z. Research on improved EGO algorithm based on Kriging surrogate model. J Eng Des 2009;4:011. [40] Jones DR, Schonlau M, Welch WJ. Efficient global optimization of expensive blackbox functions. J Global Optim 1998;13(4):455–92. [41] Bichon BJ, Eldred MS, Swiler LP, Mahadevan S, McFarland JM. Efficient global reliability analysis for nonlinear implicit performance functions. AIAA J 2008;46(10):2459–68. [42] Lv Z, Lu Z, Wang P. A new learning function for Kriging and its applications to solve reliability problems in engineering. Comput Math Appl 2015;70(5):1182–97. [43] Yang X, Liu Y, Gao Y. Unified reliability analysis by active learning Kriging model combining with random-set based Monte Carlo simulation method. Int J Numer Meth Eng 2016. https://doi.org/10.1002/nme.5255. [44] Cimellaro GP, Reinhorn A. Multidimensional performance limit state for hazard fragility functions. J Eng Mech 2010;137(1):47–60. [45] Liu XX, Wu ZY, Liang F. Multidimensional performance limit state for probabilistic seismic demand analysis. Bull Earthq Eng 2016;14(12):3389–408. [46] Gayton N, Bourinet JM, Lemaire M. CQ2RS: a new statistical approach to the response surface method for reliability analysis. Struct Saf 2003;25(1):99–121. [47] Schueremans L, Van Gemert D. Benefit of splines and neural networks in simulation based structural reliability analysis. Struct Saf 2005;27(3):246–61. [48] Echard B, Gayton N, Lemaire M. AK-MCS: an active learning reliability method combining Kriging and Monte Carlo simulation. Struct Saf 2011;33(2):145–54. [49] Guan XL, Melchers RE. Effect of response surface parameter variation on structural reliability estimates. Struct Saf 2001;23(4):429–44. [50] Mazzoni S, McKenna F, Scott MH, Fenves GL. OpenSees version 1.7.3 command language manual. Pacific Earthquake Engineering Research Center, University of California, Berkeley, CA, USA; 2006. Retrieved from http://opensees.berkeley.edu/ OpenSees/manuals/usermanual/index.html.
17