Applying the spectral stochastic finite element method in multiple-random field RC structures

Abstract

This paper uses the spectral stochastic finite element method (SSFEM) for analyzing reinforced concrete (RC) beam/slab problems. In doing so, it presents a new framework to study how the correlation length of a random field (RF) with uncertain parameters will affect modeling uncertainties and reliability evaluations. It considers: 1) different correlation lengths for uncertainty parameters, and 2) dead and live loads as well as the elasticity moduli of concrete and steel as a multi-dimensional RF in concrete structures. To show the SSFEM’s efficiency in the study of concrete structures and to evaluate the sensitivity of the correlation length effects in evaluating the reliability, two examples of RC beams and slabs have been investigated. According to the results, the RF correlation length is effective in modeling uncertainties and evaluating reliabilities; the longer the correlation length, the greater the dispersion range of the structure response and the higher the failure probability.

Graphical abstract

Keywords

uncertainty / spectral stochastic finite element method / correlation length / reliability assessment / reinforced concrete beam/slab

1 Introduction

In structural engineering, reliability assessment shows its importance when parameters involved in the analysis and design of a structural system (material properties, loading, geometry, etc.) are uncertain. This is, in fact, a well-known theoretical method of calculating the failure probability (P_f), performed through a multipurpose probabilistic integral. If f(x) and g are assumed to be, respectively, the common probability density and limit state functions (LSF) of random variables X, then P_f can be defined as follows [1–3]:

(1)P_{\mathrm{f}}=prob\left[g\left(X\right)⩽0\right]={\int }_{D_{\mathrm{f}}}^{}f\left(x\right)\mathrm{d}x,

where prob is a measure of probability and D_f is the failure domain (g(X)≤0). Since Eq. (1) is complicated for solving engineering problems with high uncertainty dimensions, various methods have been proposed for this purpose among which first-order reliability method (FORM) and second-order reliability method (SORM) are used for analysis [4,5] and Monte Carlo simulation (MCS), importance sampling (IS) and line sampling (LS) are for simulations [6–9]. However, these methods are quite costly when solving the performance function is time consuming (e.g., problems governed by partial differential equations (PDEs) and analyzed by the finite element method (FEM)).

Analysis of a structure, as an engineering system, with uncertain parameters is done by a set of stochastic partial differential equations (SPDEs), the analytical solutions of which are generally difficult; hence, many numerical approaches, developed based on the classical FEM and called “stochastic finite element method (SFEM)”, have been proposed to provide a solution. The MCS, spectral stochastic finite element method (SSFEM), Newman expansion method, stochastic response surface methodology (SRSM) and Taylor perturbation expansion method, to mention only a few, are explained in detail in Refs. [10–14], respectively.

SSFEM, that solves SPDE-based problems with less computational costs and sufficient accuracy, is a very efficient tool for measuring uncertainties in various engineering fields [15–19] and in evaluating structures’ reliabilities [20,21]. It is based on the discretization of the random field (RF) (related to the uncertainties of each parameter) using the Karhunen-Loeve (KL) expansion and estimation of response using polynomial chaos (PC). It presents a wide range of solutions for the structure’s nodal points in terms of random variables with only one analysis and, finally, yields their PDFs and CDFs by generating random samples based on their distribution functions. The P_f reliability is then easily calculable by simulation methods such as MCS [10,22].

Given that, in most studies to date, the RF correlation length has been hypothetically selected, this paper aims to study the effects of the RF correlation length on parameters’ uncertainty modeling and reliability evaluations using the SSFEM. This method was developed for reinforced concrete (RC) beams/slabs as a system with multi-dimensional RF, considering materials’ properties and dead and live loads as uncertain parameters. Such uncertainty has been ignored in SSFEM studies done so far, and some studies have investigated only the effect of correlation length in the uncertainty modeling of earth structures using the fuzzy probability method [23,24].

To achieve the goals, the rest of this paper has been organized to: describe the problem in Section 2; review the SSFEM in Section 3; develop the method for RC beams/slabs in Section 4; solve two numerical examples, examine the method and check the effects of the correlation length in uncertainty modeling and reliability evaluations in Section 5; and finally, provide conclusions in Section 6.

2 Problem description

This research considers dead and live loads and material properties as important parameters in the uncertainty modeling of concrete structures; however, to reduce calculations, the correlation length effects in uncertainty modeling and reliability evaluations are studied on dead and live loads, which are more sensitive. Fig.1 shows, for more clarity, a schematic view of loading and material properties of RC beams/slabs.

The uncertainty fluctuations of the dead and live loads are, in fact, more important for space than for properties of materials in the structures because the live load can, like a person, have variable weights ​​in different spaces, and the dead load depends on the load of such non-structural components as partition walls and flooring as well as on the weights of material structural members such as beams and slabs.

This, therefore, increases the incentive to study the dead and live loads’ correlation length sensitivity in modeling uncertainty fluctuations along the spatial domain and evaluating reliability because it is a parameter involved in RF modeling. In addition, to obtain a more accurate and realistic model of the structure, uncertainty effects of material properties have also been considered assuming the correlation length to be equal to the field dimensions, as in previous studies [11,15,21].

3 A review of the spectral stochastic finite element method

The SSFEM general framework is based on the RF discretization using the KL-expansion, response approximation using PC and FEM to find the solution [10,11].

3.1 Discretization through the KL-expansion

Consider a complete probability space (Ω, Ƒ, P) where Ω is the sample space and Ƒ is the σ-algebra on subsets Ω and P of the probability measurement. A continuous RF H(x, θ) on domain D\subset R is an operator function (H: D × Ω → R) that has a space dimension x ∈ D and a stochastic nature θ ∈ Ω. If its mean and standard deviation are µH (x) and σH (x), respectively, it can be discretized as follows using M terms from the KL-expansion [11,22]:

(2)H\left(x,\theta \right)={\mu }_{\mathrm{H}}\left(x\right)+{\sigma }_{\mathrm{H}}\left(x\right)\sum _{n=1}^M\sqrt{{\lambda }_n}{\varphi }_n\left(x\right){\xi }_n\left(\theta \right),x\in D,\theta \in \mathrm{\Omega },

where ξ_n (θ) is a set of independent random variables with mean = 0 and standard deviation = 1, λ_n and ϕ_n(x) are, respectively, the eigenvalues and eigenfunctions calculable as follows using the Fredholm quadratic integral:

(3){\int }_D^{}{Cov}_{\mathrm{H}}\left(x_1‚x_2\right){\varphi }_n\left(x_2\right)={\lambda }_n{\varphi }_n\left(x_1\right),x_1,x_2\in D,

where Cov_H is the covariance function of the RF and ϕ_n(x) should satisfy the following orthogonality properties.

(4){\int }_D^{}{\varphi }_n\left(x\right){\varphi }_m\left(x\right)\mathrm{d}{x}_1={\delta }_{nm},\{\begin{array}{c}{\delta }_{nm}=1,n=m,\\ {\delta }_{nm}=0,n\ne m,\end{array}

where δ_nm represents the Kronecker delta function.

Since the usually positive, symmetric, bounded RF covariance function (Cov_H) is unknown in engineering problems, the exponential function, a very well-known covariance function, is used in this research to specify the field properties. For one- and two-dimensional problems, it is shown in Eqs. (5) and (6), respectively [22], as follows:

(5){Cov}_{\mathrm{H}}\left(x_1‚x_2\right)=\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{\left|x_1-x_2\right|}{b}\right),

(6){Cov}_{\mathrm{H}}\left(xy\right)=\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{\left|x_1-x_2\right|}{b_x}-\frac{\left|y_1-y_2\right|}{b_y}\right),

where x and y are the dimensional coordinates of a point on the field and b is the correlation length. Hence, eigenvalues and eigenfunctions on the system domain depend on the coordinates and correlation length.

3.2 Response approximation through polynomial chaos

In analyzing a structure with SSFEM, its responses are a set of the nodes’ displacements’ stochastic vectors U(θ) the covariance function of which is unknown; the one proposed in KL-expansion cannot be used for this purpose, but U(θ) can be expressed in terms of ξ_i(θ) as a nonlinear function, the dependence of which is as follows [22,25,26]:

(7)\begin{array}{rl}U\left(\theta \right)=& {\alpha }_0{\Gamma }_0+\sum _{i_1=1}^{\mathrm{\infty }}{\alpha }_{i_1}{\Gamma }_1\left({\xi }_{i_1}\left(\theta \right)\right)\\ & +\sum _{i_1=1}^{\mathrm{\infty }}\sum _{i_2=1}^{i_1}{\alpha }_{i_1{i}_2}{\Gamma }_2\left({\xi }_{i_1}\left(\theta \right)‚{\xi }_{i_2}\left(\theta \right)\right)\\ & +\sum _{i_1=1}^{\mathrm{\infty }}\sum _{i_2=1}^{i_1}\sum _{i_3=1}^{i_2}{\alpha }_{i_1{i}_2{i}_3}{\Gamma }_3\left({\xi }_{i_1}\left(\theta \right),{\xi }_{i_2}\left(\theta \right),{\xi }_{i_3}\left(\theta \right)\right)+\cdots ,\end{array}

where {\alpha }_0,{\alpha }_{i_1},{\alpha }_{{i_{1}i}_2},...,{\alpha }_{{i_{1}i}_2\cdots i_n} are known coefficients and Γ_p is a p-order PC found as follows:

(8){\Gamma }_p\left({\xi }_{i_1}\left(\theta \right),\dots ,{\xi }_{i_n}\left(\theta \right)\right)={(-1)}^n\frac{{\mathrm{\partial }}^n}{\mathrm{\partial }{\xi }_{i_1}\left(\theta \right)\cdots {\mathrm{\partial }\xi }_{i_n}\left(\theta \right)}\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{1}{2}{\mathit{\xi }}^{\mathrm{T}}\mathit{\xi }\right),

where ξ is the vector of random variable ξ_i(θ) (\mathit{\xi }=({\xi }_{i_1}\left(\theta \right),...,{\xi }_{i_n}\left(\theta \right)). For simplicity, Eq. (7) is rewritten as follows:

(9)\begin{array}{rl}U\left(\theta \right)& =u_0{\psi }_0(\xi \left(\theta \right))+u_1{\psi }_1(\xi \left(\theta \right))+u_2{\psi }_2(\xi \left(\theta \right))+\cdots \\ & =\sum _{j=0}^{\mathrm{\infty }}{u}_j{\psi }_j(\xi \left(\theta \right)),\end{array}

where u_j is a coefficient corresponding to {\alpha }_0,{\alpha }_{i_1},{\alpha }_{{i_{1}i}_2},..., {\alpha }_{{i_{1}i}_2\cdots i_n}; and {\psi }_j\left(\theta \right) is a coefficient corresponding to orthogonal polynomial Γ_p according to ξ_i(θ) for which the following conditions should be satisfied:

(10)\begin{array}{rl}& {\psi }_0\left(\theta \right)=1,\\ & ⟨{\psi }_j\left(\theta \right)⟩=0,\\ & ⟨{\psi }_j\left(\theta \right){\psi }_k\left(\theta \right)⟩=⟨{\psi }_k^2\left(\theta \right)⟩{\delta }_{jk},\end{array}

where ⟨\cdot ⟩ is the mathematical expectation, the value of which for internal multiplication in Hilbert space is found as follows:

(11)⟨{\psi }_j\left(\theta \right){\psi }_k\left(\theta \right)⟩=\int {\psi }_j\left(\theta \right){\psi }_k\left(\theta \right)W\left(\xi \left(\theta \right)\right)\mathrm{d}\xi \left(\theta \right),

where W(ξ(θ)) is a PC-related weight function as follows:

(12)W\left(\xi \right(\theta \left)\right)=\frac{1}{\sqrt{2\text{π}}}\mathrm{e}\mathrm{x}\mathrm{p}\left(-\frac{1}{2}{\xi \left(\theta \right)}^{\mathrm{T}}\xi \left(\theta \right)\right).

To reduce computation, Eq. (9) is truncated by P terms; P is determined by P=\left(\begin{array}{c}M+p\\ p\end{array}\right)+1 where M is the number of random variables and p is the order of PC.

3.3 Formulation

In SSFEM, the general form of FE relations governing problems, where material properties and loads are both stochastic, is as follows [22]:

(13)K\left(x,\theta \right)U\left(x,\theta \right)=F\left(x,\theta \right),x\in D,\theta \in \mathrm{\Omega },

where U(x,θ), K(x,θ), and F(x,θ) are the system response, total stiffness and total applied load matrices, and the latter two are found by assembling the stiffness and element force matrices, respectively.

Considering characteristic matrix D and loading q as two separate RFs for modeling problems where material properties and loading are both stochastic. We have following relations, using the KL-expansion:

(14)D\left(x,\theta \right)={\mu }_{\mathrm{E}}\left(x\right)D_0+{\sigma }_{\mathrm{E}}\left(x\right)\sum _{n=1}^{M_1}\left(\sqrt{{\lambda }_n}{\varphi }_n\left(x\right){\xi }_n\left(\theta \right)\right)D_0,

(15)q\left(x,\theta \right)={\mu }_{\mathrm{q}}\left(x\right)+{\sigma }_{\mathrm{q}}\left(x\right)\sum _{m=1}^{M_2}\sqrt{{\lambda }_m}{\varphi }_m\left(x\right){\xi }_m\left(\theta \right).

Using Eqs. (14), (15), and FEM, each element’s stiffness matrix in 2D problems can be found as follows:

(16)k^{\left(\mathrm{e}\right)}=t\int {\int }_{D^{\left(\mathrm{e}\right)}}^{}{B}^{\mathrm{T}}D(x,\theta )B\mathrm{d}A,

where B is the transfer matrix and t is the element thickness. Substituting Eq. (14) in Eq. (16) will yield:

(17)\begin{array}{rl}{k}^{\left(\mathrm{e}\right)}\left(x,\theta \right)=& t\int {\int }_{D^{\left(\mathrm{e}\right)}}^{}{B}^{\mathrm{T}}[{\mu }_{\mathrm{E}}\left(x\right)D_0\\ & +{\sigma }_{\mathrm{E}}\left(x\right)\sum _{n=1}^{M_1}\left(\sqrt{{\lambda }_n}{\varphi }_n\left(x\right){\xi }_n\left(\theta \right)\right)D_0]B\mathrm{d}A\\ =& \underset{k_0^{\left(\mathrm{e}\right)}}{\underbrace{t\int {\int }_{D^{\left(\mathrm{e}\right)}}^{}{B}^{\mathrm{T}}{\mu }_{\mathrm{E}}\left(x\right)D_{0}B\mathrm{d}A}}\\ & +\sum _{n=1}^{M_1}{\xi }_n\left(\theta \right)\underset{k_n^{\left(\mathrm{e}\right)}}{\underbrace{t\int {\int }_{D^{\left(\mathrm{e}\right)}}^{}{\sigma }_{\mathrm{E}}\left(x\right)\sqrt{{\lambda }_n}{\varphi }_n\left(x\right)B^{\mathrm{T}}{D}_{0}B\mathrm{d}A}},\end{array}

where k_0^{\left(\mathrm{e}\right)} and k_n^{\left(\mathrm{e}\right)} are the elements’ mean and weight stiffness matrices, respectively, the assembling of which will yield the total stiffness matrix as follows:

(18)K\left(X,\theta \right)=\left[K_0\left(x\right)+\sum _{n=1}^{M_1}{K}_n\left(x\right){\xi }_n\left(\theta \right)\right].

But, the applied force matrix of each element can be found as follows:

(19)f^{\left(\mathrm{e}\right)}\left(x,\theta \right)=\int Nq\left(x,\theta \right)\mathrm{d}x,

where N is the element’s shape function. Placing Eq. (15) in Eq. (19) will yield:

(20)\begin{array}{rl}{f}^{\left(\mathrm{e}\right)}\left(x‚\theta \right)=& \int N{\mu }_{\mathrm{q}}\left(x\right)\mathrm{d}x\\ & +\int N{\sigma }_{\mathrm{q}}\left(x\right)\left(\sum _{m=1}^{M_2}\sqrt{{\lambda }_m}{\varphi }_m\left(x\right){\xi }_m\left(\theta \right)\right)\mathrm{d}x\\ =& \underset{f_0^{\left(\mathrm{e}\right)}}{\underbrace{\int N{\mu }_{\mathrm{q}}\left(x\right)\mathrm{d}x}}+\sum _{m=1}^{M_2}{\xi }_m\left(\theta \right)\left(\underset{f_m^{\left(\mathrm{e}\right)}}{\underbrace{\int N{\sigma }_{\mathrm{q}}\left(x\right)\sqrt{{\lambda }_m}{\varphi }_m\left(x\right)\mathrm{d}x}}\right)\\ =& f_0^{\left(\mathrm{e}\right)}+\sum _{m=1}^{M_2}{f}_m^{\left(\mathrm{e}\right)}{\xi }_m\left(\theta \right).\end{array}

If ξ₀(θ) = 1 for m=1:

(21)f^{\left(\mathrm{e}\right)}\left(x,\theta \right)=\sum _{m=0}^{M_2}{f_m^{\left(\mathrm{e}\right)}\xi }_m\left(\theta \right).

Assembling the elements’ force matrix will yield the system’s total force matrix as follows:

(22)F\left(x,\theta \right)=\sum _{m=0}^{M_2}{F_m\xi }_m\left(\theta \right).

Substituting Eqs. (9), (18), and (22) in Eq. (13) will yield the general form of SSFEM as follows:

\left[K_0\left(x\right)+\sum _{n=1}^{M_1}{K}_n\left(x\right){\xi }_n\left(\theta \right)\right]\sum _{j=0}^{P-1}{u}_j{\psi }_j\left[\xi \left(\theta \right)\right]=\sum _{m=0}^{M_2}{F_m\xi }_m\left(\theta \right),

(23)\begin{array}{rl}& \sum _{j=0}^{P-1}{K_0\left(x\right)u}_j{\psi }_j\left[\xi \left(\theta \right)\right]\\ & +\sum _{j=0}^{P-1}\sum _{n=1}^{M_1}{K_n\left(x\right)u}_j{{\xi }_n\left(\theta \right)\psi }_j\left[\xi \left(\theta \right)\right]=\sum _{m=0}^{M_2}{F_m\xi }_m\left(\theta \right).\end{array}

The Galerkin error minimization method is used to find the optimal solution U(θ) by first multiplying both sides of Eq. (23) by ϕ_k[ξ(θ)] and then finding their mathematical expectations as follows:

(24)\begin{array}{rl}& \sum _{j=0}^{P-1}{K_0\left(x\right)u}_j⟨{\psi }_j\left[\xi \left(\theta \right)\right]{\psi }_k\left[\xi \left(\theta \right)\right]⟩\\ & +\sum _{n=1}^{M_1}\sum _{j=0}^{P-1}{K_n\left(x\right)u}_j⟨{{\xi }_n\left(\theta \right)\psi }_j\left[\xi \left(\theta \right)\right]{\psi }_k\left[\xi \left(\theta \right)\right]⟩\\ & =\sum _{m=0}^{M_2}{F}_m⟨{\xi }_m\left(\theta \right){\psi }_k\left[\xi \left(\theta \right)\right]⟩,k=0,1,\dots ,P-1.\end{array}

Assuming:

(25)\begin{array}{rl}& {\psi }_0\left(\xi \left(\theta \right)\right)=1,\\ & ⟨{\xi }_m\left(\theta \right){\psi }_k\left(\xi \left(\theta \right)\right)⟩=D_{mk},\\ & ⟨{{\psi }_j\left(\xi \left(\theta \right)\right)\psi }_k\left(\xi \left(\theta \right)\right)⟩=C_{0jk},\\ & ⟨{\xi }_n\left(\theta \right){{\psi }_j\left(\xi \left(\theta \right)\right)\psi }_k\left(\xi \left(\theta \right)\right)⟩=C_{njk}\end{array}.

Equation (24) will take a more concise form as follows:

(26)\sum _{n=0}^{M_1}\sum _{j=0}^{P-1}{K_n\left(x\right)C_{njk}u}_j=\sum _{m=0}^{M_2}{F}_m{D}_{mk},k=0,1,\dots ,P-1.

Finally, if K_{jk}=\sum _{n=0}^{M_1}{K}_n\left(x\right)C_{njk} and F_k=\sum _{m=0}^{M_2}{F}_m{D}_{mk}, Eq. (26) can have a simplified form as follows:

(27)\sum _{j=0}^{P-1}{K_{jk}u}_j=F_k,k=0,1,...,P-1.

The matrix form of which is:

(28)\underset{K_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}}{\underbrace{\left[\begin{array}{cccc}{K}_{0,0}& K_{0,1}& \cdots & K_{0,P-1}\\ K_{1,0}& K_{1,1}& \dots & K_{1,P-1}\\ ⋮& ⋮& & ⋮\\ K_{P-1,0}& K_{P-1,1}& \dots & K_{P-1,P-1}\end{array}\right]}}\underset{U_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}}{\underbrace{\left[\begin{array}{c}{u}_0\\ u_1\\ ⋮\\ u_{P-1}\end{array}\right]}}=\underset{F_{\mathrm{t}\mathrm{o}\mathrm{t}\mathrm{a}\mathrm{l}}}{\underbrace{\left[\begin{array}{c}{F}_0\\ F_1\\ ⋮\\ F_{P-1}\end{array}\right]}}.

References [10,11,22] provide more details on the SSFEM formulation.

4 Developing spectral stochastic finite element method to solve concrete beam/slab problems with multiple RFs

According to the RF definition, multiple-RF modeling, including dead and live loads and properties of concrete and steel materials, necessitates the problem’s spatial dimension modeling. Each of these parameters is discretized as a RF using the KL-expansion, and then the effect of each is considered in each element’s stiffness and force matrices. This study has assumed 2D RC beams and slabs, and used plate and bending elements, respectively, to solve their FE equations; in RC beams, effects of the rebar’s properties and stiffness are considered as a bar element and in RC slabs as a beam element.

4.1 Concrete beam considering plate tension elements

The RC beam’s 2D FE model is shown in Fig.2 where concrete elements and steel elements are considered based on combination of plates and bar elements [27,28]. It is worth noting that for reinforcements, the RF length is calculated, in each direction, by subtracting the concrete cover (c) from the structure length.

To apply the one-time effects of the top and bottom rebars and stirrups in the total stiffness matrix, four modes of combining elements are considered conventionally as in Fig.3.

To calculate the stiffness matrix of the above elements, depending on whether there are rebars, it will suffice to update the stiffness matrix equation of each element (Eq. (17)) for each mode as follows.

For elements with top and bottom rebars with an RF domain of L-2c

(29)\begin{array}{rl}& k_{0\left\{5,7\right\},\left\{5,7\right\}}^{\left(\mathrm{e}\right)}=k_{0\left\{5,7\right\},\left\{5,7\right\}}^{(\mathrm{p}\mathrm{e})}+k_0^{(\mathrm{b})},\\ & k_{n\left\{5,7\right\},\left\{5,7\right\}}^{\left(\mathrm{e}\right)}=k_{n\left\{5,7\right\},\left\{5,7\right\}}^{(\mathrm{p}\mathrm{e})}+k_n^{(\mathrm{b})}\end{array},

or elements with stirrups with an RF domain of H-2c

(30)\begin{array}{rl}& k_{0\left\{2,8\right\},\left\{2,8\right\}}^{\left(\mathrm{e}\right)}=k_{0\left\{2,8\right\},\left\{2,8\right\}}^{(\mathrm{p}\mathrm{e})}+k_0^{(\mathrm{b})},\\ & k_{n\left\{2,8\right\},\left\{2,8\right\}}^{\left(\mathrm{e}\right)}=k_{n\left\{2,8\right\},\left\{2,8\right\}}^{(\mathrm{p}\mathrm{e})}+k_n^{(\mathrm{b})}\end{array},

where k_0^{\left(\mathrm{e}\right)} and k_n^{\left(\mathrm{e}\right)} are the mean and weight stiffness matrices of concrete beam with rebar elements, respectively. Index {5,7} indicates the degree of freedom shown in Fig.3 and such indexing is applied in Eqs. (29) and (30). k_0^{(\mathrm{p}\mathrm{e})}and k_n^{(\mathrm{p}\mathrm{e})} are, respectively, the mean and weight stiffness matrices of concrete beam without rebar elements. Since a plate element has 4 nodes each with 2 DOFs, each element’s stiffness matrix in this case is 8 × 8. Moreover, k_0^{\left(\mathrm{b}\right)} and k_n^{\left(\mathrm{b}\right)} are, respectively, the mean and weight stiffness matrices of bar elements that can be found by considering the rebar type (top, bottom, and stirrups) characteristics and the corresponding RF (since a bar element has 2 nodes each with 1 DOF, stiffness matrices are 2 × 2).

4.2 Concrete slab considering Kirchhoff-Love shell theory

Concrete slab elements are considered based on Kirchhoff-Love shell theory [29–31] and their combinations with steel elements as beam elements [27].

According to Fig.4, to apply the one-time effects of the longitudinal and transverse rebars in the total stiffness matrix, four types of elements are considered conventionally as in Fig.5.

To calculate the stiffness matrix of the above elements, depending on whether there are rebars, it will suffice to update the stiffness matrix equation of each element (Eq. (17)) for each mode as follows.

For elements with longitudinal rebars with a RF domain of L_x − 2c

(31)\begin{array}{rl}& k_{0\left\{7,9,10,12\right\},\left\{7,9,10,12\right\}}^{\left(\mathrm{e}\right)}=k_{0\left\{7,9,10,12\right\},\left\{7,9,10,12\right\}}^{\left(\mathrm{s}\mathrm{e}\right)}+k_0^{\left(\mathrm{b}\mathrm{e}\right)},\\ & k_{n\left\{7,9,10,12\right\},\left\{7,9,10,12\right\}}^{\left(\mathrm{e}\right)}=k_{n\left\{7,9,10,12\right\},\left\{7,9,10,12\right\}}^{\left(\mathrm{s}\mathrm{e}\right)}+k_n^{\left(\mathrm{b}\mathrm{e}\right)}.\end{array}

For elements with transverse rebars with an RF domain of L_y − 2c

(32)\begin{array}{rl}& k_{0\left\{1,2,10,11\right\},\left\{1,2,10,11\right\}}^{\left(\mathrm{e}\right)}=k_{0\left\{1,2,10,11\right\},\left\{1,2,10,11\right\}}^{\left(\mathrm{s}\mathrm{e}\right)}+k_0^{\left(\mathrm{b}\mathrm{e}\right)},\\ & k_{n\left\{1,2,10,11\right\},\left\{1,2,10,11\right\}}^{\left(\mathrm{e}\right)}=k_{n\left\{1,2,10,11\right\},\left\{1,2,10,11\right\}}^{\left(\mathrm{s}\mathrm{e}\right)}+k_n^{\left(\mathrm{b}\mathrm{e}\right)},\end{array}

where k_0^{\left(\mathrm{e}\right)} and k_n^{\left(\mathrm{e}\right)} are the mean and weight stiffness matrices of concrete slab with rebar elements, respectively. Index {7,9,10,12} indicates the degree of freedom which are shown in Fig.5 and such indexing is applied in Eqs. (31) and (32). k_0^{\left(\mathrm{s}\mathrm{e}\right)} and k_n^{\left(\mathrm{s}\mathrm{e}\right)} are, respectively, the mean and weight stiffness matrices of concrete slab without rebar elements. (Since a bending element has 4 nodes each with 3 DOFs, each element’s stiffness matrix in this case is 12 × 12). In addition, k_0^{\left(\mathrm{b}\mathrm{e}\right)} and k_n^{\left(\mathrm{b}\mathrm{e}\right)} are, respectively, the mean and weight stiffness matrices of beam elements that can be calculated considering the rebar type (transverse or longitudinal) features and the corresponding RF. (Since a beam element has 2 nodes each with 2 DOFs, these matrices are 4 × 4).

The equation of the external force applied to the structure can be rewritten as follows, considering dead and live loads for both problems:

(33)F\left(x,\theta \right)=\sum _{m=0}^{M_2}{{Fd}_m\xi d}_m\left(\theta \right)+\sum _{m=0}^{M_2}{{Fl}_m\xi l}_m\left(\theta \right),

where suffixes d and l relate to dead and live load parameters, respectively.

5 Numerical example

According to the points mentioned in Section 2, effort has been made in this section to use 2 practical examples of a structural ceiling’s component, including concrete beams and slabs, to apply SSFEM in related problems and study the correlation length effects of an RF dead and live load in reliability assessments done using MCS and with an LSF defined as follows:

(34)g\left(x\right)=u_{\mathrm{a}\mathrm{l}\mathrm{l}}-u,

where u is the structure’s nodal points’ displacements, u_all = L_{\mathrm{u}}/240 is the allowable displacement [32] and L_u is the span length for beams and larger span length for slabs.

In addition, a variety of methods are available to generate random samples [33–37], in which normal distribution is applied in this study to generate the random samples.

5.1 RC beam

Fig.6 shows the geometric characteristics, boundary conditions and FE meshing of an RC beam with a cover of 5 cm subjected to distributed dead and live loads (DL, LL) of, respectively, 28.5 and 13.5 kN/m. Dimensions have been so selected as to have a regular and symmetrical meshing with stirrups equal in distance to the mesh size. The beam, divided into 48 elements with 65 2-DOF-nodes, is assumed to have two Gaussian points and a concrete Poisson’s ratio of 0.2. Equivalent rebars, shown schematically, are 2ø₁₈ for the top and bottom and 2ø₈ for stirrups. The 4 uncertain parameters selected for this example are concrete and steel elasticity moduli and dead and live loads, the statistical features of which are listed in Tab.1 for their corresponding RF modeling, where the number of KL-expansion and PC terms are 5 and 3, and coefficients of variations for material properties and for loading are 0.1 and 0.2, respectively.

The number of samples for the RF modeling is 10000 and the variables’ distribution function is of the standard normal type.

Fig.7 shows the mean structure response and a 10-time greater deformation of the example analyzed with SSFEM.

To see how the correlation length affected the desired beam’s uncertainty modeling, a range varying from 0.4 to 8.2 m was analyzed with a 0.6 m step (covering the beam length) for dead and live loads (b_DL and b_LL). The PDF and CDF diagrams of the maximum deflection are shown in Fig.8–Fig.11.

According to Fig.8 and Fig.9, the mean maximum RC beam deflections for different b_DL values and b_LLvalues are almost equal and differ in deviation from the average; as b_DL values and b_LLvalues increase, the PDF diagram first becomes more elongated and responses disperse more, and then these variations converge gradually for correlation lengths > 4.1 m. As shown, the PDF diagram of the highest deflection has the highest and lowest values ​​in the peak diagram, and the lowest and highest standard deviations for correlation lengths of 0.4 and 8.2 m, respectively.

Similarly, in Fig.10 and Fig.11, turning points in CDF diagrams intersect at almost the same point for different correlation lengths, indicating that for each of them, the mean response is equal. But, CDF diagrams of the maximum deflection first become more inclined with an increase in the correlation length and then their variations converge for lengths > 4.1 m, indicating different response standard deviations for different correlation lengths.

The reliability, calculated by P_f, was also found for each b_DL−b_LL pair (Fig.12). To achieve better accuracy, the generated samples fixed 10000 samples for correlation lengths greater than 1 and 50000 for correlation lengths smaller than 1.

The P_f diagram for a b_DL equal to the concrete beam length and different b_LL values is shown in Fig.13.

As shown in Fig.12 and Fig.13, P_f of the RC beam can have a 0–0.022 variation range for different correlation lengths; for 0.4 m it is the lowest, and for 8.2 m it is the highest. The longer b_DL and b_LL are, the larger P_f is; variations are first quite considerable for increased correlation lengths, and then fall gradually until P_f values converge.

Fig.14 shows the failure probability coefficient of variation based on Eq. (35) [38].

(35)Cov=\sqrt{\frac{1-P_{\mathrm{f}}}{N_{\mathrm{s}}{P}_{\mathrm{f}}}},

where N_s is number of samples.

5.2 RC slab

Fig.15 shows an example of a concrete slab subjected to uniform dead and live loads of, respectively, 6 and 3 kN/m², with a concrete cover of 10 cm. The slab, with dimensions selected for symmetrical and regular meshing, is divided into 180 elements and has 208 3-DOF nodes. Gaussian points are 2, KL-expansion terms are 5, and PC is 3. Here too, as in the previous example, the concrete and steel elasticity moduli and dead and live loads are 4 uncertain parameters whose statistical features are presented in Tab.2.

Number of samples for the RF modeling is 10000 and the variables’ distribution function is of the standard normal type.

The desired concrete slab was analyzed with SSFEM using the mentioned data and specifications, and the mean response including deflection (displacement in z-direction) and rotations in x and y directions are shown in Fig.16–Fig.18, respectively.

This example too, like the previous one, was analyzed for a set of 0.4–8.2 m b_DL and b_LL values, and the related PDF and CDF diagrams of the maximum deflection are shown in Fig.19–Fig.22.

According to Fig.19 and Fig.20, as in the previous example, the mean maximum RC slab deflections for different b_DL and b_LL values are almost equal and differ in deviation from the average; as the values of b_DL and b_LL increase, the PDF diagram first becomes more elongated and responses disperse more, and then these variations converge gradually. As also shown in this example, the PDF diagram of the maximum deflection has the highest and lowest values ​​in the peak diagram, for correlation lengths of 0.4 and 8.2 m, respectively; it also has the lowest and highest standard deviations for correlation lengths of 0.4 and 8.2 m, respectively.

The equality of the RC slab’s mean maximum deflections are observable where CDF diagrams’ turning points (Fig.21 and Fig.22) intersect. The diagrams first become more inclined with an increase in the correlation length and then their variations converge, indicating different response standard deviations for different correlation lengths.

The reliability was also calculated for different b_DL and b_LL values (Fig.23). To achieve better accuracy in the case of RC slab, the generated samples fixed 10000 samples for correlation lengths greater than 0.4 and 50000 for correlation lengths smaller than 0.4.

The P_f diagram for a b_DL equal to the concrete beam length and different b_LL are shown in Fig.24.

As shown in Fig.23 and Fig.24, in this example too, like the previous one, P_f of the RC slab can have a 0–0.062 variation range for different correlation lengths; for 0.4 m, it is the lowest and for 8.2 m, it is the highest. The longer are the b_DL and b_LL, the more is the P_f the value of which converge over time.

Fig.25 indicates the failure probability coefficient of variation in which calculated based on Eq. (35).

6 Conclusions

In this paper, SSFEM was developed for RC beams/slabs considering multiple RFs including concrete and steel material properties and dead and live loads, and uncertainty modeling and reliability evaluation were done for different dead and live load-related RF correlation lengths (b_DL and b_LL) for RC beams/slabs. According to the results, the correlation length highly affects the dispersion of the structural response, the mean value of which is almost constant for different values of b_DL and b_LL, although the standard deviation differs. The response dispersion increases first as b_DL and b_LL get longer, then it starts reducing and gradually converging for correlation lengths more than that of the RF. Regarding the failure probability variations for different correlation lengths, the trend is the same; first it increases with an increase in the correlation length and then converges for lengths longer than that of the RF. Hence, the correlation length should be selected according to the real nature and fluctuations of the certainty parameter through such field operations as lab sampling or engineering judgment to yield more accurate and precise estimate of the structural reliability.

References

[1] Der Kiureghian A, Haukaas T, Fujimura K. Structural reliability software at the University of California, Berkeley. Structural Safety, 2006, 28( 1−2): 44–67

[2] Der Kiureghian A. Analysis of structural reliability under parameter uncertainties. Probabilistic Engineering Mechanics, 2008, 23( 4): 351–358

[3] Depina I, Le T M H, Fenton G, Eiksund G. Reliability analysis with metamodel line sampling. Structural Safety, 2016, 60 : 1–15

[4] Sakata S, Okuda K, Ikeda K. Stochastic analysis of laminated composite plate considering stochastic homogenization problem. Frontiers of Structural and Civil Engineering, 2015, 9( 2): 141–153

[5] Soltani N, Alembagheri M, Khaneghahi M H. Risk-based probabilistic thermal-stress analysis of concrete arch dams. Frontiers of Structural and Civil Engineering, 2019, 13( 5): 1007–1019

[6] Ghavidel A, Rashki M, Ghohani Arab H, Azhdary Moghaddam M. Reliability mesh convergence analysis by introducing expanded control variates. Frontiers of Structural and Civil Engineering, 2020, 14( 4): 1012–1023

[7] 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

[8] Song K, Zhang Y, Zhuang X, Yu X, Song B. Reliability-based design optimization using adaptive surrogate model and importance sampling-based modified SORA method. Engineering with Computers, 2021, 37( 2): 1295–1314

[9] Papaioannou I, Straub D. Combination line sampling for structural reliability analysis. Structural Safety, 2021, 88 : 102025

[10] Papadopoulos V, Giovanis D G. Stochastic Finite Element Methods: An Introduction. Cham: Springer, 2018, 47–70

[11] Ghanem R G, Spanos P D. Stochastic Finite Elements: A Spectral Approach. New York: Springer, 1991, 101–119

[12] Bae H R, Forster E E. Improved Neumann expansion method for stochastic finite element analysis. Journal of Aircraft, 2017, 54( 3): 967–979

[13] Sofi A, Romeo E. A unified response surface framework for the interval and stochastic finite element analysis of structures with uncertain parameters. Probabilistic Engineering Mechanics, 2018, 54 : 25–36

[14] Wu F, Yao L Y, Hu M, He Z C. A stochastic perturbation edge-based smoothed finite element method for the analysis of uncertain structural-acoustics problems with random variables. Engineering Analysis with Boundary Elements, 2017, 80 : 116–126

[15] Papadopoulos V, Kalogeris I, Giovanis D G. A spectral stochastic formulation for nonlinear framed structures. Probabilistic Engineering Mechanics, 2019, 55 : 90–101

[16] Chen N Z, Soares C G. Spectral stochastic finite element analysis for laminated composite plates. Computer Methods in Applied Mechanics and Engineering, 2008, 197( 51-52): 4830–4839

[17] Stefanou G, Papadrakakis M. Stochastic finite element analysis of shells with combined random material and geometric properties. Computer Methods in Applied Mechanics and Engineering, 2004, 193( 1−2): 139–160

[18] Kandler G, Füssl J, Eberhardsteiner J. Stochastic finite element approaches for wood-based products: theoretical framework and review of methods. Wood Science and Technology, 2015, 49( 5): 1055–1097

[19] Li K, Wu D, Gao W. Spectral stochastic isogeometric analysis for static response of FGM plate with material uncertainty. Thin-walled Structures, 2018, 132 : 504–521

[20] Zhou X Y, Gosling P D, Ullah Z, Kaczmarczyk L, Pearce C J. Stochastic multi-scale finite element based reliability analysis for laminated composite structures. Applied Mathematical Modelling, 2017, 45 : 457–473

[21] SudretBDer KiureghianA. Stochastic Finite Element Methods and Reliability: A State-of-the-Art Report. Berkeley, CA: University of California, 2000

[22] YazdaniAArabH GRashkiM. Simplified spectral stochastic finite element formulations for uncertainty quantification of engineering structures. Structures, 2020, 28: 1924–1945

[23] Schietzold F N, Schmidt A, Dannert M M, Fau A, Fleury R M, Graf W, Kaliske M, Könke C, Lahmer T, Nackenhorst U. Development of fuzzy probability based random fields for the numerical structural design. GAMM-Mitteilungen, 2019, 42( 1): e201900004

[24] Schmidt A, Henning C, Herbrandt S, Könke C, Ickstadt K, Ricken T, Lahmer T. Numerical studies of earth structure assessment via the theory of porous media using fuzzy probability based random field material descriptions. GAMM-Mitteilungen, 2019, 42( 1): e201900007

[25] Zakian P, Khaji N. A stochastic spectral finite element method for wave propagation analyses with medium uncertainties. Applied Mathematical Modelling, 2018, 63 : 84–108

[26] Wiener N. The homogeneous chaos. American Journal of Mathematics, 1938, 60( 4): 897–936

[27] ChandrupatlaT RBelegunduA D. Introduction to Finite Elements in Engineering. Upper Saddle River, NJ: Prentice Hall, 2002

[28] Nariman N A, Hamdia K, Ramadan A M, Sadaghian H. Optimum design of flexural strength and stiffness for reinforced concrete beams using machine learning. Applied Sciences, 2021, 11( 18): 8762

[29] TimoshenkoS PWoinowsky-KriegerS. Theory of Plates and Shells. New York: McGraw-hill, 1959

[30] Vu-Bac N, Duong T X, Lahmer T, Zhuang X, Sauer R A, Park H S, Rabczuk T. A NURBS-based inverse analysis for reconstruction of nonlinear deformations of thin shell structures. Computer Methods in Applied Mechanics and Engineering, 2018, 331 : 427–455

[31] Vu-Bac N, Duong T X, Lahmer T, Areias P, Sauer R A, Park H S, Rabczuk T. A NURBS-based inverse analysis of thermal expansion induced morphing of thin shells. Computer Methods in Applied Mechanics and Engineering, 2019, 350 : 480–510

[32] ACI318-08. Building Code Requirements for Structural Concrete and Commentary. Farmington Hills: American Concrete Institute, 2008

[33] Vu-Bac N, Silani M, Lahmer T, Zhuang X, Rabczuk T. A unified framework for stochastic predictions of mechanical properties of polymeric nanocomposites. Computational Materials Science, 2015, 96 : 520–535

[34] Vu-Bac N, Lahmer T, Zhuang X, Nguyen-Thoi T, Rabczuk T. A software framework for probabilistic sensitivity analysis for computationally expensive models. Advances in Engineering Software, 2016, 100 : 19–31

[35] Vu-Bac N, Zhuang X, Rabczuk T. Uncertainty quantification for mechanical properties of polyethylene based on fully atomistic model. Materials, 2019, 12( 21): 3613

[36] Liu B, Vu-Bac N, Zhuang X, Rabczuk T. Stochastic multiscale modeling of heat conductivity of polymeric clay nanocomposites. Mechanics of Materials, 2020, 142 : 103280

[37] Liu B, Vu-Bac N, Rabczuk T. A stochastic multiscale method for the prediction of the thermal conductivity of polymer nanocomposites through hybrid machine learning algorithms. Composite Structures, 2021, 273 : 114269

[38] FrangopolD M. Probability concepts in engineering: Emphasis on applications to civil and environmental engineering. Structure and Infrastructure Engineering, 2008, 4(5): 413–414

RIGHTS & PERMISSIONS

Higher Education Press 2022