Reliability mesh convergence analysis by introducing expanded control variates

Abstract

The safety evaluation of engineering systems whose performance evaluation requires finite element analysis is a challenge in reliability theory. Recently, Adjusted Control Variates Technique (ACVAT) has proposed by the authors to solve this issue. ACVAT uses the results of a finite element method (FEM) model with coarse mesh density as the control variates of the model with fine mesh and efficiently solves FEM-based reliability problems. ACVAT however does not provide any results about the reliability-based mesh convergence of the problem, which is an important tool in FEM. Mesh-refinement analysis allows checking whether the numerical solution is sufficiently accurate, even though the exact solution is unknown. In this study, by introducing expanded control variates (ECV) formulation, ACVAT is improved and the capabilities of the method are also extended for efficient reliability mesh convergence analysis of FEM-based reliability problems. In the present study, the FEM-based reliability analyses of four practical engineering problems are investigated by this method and the corresponding results are compared with accurate results obtained by analytical solutions for two problems. The results confirm that the proposed approach not only handles the mesh refinement progress with the required accuracy, but it also reduces considerably the computational cost of FEM-based reliability problems.

Keywords

finite element / reliability mesh convergence analysis / expanded control variates

Introduction

Structural reliability estimates the safety level of an engineering system by expressing the following probability integral [¹]:

(1)P_{\mathrm{f}}=Prob\left[G\left(x\right)\le \mathrm{0}\right]={\displaystyle {\int }_{G\left(x\right)\le \mathrm{0}}f\left(x\right)\mathrm{d}x},

where P_f is the failure probability, x is a vector of random variables. G(x) represents the limit state function and f(x) denotes the joint probability density function (PDF) of x.

When the performance evaluation of the system (value of G(x)) is not computable by an analytical method, finite element (FE) as the most popular and powerful numerical method is applied. Herein, a FE model with very fine mesh density is required to achieve a proper solution in reliability process [²–⁴]. Therefore, Eq. (1) reads:

(2)P_{\mathrm{f}}={\int }_{G\left(x\right)\le 0}\hspace{0.17em}f\left(x\right)\mathrm{d}x={\displaystyle {\int }_{G\left(x\right)\le \mathrm{0}}\mathbf{I}\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x},

in which G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}} is the performance of structure obtained by the finite element analysis (FEA) using a model with very fine mesh density and p\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right) is the index function of the probability integral as follows:

(3)\mathbb{I}\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)=\{\begin{array}{c}1,G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\le 0,\hfill \\ 0,G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}>0.\hfill \end{array}

Herein, proper mesh refinement (and mesh convergence analysis) plays an important role in safety evaluation of the structure. Mesh-refinement analysis allows checking whether the numerical solution is sufficiently accurate, even though the exact solution is unknown [⁵–¹¹]. In the deterministic FE problems, multiple approaches have been developed to estimate the error without excessively increasing the numerical effort [¹²–¹⁵]. For instance, mesh convergence is used as a mainstream tool to estimate discretization errors [¹⁵–¹⁷] and Grid Convergence Index (GCI) that is used incrementally as an elegant approach, particularly in Computational Fluid Dynamics (CFD) [¹⁸–²⁰]. The formal method of estimating discretization errors requires a curve of a critical result parameter, such as a specific stress or strain component [²¹], to be plotted against a range of mesh densities. If two successive runs with different mesh densities lead to the same result (or the difference is less than a specified threshold), then mesh convergence is considered to have been achieved.

In contrast to studies conducted for the deterministic mesh convergence analysis, there are no wide approaches for reliability mesh convergence analysis for FE-based reliability problems [²²–²⁴]. In these types of reliability problems, researchers often conduct a deterministic mesh convergence analysis and select a suitable FE model to evaluate G_{\mathrm{fine}\hspace{0.17em}\mathrm{mesh}}^{\mathrm{FEA}} in Eq. (2). Accordingly, an efficient reliability method may be used to estimate the system safety level. Since employing very fine mesh in reliability processes would be infeasible for many engineering problems and using a coarse mesh would also lead to improper safety evaluation, researchers often attempt to select a balanced FE model to evaluate G(x) in Eq. (2) with suitable accuracy and a reasonable computational effort.

Herein, it should be noted that although the detailed information about the deterministic mesh convergence and the error of the deterministic model are available, the effect of this error in the reliability evaluation is completely unknown. In other words, researchers generally accept an “unknown error” in reliability valuation which they have resulted from “known errors” in the deterministic performance valuations (results which obtained by the deterministic mesh convergence analysis). Therefore, they should accept the outcome as the solution of reliability analysis, since there is no approach for validation of the result. Reference [³] shows that accepting a model with 5% error in performance valuation may result in approximately 82% error for safety evaluation of a structure. Hence, this study aims to propose a new formulation based on Adjusted Control Variates Technique (ACVAT) to overcome this issue for efficient reliability mesh convergence analysis. To this aim, ACVAT and its deficiency are explained in Section 2, and then ACVAT is reformulated to improve the ability of ACVAT in Section 3. In Section 4, a convergence criterion is defined to achieve the highest accuracy with reasonable computational effort. Finally, the ability of the proposed approach is investigated in different practical engineering problems in Section 5.

Adjusted Control Variates Technique

Recently, the Control Variate Technique (CVT) has been presented as an efficient approach for solving finite element method (FEM)-based reliability problems [⁴]. The method uses a coarse mesh model (G_{\mathrm{coarse}\mathrm{mesh}}^{\mathrm{FEA}}) as control variates of the FE model with fine mesh model (G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}). In this method the failure probability integral given by Eq. (2) is rewritten as:

(4)P_f=\int \mathbf{I}\left(G_{\mathrm{coarse}\mathrm{mesh}}^{\mathrm{FEA}}\right)f\left(x\right)\mathrm{d}x+\int (\mathbf{I}\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-\mathbf{I}\left(G_{\mathrm{coarse}\mathrm{mesh}}^{\mathrm{FEA}}\right))f\left(x\right)\mathrm{d}x\mathrm{.}

The efficiency of the above equation is improved by applying linear regression and introducing three different correction factors (α) on the ACVAT [⁴]. Hence, Eq. (4) is approximated as follows:

(5)P_f=\int \mathbf{I}\left(G_{\mathrm{coarse}\mathrm{mesh}}^{\mathrm{FEA}}\right)f\left(x\right)\mathrm{d}x+\alpha ·(\int (\mathbf{I}\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-\mathbf{I}\left(G_{\mathrm{coarse}\mathrm{mesh}}^{\mathrm{FEA}}\right))f\left(x\right)\mathrm{d}x)\mathrm{.}

ACVAT however has a limitation. When the coarse and fine mesh models are selected unsuitably, this method does not provide a mesh convergence report to ensure enough accuracy. Furthermore, the accuracy of the ACVAT is highly affected by the differences between coarse and fine mesh models.

To overcome this issue, it is necessary to introduce a mesh reliability convergence with a reliable convergence criterion. To this aim, CVT is expanded in the next section.

Proposed expanded control variates (ECV) formulation

Suppose that G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{s}\mathrm{h}}^{\mathrm{F}\mathrm{E}\mathrm{A}}=G_N^{\mathrm{FEA}} is the performance of a system which obtained by a very fine mesh FE model and G_{N-\mathrm{1}}^{\mathrm{FEA}} is the performance which computed by a FE model with a coarse mesh density. By considering (G_{N-\mathrm{1}}^{\mathrm{FEA}}) as the control variates of \mathbb{I}\left(G_{N-\mathrm{1}}^{\mathrm{FEA}}\right), CVT aims to reformulate the failure probability integral as follows:

(6)P_f=\int \mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x+\int (\mathbb{I}\left(G_N^{\mathrm{FEA}}\right)-\mathbb{I}\left(G_N^{\mathrm{FEA}}\right))f\left(x\right)\mathrm{d}x\mathrm{.}

Accordingly, by employing linearly regressed CVT, Eq. (6) reads:

(7)P_f=\int (\mathbb{I}\left(G_N^{\mathrm{FEA}}\right)-\alpha ·\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))f\left(x\right)\mathrm{d}x+\alpha ·\mathrm{E}\left(\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)\mathrm{.}

In this approach, the optimum value of \alpha is determined by minimizing the variance of the estimation (mathematical expectation of Eq. (6)) with respect to \alpha [²⁵]:

\frac{\partial \mathrm{V}\mathrm{a}\mathrm{r}\left(\mathrm{E}\left(\widehat{p}\right)\right)}{\partial \alpha }=2\alpha ·\mathrm{V}\mathrm{a}\mathrm{r}\left(G_{\mathrm{c}\mathrm{o}\mathrm{a}\mathrm{r}\mathrm{s}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\left(x\right)\right)

+2\mathrm{C}\mathrm{o}\mathrm{v}(G_{\mathrm{c}\mathrm{o}\mathrm{a}\mathrm{r}\mathrm{s}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\left(x\right),(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\left(x\right)\left)\right)=0,

(8)\Rightarrow \alpha \frac{\mathrm{C}\mathrm{o}\mathrm{v}(\mathbb{I}\left(G_{\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right),\mathbb{I}\left(G_{\mathrm{c}\mathrm{o}\mathrm{a}\mathrm{r}\mathrm{s}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))}{\mathrm{V}\mathrm{a}\mathrm{r}\left(\mathbb{I}\left(G_{\mathrm{c}\mathrm{o}\mathrm{a}\mathrm{r}\mathrm{s}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)}\mathrm{.}

This study proposes to restate the last term of Eq. (7) as follows:

{\alpha }_1·\mathbb{E}\left(\mathbb{I}\left(G_{\mathrm{c}\mathrm{o}\mathrm{a}\mathrm{r}\mathrm{s}\mathrm{e}}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)={\alpha }_1({\alpha }_2·\int \mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x

+\int \mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_2·\mathbb{I}(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}f\left(x\right)\mathrm{d}x)

(9)={\alpha }_1·(\int \mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_2·\mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))f\left(x\right)\mathrm{d}x+{\alpha }_2·\mathbb{E}\left(\mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)\mathrm{.}

According to the proposed method, the CVT may be repeatedly used to break the obtained equations into more piecewise integral terms. Importance Sampling (IS) method is also employed to increase the efficiency. q(x) as an instrumental sampling PDF is introduced and the proposed ECV would be combined with the IS as follows:

P_f=\int \left(\mathbb{I}\left(G_N^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_1·\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)\frac{f\left(x\right)}{q\left(x\right)}q\left(x\right)\mathrm{d}x+{\alpha }_1·\int (\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_2·\mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))\frac{f\left(x\right)}{q\left(x\right)}q\left(x\right)\mathrm{d}x+{\alpha }_1·{\alpha }_2·\int (\mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_3·\mathbb{I}\left(G_{N-3}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))\frac{f\left(x\right)}{q\left(x\right)}q\left(x\right)\mathrm{d}x+\mathrm{\cdots }+{\alpha }_1·{\alpha }_2·\mathrm{\cdots }·{\alpha }_{N-2}·\int (\mathbb{I}\left(G_2^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_{N-1}·\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))\frac{f\left(x\right)}{q\left(x\right)}q\left(x\right)\mathrm{d}x

(10)+{\alpha }_1·{\alpha }_2·\mathrm{\cdots }·{\alpha }_{N-1}·{\alpha }_N·\int \mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x\mathrm{.}

By this approach, only a portion of estimations is engaged with a very fine mesh FE model and some coarse mesh FE models to reduce the computational cost. The outcome is a high efficient formulation for reliability analysis of problems involving FE. A Simplified ECV formulation can be defined as follow (based on Eqs. (8) and (10)):

(11)P_f={\alpha }_1·{\alpha }_2·\mathrm{\cdots }·{\alpha }_{N-1}·{\alpha }_N·\int \mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x�{\alpha }_1·{\alpha }_2·\mathrm{\cdots }·{\alpha }_{N-1}·{\alpha }_N·\mathbb{E}\left(\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)�\prod _{i=1}^N{\alpha }_i·\mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)\mathrm{.}

This equation means that it is enough to estimate the failure probability using a coarse mesh model and regulate it by means of correction factors {\displaystyle \prod _{i=\mathrm{1}}^N{\alpha }_i}.

Reliability mesh convergence analysis based on ECV

It was shown that ECV introduces chain-like control variates for the index function of original failure probability integral. In the Eq. (10), the ECV terms originate from \mathbb{E}\left(\mathbb{I}\left(G_N^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right) (involving FE model with a very fine mesh) and are terminated at \mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right) (involving FE model with a coarse mesh model). In the ECV, it was assumed that the model #N was the most suitable FE model for proper reliability analysis of system. In the simplified form, the proposed formulation regulated as P_{\mathrm{f}}={\displaystyle \prod _{i=1}^N}{\alpha }_i·\mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right) that just uses correction factor α to regulate the failure probability obtained by the coarse mesh model.

To propose an efficient method for reliability mesh convergence, the proposed strategy is contrariwise estimation of the ECV terms of Eq. (10), beginning from \mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right), continually to meet the original CVT equation as \int (\mathbb{I}\left(G_N^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_1·\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x) as follows:

P_f={\alpha }_1·{\alpha }_2·\cdots ·{\alpha }_{N-1}·{\alpha }_N·\mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)+{\alpha }_1·{\alpha }_2·\cdots ·{\alpha }_{N-2}·\int (\mathbb{I}\left(G_2^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_{N-1}·\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))\frac{f\left(x\right)}{q\left(x\right)}\mathrm{d}x+\cdots +{\alpha }_1·{\alpha }_2·\int (\mathbb{I}\left(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_3·\mathbb{I}\left(G_{N-3}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right))\frac{f\left(x\right)}{q\left(x\right)}\mathrm{d}x+{\alpha }_1·\int (\mathbb{I}\left(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_2·\mathbb{I}(G_{N-2}^{\mathrm{F}\mathrm{E}\mathrm{A}}\left)\right)\frac{f\left(x\right)}{q\left(x\right)}\mathrm{d}x

(12)+\int (\mathbb{I}\left(G_N^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)-{\alpha }_1·\mathbb{I}(G_{N-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\left)\right)\frac{f\left(x\right)}{q\left(x\right)}\mathrm{d}x\mathrm{.}

Equation (12) would present the following steps for simplified ECV:

To apply the proposed idea, the required steps to perform reliability mesh converge analysis are presented in the following subsections.

Initial failure probability using a coarse mesh model

In the proposed approach, the reliability analysis begins by employing the model #1 with a coarse mesh model. For the proposed mesh density, the reliability estimation would not be time-consuming and the crude MCS may be used to approximate the failure probability for high-dimensional problems:

(14)p_{\mathrm{f}}^1=p_{\mathrm{f}}^{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯1}=\int \mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)f\left(x\right)\mathrm{d}x=\mathbb{E}\left(\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\right)=\frac{{\displaystyle \mathrm{\sum }_{i=1}^{N_{\mathrm{M}\mathrm{C}\mathrm{S}}}}\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\left(i\right)\right)}{N_{\mathrm{M}\mathrm{C}\mathrm{S}}},

or alternatively the Weighted Average Simulation Method (WASM) may be used for small-size problems [²⁶,²⁷]:

(15)p_{\mathrm{f}}^{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯1}=\frac{{\displaystyle }\mathrm{\sum }_{i=1}^{N_{\mathrm{W}\mathrm{A}\mathrm{S}\mathrm{M}}}[{\displaystyle \prod _{i=\mathrm{1}}^{D}PDF\left(i\right)}·\mathbb{I}\left(G_1^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)\left(i\right)}{\mathrm{\sum }_{i=1}^{N_{\mathrm{W}\mathrm{A}\mathrm{S}\mathrm{M}}}\left({\displaystyle \prod _{i=\mathrm{1}}^{D}PDF\left(i\right)}\right)}

In these equations, N_{\mathrm{MCS}} and N_{\mathrm{WASM}} are the number of samples employed in MCS and WASM for reliability analysis, respectively, and D is the dimension of reliability problem. Other simulation methods (i.e., subset simulation [²⁸]) may be also used for reliability analysis in this step. Accordingly, Most Probable Point (MPP) may be simply determined by employing the definition presented in Ref. [²⁶] for simulation methods as the point with maximum weight in failure region. The location of MPP (L_MPP) may be simply determined as follows:

(16)L_{\mathrm{MPP}}=\mathrm{Find}\hspace{0.17em}\max \left({\displaystyle \prod _{i=\mathrm{1}}^{D}PDF\left(i\right)}\right)\mathrm{.}

Once the failure probability and MPP are estimated for model #1, a repetitive cycle should be performed until the method meets the convergence criterion which is presented in Section 4.2.

Convergence cycle

Suppose model #i−1 and P_f^{i-\mathrm{1}} are the employed mesh model and the failure probability approximated in the former step, respectively. A new failure probability (P_f^i) based on the newly refined FE model (model #i), may be estimated as follows:

(17)P_f^i={\alpha }_{i-\mathrm{1}}·P_f^{i-\mathrm{1}},

where

(18){\alpha }_{i-1}=\frac{{\mathbb{E}}_q\left(\mathbb{I}\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right){\displaystyle \frac{f\left(x\right)}{q\left(x\right)}}\right)}{{\mathbb{E}}_q\left(\mathbb{I}\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right){\displaystyle \frac{f\left(x\right)}{q\left(x\right)}}\right)}\mathrm{.}

In this equation, q(x) is the instrumental PDF with the mean of MPP determined in the previous step by model #i−1. The proposed step should be repetitively performed until the method meets model #N that corresponds to the best-refined mesh model for reliability analysis. To find the qualified FE model, a new reliability convergence criterion is proposed in the following subsection.

Proposed reliability mesh convergence criterion

At the initial steps of mesh convergence analysis, the differences between performance valuations of coarse and refined mesh densities may be calculated rationally large values. However, as the number of analyzed models increases and more refined models are used in analyses, the differences in performance valuation of the two models become smaller. Supposing two employed FE models #i and i-\mathrm{1} presented high correlation (see Eq. (8)), one may reformulate ECV as:

p_{\mathrm{f}}^i={\alpha }_{i-1}·P_{\mathrm{f}}^{i-1}+(\frac{1}{N}\mathrm{\sum }_{k=1}^N\mathbb{I}\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)·\frac{f\left(x\right)}{q\left(x\right)}

(19)-\frac{{\alpha }_{i-1}}{N}\mathrm{\sum }_{k=1}^N\mathbb{I}\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i-1}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)·\frac{f\left(x\right)}{q\left(x\right)},

(20)P_{\mathrm{f}}^i-{\alpha }_{i-1}·P_{\mathrm{f}}^{i-1}=(1-{\alpha }_{i-1})·\left(\frac{1}{N}\mathrm{\sum }_{k=1}^N\mathbb{I}\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)·\frac{f\left(x\right)}{q\left(x\right)}\right)\mathrm{.}

Accordingly, using correction factor, an scaled- convergence criterion may derive from Eq. (20) as follows:

(21)|1-{\alpha }_{i-1}|=\frac{P_{\mathrm{f}}^i-{\alpha }_{i-1}·P_{\mathrm{f}}^{i-1}}{{\displaystyle \frac{1}{N}}g\left(G_{\mathrm{m}\mathrm{o}\mathrm{d}\mathrm{e}\mathrm{l}♯i}^{\mathrm{F}\mathrm{E}\mathrm{A}}\right)·{\displaystyle \frac{f\left(x\right)}{q\left(x\right)}}}\le {ϵ}_{\alpha }\mathrm{.}

A common deterministic mesh convergence evaluates the accuracy of FE model for a single sample (i.e., at the mean point of random variables), while the proposed criterion evaluates the accuracy of FE model using the complete results of simulation process (see Eqs. (8) and (18)).

For easy comprehension, Fig. 1 shows the proposed algorithm.

Engineering examples

The reliability-based FE mesh convergence analysis of four engineering problems are investigated. Reliability analysis in the first step (model #1) is performed by WASM which provides failure probability and MPP with suitable accuracy.

Clamped non-prismatic beam

In this example, reliability of a clamped beam with the rectangular section is investigated (see Fig. 2), where length (L), width (B), height of start section (d_s), height of end section (d_e), the applied load intensity (q), and modulus of elasticity (E) have been selected as random variables and the high of the section in the span of the beam has been considered linear. Table 1 presents the statistical parameters of the mentioned structure and the performance function of the problem is defined as below:

(22)G=y_{\mathrm{a}}-y,

in which y_{\mathrm{a}} is the maximum allowable deflection and y is the maximum deflection that can be calculated by analytical solution (y_{\mathrm{analytical}}) or FEA (y_{\mathrm{FEA}}).

To perform FEA of this structure, the beam was divided into the same lengths and geometrical properties of the section for each element have been considered as the mean of start and end of elements. Double integral method has been applied as the theoretical solution of this problem to calculate deflection of the beam [²⁹].

Figure 3 shows the deterministic FE mesh convergence of the structure at the mean value of random variables. For safety evaluation, reliability mesh convergence analysis was performed using the proposed approach in which the FE model with two elements was selected as the first coarse mesh model (model #1). In refining purpose, one element was added to the coarse model in each step in order to increase the FE accuracy. Figure 4 presents the reliability mesh convergence of the problem. The figure shows that the reliability index (β) of the problem converged to the \beta =\mathrm{3.98} that it is highly in agreement with the reliability results which obtained by theoretical solution. It should be noted that only 100 samples employed in each cycle; therefore, the reliability analysis is conducted very efficiently.

Figure 5 shows the value of the proposed correction factor in each step (each number step is related to (n−1)th and nth mesh density models). This figure indicates that the proposed convergence criterion in Eq. (21) is satisfied in Step 8 (related to 9 elements as coarse mesh model and 10 elements as fine mesh model in the proposed approach) that is in good agreement with the result of theoretical solution and confirms the efficiency of the proposed convergence criterion.

Figure 6 also shows the normal plot of performance function G for different mesh densities. In this figure, among 14 FE models presented in Fig. 6, four FE models have been selected (number of elements: 2, 6, 10, and 14). This figure indicates that models with 10 and 14 numbers of elements are in a good agreement. In addition, this figure shows that the first step of the proposed approach requiring huge samples for a low-cost coarse mesh model and other mesh refinement steps require very small size for proper safety evaluation of structure. In the last FE models, the number of failure samples (G\le \mathrm{0} ) increased to prevent the wasted samples that are generally time-consuming FE analyses. For instance, FEA-based reliability evaluation of this structure takes 5000 separate FEA for fine mesh density to calculate failure probability using WASM as an efficient reliability method (definitely, reliability analysis using crude Monte Carlo Simulation requires huge number of FEA), while the proposed approach saves computational effort by only 100 separate FEAs for each FE model.

Non-prismatic column

In this example, reliability analysis of a non-prismatic column is investigated (see Fig. 7). The performance function for the problem is defined as follows:

(23)G=P_{\mathrm{C}}-P,

where P is the applied load and P_{\mathrm{C}} is the critical load.

In this example, flange thickness (t_f), web thickness (t_w), web height (h) at the start (i) and end point (j), the length of the column (L), modulus of elasticity (E), and the applied load (P) are considered random variables. The variation of web high is considered linear in the length.

Table 2 shows the statistical parameters of this problem. It is worthy to be noted that the proposed approach is not under influence of distribution type because the proposed approach is a simulation based technique and it is possible to transform any distributions to the standard normal space [¹,²⁵–²⁸]. In addition, it was assumed that the elastic buckling occurred in the structure for calculating the critical load using FE. Bifurcation concept was employed to calculate the critical load. Geometrical stiffness of elements were calculated based on equations presented in Ref. [³⁰]. The critical load was also calculated by solving specific value problem and total determinant of stiffness matrix equal to zero. In FEA, the non-prismatic column is divided into the same length, and geometrical properties of each element are considered as the mean of start and end points. In addition, FEA of this structure was performed at the mean points which are presented in Table 2.

Different numbers of elements including 2, 4, 6, 8, 10, and 12 are used for mesh-refinement purpose to evaluate the structure safety using the proposed method.

Figure 8 shows the accuracy and consuming time by each FE model to achieve the proper critical load (P_{\mathrm{C}}) at the mean point of random variables. This figure confirms that the finer meshes need more computational time and computational effort increased, dramatically. This finding is a clear explanation for high efficiency of the proposed approach that requires very few samples during the refinement process.

Figures 9 and 10 show the convergence of reliability index for different mesh densities and the convergence of the defined correction factor, respectively. Figure 9 shows that the FE model with eight elements is enough to be converged appropriately, and it is also proved by Fig. 10.

Single edge notch test specimen

The reliability analysis of a rectangular plate with 60 mm ×30 mm dimensions and the single edge crack by 15 mm is investigated. Figure 11 shows the specimen. For this problem, an analytical formulation for determination of the stress intensity factor is presented as [³¹]:

K_{\mathrm{I}}=t\sqrt{\mathrm{\pi }a}\sqrt{{\displaystyle \frac{\mathrm{2}B}{\mathrm{\pi }a}}\tan \left({\displaystyle \frac{\mathrm{\pi }a}{\mathrm{2}B}}\right)}

(24)·{\displaystyle \frac{\mathrm{0.752}+\mathrm{2.02}\left({\displaystyle \frac{a}{B}}\right)+\mathrm{0.37}{\left(\mathrm{1}-\sin \left({\displaystyle \frac{\mathrm{\pi }a}{\mathrm{2}B}}\right)\right)}^{\mathrm{3}}}{\cos \left({\displaystyle \frac{\mathrm{\pi }a}{\mathrm{2}B}}\right)}}\mathrm{.}

The FEA of the problem to obtain stress intensity was conducted in different mesh densities (5 × 10, 25 × 50, 50 × 100, 75 × 150, and 150 × 300) and then compared with the analytical solution presented in Eq. (24).

For reliability evaluation of this problem, description of basic random variables is presented in Table 3 and the limit state function of this problem may be presented as:

(25)G=K_{\mathrm{IC}}-K_{\mathrm{I}}.

In the above formula, K_{\mathrm{IC}} is the critical stress intensity factor and K_{\mathrm{I}} is the existent stress intensity factor.

Table 3 presents the convergence of FE models at mean point, and Fig. 12 shows the relevant scaled CPU time. This figure shows that employing very fine mesh model for reliability analysis would be highly time-consuming. For this problem, the required time to solve problem by the fine mesh model (150 × 300 element) is approximately 3650 times more than the required time for coarse mesh model (5 × 10). The problem is solved by the proposed approach which its results are shown in Figs. 13 and 14. These figures indicate that employing the proposed approach would lead to proper approximation of failure probability with suitable mesh convergence criterion. Accurate reliability analysis by common approach required huge simulations with very fine mesh model (150 × 300 elements), while solving this problem by the proposed approach required evaluation of only 100 FEA with the very fine model (150 × 300 elements), which is confirming high efficiency of the method for solving FE-based reliability problems.

The nonlinear simply supported beam

In this example, reliability of a simply supported beam with nonlinear material as shown in Fig. 15 is investigated. A hypoelastic material with the following stress-strain formulation [³²] is considered as the material of the beam:

(26){\displaystyle \frac{{\sigma }_{\mathrm{e}}}{{\sigma }_{\mathrm{0}}}}=\{\begin{array}{ll}\sqrt{{\displaystyle \frac{\mathrm{1}+n^{\mathrm{2}}}{{\left(n-\mathrm{1}\right)}^{\mathrm{2}}}}-{\left({\displaystyle \frac{n}{n-\mathrm{1}}}-{\displaystyle \frac{{ϵ}_{\mathrm{e}}}{{ϵ}_{\mathrm{0}}}}\right)}^{\mathrm{2}}}-{\displaystyle \frac{\mathrm{1}}{n-\mathrm{1}}},\hfill & {ϵ}_{\mathrm{e}}\le {ϵ}_{\mathrm{0}},\hfill \\ \hfill \\ {\left({\displaystyle \frac{{ϵ}_{\mathrm{e}}}{{ϵ}_{\mathrm{0}}}}\right)}^{\mathrm{1}/n},\hfill & {ϵ}_{\mathrm{e}}>{ϵ}_{\mathrm{0}}\mathrm{.}\hfill \end{array}

ƒSpan length (L), width (W), the distributed load (q) and parameters related to the stress-strain curve ({\sigma }_{\mathrm{0}}, {ϵ}_{\mathrm{0}}, n) are considered as random variables and the performance function of the problem is defined as follow for reliability analysis:

(27)G=u_{\mathrm{a}}-u,

in which u_{\mathrm{a}} is the allowable deflection (it is considered 13 cm in this study) and is the maximum deflection that is obtained by FEA. Description of random variables are presented in Table 4.

Six different mesh densities have been selected for FEA and results are shown in Fig. 16. It is shown that the coarse and the finest mesh models (30 × 5 and 66 × 11 elements, respectively) have 9.74% error in results while employing fine mesh model in FEA takes 5.4 times more CPU time in comparison with the coarse mesh model.

Figure 17 presents the reliability indices which are obtained by the proposed approach. The reliability index of problem from b = 3.27 turned to b = 2.91 by using the proposed approach. In this example, accuracy of the result and the saved computational cost confirms the capability of the proposed method to handle problems with nonlinear material. Figure 18 shows the convergence of the proposed correction to one value.

Conclusions

Reliability-based mesh convergence analysis was an overlooked topic in FE-based reliability problems. It was often impractical for FE problems, since for each FE model, a separate reliability analysis was required. Accordingly, it was regular to accept “unknown errors” in reliability valuation that resulted from “known errors” distinguished by the deterministic mesh convergence analysis. In this study, ECV is introduced for an efficient reliability mesh convergence approach using inversely computing the proposed ECV terms. Based on the proposed formulation, only a portion of reliability estimation is engaged with very fine mesh models in FEA and the main portion of computations required employing models with coarse meshes. In addition, the proposed strategy provides a step-by-step mesh convergence history of FE models in reliability problems to reduce the computational cost. It becomes more desirable by recording the reliability results at each step, instead of performing a separate reliability analysis, obtained from upgrading former models results. Four practical engineering problems solved by the proposed approach. For all solved examples, it was demonstrated that the proposed approach not only provided a useful mesh convergence for the problem, but also considerably reduced the computational cost of reliability problems.

References

[1] Nowak A S, Collins K R. Reliability of Structures. New York: McGraw-Hill, 2000

[2] Gallimard L. Error bounds for the reliability index in finite element reliability analysis. International Journal for Numerical Methods in Engineering, 2011, 87(8): 781–794

[3] Ghavidel A, Rashki M, Mousavi S R. The effect of FEM mesh density on the failure probability analysis of structures. KSCE Journal of Civil Engineering, 2018, 22(7): 2371–2383

[4] Rashki M, Ghavidel A, Ghohani Arab H, Mousavi S R. Low-cost finite element method-based reliability analysis using adjusted control variate technique. Structural Safety, 2018, 75: 133–142

[5] Javed A, Djijdeli K, Xing J T. A coupled meshfree-mesh-based solution scheme on hybrid grid for flow-induced vibrations. Acta Mechanica, 2016, 227(8): 2245–2274

[6] Li Z, Wood R. Accuracy verification of a 2D adaptive mesh refinement method for incompressible or steady flow. Journal of Computational and Applied Mathematics, 2017, 318: 259–265

[7] Morgan N R, Waltz J I. 3D level set methods for evolving fronts on tetrahedral meshes with adaptive mesh refinement. Journal of Computational Physics, 2017, 336: 492–512

[8] Anitescu C, Hossain M N, Rabczuk T. Recovery-based error estimation and adaptivity using high-order splines over hierarchical T-meshes. Computer Methods in Applied Mechanics and Engineering, 2018, 328: 638–662

[9] Hennig P, Kästner M, Morgenstern P, Peterseim D. Adaptive mesh refinement strategies in isogeometric analysis—A computational comparison. Computer Methods in Applied Mechanics and Engineering, 2017, 316: 424–448

[10] Ray N, Grindeanu I, Zhao X, Mahadevan V, Jiao X. Array-based, parallel hierarchical mesh refinement algorithms for unstructured meshes. Computer Aided Design, 2017, 85: 68–82

[11] Bespalov A, Haberl A, Praetorius D. Adaptive FEM with coarse initial mesh guarantees optimal convergence rates for compactly perturbed elliptic problems. Computer Methods in Applied Mechanics and Engineering, 2017, 317: 318–340

[12] Rabczuk T, Ren H, Zhuang X. A nonlocal operator method for partial differential equations with application to electromagnetic waveguide problem. Computers, Materials & Continua, 2019, 59(1): 31–55

[13] Anitescu C, Atroshchenko E, Alajlan N, Rabczuk T. Artificial neural network methods for the solution of second order boundary value problems. Computers, Materials & Continua, 2019, 59(1): 345–359

[14] Guo H, Zhuang X, Rabczuk T. A deep collocation method for the bending analysis of Kirchhoff plate. Computers, Materials & Continua, 2019, 59(2): 433–456

[15] Weiber S. Residual based error estimate and quasi-interpolation on polygonal meshes for high order BEM-based FEM. Computers & Mathematics with Applications (Oxford, England), 2017, 73(2): 187–202

[16] Wu J, Zheng H. Uniform convergence of multigrid methods for adaptive meshes. Applied Numerical Mathematics, 2017, 113: 109–123

[17] Liu G R. The smoothed finite element method (S-FEM): A framework for the design of numerical models for desired solutions. Frontiers of Structural and Civil Engineering, 2019, 13(2): 456–477

[18] Brambilla P, Brivio P, Guardone A, Romanelli G. Grid convergence assessment for adaptive grid simulations of normal drop impacts onto liquid films in axi-symmetric and three-dimensional geometries. Applied Mathematics and Computation, 2015, 267: 487–497

[19] Lockard D P. In search of grid converged solutions. Procedia Engineering, 2010, 6: 224–233

[20] Płaszewski P, Banaś K. Performance analysis of iterative solvers of linear equations for Hp-adaptive finite element method. Procedia Computer Science, 2013, 18: 1584–1593

[21] Jones A C, Wilcox R K. Finite element analysis of the spine: Towards a framework of verification, validation and sensitivity analysis. Medical Engineering & Physics, 2008, 30(10): 1287–1304

[22] Hamdia K M, Ghasemi H, Zhuang X, Alajlan N, Rabczuk T. Sensitivity and uncertainty analysis for flexoelectric nanostructures. Computer Methods in Applied Mechanics and Engineering, 2018, 337: 95–109

[23] Liu P L, Liu K G. Selection of random field mesh in finite element reliability analysis. Journal of Engineering Mechanics, 1993, 119(4): 667–680

[24] Manjuprasad M, Manohar C S. Adaptive random field mesh refinements in stochastic finite element reliability analysis of structures. Computer Modeling in Engineering & Sciences, 2007, 19(1): 23–54

[25] Tracey B, Wolper D. Reducing the error of Monte Carlo algorithms by learning control variates. In: The 29th Conference on Neural Information Processing Systems (NIPS). Barcelona, 2016

[26] Rashki M, Miri M, Azhdary Moghaddam M. A new efficient simulation method to approximate the probability of failure and most probable point. Structural Safety, 2012, 39: 22–29

[27] Arab H G, Rashki M, Rostamian M, Ghavidel A, Shahraki H. Refined first-order reliability method using cross-entropy optimization method. Engineering with Computers, 2018, 35: 1–13

[28] Green D K E. Efficient Markov Chain Monte Carlo for combined subset simulation and nonlinear finite element analysis. Computer Methods in Applied Mechanics and Engineering, 2017, 313: 337–361

[29] Megson T H G. Structural and Stress Analysis. 2nd ed. Amsterdam/London: Elsevier Butterworth-Heinemann, 2005

[30] Clough R W, Penzien J. Dynamics of Structures. 2nd ed. New York: McGraw-Hill, 1993

[31] Tada H, Paris P C, Irwin G R. The Stress Analysis of Cracks Handbook. New York: ASME press, 2000

[32] Bower A F. Applied Mechanics of Solids. Boca Raton: CRC Press, Taylor & Francis, 2009

RIGHTS & PERMISSIONS

Higher Education Press