An adaptive trivariate dimension-reduction method for statistical moments assessment and reliability analysis

An adaptive trivariate dimension-reduction method for statistical moments assessment and reliability analysis

Journal Pre-proof

An adaptive trivariate dimension-reduction method for statistical moments assessment and reliability analysis Jun Xu, Lijuan Zhou PII: DOI: Reference:

S0307-904X(20)30076-7 https://doi.org/10.1016/j.apm.2020.01.065 APM 13278

To appear in:

Applied Mathematical Modelling

Received date: Revised date: Accepted date:

25 July 2019 20 January 2020 30 January 2020

Please cite this article as: Jun Xu, Lijuan Zhou, An adaptive trivariate dimension-reduction method for statistical moments assessment and reliability analysis, Applied Mathematical Modelling (2020), doi: https://doi.org/10.1016/j.apm.2020.01.065

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier Inc.

Highlights • An adaptive trivariate dimension-reduction method is proposed for statistical moments assessment. • Delineating the existence of trivariate and bivariate cross terms is performed. • High-order unscented transformation is employed for numerical integrations. • Numbers of low-dimensional integrals and points for integrations are reduced. • Numerical examples demonstrate the efficacy of the proposed method.

1

An adaptive trivariate dimension-reduction method for statistical moments assessment and reliability analysis Jun Xua,b,∗, Lijuan Zhoua b Hunan

a College of Civil Engineering, Hunan University, Changsha, 410082, PR China. Provincial Key Lab on Damage Diagnosis for Engineering Structures, Hunan University, Changsha, 410082, PR China.

Abstract An adaptive trivariate dimension-reduction method is proposed for statistical moments evaluation and reliability analysis in this paper. First, the raw moments of the performance function can be estimated by means of the trivariate dimension-reduction method, where the trivariate, bivariate and univariate Gaussian-weighted integrals are involved. Since the trivariate and bivariate integrals control the efficiency and accuracy, delineating the existence of bivariate and trivariate cross terms is performed, which could significantly reduce the numbers of trivariate and bivariate integrals to be evaluated. When the cross terms exist, the trivariate and bivariate integrals are numerically evaluated directly by the high-order unscented transformation, where the involved free parameters are provided. When the cross terms don’t exist, the trivariate and bivariate integrals can be further decomposed to be the lower-dimensional integrals, where the high-order unscented transformation is again adopted for numerical integrations. In that regard, the first-four central moments can be computed accordingly and the performance function’s probability density function can be reconstructed by fitting the shifted generalized lognormal distribution model based on the first-four central moments. Then, the failure probability can be computed by a one-dimensional integral over the performance function’s probability density function in the failure domain. Three numerical examples, including both the explicit and implicit performance functions, are investigated, to demonstrate the efficacy of the proposed method for both the statistical moments assessment and reliability analysis. Keywords: Reliability analysis, Statistical moments, Trivariate dimension-reduction method, Cross terms, High-order unscented transformation

1. Introduction As the recognition and emphasis of probability-based design for engineering structures, structural reliability analysis plays an increasingly important role in the field of structural engineering [1, 2]. The key problem of structural reliability analysis is to compute the failure probability, denoted as pf . The performance function of the structure, or called the limit state function, can be written as: Z = G (X) (1) where X = [X1 , X2 , . . . , Xn ] is a vector containing a total of n basic random variables, G represents the deterministic function, and Z < 0 indicates the failure of the structure. Usually, the failure probability pf could be expressed as Z pf = Pr [G (X) ≤ 0] = fX (x)dx (2) G(X)≤0

where Pr denotes the probability and fX (x) represents the joint probability density functions (PDF) of input random vector X. Despite the simplicity of Eq. (2) for assessing the failure probability, great difficulties still arise for practical ∗ Corresponding

author Email addresses: [email protected] (Jun Xu), [email protected] (Lijuan Zhou)

Preprint submitted to Elservier

January 20, 2020

Notations α3Z

skewness of the performance function

Jp

the value represents the existence or absence of the cross terms

α4Z

kurtosis of the performance function

β

reliability index

N

δ

limit used to determine the cross terms exist or not

total number of deterministic model evaluations

P (·)

γ

free parameters in HUT

Γ (·)

gamma function

µZ

mean of the performance function

Pr

probability

ρ

shape parameters for SGLD

σZ

standard deviation of the performance function

τ

scale parameters for SGLD

ε

shape parameters for SGLD

A0

weight of the first types of sigma points in HUT

A1

weight of the second types of sigma points in HUT

A2

weight of the third types of sigma points in HUT

b

location parameters for SGLD

di

numbers of quadrature points dimensional numerical integrations

deterministic function with respect to Y
P yj1 , yj2 , . . . , yjs−i , 0 a summation of terms that contain at most s − i variables

pf

failure probability

s1

coodinate of the second types of sigma points in HUT

s2

coodinate of the third types of sigma points in HUT

W

auxiliary random variable

Xi (i = 1, 2, . . . , n) basic random variables

in

Yi (i = 1, 2, . . . , n) standard normal random variables Z

performance function

0

zero vector

u

vector containing ρ and ε

X

vector containing n basic random variables

Y

vector containing n independent standard normal random variables

FZ (z) cumulative distribution function of Z

i-

fY (y) the joint probability density functions of Y

fX (x) the joint probability density functions of X

fYj1 (yj1 ) the probability density function of Yj1

G deterministic function  k  I P (yj1 , yj2 , yj3 , 0) three-dimensional Gaussianweighted integrals  k  I P (yj1 , yj2 , 0) two-dimensional Gaussianweighted integrals   I P k (yj1 ), 0 one-dimensional Gaussian-weighted integrals  k I Z kth-order raw moments of Z

fZ (z) the probability density function of Z T −1 (·) inverse Nataf transformation ω1,l

weight for the one-dimensional numerical integrations

ω2,l

weight for the two-dimensional numerical integrations

ω3,l

weight for the three-dimensional numerical integrations

cases since the explicit expression of the performance function is generally not available, especially when the nonlinear behavior needs to be addressed [3]. Various approximation techniques have been developed to circumvent this difficulty for calculating the failure probability. The surrogate models are now quite popular in reliability community, where the response surface method [4–6], kriging method [7–9], support vector machine [10, 11], and polynomial chaos expansion [12–15] have been extensively investigated. Nevertheless, quantifying the errors of surrogate models for structural reliability analysis could be still an open challenge. The first-order and second-order reliability method 2

(FORM/SORM) [16, 17] have been gradually developed to be mature for structural reliability analysis. However, when the nonlinearity of the performance function around the design point is relatively high, the errors of using FORM/SORM for structural reliability analysis could be relatively significant. Besides, since the gradient information and iterations are required in FORM/SORM, these intrusive methods could be computationally inefficient for complex structures since obtaining the gradient information of a complex structure is not an easy task. Alternatively, the method of moments [18], which dose not require the gradient information and iterations, is a non-intrusive method, and mainly includes the following two steps for structural reliability analysis. First, the first-four central moments of the performance function are computed. Then, the performance function’s PDF is recovered based on the first-four central moments, where some parametric PDF models could be employed, e.g., the Pearson’s system. In that regard, the failure probability could be straightforwardly evaluated by integrating the performance function’ PDF over the failure domain. Since the evaluation of the first-four central moments involves repeatedly deterministic model evaluations, the computational effort of the method of moments is mainly determined by the numerical evaluation of the first-four central moments, where the computational time of deriving the performance function’s PDF could be negligible. In this regard, the approaches for numerical evaluation of the first-four central moments are of paramount importance to the accuracy and efficiency for structural reliability analysis based on the method of moments. Usually, the statistical moments of the performance function can be numerically evaluated by two kinds of methods. The first kind is the so-called point estimate method, where the statistical moments are evaluated by the weighted summation of values of the performance function at a finite number of points of input random variables. The earliest point estimate methods include the Rosenblueth’s method [19, 20], Harr’s method[21], Chang’s method[22], etc.. However, these methods may not give sufficient precision for statistical moments assessment of nonlinear systems. In fact, the original random variables can be converted into the independent standard normal random variables through the non-normal to normal transformation [23] e.g., the inverse Rosenblatt transformation or inverse Nataf transformation. In that regard, the Gaussian-weighted integrals (GWIs) are of great concern for statistical moments assessment of the performance function. In the last decade, many numerical methods have been developed for the GWIs, which could be also regarded as the point estimate method, e.g., the random configuration method [24], sparse grid numerical integration [25], high-order unscented transformation (HUT) [26], 5th-degree cubature formula [27], sparse grid integration method of spreading [28], etc.. Nevertheless, the tradeoff of accuracy and efficiency for statistical moments assessment may not be guaranteed by using those methods above, especially for the high-order central moments assessment. The second kind is the so-called dimension-reduction method (DRM) [29–31], which decomposes a multidimensional moment integral by the summation of a series of low-dimensional integrals. The univariate dimensionreduction method (UDRM) [29] is very efficient since the computational effort only grows linearly with the dimension of the problem. The multiplicative UDRM (M-UDRM) is also developed to further enhance the accuracy for fractional moments assessment in the maximum entropy method [32], which has been further applied to complex structural models where the deterministic analysis could be time-consuming [33, 34]. However, the UDRM and M-UDRM may still provide poor accuracy for estimating the high-order central moments, especially for nonlinear problems [35], as will be shown in Example 1. On the other hand, the multivariate DRM, e.g., the bivariate or the trivariate DRM is more accurate and theoretically appropriate. However, a multivariate DRM may not be able to appropriately keep the efficiency [36]. Therefore, it is necessary to reduce the computational effort of a multivariate DRM without losing accuracy. A criterion to judge the cross terms among variables in bivariate DRM (BDRM) is proposed to reduce unnecessary two-dimensional integrals [36] to enhance the computational efficiency of BDRM for statistical moments assessment. However, for nonlinear problems, the accuracy of BDRM for estimating the high-order central moments may not be sufficient. The similar criterion for determining cross terms is further applied in trivariate DRM (TDRM) with high accuracy. Unfortunately, this method can be only applied for very low-dimensional cases, e.g., the number of random variables is less than 3 [37]. In that regard, further improving the efficiency of TDRM without losing accuracy could be a plausible way to generally evaluate the central moments of the performance function. The objective of this paper is to develop an adaptive TDRM (ATDRM) for estimating the first-four central moments and structural reliability with the tradeoff of accuracy and efficiency. The rest of the paper is organized as follows. In Section 2, the problem formulation is presented, where the raw moments of the performance function is estimated by TDRM. Section 3 devotes to establishing the ATDRM with delineating the existence of cross terms and high-order unscented transformation. In that regard, the first-four central moments can be estimated by ATDRM with accuracy and efficiency. The shifted generalized lognormal distribution (SGLD) model is fitted based on the first-four central moments to recover the probability density function (PDF) of the performance function in Section 4, where the failure 3

probability can be evaluated. In Section 5, numerical examples are investigated to validate the accuracy and efficiency of the proposed method. The final section contains the concluding remarks. 2. Problem formulation The first-four central moments of the performance function, denoted as µZ , σZ , α3Z , α4Z , can be expressed as [35] µZ = I [Z] (3)   2 σZ = I Z 2 − µZ 2

(4)

      3 α3Z = I Z 3 − 3I Z 2 + 2µZ 2 σZ

(5)

        4 α4Z = I Z 4 − 4I Z 3 µZ + 6I Z 2 µZ 2 − 3µZ 4 σZ

where

    I Z k = I Gk (X) =

Z

(6)

+∞

−∞

Gk (X) · fX (x) dx, k = 1, 2, 3, 4

(7)

denote the raw moments of Z and the symbol I denotes the integration for short. Further, Eq. (7) can be transformed to be the GWIs such that   
 I Z k = I Gk T −1 (y) Z +∞
= Gk T −1 (y) · fY (y) dy

(8)

−∞ +∞

=

Z

−∞

P k (y) · fY (y) dy, k = 1, 2, 3, 4

where, T −1 (·) represents the inverse Rosenblatt transformation or inverse Nataf transformation ,Y = [Y1 , Y2 , . . . , Yn ]  T  n Q 1 is a vector containing n independent standard normal random variables, fY (y) = fYi (yi ) = √2π − y2y n exp ( ) i=1 is the joint PDF of Y, and P (·) is a deterministic function with respect to Y. Since the direct integration of the multi-dimensional GWIs in Eq. (8) could not be an easy task, the dimension-reduction method (DRM) can be employed to alleviate the difficulty for central moments assessment of the performance function. The generalized DRM for decomposing the response function that contains n random variables can be expressed as P (y) =

s X i=0

(−1)

i



n−s+i−1 i



×

X

j1

P yj1 , yj2 , . . . , yjs−i , 0

(9)


where s is the number of variables in DRM, P yj1 , yj2 , . . . , yjs−i , 0 denotes a summation of terms that contain at most s − i variables, which can be expressed as

P yj1 , yj2 , . . . , yjs−i 0 = P 0, . . . , 0, yj1 , 0, . . . , 0, yj2 , 0, . . . , 0, yjs−i , 0, . . . , 0 (10) When s = 3 , the trivariate dimension-reduction method (TDRM) could be established, and Eq. (9) could be reduced as [30] P (y) =

X

j1
P (yj1 , yj2 , yj3 , 0) − (n − 3)

(n − 3)(n − 2)(n − 1) − P (0) 6

X

P (yj1 , yj2 , 0) +

j1
4

n (n − 3)(n − 2) X P (yj1 , 0) 2 j =1 1

(11)

where 0 is the zero vector . Substituting Eq. (11) into Eq. (8), we have X X       I P k (yj1 , yj2 , yj3 , 0) − (n − 3) I Zk ≈ I P k (yj1 , yj2 , 0) j1
j1
n  (n − 3) (n − 2) (n − 1) k (n − 3) (n − 2) X  k + I P (yj1 ) , 0 − P (0) , k = 1, 2, 3, 4 2 6 j =1

(12)

1



     in which I P k (yj1 , yj2 , yj3 , 0) , I P k (yj1 , yj2 , 0) and I P k (yj1 ), 0 denote the three-, two- and one-dimensional GWIs respectively, i.e., Z +∞ Z +∞ Z +∞   I P k (yj1 , yj2 , yj3 ), 0 = P k (yj1 , yj2 , yj3 , 0) · fYj1 (yj1 ) · fYj2 (yj2 ) · fYj3 (yj3 ) dyj1 dyj2 dyj3 −∞

−∞

  I P k (yj1 , yj2 , 0) =

−∞

Z

(13)

+∞ Z

−∞

+∞

−∞

  I P k (yj1 , 0) =

P k (yj1 , yj2 , 0) · fYj1 (yj1 ) · fYj2 (yj2 ) dyj1 dyj2

Z

(14)

+∞

−∞

P k (yj1 , 0) · fYj1 (yj1 ) dyj1

(15)

It can be seen from Eqs. (12)-(15), the TDRM can additively decompose an n-dimensional GWI into a total of Cn3 = n(n−1)(n−2) three-dimensional GWIs , Cn2 = n(n−1) two-dimensional GWIs and n one-dimensional GWIs, 6 2 respectively. Numerical integrations are usually employed to approximate the GWIs above. If the numbers of integration points required for the three-, two- and one-dimensional numerical integrations are denoted as d3 , d2 and d1 , we have d3  k  X I P (yj1 , yj2 , yj3 ) , 0 = ω3,l P k (yj1 ,l , yj2 ,l , yj3 ,l , 0)

(16)

l=1

d2   X I P k (yj1 , yj2 , 0) = ω2,l P k (yj1 ,l , yj2 ,l , 0)

(17)

l=1

d1   X I P k (yj1 , 0) = ω1,l P k (yj1 ,l , 0)

(18)

l=1

where ω3,l , ω2,l , ω1,l denote the weights for the three-, two- and one-dimensional numerical integrations; then the total number of deterministic model evaluations for assessing the raw moments (Eq. (12)) of the performance function, denoted as N , can be expressed as N=

n (n − 1) (n − 2) n(n − 1) · d3 + · d2 + n · d1 + 1 6 2

(19)

In general, the Gauss-Hermite quadratures (GHQs) are often used to numerically solve the low-dimensional GWIs (1 ≤ s ≤ 3, s ∈ N) involved in Eqs. (16)-(18), where the numbers of quadrature points are d1 = d, d2 = d2 and d3 = d3 . Since the model evaluation on the origin point 0 is repeatedly required in all these low-dimensional numerical integrations (1 ≤ s ≤ 3, s ∈ N), the total number of model evaluations can be rewritten as N=

n (n − 1) (n − 2) n(n − 1) 3 2 · (d − 1) + · (d − 1) + n · (d − 1) + 1 6 2

(20)

Originally, the seven-point GHQ (d = 7), denoted as GHQ7, is usually adopted in TDRM and the total number of points required in TDRM is N = 36n3 − 90n2 + 60n + 1 (21) 5

It is seen that the number of points in the original TDRM cubically grows with the dimension, which indicates that the prohibitively large computational effort may be required for statistical moments assessment by TDRM. It is also noted that the accuracy for calculating the statistical moment of the performance function with TDRM is mainly controlled by the accuracy of three- and two-dimensional numerical integrations. Besides, the required numbers of three- and two-dimensional numerical integrations are much larger than that of one-dimensional numerical integration. In that regard, even if the error of each three- and two-dimensional GWI is small, the cumulative error may be quite large in TDRM for statistical moments assessment. As a matter of fact, the computational effort of TDRM for statistical moments assessment of the performance function mainly depends on two aspects: (1) the numbers of three- and two-dimensional GWIs (see Eqs. (19) and (20)); (2). the numbers of integration points for three- and two-dimensional numerical integrations. If the numbers of three- and two-dimensional GWIs and the numbers of points required for three- and two-dimensional numerical integrations can be both decreased, the total computational burden of TDRM could be significantly reduced. In the subsequent section, an adaptive TDRM (ATDRM) will be established based on these two aspects for statistical moments assessment of the performance function, which can significantly reduce the computational time without losing accuracy. 3. An adaptive trivariate dimension-reduction method for statistical moments assessment From Eq. (11), one can note that the response function can be decomposed into several low-dimensional functions(1 ≤ s ≤ 3, s ∈ N), where the bivariate and trivariate cross terms could be involved. However, not all cross terms are necessary for the contribution to statistical moments assessment. In other words, some bivariate and trivariate cross terms in Eq. (11) could be unnecessary [36], indicating the numbers of three- and two-dimensional GWIs could be reduced when some bivariate and trivariate cross terms are detected as unnecessary. In that regard, the computational effort of TDRM for statistical moments assessment could be also reduced since a small number of three- and twodimensional numerical integrations need to be performed. Hence, the existence of bivariate and trivariate cross terms is first delineated in the performance function. 3.1. Delineating the existence of cross terms According to Ref. [36], the criterion for delineating the existence of cross term of Yj1 and Yj2 is illustrated. To further ensure the non-existence of cross terms, the criterion must be unanimously satisfied by two different sets of samples. The practical procedure for delineating the existence of cross term of Yj1 and Yj2 is explained as follows: Step 1: Select an appropriate reference point of yc = (yj1 ,c, yj2 ,c , 0), and evaluate P (yc ).  (1)



Step 2: Determine another three points, namely yj1 =  (1) (1) yj1 , yj2 , 0 . Then, calculate the following equation ∆P (1)

(1)

(1)

yj1 , yj2 ,c , 0 , yj2 =

(1)

(1)

yj1 ,c , yj2 , 0 , and yj1 j2 =

h   i h    i (1) (1) (1) P yj1 − P (yc ) − P yj1 j2 − P yj2 n h   i h    i o = (1) (1) (1) min P yj1 − P (yc ) , P yj1 j2 − P yj2

(22)

If ∆P (1) > δ, there must exist cross term of Yj1 and Yj2 ; otherwise, turn to Step 3. Herein, δ = 0.01 is specifically considered.       (2) (2) (2) (2) (2) (2) (2) Step 3: Select yj1 = yj1 , yj2 ,c , 0 , yj2 = yj1 ,c , yj2 , 0 , yj1 j2 = yj1 , yj2 , 0 , and then similarly compute h   i h    i (2) (2) (2) P yj1 − P (yc ) − P yj1 j2 − P yj2 (2) n h   i h    i o ∆P = (23) (2) (2) (2) min P yj1 − P (yc ) , P yj1 j2 − P yj2 If ∆P (2) > δ, there must exist cross term of Yj1 and Yj2 ; otherwise, the cross term of Yj1 and Yj2 do not exist. The rigorous proof of delineating the existence of cross terms can be found in Ref. [36]

6

According to the existence of cross terms of Yj1 and Yj2 , we define an indication function for the bivariate cross term of Yj1 and Yj2 as  1, if the cross terms of Yj1 and Yj2 exist Jp(2) (Yj1 , Yj2 ) = (24) 0, otherwise Then, define the indication function for the trivariate cross term of Yj1 , Yj2 , Yj3 , j1 < j2 < j3 , j1 , j2 , j3 ∈ {1, 2, . . . , n},namely Jp(3) (Yj1 , Yj2 , Yj3 ) = Jp(2) (Yj1 , Yj2 ) · Jp(2) (Yj1 , Yj3 ) · Jp(2) (Yj2 , Yj3 ) and we have Jp(3)

(Yj1 , Yj2 , Yj3 ) =



1, if the cross terms among Yj1 , Yj2 , Yj3 exist 0, otherwise

(25)

(26)

3.2. Adaptive trivariate dimension-reduction method (ATDRM) Based on delineating the existence of cross terms (bivariate and trivariate cross terms) above, an ATDRM can be established for statistical moments assessment of the performance function with the balance between accuracy and efficiency. As noted in Eq. (12), the univariate, bivariate and trivariate GWIs are involved in TDRM. Since the computational time of the univariate GWIs just grows linearly with the dimension, the bivariate and trivariate GWIs actually control the efficiency and accuracy of TDRM. If the numbers of the bivariate and trivariate GWIs can be reduced, the total computational effort of TDRM may be greatly eased. To this end, delineating the existence of cross terms (bivariate and trivariate terms), which indicates the bivariate and trivariate GWIs could be unnecessary if the cross terms do not exist, could be employed to alleviate the computational burden of TDRM. (2) For a bivariate GWI, if Jp (Yj1 , Yj2 ) = 1, which means the bivariate cross term of Yj1 and Yj2 exists, then Eq. (2) (17) is employed for the integration (two-dimensional numerical integration). However, when Jp (Yj1 , Yj2 ) = 0, which reveals the bivariate cross term of Yj1 and Yj2 does not exist, the bivariate function P (yj1 , yj2 , 0) could be further simplified by the univariate dimension-reduction such that (Eq. (9) when s = 1) P (yj1 , yj2 , 0) = P (yj1 , 0) + P (yj2 , 0) − P (0)

(27)

h i   k I P k (yj1 , yj2 , 0) = I (P (yj1 , 0) + P (yj2 , 0) − P (0)) X k! q q q = I [(P (yj1 , 0)) 1 ] · I [(P (yj2 , 0)) 2 ] · I [(−P (0)) 3 ] q1 !q2 !q3 !

(28)

Then, we have

q1 +q2 +q3 =k

which only involves one-dimensional GWIs. (3) For a trivariate GWI, if Jp (Yj1 , Yj2 , Yj3 ) = 1 , which means the trivariate cross term among Yj1 , Yj2 , Yj3 exists, the three-dimensional numerical integration could be employed for this purpose, i.e., Eq. (16). On the other hand, if (3) (2) (2) (2) Jp (Yj1 , Yj2 , Yj3 ) = 0, which means Jp (Yj1 , Yi2 ) · Jp (Yj1 , Yj3 ) · Jp (Yi2 , Yi3 ) = 0, at least one of the indication (2) (2) (2) functions Jp (Yj1 , Yi2 ), Jp (Yj1 , Yj3 ), Jp (Yi2 , Yi3 ) equals to zero. Then, according to Eqs. (8)-(9), we can expand k the trivariate function P (yj1 , yj2 , yj3 , 0) by using the BDRM such that P k (yj1 , yj2 , yj3 , 0) ≈

X

j1
P k (yj1 , yj2 , 0) − (n − 2)

n X

j1 =1

P k (yj1 , 0) +

(n − 2) (n − 1) k P (0) 2

(29)

where n = 3, Eq. (29) can be further expressed as P k (yj1 , yj2 , yj3 , 0) ≈P k (yj1 , yj2 , 0) + P k (yj1 , yj3 , 0) + P k (yj2 , yj3 , 0)

− P k (yj1 , 0) − P k (yj2 , 0) − P k (yj3 , 0) + P k (0) 7

(30)

  In this regard, the trivariate GWI, i.e., I P k (yj1 , yj2 , yj3 0) can be estimated by         I P k (yj1 , yj2 , yj3 , 0) =I P k (yj1 , yj2 , 0) + I P k (yj1 , yj3 , 0) + I P k (yj2 , yj3 , 0)       − I P k (yj1 , 0) − I P k (yj2 , 0) − I P k (yj3 , 0) + P k (0) , k = 1, 2, 3, 4

(31)

It can be seen that there only exist two- and one-dimensional GWIs, in which the corresponding numerical integrations can be employed. Then, the ATDRM can be established for statistical moments assessment such that (1). First, evaluate the one-dimensional GWIs involved in Eq. (12) by Eq. (18); (2) (2). Delineate the existence of the bivariate cross terms. When Jp (Yj1 , Yj2 ) = 1, the two-dimensional numerical integration (Eq. (17)) is employed to evaluate the bivariate GWIs involved in TDRM, otherwise, use Eq. (28) to approximate those integrals, where only the one-dimensional numerical integration can be applied. (3) (3). Delineate the existence of the trivariate cross terms. When Jp (Yj1 , Yj2 , Yj3 ) = 1, the three-dimensional numerical integration (Eq. (16)) is employed to evaluate the trivariate GWIs involved in TDRM, otherwise, use Eq. (31) to approximate those integrals, where only the two-dimensional and one-dimensional GWIs can be applied. Since the bivariate cross terms have been previously delineated in Step (2), the results can be directly utilized in Eq. (31) to get the trivariate GWIs when the trivariate cross terms do not exist. The difference between of the proposed method and TDRM in Ref. [37] is the way for calculating a trivariate GWI (3) without existence of cross term. For the TDRM in Ref. [37], if Jp (Yj1 , Yj2 , Yj3 ) = 0, at least one of the indication (2) function equals to zero. Assume that Jp (Yj1 , Yj2 ) = 0, Eq. (13) can be further simplified as [37]   I P k (yj1 , yj2 , yj3 , 0) =

X

q1 +q2 +q3 =k

k! q q q I [(P (yj2 , yj3 , 0)) 1 ] · I [(P (yj1 , yj3 , 0)) 2 ] · I [(−P (yj3 , 0)) 3 ] q1 !q2 !q3 ! (32)

Then, the GHQ is then used to compute the two- and one-dimensional numerical integrations involved in Eq. (32). Herein, a simple example is employed to illustrate why the TDRM in Ref. [37] is not accurate for statistical moments assessment. We consider a simple performance function with an explicit expression such that Z = P (Y) = (2) Y1 · Y2 + Y3 , where Yi (i = 1, 2, 3) are independent standard normal variables. It’s easy to know Jp (Y1 , Y2 ) 6= 0, (2) (2) (3) Jp (Y1 , Y3 ) = 0, Jp (Y2 , Y3 ) = 0 and Jp (Y1 , Y2 , Y3 ) = 0. When k = 2, according to Eq. (32), the result of 2 I P (y1 , y2 , y3 ) in Ref. [37] could be expressed as         I P 2 (y1 , y2 , y3 ) =I P 2 (0, y2 , y3 ,) + I P 2 (y1 , 0, y3 ) + I P 2 (0, 0, y3 ) + 2I [P (0, y2 , y3 )] · I [P (y1 , 0, y3 )] − 2I [P (0, y2 , y3 )] · I [P (0, 0, y3 )]

− 2I [P (y1 , 0, y3 )] · I [P (0, 0, y3 )]   =3I Y32 − 2I 2 [Y3 ] = 3

where I[Y32 ] = 1 and I 2 [Y3 ] = 0.   However, the exact result of I P 2 (y1 , y2 , y3 ) gives   I P k (y1 , y2 , y3 ) = Ih[(Y1 · Y2 +i Y3 ) · (Y1 · Y2 + Y3 )] h i 2 2 = I (Y1 · Y2 ) + 2I [Y1 · Y2 · Y3 ] + I (Y3 ) i h i h 2 2 = I (Y1 · Y2 ) + I (Y3 ) =2

(33)

(34)

i h 2 where I (Y1 · Y2 ) = 1 and I [Y1 · Y2 · Y3 ] = 0. By comparing Eq. (33) and Eq. (34), we can be find that they are not identical with each other, which demonstrates Eq. (33) could be incorrect. If Eq. (33) is incoperated into the ATDRM, it could be anticipated that inaccurate results for statistical moments could be computed. On the contrary, if the proposed method is employed, we can easily calculate   the result of I P 2 (y1 , y2 , y3 ) from Eq. (31) such that h i h i   2 2 I P k (y1 , y2 , y3 ) = I (Y1 · Y2 ) + I (Y3 ) = 2 (35) 8

which is identical with the exact result. In that regard, when Eq. (31) is employed in the proposed method, accurate statistical moments could be computed. Compared with the original TDRM, the numbers of three- and two-dimensional numerical integrations could be significantly reduced if the cross terms do not exist. The worst case is that the cross terms always exist, which means (3) (2) Jp (Yj1 , Yj2 , Yj3 ) = 1 and Jp (Yj1 , Yj2 ) = 1 always hold. However, the worst case actually dose not occur in practical cases, where delineating the existence of cross terms in TDRM adaptively is of great necessity. In that regard, the computational efficiency of ATDRM could be significantly improved compared with the original one since only the lower-dimensional GWIs are involved. To further reduce the computational effort, efficient numerical integrations for the three- and two-dimensional GWIs involved in ATDRM need to be employed. As mentioned, if the GHQs are used to evaluate these numerical integrations, where d1 = d, d2 = d2 and d3 = d3 , the computational burden still exists especially when the number of trivariate cross terms is relatively large. Alternatively, an efficient scheme needs to be employed to substitute the GHQ for numerical integrations involved in ATDRM to further enhance the computational efficiency without losing accuracy. In this regard, a high-order unscented transformation will be introduced to fulfill this aim. 3.3. High-order unscented transformation Three- and two-dimensional GWIs involved in ATDRM can be calculated by using the high-order unscented transformation (HUT) [38]. The central idea of HUT could be divided into three parts: first, select some sigma points and weights to match the probability information of input random variables; second, obtain the corresponding sample points of output random variable by means of the nonlinear transformation of sigma points; finally, the statistical moments of output random variable can be estimated by the weights of input variables and sample points of output random variable [39]. Three types of symmetrically located sigma points and their weights are involved in HUT to match the first-four principal moments at least for the independent standard normal vector. The first type of sigma point is located at the origin and the weight is denoted as A0 . The second type of sigma points are symmetrically distributed on the axes with the distance s1 from the origin and the weight of each point is denoted as A1 . In that regard, a total of 2n sigma points are involved in the second type. The third type of sigma points with the weight A2 are located on (0, · · · , ±s2 , · · · , ±s2 , · · · , 0), where only two elements are nonzero. It is noted that a total of five parameters, i.e., A0 , A1 , A2 , s1 and s2 , need to be specified. To accurately match the mean, standard deviation, skewness and kurtosis of the independent standard normal random vector, the following conditions should be satisfied:  A0 + 2nA1 + 2n(n − 1)A2 − 1 = 0    2A1 s21 + 4(n − 1)A2 s22 − 1 = 0 (36) 2A1 s41 + 4(n − 1)A2 s42 − 3 = 0    4A2 s42 − 3 = 0

Since only four equations are available, the sigma points and their weights cannot be solved from these equations, 1 where five parameters are involved. If we assume A2 = (n+γ) 2 , where γ is a free parameter, then the solution to Eq. (36) can be derived such that [38]   −2n2 +(4−2n)γ 2 +(4γ+4)n q    A0 =  s1 = (4−n)(n+γ) (n+γ)2 (4−n) (γ+2−n)2 q (γ+2−n) , (37) A1 = 2(n+γ) 2 (4−n)   s = (n+γ)  2  A2 = 1 2 2 (n+γ)

It is noted that the HUT possesses a total of N = 4Cn2 + 2n + 1 = 2n2 + 1 points for numerical integration. Specifically, when n = 2 and n = 3, the HUT contains 9 and 19 sigma points, respectively, which are far less than the GHQ points with the same degree of accuracy. According to computational experiences, γ = 7.2 is recommended in HUT to evaluate the two- and one-dimensional GWIs and and γ = 5.8 is adopted for evaluating the three-dimensional GWIs in this paper. The reason of using these values for γ is that the HUT can practically achieve the best accuracy for the low-dimensional Gaussian-weighted numerical integrations. Since only the low-dimensional GWIs are involved in ATDRM, the recommendations of these values for γ, which ensure the accuracy of three- and two-dimensional Gaussian-weighted numerical integrations, could also result in the accuracy of ATDRM for statistical moments 9

assessment in all cases. In that regard, the points and weights can be utilized in ATDRM to further enhance the computational efficiency for statistical moments assessment. Then, the first-four central moments can be evaluated by using the ATDRM (Eqs. (3)-(6)), accordingly. 3.4. Computational effort of the proposed ATDRM with HUT To recycle the points for statistical moments assessment, the points generated by HUT in Table 1 are also utilized for delineating the existence of cross terms in ATDRM. It could be observed that the origin point 0 is repeatedly used in both judging the existence of cross terms and evaluating the three-, two- and one-dimensional GWIs, where only one time of deterministic model evaluation at the origin point actually needs to be performed. Suppose the numbers of bivariate and trivariate cross terms are l2 and l3 , respectively, the total number of points in the proposed ATDRM with HUT gives N = Cn2 · (7 − 1) + n · (3 − 1) + l2 · (9 − 5) + l3 · (19 − 1) + 1 (38) = 18l3 + 4l2 + 3n2 − n + 1

In Eq. (38), Cn2 means the number of bivariate terms, where the cross terms need to be delineated by the seven HUT points in Table 1, and thus the total number of points required for judging the existence of cross terms is Cn2 · (7 − 1) excluding the origin point; the number of integration points for computing all one-dimensional GWIs with HUT is n · (3 − 1) without the origin point; the two-dimensional GWIs should be evaluated by HUT with 9 points, where 5 of these points could be recycled from judging the existence of cross terms in Table 1 and thus the number of integration points for all two-dimensional GWIs is determined as l2 · (9 − 5); the three-dimensional GWIs could be evaluated by HUT with 19 points, where the number of required integration points for all three-dimensional GWIs is l3 · (19 − 1) excluding the origin point; and the term “1” represents the 1 time of model evaluation over the origin point, which is repeatedly applied. On the other hand, if the seven-point GHQ (GHQ7) is employed in ATDRM, denoted as ATDRM with GHQ7, the numbers of integration points for calculating the one-, two- and three-dimensional GWIs are 7, 49 and 343, respectively. It should be mentioned that delineating the existence of cross terms in ATDRM with GHQ7 is also performed based on GHQ7, where the seven points are also listed in Table 1 and can be all recycled. Similarly, the number of ATDRM with GHQ7 could be expressed as N = Cn2 · (7 − 1) + n · (7 − 1) + l2 · (49 − 7) + l3 · (343 − 7) + 1

(39)

= 336l3 + 42l2 + 3n2 + 3n + 1

In the worst case, all cross terms exist, which means l2 = Cn2 and l3 = Cn3 , the total number of points in ATDRM with HUT is N = 3n3 − 4n2 + 3n + 1 and the total number of points for the ATDRM with GHQ7 is as large as N = 56n3 − 144n2 + 94n + 1. When all the cross terms don’t exist, i.e., l2 = 0 and l3 = 0, the total numbers of points in ATDRM with HUT and ATDRM with GHQ7 could be reduced to be N = 3n2 − n + 1 and N = 3n2 + 3n + 1. As mentioned in Eq. (21), the number of points required in the original TDRM with GHQ7 is N = 36n3 − 90n2 + 60n + 1. Since only a few of cross terms, especially for the trivariate cross terms, exist in practical cases, it could be known that the number of points required in the proposed ATDRM with HUT is much smaller than those of the ATDRM with GHQ7 and original TDRM with GHQ7. Table 1: Points for delineating the existence of cross terms of Yj1 and Yj2

Method

l

1

2

3

4

5

6

7

HUT

yj1 ,l yj2 ,l

0.00000 0.00000

1.59861 0.00000

0.00000 1.59861

1.59861 1.59861

- 1.59861 0.00000

0.00000 - 1.59861

- 1.59861 - 1.59861

GHQ7

yj1 ,l yj2 ,l

0.00000 0.00000

3.75044 0.00000

0.00000 3.75044

3.75044 3.75044

-3.75044 0.00000

0.00000 -3.75044

-3.75044 -3.75044

10

4. Evaluation of the probability density function of the performance function Since the first-four central moments are available by using the proposed ATDRM with HUT, to reconstruct the entire distribution range of PDF fZ (z) for reliability analysis, a flexible PDF model, called the shifted generalized lognormal distribution (SGLD) [40], where only four parameters need to be specified, is employed. The reasons of choosing the SGLD model to represent the unknown PDF for structural reliability analysis are as follows: [40] (1). the SGLD model, which has a rich flexibility in shape, is able to assimilate a plenty of unimodal distributions; (2). the SGLD model can reduce to be many commonly used distributions; (3). the involved parameters of the SGLD model can be easily computed from the statistical moments. The PDF represented by the SGLD model can be expressed as    ε  z − b α 1 (40) fZ (z) = exp − ε ln , b < z < +∞ z−b ερ τ   where α = 1/ 2ε1/ε ρΓ(1 + 1/ε) , Γ (·) is the gamma function, b and τ are the parameters for the location and scale of the SGLD model, and ρ and ε are the parameters for the shape of the SGLD model. The corresponding cumulative distribution function (CDF) can be written as : 
ε    ρ1 ln z−b τ 1 1 z−b 1  , b < z < +∞ (41) −1 g , FZ (z) = + sgn 2 2 τ ε ε

Rx where g(v, x) = 0 tv−1 e−t dt/Γ(v) . To solve the four undetermined parameters, i.e., b, τ , ρ and ε, an auxiliary random variable W = (Z − b)/τ is introduced. Then, the raw moments of W can be expressed as [41]   I Wk =

  +∞ 2m X (kρ) 2m + 1 1 ε2m/ε Γ Γ(1/ε) m=0 (2m)! ε

(42)

Once ρ and ε are known, the raw moments of W can be easily computed, and the corresponding first-four central moments can be computed accordingly, denoted as µW , σW , α3W and α4W . According to the probability theory, the following relations hold for the first-four central moments of Z and W :  µ = µZτ−b    W σW = στZ (43) α = α3Z    3W α4W = α4Z The scale and location parameters can be computed as τ = σZ /σW and b = µZ − τ µW . In that regard, the shape parameters ε and ρ can be solved by  α3W (u) − α3Z = 0 (44) α4W (u) − α4Z = 0 T

where u = [ρ, ε] . It should be pointed out that the absolute value of α3Z should be considered to compute the parameters. When α3Z < 0, the symbol z needs to be replaced by 2µZ − z on the right-hand side of Eq. (40) [40]. To promptly compute the solution to Eq. (44), the initial values of ρ and ε are suggested as [42]: s     − |α3Z | 2 1 −1 sinh , ε0 = 2 (45) ρ0 = ln 1 + 4sinh 3 2 Once the PDF of the performance function is recovered, the failure probability can be directly computed by a integral over the failure domain. Then, the reliability index β can be computed by 1 − pf = Φ (β), where Φ denotes the CDF of a standard normal random variable. 11

5. Computational procedure In summary, the statistical moments and reliability can be evaluated by using the proposed ATDRM with HUT, where the step-by-step procedure is explained as follows: Step 1: Delineate the existence of cross terms (bivariate and trivariate) by the HUT points in Table 1 and Eqs. (24)-(26); Step 2: Calculate the univariate, bivariate and trivariate GWIs involved in TDRM for raw moments assessment. First, calculate the univariate GWIs by Eq. (18) with the one-dimensional HUT. Second, if the bivariate cross terms exist, Eq. (17) with the two-dimensional HUT is adopted for evaluating the two-dimensional GWIs; otherwise, Eq. (28) is adopted to approximate the bivariate GWIs involved in TDRM and the one-dimensional GWIs in Eq. (28) can be also assessed by HUT. Then, compute the trivariate GWIs involved in TDRM. If the trivariate cross terms exist, the trivariate GWIs in TDRM (Eq. (16)) are assessed by using the three-dimensional HUT; otherwise, Eq. (31) is adopted to approximate the trivariate GWIs, where the two- and one-dimensional GWIs in Eq. (31) have been already computed. Step 3: Compute the first-four central moments from the raw moments evaluated by the proposed ATDRM with HUT. It is noted that all the three-, two- and one-dimensional GWIs have been computed from Step 2, and then Eq.(12) can be used to calculate the first-four raw moments. Hence, the first-four central moments can be transformed from the raw moments by Eqs. (3)-(6). Step 4: Reconstruct the performance function’s PDF based on the first-four central moments by fitting the SGLD model. First, the first-four central moments of W are calculated according to the initial values of ρ and ε, as shown in Fig.1. Second, according to the relations between the first-four central moments of W and Z, the shape parameters ρ and ε can be computed by solving the nonliear equations (Eq. (44)) via the Levenberg-Marquardt algorithm in Matlab. The scale and location parameters τ and b can be easily determined by the means and standard deviations of W and Z (Eq. (43)). Finally, the performance function’s PDF can be recovered by substituting those parameters into the SGLD model (Eq. (40)). Step 5: Calculate the failure probability and reliability index based on the recovered PDF of the performance function. 6. Numerical Examples Three numerical examples are studied in this section to verify the efficacy of the proposed method for statistical moments assessment and reliability analysis. The pertinent MCS is employed to provide the reference results for the moments, PDF and failure probability for comparisons. Besides, the results of statistical moments by some other approaches are also presented to demonstrate the accuracy and efficency of the proposed method. 6.1. Example 1 The moment bearing capacity of a reinforced concrete section, shown in Fig. 2, is first considered [35], where the performance function could be expressed as: Z = G (X) = Afy d −

A2 f y 2 K −M fc b

(46)

where A is the area of reinforcement; fy is the yielding stress of reinforcement; d is the effective depth of reinforcement; K is the stress–strain factor of concrete; fc is the compressive strength of concrete; b is the width of section; and M is the applied moment. These parameters above are all considered as independent random variables, where the statistical information is provided in Table 2[35]. Hereafter, the abbreviation “C.O.V.” denotes the coefficient of variation for short in Tables. Table 3 presents the first-four central moments of the performance function evaluated by different approaches for comparisons. The reference results provided by MCS (107 runs) are 180.3729, 67.8082,0.5922 and 3.7323, respectively. The results computed from the UDRM, M-UDRM, BDRM with the five-point GHQ (denoted as BDRM with GHQ5), BDRM with the three-point GHQ (denoted as BDRM with GHQ3) and BDRM with HUT are all taken from Ref. [35] (Section 5.3, Table 6). It is noted that the UDRM and M-UDRM gives spurious results for 12

start Input  0 ,  0

k 1 m 1 E  m  1  E  m  2 

1 k   2m  1   2m      1/    2m !    2m

2 m1

 2  m  1  1  1 k   2 m1      1/    2  m  1 !   

E  m  2   E  m  1  106 E  m  1

k  k 1

m  m 1

Yes

No

I W k     E  m  1  Yes

k4 No

Transform raw moments into central moments by Eqs. (7)-(10) Output W ,  W ,3W ,4W end Figure 1: Flow-chart of computing the first-four central moments of W Table 2: Random variables for Example 1

Variable

Distribution

Mean

C.O.V

A fy d K fc b M

Lognormal Lognormal Lognormal Normal Weibull Normal Lognormal

1260.00 mm 300.00 N/mm2 770.00 mm 0.35 25.00 N/mm2 200.00 mm 100.00 kN·m 2

0.15 0.15 0.10 0.25 0.10 0.10 0.20

the high-order moments, where the relative errors of skewness are up to 68.0528% and 33.1506%, respectively. The BDRM with GHQ5 and BDRM with GHQ3 improve the calculation accuracy, however, the relative errors of kurtosis are still quite large, says 7.0525% and 9.4880%. The BDRM with HUT can utilize a small number of samples to achieve fairly accurate results of statistical moments. Nevertheless, the maximum error of moments by this method is still much larger than that of the proposed method. Besides, the results computed from the original TDRM, TDRM in Ref. [37] and ATDRM with GHQ7 are also computed and listed in Table 3. It is seen that although the original TDRM gives very accurate results for the first-four central moments, a total of 8359 times of deterministic model evaluations 13

fc

d

Uncertainty shape (modeled by factor K)

fy

b a) Cross section

b) Ultimate stress distribution

Figure 2: Ultimate state for reinforced concrete section

need to be carried out. The TDRM in Ref. [37] can significantly reduce the number of points, however, it may also render totally incorrect results for statistical moments assessment, where the relative error of standard deviation is even up to 166.55%. Further, the ATDRM with GHQ7 also provides quite accurate results for the statistical moments while the computational effort is still very large for this 7-dimensional problem, where 4369 times of repeatedly deterministic model evaluations need to be implemented. The proposed ATDRM with HUT can achieve a good tradeoff of accuracy and efficiency for statistical moments assessment, where the maximum relative error is just 1.31% with only using 387 times of repeatedly deterministic model evaluations, which is far less than that of the ATDRM with GHQ7 and original TDRM. Table 3: Comparisons of statistical moments for Example 1

Method

N

UDRM

29

µZ (R.E.)

σZ (R.E.)

α3Z (R.E.)

α4Z (R.E.)

Reference

180.3610 66.4259 0.1892 3.0907 [35](Sec.5.3, Tab. (0.0066) (2.0385) (68.0528) (17.1910) M-UDRM 29 180.3641 66.7892 0.3959 3.3057 [35](Sec.5.3, Tab. (0.0049) (1.5028) (33.1506) (11.4305) BDRM with GHQ5 533 180.3772 67.7819 0.5755 3.4691 [35](Sec.5.3, Tab. (0.0024) (0.0388) (2.8244) (7.0525) BDRM with GHQ3 183 180.3775 67.7766 0.5714 3.3782 [35](Sec.5.3, Tab. (0.0025) (0.0466) (3.5167) (9.4880) BDRM with HUT 183 180.3756 67.7963 0.6070 3.9263 [35](Sec.5.3, Tab. (0.0189) (0.0158) (2.4992) (5.1978) Original TDRM 8359 180.3765 67.7946 0.5917 3.7337 [30] (0.0001) (0.0038) (0.2080) (0.0508) TDRM in Ref. [37] 4369 180.3765 180.6699 0.0132 0.7080 [37] (0.0181) (166.5507) (97.7577) (81.0065) ATDRM with GHQ7 4369 180.3765 67.7965 0.5967 3.7816 —(0.0181) (0.0233) (1.0515) (1.4423) ATDRM with HUT 387 180.3764 67.7928 0.5983 3.7582 —(0.0181) (0.0179) (1.3146) (0.8150) Note:N =number of repeatedly deterministic model evaluations; R.E. = relative error (%)

6) 6) 6) 6) 6)

Once the statistical moments evaluated by the proposed ATDRM with HUT are available, the performance function’s PDF can be reconstructed by using the SGLD model. Fig. 3 pictures the PDF and CDF in logarithmic scale by the proposed method, where the distribution histogram and CDF curve by MCS (107 trails) are also provided for comparisons. It is easy to observe that PDF by the proposed method is quite close to the histogram 14

by MCS. Besides, the CDF in logarithmic scale by the proposed method also agrees well with that by MCS, which demonstrates that the distribution tail can be also accurately recovered by the proposed method. Then, one can straightforwardly compute the failure probability by the proposed method by integrating the PDF over the failure domain. Or, the value of CDF at the abscissa 0 gives the failure probability. The proposed method and MCS yield the failure probabilities as 2.10 × 10−4 and 3.68 × 10−4 , which further demonstrates the proposed ATDRM with HUT can compute fairly accurate result with reasonable computational effort for structural reliability assessment. 7

10-3

100 MCS Proposed method

6

MCS Proposed method

10-2

4

CDF

PDF

5

3

10-4

2 1 0 -200

0

200

400

600

10-6 -200

800

0

200

z (kN m) (a) PDF

400

600

800

z (kN m) (b) CDF

Figure 3: Recovered PDF and CDF in Example 1.

6.2. Example 2 In this example, an annular section column under axial compressive load is investigated, which is shown in Fig. 4[43]. It is known that the buckling failure could occur since this column is quite slender, where the performance function is given as Z = G (X) =

i π2 E π h 4 4 · (D + T ) − D −P L2 64

(47)

where E is the elastic modulus of material; D is the internal diameter; T is the thickness; L is the height of column and P is the axial load. The distribution information of basic random variables is provided in Table 4. Table 4: Random variables for Example 2

Variable

Distribution

Mean

Standard deviation

E D T L P

Lognormal Normal Normal Normal Lognormal

203.00 GPa 30.00 mm 6.00 mm 2500.00 mm 6000.00 N

40.00 GPa 1.00 mm 0.40 mm 90.00 mm 400.00 N

Like Example 1, the first-four central moments of the performance function computed from different approaches are shown in Table 5. The MCS (107 traills) is also employed to provide the “exact” results for comparisons, which are 8270.6973, 3457.0031, 0.7095 and 3.9214, respectively. In this example, the results computed from the UDRM, BDRM and the numerical integration based BDRM (NI-BDRM) are taken from Ref. [43] (Section 4.2, Table 9). It is clear that the UDRM still cannot provide accuate high-order statistical moments although it is highly efficient. Although the BDRM uses a large number of samples, say 442, the relative error of kurtosis is still quite large, say 15

P

D

A

A T

A A

Figure 4: An annular section column.

7.49%. The NI-BDRM can yield better results than the BDRM with the same computational effort. Nevertheless, it is seen the proposed ATDRM with HUT can achieve the comparable results with much smaller sample size. Similarly, the original TDRM gives accurate results with a large amount of computational effort and the TDRM in Ref. [37] could not yield accurate results, where the relative error of standard deviation is as large as 125.7830%. Though the ATDRM with GHQ7 can significantly improve the accuracy for statistical moments assessment, the computational effort of this method is still unacceptable for practical applications, where 1687 times of model evaluations are necessary. It is again noted that the proposed ATDRM with HUT can utilize a small number of samples (only 167) to compute the very accurate statistical moments, where the maximum relative error of moments is just 2.32%. Table 5: Comparisons of statistical moments for Example 2

Method UDRM

N 31

BDRM

391

NI-BDRM

391

Original TDRM

2551

TDRM in Ref. [37]

1687

ATDRM with GHQ7

1687

ATDRM with HUT

167

µZ (R.E.) 8270.1952 (0.0061) 8270.4415 (0.0031) 8270.4415 (0.0031) 8270.4403 (0.0031) 8270.4403 (0.0031) 8270.4403 (0.0031) 8270.4438 (0.0031)

σZ (R.E.) 3407.9319 (1.4195) 3457.4557 (0.0131) 3457.4557 (0.0131) 3457.9035 (0.0260) 5194.5198 (50.2608) 3457.9123 (0.0263) 3457.8933 (0.0258)

α3Z (R.E.) 0.3094 (56.3888) 0.6921 (2.4521) 0.6998 (1.354) 0.7085 (0.1397) 0.1521 (78.5618) 0.709 (0.0698) 0.7058 (0.5087)

α4Z (R.E.) 3.2386 (17.4116) 3.6276 (7.4915) 3.8537 (1.7255) 3.9100 (0.2902) 1.6988 (56.6794) 3.9147 (0.1703) 3.8303 (2.3224)

Reference [43] (Sec. 4.2, Tab. 9) [43] (Sec. 4.2, Tab. 9) [43] (Sec. 4.2, Tab. 9) [30] [37] ——-

Then, the SGLD model is employed to rebuild the PDF of the performance function, where the first-four central moments evaluated by the proposed method are served as constraints. The PDF and CDF in logarithmic scale by the proposed method are plotted and compared with those by MCS (106 runs) in Fig. 5. Again, it is clear that both the PDF and CDF in logarithmic scale by the proposed method still have very good accordance with those by MCS. The failure probabilities calculated by the proposed method and MCS are 2.1217 × 10−4 and 3.6150 × 10−4 , respectively. In this regard, the effectiveness of the proposed method for statistical moments assessment and reliability analysis is validated through the analyses and comparisons above. 16

1.4

10-4

100 MCS Proposed method

1.2

MCS Proposed method

10-2

0.8

CDF

PDF

1

0.6

10-4

0.4 0.2 0 -0.5

0

0.5

1

1.5

2

10-6 -0.5

2.5

0

0.5

104

z (N) (a) PDF

1

z (N) (b) CDF

1.5

2

2.5 104

Figure 5: Recovered PDF and CDF in Example 2.

6.3. Example 3 The last example considers a practical three-bay seven-storey planar frame subjected to lateral loads, which could exhibits nonlinear behaviors and is shown in Fig. 6. The implicit performance function can be expressed as 6m

2.7m

6m



F6 3.6m

fc fu

3.6m

F5 u

c

F4



3.6m

Concrete 01

F3 3.6m

 b  Es

fy

F2 3.6m

Es



F1 4.2m

fy

Steel 01 3"25 550mm

4"20 550mm

4"20

3"22

4"20 550mm

300mm

Column

Beam

Figure 6: A nonlinear three-bay seven-storey frame structure

Z = G (X) = ub − |utop (X)|

(48)

where utop denotes the interstory drift of the top floor and ub is the deterministic threshold. Statistical information of basic independent random variables are given in Table 6, where the elastic modulus of rebar E is a nomarlly 17

distributed random variable and the external forces F1 − F6 are all lognormally distributed. Nonlinear finite-element analyses of the structure are carried out for structural reliability analysis by using the Opensees software. The uniaxial Concrete01 and Steel01 are adopted to model the nonlinear constitutive laws of concrete and steel bar, respectively. For confined concrete, we have the compressive strength fcc = 35 MPa, the strain at maximum strength εcc = 0.005, the crushing strength fcu = 25 MPa and the strain at crushing strength εcu = 0.02. For unconfined concrete, we have the compressive strength fc = 27 MPa, the strain at maximum strength εc = 0.002, the crushing strength fu = 0 MPa and the strain at crushing strength εu = 0.006. For steel bar, the yielding strength fy is 400 MPa, the initial elastic modulus Es is 200 GPa, and the strain-hardening ratio b is 0.007. R∞ Actually, we can define Z1 = |utop (X)| and the failure probability can be evaluated by pf = ub fZ1 (z1 )dz1 . In that regard, the statistical moments of Z1 need to be first evaluated. Table 6: Statistical characteristics of random variables for Example 4

Variable

E(GPa)

F1 (kN)

F2 (kN)

F3 (kN)

F4 (kN)

F5 (kN)

F6 (kN)

Mean C.O.V.

200.00 0.10

10.00 0.25

20.00 0.25

30.00 0.25

40.00 0.25

50.00 0.25

60.00 0.25

Table 7 shows the comparions of the statistical moments computed by different approaches. Since the timeconsuming nonlinear finite-element analyses of the structure are involved, the computational times of these methods are also reported in the table, where a computer with Inter Core i7-7700K CPU at 4.20 GHz and 8 GB of RAM is used. The results by MCS (106 trails) are served as the references, where the first-four central moments are 0.0573, 0.0090, 0.5853 and 3.6561, respectively. Besides, the total time comsumed by MCS is 1.2416 ×105 s (34.39 hours). Remarkably, the proposed ATDRM with HUT can still provide very accurate results for the statistical moments within a low computational effort, where the maximum relative error of moments is just 2.20% and only 161 times of repeatedly deterministic model evaluations need to be implemented. The UDRM still gives spurious results of the high-order central moments although it requires the smallest computational time. The BDRM and ATDRM with GHQ7 cannot provide precise results, where the maximum relative error of kurtosis is up to 9%. Although the TDRM enables to calculate very accurate statistical moments, the computational effort of this method is far larger than that of the proposed method. It can be seen that the proposed method only costs about 20s, while the required computational time of TDRM is about 919s when nonlinear finite-element analyses are involved. Besides, the TDRM in Ref. [37] is still unable to provide satisfactory results for the statistical moments assessment. Table 7: Comparisons of statistical moments for Example 3

Method

N

Time (sec)

µZ (R.E.)

σZ (R.E.)

α3Z (R.E.)

α4Z (R.E.)

Reference

UDRM

43

1.2745

799

85.0182

Original TDRM

8359

918.6564

TDRM in Ref. [37]

379

58.5264

ATDRM with GHQ7

379

41.1128

ATDRM with HUT

161

19.8482

0.0089 (1.7071) 0.0091 (0.3088) 0.0091 (0.3427) 0.0392 (332.9804) 0.0091 (0.2801) 0.0091 (0.3342)

0.3843 (34.2195) 0.5619 (3.9904) 0.5873 (0.3423) -0.0744 (112.7184) 0.5603 (4.2730) 0.5858 (0.0866)

3.2972 (9.7705) 3.3026 (9.6683) 3.6645 (0.2291) 0.2768 (92.4293) 3.3029 (9.6615) 3.5757 (2.1998)

[29]

BDRM

0.0572 (0.0231) 0.0572 (0.0574) 0.0572 (0.0582) 0.0572 (0.0591) 0.0572 (0.0591) 0.0572 (0.0590)

[30] [30] [37] ——-

Meanwhile, the SGLD model is fitted to recover the performance function’s PDF based on the evaluated statistical moments by the proposed method. Again, it is clear that the PDF and probability of exceedance (POE=1-CDF) in logarithmic scale by the proposed method still agree very well with those by MCS(106 runs), as shown in Fig. 7. When 18

the threshold ub is set to be 0.1 m, the failure probability calculated from the proposed method is 3.40 × 10−4 , which is fairly close to that computed from MCS, say 2.33 × 10−4 . 100

50 MCS Proposed method

MCS Proposed method -1

10

30

10-2

PDF

POE

40

20

10-3

10

10-4

0 0

0.02

0.04

0.06

0.08

0.1

0.12

0

0.02

0.04

z1 (m) (a) PDF

0.06

0.08

0.1

0.12

z1 (m) (b) POE

Figure 7: Recovered PDF and POE of Z1 in Example 4.

Further, some other state-of-the art methods, e.g., importance sampling and subset sampling are also employed to evaluate the failure probabilities for comparisons. The failure probabilities are listed in Table 8. It is clear that the importance sampling and subset simulation method can give fairly accurate failure probabilities, however, the computational efforts of these methods are also quite large, where 1265 and 3700 times of model evaluations need to be performed. Besides, the variation of the failure probability by subset simulation method is somewhat large, where the C.O.V. is 0.29. Table 8: Comparisons of failure probabilities for Example 4

Method

N

Failure probability

ATDRM with HUT Importance sampling Subset simulation

161 1265 3700

−4

3.40 × 10 3.10 × 10−4 3.04 × 10−4

C.O.V. – 0.10 0.29

7. Concluding remarks This paper proposes an adaptive trivariate dimension-reduction method (TDRM) with the high-order unscented transformation (HUT) for statistical moments assessment and structural reliability analysis. The statistical moments of the performance function are first approximated by using the TDRM, where trivariate, bivariate and univariate Gaussian-weighted integrals (GWIs) need to be evaluated. First, delineating the existence of cross terms is performed with the purpose of reducing the numbers of trivariate and bivariate GWIs in TDRM. When the cross terms exist, the three- and two-dimensional numerical integrations based on HUT are employed to compute the trivariate and bivariate GWIs involved in TDRM, otherwise, the trivariate and bivariate GWIs can be further decomposed as the lower-dimensional GWIs and the HUT is again employed for the lower-dimensional numerical integrations. For univariate GWIs in TDRM, the HUT can be directly employed to provide the one-dimensional numerical integrations. In this regard, an adaptive TDRM (ATDRM) with HUT can be established for the first-four central moments assessment of the performance function. Once the first-four central moments are computed, the performance function’s probability density function (PDF) can be reconstructed by the shifted generalized lognormal distribution (SGLD) and the failure probability and reliability index can be easy computed accordingly. Three numerical examples are investigated to illustrate the efficiency and accuracy of the proposed method and the following conclusions could be arrived: 19

(1). The proposed method first reduces the number of three- and two-dimensional GWIs by judging the cross terms, and then uses the HUT in a simple and effective way for numerical integrations. In that regard, the proposed ATDRM with HUT could keep the tradeoff of accuracy and efficiency for statistical moments assessment and reliability analysis. (2). The bivariate dimension-reduction method (BDRM) may not always result in fairly accurate results for statistical moments assessment. Besides, the computational effort of BDRM is still much larger than that of the proposed ATDRM with HUT. The BDRM with HUT could be more efficient than the proposed ATDRM with HUT for statistical moments assessment, however the accuracy of BDRM with HUT may not be competitive with that of the proposed one. (3). Although the original TDRM is always very accurate for statistical moments assessment, the prohibitively large computational burden makes this method intractable for practical applications. Spurious results could be computed if the TDRM in Ref. [37] is employed for statistical moments assessment. (4). If the ATDRM with the seven-point Gauss-Hermite quadrature (GHQ7) is employed, although the accuracy of statistical moments assessment could be satisfactory, the computational efficiency cannot be ensured since a large number of repeatedly deterministic model evaluations still need to be performed. (5). The univariate dimension-reduction method (UDRM) and multiplicative UDRM (M-UDRM) may not be capable of evaluating the high-order statistical moments of nonlinear performance functions with sufficient accuracy. In summary, the proposed ATDRM with HUT is an effective way to compute the statistical moments with accuracy and efficiency for low- and moderate- dimensional problems. However, the proposed ATDRM with HUT still has some shortcomings. It is seen that the proposed ATDRM could also encounter the “curse of dimensionality”, where the required number of model evaluations still quadratically grows with the dimension. That means the computational effort of this method could be prohibitively large for high-dimensional problems, e.g., n > 20. Besides, if the performance function involves too many trivariate terms, the computational burden also exists since judging the existence of trivariate terms and the numerical computations of all three-dimensional GWIs cost considerably large computational times. Further investigations will deal with these issues. Acknowledgments: The research reported in this paper is supported by the National Natural Science Foundation of China (No. 51978253) and the Fundamental Research Funds for the Central Universities (No. 531118040110). The anonymous reviewers and editors are highly appreciated for their constructive comments for improving the paper. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13]

O. Ditlevsen, H. O. Madsen, Structural Reliability Methods, volume 178, Wiley New York, 1996. P. Thoft-Cristensen, M. J. Baker, Structural Reliability Theory and Its Applications, Springer Science & Business Media, 2012. J. Li, J. Chen, Stochastic Dynamics of Structures, John Wiley & Sons, 2009. N. Roussouly, F. Petitjean, M. Salaun, A new adaptive response surface method for reliability analysis, Probabilistic Engineering Mechanics 32 (2013) 103–115. S. Goswami, S. Ghosh, S. Chakraborty, Reliability analysis of structures by iterative improved response surface method, Structural Safety 60 (2016) 56–66. A. Hadidi, B. F. Azar, A. Rafiee, Efficient response surface method for high-dimensional structural reliability analysis, Structural Safety 68 (2017) 15–27. I. Kaymaz, Application of kriging method to structural reliability problems, Structural Safety 27 (2005) 133–151. W. Jian, S. Zhili, Y. Qiang, L. Rui, Two accuracy measures of the kriging model for structural reliability analysis, Reliability Engineering & System Safety 167 (2017) 494–505. B. Echard, N. Gayton, M. Lemaire, Ak-mcs: an active learning reliability method combining kriging and monte carlo simulation, Structural Safety 33 (2011) 145–154. J.-M. Bourinet, F. Deheeger, M. Lemaire, Assessing small failure probabilities by combined subset simulation and support vector machines, Structural Safety 33 (2011) 343–353. H. Dai, B. Zhang, W. Wang, A multiwavelet support vector regression method for efficient reliability assessment, Reliability Engineering & System Safety 136 (2015) 132–139. G. Blatman, B. Sudret, Adaptive sparse polynomial chaos expansion based on least angle regression, Journal of Computational Physics 230 (2011) 2345–2367. G. Blatman, B. Sudret, An adaptive algorithm to build up sparse polynomial chaos expansions for stochastic finite element analysis, Probabilistic Engineering Mechanics 25 (2010) 183–197.

20

[14] J. Xu, F. Kong, A cubature collocation based sparse polynomial chaos expansion for efficient structural reliability analysis, Structural Safety 74 (2018) 24–31. [15] J. Xu, D. Wang, Structural reliability analysis based on polynomial chaos, voronoi cells and dimension reduction technique, Reliability Engineering & System Safety 185 (2019) 329–340. [16] X. Du, Unified uncertainty analysis by the first order reliability method, Journal of Mechanical Design 130 (2008) 091401. [17] A. D. Kiureghian, M. D. Stefano, Efficient algorithm for second-order reliability analysis, Journal of Engineering Mechanics 117 (1991) 2904–2923. [18] Y.-G. Zhao, T. Ono, Moment methods for structural reliability, Structural safety 23 (2001) 47–75. [19] E. Rosenblueth, Point estimates for probability moments, Proceedings of the National Academy of Sciences 72 (1975) 3812–3814. [20] E. Rosenblueth, Two-point estimates in probabilities, Applied Mathematical Modelling 5 (1981) 329–335. [21] M. E. Harr, Probabilistic estimates for multivariate analyses, Applied Mathematical Modelling 13 (1989) 313–318. [22] C.-H. Chang, J.-C. Yang, Y.-K. Tung, Uncertainty analysis by point estimate methods incorporating marginal distributions, Journal of Hydraulic Engineering 123 (1997) 244–250. [23] H. Li, Z. L¨u, X. Yuan, Nataf transformation based point estimate method, Chinese Science Bulletin 53 (2008) 2586. [24] P. Bressolette, M. Fogli, C. Chauvi`ere, A stochastic collocation method for large classes of mechanical problems with uncertain parameters, Probabilistic Engineering Mechanics 25 (2010) 255–270. [25] F. Xiong, Y. Xiong, S. Greene, W. Chen, S. Yang, A new sparse grid based method for uncertainty propagation, in: ASME 2009 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, American Society of Mechanical Engineers, 2009, pp. 1205–1215. [26] S. Xiao, Z. Lu, Reliability analysis by combining higher-order unscented transformation and fourth-moment method, ASCE-ASME Journal of Risk and Uncertainty in Engineering Systems, Part A: Civil Engineering 4 (2017) 04017034. [27] J. Xu, Z.-H. Lu, Evaluation of moments of performance functions based on efficient cubature formulation, Journal of Engineering Mechanics 143 (2017) 06017007. [28] X. Jia, C. Jiang, C. Fu, B. Ni, C. Wang, M. Ping, Uncertainty propagation analysis by an extended sparse grid technique, Frontiers of Mechanical Engineering 14 (2019) 33–46. [29] S. Rahman, H. Xu, A univariate dimension-reduction method for multi-dimensional integration in stochastic mechanics, Probabilistic Engineering Mechanics 19 (2004) 393–408. [30] H. Xu, S. Rahman, A generalized dimension-reduction method for multidimensional integration in stochastic mechanics, International Journal for Numerical Methods in Engineering 61 (2004) 1992–2019. [31] G. Li, K. Zhang, A combined reliability analysis approach with dimension reduction method and maximum entropy method, Structural and Multidisciplinary Optimization 43 (2011) 121–134. [32] X. Zhang, M. D. Pandey, Structural reliability analysis based on the concepts of entropy, fractional moment and dimensional reduction method, Structural Safety 43 (2013) 28–40. [33] G. P. Balomenos, M. D. Pandey, Finite element reliability and sensitivity analysis of structures using the multiplicative dimensional reduction method, Structure and Infrastructure Engineering 12 (2016) 1553–1565. [34] G. P. Balomenos, A. S. Genikomsou, M. A. Polak, M. D. Pandey, Efficient method for probabilistic finite element analysis with application to reinforced concrete slabs, Engineering Structures 103 (2015) 85–101. [35] J. Xu, C. Dang, A new bivariate dimension reduction method for efficient structural reliability analysis, Mechanical Systems and Signal Processing 115 (2019) 281–300. [36] W. Fan, J. Wei, A. H. Ang, Z. Li, Adaptive estimation of statistical moments of the responses of random systems, Probabilistic Engineering Mechanics 43 (2016) 50–67. [37] Q. Zhou, Z. Li, W. Fan, A. H. Ang, R. Liu, System reliability assessment of deteriorating structures subjected to time-invariant loads based on improved moment method, Structural Safety 68 (2017) 54–64. [38] Y.-G. Zhang, Y.-L. Huang, Z.-M. Wu, N. Li, A high order unscented kalman filtering method, Acta Automatica Sinica 40 (2014) 838–848(in Chinese). [39] S. J. Julier, J. K. Uhlmann, Consistent debiased method for converting between polar and cartesian coordinate systems, in: Acquisition, Tracking, and Pointing XI, volume 3086, International Society for Optics and Photonics, 1997, pp. 110–121. [40] Y. Low, A new distribution for fitting four moments and its applications to reliability analysis, Structural Safety 42 (2013) 12–25. [41] A. Pollastri, A. Brunazzo, Proposta di una nuova distribuzione: la lognormale generalizzata (1986). [42] C. Dang, J. Xu, Novel algorithm for reconstruction of a distribution by fitting its first-four statistical moments, Applied Mathematical Modelling 71 (2019) 505 – 524. [43] W. Fan, R. Liu, A. H. Ang, Z. Li, A new point estimation method for statistical moments based on dimension-reduction method and direct numerical integration, Applied Mathematical Modelling 62 (2018) 664–679.

21