# Reliability Sensitivity Analysis Method for Mechanical Components

Authors: Yan-Fang Zhang; Yan-Lin Zhang
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

Based on the univariate dimension-reduction method (UDRM), Edgeworth series, and sensitivity analysis, a new method for reliability sensitivity analysis of mechanical components is proposed. The univariate dimension-reduction method is applied to calculate the response origin moments and their sensitivity with respect to distribution parameters (e.g., mean and standard deviation) of fundamental input random variables. Edgeworth series is used to estimate failure probability of mechanical components by using first few response central moments. The analytic formula of reliability sensitivity can be derived by calculating partial derivative of the failure probabilityPf with respect to distribution parameters of basic random variables. The nonnormal random parameters need not to be transformed into equivalent normal ones. Three numerical examples are employed to illustrate the accuracy and efficiency of the proposed method by comparing the failure probability and reliability sensitivity results obtained by the proposed method with those obtained by Monte Carlo simulation (MCS).


---

## Body

## 1. Introduction

The responses of mechanical components or engineering structures are often random due to random inputs of them including loads, material properties, and geometry. Many reliability analysis methods have been developed to calculate the failure probabilities of these random structures [1–4]. In reliability analysis, reliability sensitivity is defined as the partial derivative of the failure probability with respect to the distribution parameters (e.g., mean and standard deviation) of fundamental random input variables. Reliability sensitivity provides information about the importance of each input random variable to a structure’s failure probability. Reliability sensitivity analysis methods can be divided into the numerical simulation method and approximate analytic calculation method.The numerical simulation method can be divided into many different methods due to different sampling methods including importance sampling, direction sampling, line sampling, subset simulation, and low-discrepancy sampling [5–10]. These different simulation methods are all based on Monte Carlo simulation but have different sampling methods. Although the program of these methods is not complicated, lots of samples are needed to get precise results because the structural failure probability is often very small and the computational efficiency is therefore reduced.Analytic methods for reliability sensitivity analysis are always based on analytic reliability methods. These reliability analysis methods are the mean-value first-order reliability method (MVFORM)/mean-value second-order reliability method (MVSORM) [11, 12], JC method, mean-value first-order saddlepoint approximation (MVFOSA) [13], and moment method. FORM/SORM needs to expand the performance function into first-order or second-order Taylor series at the most likely failure point. The expansion is widely used because of its relatively moderate precision and efficiency. The corresponding reliability sensitivity calculation methods have been derived from these reliability methods [14–23]. Analytic methods for reliability sensitivity analysis always have higher calculation efficiency than MCS.In this study, a feasible method to compute the sensitivity of failure probability with respect to distribution parameters of basic random input variables is proposed. The method can be applied to estimate failure probability and reliability sensitivity based on Edgeworth series. The proposed reliability sensitivity analysis method includes the following steps. (1) Choose random input variables and determine their distributions and distribution parameters. (2) Calculate the first four response origin moments and their sensitivity with respect to distribution parameters of random input variables by UDRM. (3) Calculate failure probability by the reliability method based on Edgeworth series. (4) Calculate partial derivatives of failure probability with respect to response central moments. (5) Calculate partial derivatives of response central moments with respect to response origin moments. (6) Calculate reliability sensitivity defined by partial derivatives of failure probability with respect to distribution parameters of input random variables.


## 2. Response Origin Moments and Response Sensitivity

### 2.1. Response Origin Moments by UDRM

Structures subject to random input vectorX=X1,X2,Λ,XNT∈RN, which characterizes uncertainty in loads, material properties, and geometry. Let YX represent a response of interest that depends on independent random variables X=X1,X2,Λ,XNT; then, the lth order origin moment of YX can be written as follows:(1)ml=EYlX=∫RNYlXfXxdx,where mlis the lth order origin moment of YX, fXxis the joint probability density function of X, and E⋅is the expectation operator. According to the univariate dimension reduction method proposed by Rahman and Xu [24–26], YX can be approximately written as follows:(2)YX≈Y^X=Y^X1,X2,Λ,XN=∑j=1NYμ1,Λ,μj−1,xj,μj+1,Λ,μN−N−1Yμ1,Λ,μN,where Yμ1,Λ,μj−1,xj,μj+1,ΛμN is the response which depends on the jth random input Xj, Yμ1,Λ,μj−1,xj,μj+1,ΛμN is the response where X=X1,X2,Λ,XNT, and μj=EXj is the mean of Xj.Then, thelth order origin moment ml of YX can be written as follows:(3)ml≅EY^lX=E∑j=1NYμ1,Λ,μj−1,Xj,μj+1,Λ,μN−N−1Yμ1,Λ,μNl.Applying the binomial formula on the right-hand side of equation (3), the lth response origin moment ml can be written as follows:(4)ml≅∑i=0lCliE∑j=1NYμ1,Λ,μj−1,Xj,μj+1,Λ,μNi−N−1Yμ1,Λ,μNl−i,where C⋅ is the combination operator.Define(5)Sji=E∑i=1jYμ1,Λ,μj−1,Xj,μj+1,Λ,μNi,j=1,…N,i=1,…,l.Equation (5) can be expressed using the recursive formula as follows [23]:(6)S1i=EYiX1,μ2,Λ,μN,i=1,…,N,S2i=∑k=0iCikS1iEYi−kμ1,X2,Λ,μN,i=1,…,N,Sji=∑k=0iCikSj−1iEYi−kμ1,Λ,μj−1,Xj,μj+1,Λ,μN,i=1,…,N,SNi=∑k=0iCikSN−1iEYi−kμ1,Λ,μN−1,XN,i=1,…,N.Then, thelth origin moment of YX can be written as follows:(7)ml≅EY^lX=∑i=0lCliSNi−N−1Yμ1,Λ,μN−1,μNl−i.From equations (6) and (7), the lth response origin moment ml can be calculated out by using the origin moments of univariate functions Yμ1,Λ,μj−1,Xj,μj+1,Λ,μNj=1,2,Λ,N. The qth origin moment of univariate function Yμ1,Λ,μj−1,Xj,μj+1,Λ,μNj=1,2,Λ,N can be written as follows:(8)mYjq=EYqμ1,Λ,μj−1,Xj,μj+1,Λ,μN=∫Yqμ1,Λ,μj−1,xj,μj+1,Λ,μNfXjxjdxj,where fXjxj is the probability density function of random variable Xj.Theqth origin moment mYjq in equation (8) can be calculated out by the numerical integration method and can be written as follows:(9)mYjq=∫Yqμ1,Λ,μj−1,xj,μj+1,Λ,μNfxjdxj=∑k=0nWjkYqμ1,Λ,μj−1,xjk,μj+1,Λ,μN,where Wjk is the kth integral weight, xjk is the kth integration point, and n is the number of integration points.Ifxj has different probability density (e.g., normal, lognormal, and Weibull), the corresponding integral weights and integration points in equation (9) are different. The numerical integration points and the corresponding orthogonal polynomial for different probability densities are listed in Table 1 [27].Table 1
Orthogonal polynomial and Gauss integration.
DistributionOrthogonal polynomialIntegral formulaNormalHermite1/π∑k=0nWjkYq2σtjk+μLognormalHermite1/π∑k=0nWjkYqe2σ1tjk+μ1UniformLegendre1/2∑k=0nWlkYqb−a/2tjk+a+b/2ExponentialLaguerre1/λ∑k=0nWjkYqtjk/λWeibullLaguerre∑k=0nWjkYqαtj1/β


### 2.2. Sensitivity of Response Moment with Respect to Distribution Parameters of Inputs

Calculate the partial derivative on both side of equation (1) with respect to distribution parameters Tji of random input variables X=X1,X2,Λ,XNT. The partial derivative can be written as follows:(10)∂ml∂Tji=∂EYlX∂Tji=∂∫RNYlXfXdx∂Tji=∫RNYlX1fjx⋅∂fjx∂TjifXdx=EYlX1fjx⋅∂fjx∂Tji≈E∑j=1NY^jXj−N−1Yμ1,Λ,μNl1fjx⋅∂fjx∂Tji,j=1,2,…,N,i=1,2,where Tjij=1,2,…,N,i=1,2, is the ith distribution parameter of the jth random input variable (i=1, represents the mean; i=2 represents standard deviation), fjx is the probability density function of the jth random variable Xj, and Y^jXj=Yμ1,Λ,μj−1,Xj,μj+1,Λ,μN.Applying the binomial formula on the right-hand side of equation (10), ∂ml/∂Tji can be written as follows:(11)∂ml∂Tji≈∑z=0lClzE∑j=1NY^jXjz1fjx⋅∂fjx∂Tji−N−1Yμ1,Λ,μNl−z,where ∂ml/∂Tji can be calculated out by equation (11) and recursive formula similar to equations (5) and (6).Suppose(12)kji=1fjx⋅∂fjx∂Tji,j=1,2,…,N,i=1,2,where kji is the kernel function of the ith distribution parameter Tji of jth the random input variable.The kernel function in equation (12) of Xj with different probability density functions can be obtained as follows: here, xj is denoted by x.IfX is a normal distributed variable, its probability density function can be written as follows:(13)fx=12πσe−x−μ2/2σ2,where μ is the mean of X, and σ is the standard deviation of X.The kernel function for normal variableX with respect to its distribution parameters (mean and standard deviation) can be derived from equation (12) directly and can be written as follows:(14)kμ=x−μσ2,kσ=1σx−μσ2−1,where kμ and kσ are the kernel functions with respect to mean and standard deviation, respectively. μ and σ are the mean and standard deviation of normal variable X, respectively.IfX is a lognormal distributed variable, its probability density function can be written as follows:(15)fx=12πσ1xe−lnx−μ12/2σ12,where μ1 is the mean of lnX, and σ1 is the standard deviation of lnX.The kernel function of lognormal variableX with respect to μ1 and σ1 can be derived from equation (12) directly and can be written as follows:(16)kμ1=lnx−μ1σ12,kσ1=lnx−μ12−σ12σ13,where kμ1 and kσ1 are the kernel functions with respect to μ1 and σ1, respectively.Distribution parametersμ1 and σ1 of the lognormal distributed variable x can be written as follows:(17)μ1=lnμ2σ2+μ2,σ1=lnσ2μ2+1,where μ and σ are the mean and standard deviation of the lognormal distributed variable x, respectively.Calculate partial derivative ofμ1 and σ1 with respect to μ and σ, respectively, according to equation (16); then,(18)∂μ1∂μ=2σ2+μ2μσ2+μ2,∂μ1∂σ=−2σσ2+μ2,∂σ1∂μ=μlnσ2/μ2+1⋅μ2+1,∂σ1∂σ=1σlnσ2/μ2+1.Jacobi matrixJ can be written as follows:(19)J=∂μ1∂μ∂σ1∂μ∂μ1∂σ∂σ1∂σ.The kernel function of the lognormal distributed variableX with respect to its mean and standard deviation can be written as follows:(20)kμkσ=∂μ1∂μkμ1+∂σ1∂μkσ1∂μ1∂σkμ1+∂σ1∂σkσ1=∂μ1∂μ∂σ1∂μ∂μ1∂σ∂σ1∂σ⋅kμ1kσ1=J⋅kμ1kσ1,where kμ and kσ are the kernel functions with respect to mean and standard deviation of the lognormal distributed variable x, respectively.Ifx is two parameter Weibull distributed variable, its probability density function can be written as follows:(21)fx=βαxαβ−1e−x/αβ,where β is the shape parameter, and α is the scale parameter.The kernel function of the Weibull distributed variableX with respect to scale parameter α and shape parameter β can be derived from equation (12) directly, which can be written as follows:(22)kα=βαxαβ−1,kβ=1β+1−xββlnxα.The mean and standard deviationμ and σof X can be written as follows:(23)μ=α⋅Γ1β+1,σ=α⋅Γ2β+1−Γ1β+12,where μ and σ are the mean and standard deviation of X, respectively, and Γ⋅ is the gamma function.Calculate partial derivative ofμ and σ with respect to α and β, respectively, from equation (21); then,(24)∂μ∂α=Γ1β+1,∂σ∂α=σα,∂μ∂β=−αΓ1/β+1Ψ1/β+1β2,(25)∂σ∂β=α2Γ21/β+1Ψ1/β+1−Γ2/β+1Ψ2/β+1β2σ,where Ψ⋅ is the psi function.Jacobi matrixJ can be written as follows: J=∂μ/∂α∂σ/∂α∂μ/∂β∂σ/∂β; then, the inverse matrix of J can be written as J−1=∂α/∂μ∂β/∂μ∂α/∂σ∂β/∂σ.The kernel function of the Weibull distributed variableX with respect to its mean and standard deviation can be written as follows:(26)kμkσ=∂α∂μkα+∂β∂μkβ∂α∂σkα+∂β∂σkβ=∂α∂μ∂β∂μ∂α∂σ∂β∂σ⋅kαkβ=J−1⋅kαkβ,where kα and kβ are the kernel functions with respect to scale parameter and shape parameter, respectively. α and β are the scale parameter and shape parameter of the two parameter Weibull distributed X, respectively.


## 2.1. Response Origin Moments by UDRM

Structures subject to random input vectorX=X1,X2,Λ,XNT∈RN, which characterizes uncertainty in loads, material properties, and geometry. Let YX represent a response of interest that depends on independent random variables X=X1,X2,Λ,XNT; then, the lth order origin moment of YX can be written as follows:(1)ml=EYlX=∫RNYlXfXxdx,where mlis the lth order origin moment of YX, fXxis the joint probability density function of X, and E⋅is the expectation operator. According to the univariate dimension reduction method proposed by Rahman and Xu [24–26], YX can be approximately written as follows:(2)YX≈Y^X=Y^X1,X2,Λ,XN=∑j=1NYμ1,Λ,μj−1,xj,μj+1,Λ,μN−N−1Yμ1,Λ,μN,where Yμ1,Λ,μj−1,xj,μj+1,ΛμN is the response which depends on the jth random input Xj, Yμ1,Λ,μj−1,xj,μj+1,ΛμN is the response where X=X1,X2,Λ,XNT, and μj=EXj is the mean of Xj.Then, thelth order origin moment ml of YX can be written as follows:(3)ml≅EY^lX=E∑j=1NYμ1,Λ,μj−1,Xj,μj+1,Λ,μN−N−1Yμ1,Λ,μNl.Applying the binomial formula on the right-hand side of equation (3), the lth response origin moment ml can be written as follows:(4)ml≅∑i=0lCliE∑j=1NYμ1,Λ,μj−1,Xj,μj+1,Λ,μNi−N−1Yμ1,Λ,μNl−i,where C⋅ is the combination operator.Define(5)Sji=E∑i=1jYμ1,Λ,μj−1,Xj,μj+1,Λ,μNi,j=1,…N,i=1,…,l.Equation (5) can be expressed using the recursive formula as follows [23]:(6)S1i=EYiX1,μ2,Λ,μN,i=1,…,N,S2i=∑k=0iCikS1iEYi−kμ1,X2,Λ,μN,i=1,…,N,Sji=∑k=0iCikSj−1iEYi−kμ1,Λ,μj−1,Xj,μj+1,Λ,μN,i=1,…,N,SNi=∑k=0iCikSN−1iEYi−kμ1,Λ,μN−1,XN,i=1,…,N.Then, thelth origin moment of YX can be written as follows:(7)ml≅EY^lX=∑i=0lCliSNi−N−1Yμ1,Λ,μN−1,μNl−i.From equations (6) and (7), the lth response origin moment ml can be calculated out by using the origin moments of univariate functions Yμ1,Λ,μj−1,Xj,μj+1,Λ,μNj=1,2,Λ,N. The qth origin moment of univariate function Yμ1,Λ,μj−1,Xj,μj+1,Λ,μNj=1,2,Λ,N can be written as follows:(8)mYjq=EYqμ1,Λ,μj−1,Xj,μj+1,Λ,μN=∫Yqμ1,Λ,μj−1,xj,μj+1,Λ,μNfXjxjdxj,where fXjxj is the probability density function of random variable Xj.Theqth origin moment mYjq in equation (8) can be calculated out by the numerical integration method and can be written as follows:(9)mYjq=∫Yqμ1,Λ,μj−1,xj,μj+1,Λ,μNfxjdxj=∑k=0nWjkYqμ1,Λ,μj−1,xjk,μj+1,Λ,μN,where Wjk is the kth integral weight, xjk is the kth integration point, and n is the number of integration points.Ifxj has different probability density (e.g., normal, lognormal, and Weibull), the corresponding integral weights and integration points in equation (9) are different. The numerical integration points and the corresponding orthogonal polynomial for different probability densities are listed in Table 1 [27].Table 1
Orthogonal polynomial and Gauss integration.
DistributionOrthogonal polynomialIntegral formulaNormalHermite1/π∑k=0nWjkYq2σtjk+μLognormalHermite1/π∑k=0nWjkYqe2σ1tjk+μ1UniformLegendre1/2∑k=0nWlkYqb−a/2tjk+a+b/2ExponentialLaguerre1/λ∑k=0nWjkYqtjk/λWeibullLaguerre∑k=0nWjkYqαtj1/β


## 2.2. Sensitivity of Response Moment with Respect to Distribution Parameters of Inputs

Calculate the partial derivative on both side of equation (1) with respect to distribution parameters Tji of random input variables X=X1,X2,Λ,XNT. The partial derivative can be written as follows:(10)∂ml∂Tji=∂EYlX∂Tji=∂∫RNYlXfXdx∂Tji=∫RNYlX1fjx⋅∂fjx∂TjifXdx=EYlX1fjx⋅∂fjx∂Tji≈E∑j=1NY^jXj−N−1Yμ1,Λ,μNl1fjx⋅∂fjx∂Tji,j=1,2,…,N,i=1,2,where Tjij=1,2,…,N,i=1,2, is the ith distribution parameter of the jth random input variable (i=1, represents the mean; i=2 represents standard deviation), fjx is the probability density function of the jth random variable Xj, and Y^jXj=Yμ1,Λ,μj−1,Xj,μj+1,Λ,μN.Applying the binomial formula on the right-hand side of equation (10), ∂ml/∂Tji can be written as follows:(11)∂ml∂Tji≈∑z=0lClzE∑j=1NY^jXjz1fjx⋅∂fjx∂Tji−N−1Yμ1,Λ,μNl−z,where ∂ml/∂Tji can be calculated out by equation (11) and recursive formula similar to equations (5) and (6).Suppose(12)kji=1fjx⋅∂fjx∂Tji,j=1,2,…,N,i=1,2,where kji is the kernel function of the ith distribution parameter Tji of jth the random input variable.The kernel function in equation (12) of Xj with different probability density functions can be obtained as follows: here, xj is denoted by x.IfX is a normal distributed variable, its probability density function can be written as follows:(13)fx=12πσe−x−μ2/2σ2,where μ is the mean of X, and σ is the standard deviation of X.The kernel function for normal variableX with respect to its distribution parameters (mean and standard deviation) can be derived from equation (12) directly and can be written as follows:(14)kμ=x−μσ2,kσ=1σx−μσ2−1,where kμ and kσ are the kernel functions with respect to mean and standard deviation, respectively. μ and σ are the mean and standard deviation of normal variable X, respectively.IfX is a lognormal distributed variable, its probability density function can be written as follows:(15)fx=12πσ1xe−lnx−μ12/2σ12,where μ1 is the mean of lnX, and σ1 is the standard deviation of lnX.The kernel function of lognormal variableX with respect to μ1 and σ1 can be derived from equation (12) directly and can be written as follows:(16)kμ1=lnx−μ1σ12,kσ1=lnx−μ12−σ12σ13,where kμ1 and kσ1 are the kernel functions with respect to μ1 and σ1, respectively.Distribution parametersμ1 and σ1 of the lognormal distributed variable x can be written as follows:(17)μ1=lnμ2σ2+μ2,σ1=lnσ2μ2+1,where μ and σ are the mean and standard deviation of the lognormal distributed variable x, respectively.Calculate partial derivative ofμ1 and σ1 with respect to μ and σ, respectively, according to equation (16); then,(18)∂μ1∂μ=2σ2+μ2μσ2+μ2,∂μ1∂σ=−2σσ2+μ2,∂σ1∂μ=μlnσ2/μ2+1⋅μ2+1,∂σ1∂σ=1σlnσ2/μ2+1.Jacobi matrixJ can be written as follows:(19)J=∂μ1∂μ∂σ1∂μ∂μ1∂σ∂σ1∂σ.The kernel function of the lognormal distributed variableX with respect to its mean and standard deviation can be written as follows:(20)kμkσ=∂μ1∂μkμ1+∂σ1∂μkσ1∂μ1∂σkμ1+∂σ1∂σkσ1=∂μ1∂μ∂σ1∂μ∂μ1∂σ∂σ1∂σ⋅kμ1kσ1=J⋅kμ1kσ1,where kμ and kσ are the kernel functions with respect to mean and standard deviation of the lognormal distributed variable x, respectively.Ifx is two parameter Weibull distributed variable, its probability density function can be written as follows:(21)fx=βαxαβ−1e−x/αβ,where β is the shape parameter, and α is the scale parameter.The kernel function of the Weibull distributed variableX with respect to scale parameter α and shape parameter β can be derived from equation (12) directly, which can be written as follows:(22)kα=βαxαβ−1,kβ=1β+1−xββlnxα.The mean and standard deviationμ and σof X can be written as follows:(23)μ=α⋅Γ1β+1,σ=α⋅Γ2β+1−Γ1β+12,where μ and σ are the mean and standard deviation of X, respectively, and Γ⋅ is the gamma function.Calculate partial derivative ofμ and σ with respect to α and β, respectively, from equation (21); then,(24)∂μ∂α=Γ1β+1,∂σ∂α=σα,∂μ∂β=−αΓ1/β+1Ψ1/β+1β2,(25)∂σ∂β=α2Γ21/β+1Ψ1/β+1−Γ2/β+1Ψ2/β+1β2σ,where Ψ⋅ is the psi function.Jacobi matrixJ can be written as follows: J=∂μ/∂α∂σ/∂α∂μ/∂β∂σ/∂β; then, the inverse matrix of J can be written as J−1=∂α/∂μ∂β/∂μ∂α/∂σ∂β/∂σ.The kernel function of the Weibull distributed variableX with respect to its mean and standard deviation can be written as follows:(26)kμkσ=∂α∂μkα+∂β∂μkβ∂α∂σkα+∂β∂σkβ=∂α∂μ∂β∂μ∂α∂σ∂β∂σ⋅kαkβ=J−1⋅kαkβ,where kα and kβ are the kernel functions with respect to scale parameter and shape parameter, respectively. α and β are the scale parameter and shape parameter of the two parameter Weibull distributed X, respectively.


## 3. Reliability and Reliability Sensitivity Based on Edgeworth Series

### 3.1. Reliability Based on Edgeworth Series

The failure probability of a mechanical component or a random structure can be calculated out by the following multidimensional integral:(27)Pf=PgX<0=∫gX<0fXxdx,where Pf is the failure probability of a random structure, X=X1,X2,Λ,XNT∈RN represents the N-dimensional random input variables of the random structure, fXx is the joint probability density function of random input variable X, gX is the performance function, gX<0 represents the failure domain, gX>0 represents the safety domain, and gX=0 represents the limit state.In practice, it is difficult to obtain the analytical solution of the multidimensional integral in equation (25) because of the nonlinear integration boundary gX=0 and high dimension of X. There are two classes of approaches available for estimating the failure probability Pf in equation (25), the analytical method and numerical simulation method. The numerical simulation method is the Monte Carlo simulation as known. Analytical methods (e.g., first-order reliability method (FORM)/second-order reliability method (SORM)), higher-order moment method, and saddlepoint approximation) can always be applied to calculate failure probability Pf because of high computational efficiency and accuracy of them.Suppose the first four central moments of performance functiongX are μg,vg,θg, and ηg. The failure probability Pf in equation (27) can be approximately estimated by Edgeworth series which can be written as follows:(28)Pfζ=Φζ−ϕζ16θgσg3H2ζ+124ηgσg4−3H3ζ+172θgσg32H5ζ,where Φ⋅ is the standard normal distribution function, φ⋅ is the standard normal probability density, μg,vg,θg, and ηg are the mean, deviation, the third central moment, and the fourth central moment of performance function gX, σg is the standard deviation of gX, ς=−μg/σh is the negative reliability index, and Hj⋅ is the jth Hermite polynomial and can be written as follows:(29)Hj+1ζ=2ζHjζ−2jHj−1ζ,H0ζ=1,H1ζ=2ζ.


### 3.2. Reliability Sensitivity Based on Edgeworth Series

Partial derivative of failure probabilityPf in equation (28) with respect to the first four central moments μg,vg,θg, and ηg of performance function gX can be written as follows:(30)∂Pf∂μg=∂Pf∂ζ∂ζ∂μg,(31)∂Pf∂νg=∂Pf∂ζ∂ζ∂σg+∂Pf∂σg∂σg∂νg,(32)∂Pf∂θg=−ϕζ16⋅σg3H2ζ+136θgσg6H5ζ,∂Pf∂ηg=−ϕζ124ηgσg4H3ζ,where(33)∂Pf∂ζ=ϕζ+ζϕζ16θgσg3H2ζ+124ηgσg4−3H3ζ+172θgσg32H5ζ−ϕζ13θgσg3H1ζ+18ηgσg4−3H2ζ+572θgσg32H4ζ,∂ζ∂μg=−1σg,∂ζ∂σg=μgσg2,∂σg∂νg=12σg,∂Pf∂σg=ϕζ12θgσg4H2ζ+16ηgσg5H3ζ+112θgσg7/22H5ζ.The first four central momentsμg,vg,θg, and ηg of performance function gX can be written as follows:(34)μg=m1,νg=m2−m12,θg=m3−3m1m2+2m13,ηg=m4−4m1m3+6m12m2−3m14,where mll=1,2,3,4 are the first four origin moments of gX=0 which can be estimated by equation (7).From equation (34), the partial derivative of central moments μg,vg,θg, and ηg of performance function gX with respect to the distribution parameters Tjij=1,2,…,N,i=1,2 (i=1 represents mean; i=2 represents standard deviation) of random inputs X=X1,X2,…,XNT can be written as follows:(35)∂μg∂Tji=∂m1∂Tji,∂νg∂Tji=∂m2∂Tji−2m1∂m1∂Tji,∂θg∂Tji=∂m3∂Tji−3∂m1∂Tjim2+∂m2∂Tjim1+6m12∂m2∂Tji,∂ηg∂Tji=∂m4∂Tji−4∂m1∂Tjim3+∂m3∂Tjim1+62m1m2∂m1∂Tji+m12∂m2∂Tji−12m13∂m1∂Tji,i=1,2,j=1,2,…,N,where ∂ml/∂Tjil=1,2,3,4 can be calculated out by equation (11).The failure probability sensitivity defined as partial derivative of failure probabilityPf with respect to distribution parameters Tji of random inputs X=X1,X2,…,XNT can be written as follows:(36)∂Pf∂Tji=∂Pf∂μg∂μg∂Tji+∂Pf∂νg∂νg∂Tji+∂Pf∂θg∂θg∂Tji+∂Pf∂ηg∂ηg∂Tji,j=1,2,…,N,i=1,2,.From equation (36), ∂Pj/∂Tji can be calculated out according to equations (28)–(31) and equation (35). The flow process diagram of the computational procedure is shown in Figure 1.Figure 1
Flow process diagram for computations.


## 3.1. Reliability Based on Edgeworth Series

The failure probability of a mechanical component or a random structure can be calculated out by the following multidimensional integral:(27)Pf=PgX<0=∫gX<0fXxdx,where Pf is the failure probability of a random structure, X=X1,X2,Λ,XNT∈RN represents the N-dimensional random input variables of the random structure, fXx is the joint probability density function of random input variable X, gX is the performance function, gX<0 represents the failure domain, gX>0 represents the safety domain, and gX=0 represents the limit state.In practice, it is difficult to obtain the analytical solution of the multidimensional integral in equation (25) because of the nonlinear integration boundary gX=0 and high dimension of X. There are two classes of approaches available for estimating the failure probability Pf in equation (25), the analytical method and numerical simulation method. The numerical simulation method is the Monte Carlo simulation as known. Analytical methods (e.g., first-order reliability method (FORM)/second-order reliability method (SORM)), higher-order moment method, and saddlepoint approximation) can always be applied to calculate failure probability Pf because of high computational efficiency and accuracy of them.Suppose the first four central moments of performance functiongX are μg,vg,θg, and ηg. The failure probability Pf in equation (27) can be approximately estimated by Edgeworth series which can be written as follows:(28)Pfζ=Φζ−ϕζ16θgσg3H2ζ+124ηgσg4−3H3ζ+172θgσg32H5ζ,where Φ⋅ is the standard normal distribution function, φ⋅ is the standard normal probability density, μg,vg,θg, and ηg are the mean, deviation, the third central moment, and the fourth central moment of performance function gX, σg is the standard deviation of gX, ς=−μg/σh is the negative reliability index, and Hj⋅ is the jth Hermite polynomial and can be written as follows:(29)Hj+1ζ=2ζHjζ−2jHj−1ζ,H0ζ=1,H1ζ=2ζ.


## 3.2. Reliability Sensitivity Based on Edgeworth Series

Partial derivative of failure probabilityPf in equation (28) with respect to the first four central moments μg,vg,θg, and ηg of performance function gX can be written as follows:(30)∂Pf∂μg=∂Pf∂ζ∂ζ∂μg,(31)∂Pf∂νg=∂Pf∂ζ∂ζ∂σg+∂Pf∂σg∂σg∂νg,(32)∂Pf∂θg=−ϕζ16⋅σg3H2ζ+136θgσg6H5ζ,∂Pf∂ηg=−ϕζ124ηgσg4H3ζ,where(33)∂Pf∂ζ=ϕζ+ζϕζ16θgσg3H2ζ+124ηgσg4−3H3ζ+172θgσg32H5ζ−ϕζ13θgσg3H1ζ+18ηgσg4−3H2ζ+572θgσg32H4ζ,∂ζ∂μg=−1σg,∂ζ∂σg=μgσg2,∂σg∂νg=12σg,∂Pf∂σg=ϕζ12θgσg4H2ζ+16ηgσg5H3ζ+112θgσg7/22H5ζ.The first four central momentsμg,vg,θg, and ηg of performance function gX can be written as follows:(34)μg=m1,νg=m2−m12,θg=m3−3m1m2+2m13,ηg=m4−4m1m3+6m12m2−3m14,where mll=1,2,3,4 are the first four origin moments of gX=0 which can be estimated by equation (7).From equation (34), the partial derivative of central moments μg,vg,θg, and ηg of performance function gX with respect to the distribution parameters Tjij=1,2,…,N,i=1,2 (i=1 represents mean; i=2 represents standard deviation) of random inputs X=X1,X2,…,XNT can be written as follows:(35)∂μg∂Tji=∂m1∂Tji,∂νg∂Tji=∂m2∂Tji−2m1∂m1∂Tji,∂θg∂Tji=∂m3∂Tji−3∂m1∂Tjim2+∂m2∂Tjim1+6m12∂m2∂Tji,∂ηg∂Tji=∂m4∂Tji−4∂m1∂Tjim3+∂m3∂Tjim1+62m1m2∂m1∂Tji+m12∂m2∂Tji−12m13∂m1∂Tji,i=1,2,j=1,2,…,N,where ∂ml/∂Tjil=1,2,3,4 can be calculated out by equation (11).The failure probability sensitivity defined as partial derivative of failure probabilityPf with respect to distribution parameters Tji of random inputs X=X1,X2,…,XNT can be written as follows:(36)∂Pf∂Tji=∂Pf∂μg∂μg∂Tji+∂Pf∂νg∂νg∂Tji+∂Pf∂θg∂θg∂Tji+∂Pf∂ηg∂ηg∂Tji,j=1,2,…,N,i=1,2,.From equation (36), ∂Pj/∂Tji can be calculated out according to equations (28)–(31) and equation (35). The flow process diagram of the computational procedure is shown in Figure 1.Figure 1
Flow process diagram for computations.


## 4. Numerical Examples

### 4.1. Example 1: Nonlinear Performance Function with Normal Variables

This example considers a nonlinear performance function which is the stress limit state function of a multileaf spring written as follows:(37)gX=r−3Pl7bh2,where basic random variables X=r,p,l,b,h are the independent and identically distributed Gaussian random variables, r is the material strength of the leaf spring, p is the load, and b, h, and l are the geometric dimension width, thickness, and span of the leaf spring, respectively. The mean and standard deviation of random variables X=r,p,l,b,h are listed in Table 2 [22].Table 2
Statistic moments and distributions of input variables of Example 1.
DistributionNormalNormalNormalNormalNormalVariablesr (MPa)P (N)l (mm)b (mm)h (mm)Mean61455371100706.5Standard deviation45.8276.855.50.350.0325The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table3 for comparison. The MCS results can be considered as accurate results and N is the number of samples.Table 3
Failure probability and sensitivity results of Example 1.
Proposed methodMCS (N = 5 × 106)Pf3.6904 × 10−43.7501 × 10−4∂Pf/∂μr−2.6160 × 10−5−2.6826 × 10−5∂Pf/∂σr7.9061 × 10−57.8365 × 10−5∂Pf/∂μh−3.5000 × 10−3−3.9000 × 10−3∂Pf/∂σh1.1000 × 10−39.4906 × 10−4∂Pf/∂μb−1.6736 × 10−4−1.5067 × 10−4∂Pf/∂σb2.5960 × 10−52.7601 × 10−5∂Pf/∂μl1.0535 × 10−51.1517 × 10−5∂Pf/∂σl1.5250 × 10−61.8140 × 10−6∂Pf/∂μP2.0798 × 10−62.1164 × 10−6∂Pf/∂σP3.0267 × 10−63.0777 × 10−6In Table3, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 1 indicates that the proposed method is suitable for problems with nonlinear performance functions.


### 4.2. Example 2: Linear Performance Function with Six Lognormal Variables

This example considers a linear performance function, a plastic collapse mechanism of a one-bay frame which has been used as example 5 by Der et al. [28] and Zhao and Ono [29], written as follows:(38)gX=x1+2x2+2x3+x4−5x5−5x6,where basic random variables X=x1,x2,x3,x4,x5,x6 are the independent and identically distributed lognormal random variables. The mean and standard deviation of random variables X=x1,x2,x3,x4,x5,x6 are listed in Table 4.Table 4
Statistic moments and distributions of input variables of Example 2.
DistributionLognormalLognormalLognormalLognormalLognormalLognormalVariablesx1x2x3x4x5x6Mean1201201201205040Standard deviation121212121512The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table5 for comparison.Table 5
Failure probability and sensitivity results of Example 2.
Proposed methodMCS (N = 5 × 106)Pf1.3600 × 10−21.204 × 10−2∂Pf/∂μ1−2.2290 × 10−4−1.942 × 10−4∂Pf/∂σ13.4596 × 10−54.104 × 10−5∂Pf/∂μ2−4.4554 × 10−4−4.068 × 10−4∂Pf/∂σ21.3341 × 10−41.610 × 10−4∂Pf/∂μ3−4.4554 × 10−4−4.072 × 10−4∂Pf/∂σ31.3341 × 10−41.429 × 10−4∂Pf/∂μ4−2.2290 × 10−4−2.342 × 10−4∂Pf/∂σ43.4596 × 10−57.114 × 10−5∂Pf/∂μ59.5913 × 10−49.157 × 10−4∂Pf/∂σ52.9000 × 10−32.482 × 10−3∂Pf/∂μ61.0000 × 10−39.500 × 10−4∂Pf/∂σ61.9000 × 10−31.732 × 10−3In Table5, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 2 indicates that the proposed method is suitable for problems with nonnormal random input variables.


### 4.3. Example 3: Nonlinear Performance Function with Nonnormal Variables

This example considers a nonlinear performance function which describes the displacement response of aI section cantilever beam shown in Figure 2 written as follows:(39)gX=d−4Fl3Ebh3−bh−2t3+ah−2t3,where basic random variables X=a,t,b,h,l,F,d,E are the independent random variables. Here, a is the thickness of web, t is the thickness of flange, b is the width of I-beam section, h is the height of I-beam section, l is the length of cantilever beam, f is the load force, d is the allowable maximum deformation, and E is the elastic modulus. The distributions and distribution parameters (mean and standard deviation) of these input random variables are listed in Table 6.Figure 2
Section cantilever beam.Table 6
Statistic moments and distributions of input variables of Example 2.
DistributionNormalNormalNormalNormalNormalWeibullLognormalLognormalVariablesa (mm)t (mm)b (mm)h (mm)l (mm)F (N)d (mm)E (MPa)Mean10148820520005.187 × 104272.03 × 105Standard deviation0.050.070.441.02510213.8650.1351.015 × 103The scale parameter and shape parameter of random variableF are α = 5.197 × 104 and β = 310.36, respectively.The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table7 for comparison.Table 7
Failure probability and sensitivity results of Example 3.
Proposed methodMCS (N = 5 × 106)Pf9.7363 × 10−42.8000 × 10−4∂Pf/∂μa−8.2470 × 10−4−5.0299 × 10−4∂Pf/∂σa1.0579 × 10−41.1908 × 10−4∂Pf/∂μt−2.2000 × 10−3−2.0500 × 10−3∂Pf/∂σt1.1000 × 10−31.4000 × 10−3∂Pf/∂μb−4.4905 × 10−4−5.8711 × 10−4∂Pf/∂σb3.2129 × 10−41.8541 × 10−4∂Pf/∂μh−5.5812 × 10−4−6.4358 × 10−4∂Pf/∂σh1.1000 × 10−31.5000 × 10−3∂Pf/∂μl9.3914 × 10−58.9505 × 10−5∂Pf/∂σl1.9045 × 10−41.9436 × 10−4∂Pf/∂μF5.5696 × 10−61.5504 × 10−6∂Pf/∂σF8.4887 × 10−58.8620 × 10−7∂Pf/∂μd−1.9000 × 10−3−2.1000 × 10−3∂Pf/∂σd1.6000 × 10−31.5000 × 10−3∂Pf/∂μE−2.2514 × 10−7−3.6546 × 10−7∂Pf/∂σE1.8949 × 10−72.5002 × 10−7In Table7, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 3 indicates that the proposed method is suitable for problems with nonlinear performance functions with nonnormal random input variables. The nonnormal variables need not to be transformed into equivalent normal ones. But the coefficient of variation of the input random variables should be small (≤0.1 for Weibull distribution variables).


## 4.1. Example 1: Nonlinear Performance Function with Normal Variables

This example considers a nonlinear performance function which is the stress limit state function of a multileaf spring written as follows:(37)gX=r−3Pl7bh2,where basic random variables X=r,p,l,b,h are the independent and identically distributed Gaussian random variables, r is the material strength of the leaf spring, p is the load, and b, h, and l are the geometric dimension width, thickness, and span of the leaf spring, respectively. The mean and standard deviation of random variables X=r,p,l,b,h are listed in Table 2 [22].Table 2
Statistic moments and distributions of input variables of Example 1.
DistributionNormalNormalNormalNormalNormalVariablesr (MPa)P (N)l (mm)b (mm)h (mm)Mean61455371100706.5Standard deviation45.8276.855.50.350.0325The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table3 for comparison. The MCS results can be considered as accurate results and N is the number of samples.Table 3
Failure probability and sensitivity results of Example 1.
Proposed methodMCS (N = 5 × 106)Pf3.6904 × 10−43.7501 × 10−4∂Pf/∂μr−2.6160 × 10−5−2.6826 × 10−5∂Pf/∂σr7.9061 × 10−57.8365 × 10−5∂Pf/∂μh−3.5000 × 10−3−3.9000 × 10−3∂Pf/∂σh1.1000 × 10−39.4906 × 10−4∂Pf/∂μb−1.6736 × 10−4−1.5067 × 10−4∂Pf/∂σb2.5960 × 10−52.7601 × 10−5∂Pf/∂μl1.0535 × 10−51.1517 × 10−5∂Pf/∂σl1.5250 × 10−61.8140 × 10−6∂Pf/∂μP2.0798 × 10−62.1164 × 10−6∂Pf/∂σP3.0267 × 10−63.0777 × 10−6In Table3, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 1 indicates that the proposed method is suitable for problems with nonlinear performance functions.


## 4.2. Example 2: Linear Performance Function with Six Lognormal Variables

This example considers a linear performance function, a plastic collapse mechanism of a one-bay frame which has been used as example 5 by Der et al. [28] and Zhao and Ono [29], written as follows:(38)gX=x1+2x2+2x3+x4−5x5−5x6,where basic random variables X=x1,x2,x3,x4,x5,x6 are the independent and identically distributed lognormal random variables. The mean and standard deviation of random variables X=x1,x2,x3,x4,x5,x6 are listed in Table 4.Table 4
Statistic moments and distributions of input variables of Example 2.
DistributionLognormalLognormalLognormalLognormalLognormalLognormalVariablesx1x2x3x4x5x6Mean1201201201205040Standard deviation121212121512The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table5 for comparison.Table 5
Failure probability and sensitivity results of Example 2.
Proposed methodMCS (N = 5 × 106)Pf1.3600 × 10−21.204 × 10−2∂Pf/∂μ1−2.2290 × 10−4−1.942 × 10−4∂Pf/∂σ13.4596 × 10−54.104 × 10−5∂Pf/∂μ2−4.4554 × 10−4−4.068 × 10−4∂Pf/∂σ21.3341 × 10−41.610 × 10−4∂Pf/∂μ3−4.4554 × 10−4−4.072 × 10−4∂Pf/∂σ31.3341 × 10−41.429 × 10−4∂Pf/∂μ4−2.2290 × 10−4−2.342 × 10−4∂Pf/∂σ43.4596 × 10−57.114 × 10−5∂Pf/∂μ59.5913 × 10−49.157 × 10−4∂Pf/∂σ52.9000 × 10−32.482 × 10−3∂Pf/∂μ61.0000 × 10−39.500 × 10−4∂Pf/∂σ61.9000 × 10−31.732 × 10−3In Table5, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 2 indicates that the proposed method is suitable for problems with nonnormal random input variables.


## 4.3. Example 3: Nonlinear Performance Function with Nonnormal Variables

This example considers a nonlinear performance function which describes the displacement response of aI section cantilever beam shown in Figure 2 written as follows:(39)gX=d−4Fl3Ebh3−bh−2t3+ah−2t3,where basic random variables X=a,t,b,h,l,F,d,E are the independent random variables. Here, a is the thickness of web, t is the thickness of flange, b is the width of I-beam section, h is the height of I-beam section, l is the length of cantilever beam, f is the load force, d is the allowable maximum deformation, and E is the elastic modulus. The distributions and distribution parameters (mean and standard deviation) of these input random variables are listed in Table 6.Figure 2
Section cantilever beam.Table 6
Statistic moments and distributions of input variables of Example 2.
DistributionNormalNormalNormalNormalNormalWeibullLognormalLognormalVariablesa (mm)t (mm)b (mm)h (mm)l (mm)F (N)d (mm)E (MPa)Mean10148820520005.187 × 104272.03 × 105Standard deviation0.050.070.441.02510213.8650.1351.015 × 103The scale parameter and shape parameter of random variableF are α = 5.197 × 104 and β = 310.36, respectively.The proposed method and MCS are applied to analyze the failure probability and reliability sensitivity, respectively, and their results are listed in Table7 for comparison.Table 7
Failure probability and sensitivity results of Example 3.
Proposed methodMCS (N = 5 × 106)Pf9.7363 × 10−42.8000 × 10−4∂Pf/∂μa−8.2470 × 10−4−5.0299 × 10−4∂Pf/∂σa1.0579 × 10−41.1908 × 10−4∂Pf/∂μt−2.2000 × 10−3−2.0500 × 10−3∂Pf/∂σt1.1000 × 10−31.4000 × 10−3∂Pf/∂μb−4.4905 × 10−4−5.8711 × 10−4∂Pf/∂σb3.2129 × 10−41.8541 × 10−4∂Pf/∂μh−5.5812 × 10−4−6.4358 × 10−4∂Pf/∂σh1.1000 × 10−31.5000 × 10−3∂Pf/∂μl9.3914 × 10−58.9505 × 10−5∂Pf/∂σl1.9045 × 10−41.9436 × 10−4∂Pf/∂μF5.5696 × 10−61.5504 × 10−6∂Pf/∂σF8.4887 × 10−58.8620 × 10−7∂Pf/∂μd−1.9000 × 10−3−2.1000 × 10−3∂Pf/∂σd1.6000 × 10−31.5000 × 10−3∂Pf/∂μE−2.2514 × 10−7−3.6546 × 10−7∂Pf/∂σE1.8949 × 10−72.5002 × 10−7In Table7, the results obtained by the proposed method have small errors compared with those obtained by MCS. Example 3 indicates that the proposed method is suitable for problems with nonlinear performance functions with nonnormal random input variables. The nonnormal variables need not to be transformed into equivalent normal ones. But the coefficient of variation of the input random variables should be small (≤0.1 for Weibull distribution variables).


## 5. Conclusions

The univariate dimension-reduction method, reliability analysis based on Edgeworth series, and reliability sensitivity analysis are employed to present a computational procedure for estimating failure probability sensitivity of random engineering structures. Three numerical examples are employed to illustrate the performance of the proposed method. To illustrate the feasibility of the proposed method, failure probability and failure probability sensitivity results obtained by the new method are compared with those obtained by direct MCS. It is concluded as follows: The proposed method is suitable for models of which the input random variables have small coefficient of variation (≤0.1). The nonnormal input variables need not to be transformed into equivalent normal ones, so the new method can be conveniently applied to solve problems with nonnormal input variables.



---

*Source: 7867003-2021-07-16.xml*