# A Stochastic Wavelet Finite Element Method for 1D and 2D Structures Analysis

Authors: Xingwu Zhang; Xuefeng Chen; Zhibo Yang; Bing Li; Zhengjia He
Journal: Shock and Vibration
 (2014)

---

## Abstract

A stochastic finite element method  based on B-spline wavelet on the interval (BSWI-SFEM) is presented for static analysis of 1D and 2D structures in this paper. Instead of conventional polynomial interpolation, the scaling functions of BSWI are employed to construct the displacement field. By means of virtual work principle and BSWI, the wavelet finite elements of beam, plate, and plane rigid frame are obtained. Combining the Monte Carlo method and the constructed BSWI elements together, the BSWI-SFEM is formulated. The constructed BSWI-SFEM can deal with the problems of structural response uncertainty caused by the variability of the material properties, static load amplitudes, and so on. Taking the widely used Timoshenko beam,  the Mindlin plate, and the plane rigid frame as examples, numerical results have demonstrated that the proposed method can give a higher accuracy and a better constringency than the conventional stochastic finite element methods.


---

## Body

## 1. Introduction

The finite element method (FEM) is well accepted now in numerous industrial areas [1–4]. However, the validity of the results obtained by FEM can be drastically limited by the uncertainties of the parameters introduced in the model. Many factors may cause the uncertainties, such as a consequence of ignorance and spatial variability [5, 6].On the other hand, the reliability design and analysis for engineering structures also require a kind of FE model which could deal with the variability and uncertainty. For example, a slight change in design variables may result in a poignant change of structural properties. As a powerful tool in computational stochastic mechanics, the stochastic finite element method (SFEM) appeared.SFEM is an extension of the classical FEM to stochastic framework. The most important factor of SFEM is the uncertainty modeling. Over the past several decades, many methods were combined with FEM to approximate the uncertainties mentioned above. Chakraborty and Dey investigated the stochastic finite element simulation based on the Neumann expansion of uncertain structures subjected to earthquake [7]. A method named first-order, second-moment (FOSM) method was proposed by Cornell [8]; in fact it is a special case of perturbation method [9, 10]. Valdebenito et al. studied the application of first-order expansion considering intervening variables for estimating second-order statistics for stochastic finite element analysis [11]. The Karhunen-Loeve [12] expansion, which can be seen as a kind of the orthogonal series expansion, is also an important representation of SFEM fields. Sakamoto and Ghanem [13] proposed a method based on polynomial chaos expansion for presenting the variability. Proppe implemented multiresolution analysis for stochastic finite element problems with Karhunen-Loeve expansion [14]. Based on Karhunen-Loeve series expansion and hierarchical approach, Yang et al. applied the hierarchical stochastic finite element method for structural analysis [15]. Another method is Monte Carlo simulation (MCS) [16]; the name “Monte Carlo” was coined in the 1940s by scientists working in Los Alamos to designate a class of numerical methods based on the use of random numbers. Compared with other technologies, MCS is more direct and widely used. Taking a certain amount of samples as the priori testing data, MCS can give an accurate approximation for the stochastic field [17]. Besides, Diazdelao and Adhikari explored a method to reduce the computation cost of stochastic finite element codes based on the Gaussian process emulation [18]. Han and Bang performed the probabilistic optimal safety assessment of steel cable-stayed bridges using SFEM and expected life cycle cost concept [19]. Lang et al. investigated the prediction of heat transfer in composite materials with uncertain geometry based on the combination of the extended finite element method and the spectral stochastic finite element method [20]. Augustin and Rentrop constructed the stochastic Galerkin techniques for the random ordinary differential equations and demonstrated the advantages through numerical examples [21]. Shang and Yun proposed a stochastic finite element within a general purpose finite element analysis program ABAQUS to simulate the probabilistic structural response of stochastic materials [22].Another important factor of SFEM is the properties of FEM element which influence the presentation directly. Monte Carlo simulation is widely used to solve complex mechanics and physics problems, particularly those involving multiple independent stochastic variables. The essential characteristic of Monte Carlo method is the use of random numbers and random variables. A random variable is a quantity that results from a repeatable process and whose actual values cannot be predicted with certainty [23]. The occurrence frequency of a sample value would be deemed to be the approximation to its probability. With the famous “law of large number,” the correctness of Monte Carlo simulation can be understood easily. Usually, a sample size N of the order 30–100 could already provide estimates with acceptable accuracy [24]. On the other hand, the law makes the exactness of Monte Carlo simulation heavily depend on the sampling counter. With the growing number of sampling counters, the calculation efficiency seems to be particularly important. That is to say, the realization of MCS should be based not only on an effective computer system but also on an efficient algorithm. Therefore, a kind of new stochastic finite element method based on B-spline wavelet on the interval and MCS (BSWI-SFEM) is constructed to deal with the current problems in SFEM analysis.It is necessary to recommend wavelet analysis and BSWI briefly before we introduce BSWI-FEM. Wavelet analysis is a new method developed in recent years [25]. The wavelet method can be viewed as a method in which the approximating function is defined by a multiresolution technique based on scaling or wavelet functions, similar to those used in signal and image processing. With its desirable characteristics, such as multiresolution properties and various basis functions for structural analysis, wavelet method is well argued by many researchers not only in numerical analysis domains [26–28] but also in structural analysis fields [29–31]. As a kind of wavelet, BSWI basis has the good characteristics of compact support, smoothness, and symmetry besides the multiresolution analysis. Moreover, it has the explicit expression, which will not lead to any trouble for differentiation and integration. So it is the best one among all existing wavelets in approximation of numerical calculation [32]. Based on the good characteristics of BSWI, several investigators have tried to combine it with FEM to solve practical problems [33–37]. However, the above-mentioned methods, which are all constructed in a determined wavelet space, can hardly make a satisfying approximation for the influence caused by structural uncertainties.In this paper, a new BSWI-SFEM is presented. Combining BSWI element with Monte Carlo method, a BSWI-SFEM formulation is constructed to deal with the static and vibration analysis of 1D and 2D structures. Taking the Timoshenko beam, the Mindlin plate, and the plane rigid frame as examples, numerical results indicate that the proposed method has a higher accuracy and a better constringency than the conventional SFEM.


## 2. The Construction of BSWI Element

### 2.1. One-Dimensional BSWI on the Interval[0,1]

Goswami et al. constructed the scaling functions of BSWI and gave themth order wavelet formulations at scale j=0 in 1995 [39]. To any scale j, the mth order BSWI scaling functions must satisfy the following condition [39] in order to have at least one inner wavelet on the interval [0,1]:
(1)2j≥2m-1.Thenj scale mth order BSWI (BSWImj) scaling functions ϕm,kj(ξ) and the corresponding wavelets ψm,kj(ξ) can be evaluated by the following formulas:(2)ϕm,kj(ξ)={ϕm,kl(2j-lξ),k=-m+1,…,-1(0boundaryscalingfunctions)ϕm,2j-m-kl(1-2j-lξ),k=2j-m+1,…,2j-1(1boundaryscalingfunctions)ϕm,0l(2j-lξ-2-lk),k=0,…,2j-m(innerscalingfunctions)ψm,kj(ξ)={ψm,kl(2j-lξ),k=-m+1,…,-1(0boundarywaveletfunctions)ψm,2j-m-kl(1-2j-lξ),k=2j-m+1,…,2j-1(1boundarywaveletfunctions)ψm,0l(2j-lξ-2-lk),k=0,…,2j-m(innerwaveletfunctions).Letj0 be the scale which satisfies (1); then let l=0 in (2), for any j>j0; according to the 0 scale mth order scaling functions and wavelets given in [39], we can get the j scale mth order scaling functions and wavelets. There are m-1 boundary scaling functions and wavelets, 2j-m+1 inner scaling functions, and 2j-2m+2 inner wavelets. The scaling functions and the wavelets on the interval [0,1] can be written in vector form as follows:
(3)Φ={ϕm,-m+1j(ξ)ϕm,-m+2j(ξ)⋯ϕm,2j-1j(ξ)},(4)Ψ={ψm,-m+1j(ξ)ψm,-m+2j(ξ)⋯ψm,2j-1j(ξ)}.
The corresponding eleven scaling functions and eleven wavelets of BSWI43 are given in Figure 1.BSWI43 on the interval [0,1]: (a) scaling functions and (b) wavelets.
(a)
(b)


### 2.2. Two-Dimensional Tensor Product BSWI on the Interval[0,1]

By means of tensor product, 2D BSWI can be easily extended from 1D BSWI. A separable wavelet semiorthonormal basis at scalej of L2(R2) is constructed with tensor products of the multiresolution approximation spaces Vj1 and Vj2, so the tensor product subspace can be expressed as Fj=Vj1⊗Vj1, and the scaling functions can be expressed as
(5)Φ=Φ1⊗Φ2,
where Φ1={ϕm,-m+1j(ξ)ϕm,-m+2j(ξ)⋯ϕm,2j-1j(ξ)}, indicating one row vector combined by the mth order scaling functions at scale j. Φ2={ϕm,-m+1j(η)ϕm,-m+2j(η)⋯ϕm,2j-1j(η)}, indicating another row vector combined by the mth order scaling functions at scale j. The wavelets functions are
(6)Ψ1=Φ1⊗ψ2Ψ2=ψ1⊗Φ2Ψ3=ψ1⊗ψ2,
where ψ1 and ψ2 are wavelets of BSWImj in (4). Figure 2(a) shows all the scaling functions of 2D tensor products BSWI, and the wavelets of tensor products BSWI are shown in Figures 2(b), 2(c), and 2(d).Two-dimensional tensor product BSWI: (a) scaling functionsΦ=Φ1⊗Φ2, (b) wavelets Ψ1=Φ1⊗ψ2, (c) wavelets Ψ2=ψ1⊗Φ2, and (d) wavelets Ψ3=ψ1⊗ψ2.
(a)
(b)
(c)
(d)


### 2.3. The BSWI Timoshenko Beam Element

Considering the influence of shear deformation, the generalized function of potential energy of the Timoshenko beam element can be given as [32]
(7)Πp(w,θ)=∫0LEI2(-dθdx)2dx+∫0LGA2k(dwdx-θ)2dx-∫0Lf(x)wdx-∑iPiw(xi)+∑iMiθ(xi),
where EI is the bending rigidity, θ is cross section rotation, G is shear modulus, A is the area of cross section, Pi is the lump force, Mi is the lump bending moment, and k is the shear deformation coefficient. According to the Timoshenko beam theory, w(ξ) and θ(ξ) can be independently interpolated by BSWI scaling functions as
(8)w(ξ)=ΦTewe,θ(ξ)=ΦTeθe.Translating (7) into interval [0,1] and taking (8) into (7), let δΠp=0, and the element solving equation of bending problem can be obtained as follows:
(9)[Ke,1Ke,2Ke,3Ke,4][weθe]=[Pwe0]+[PwieMθle],
where
(10)Pwe=le∫01(Te)TΦTf(ξ)dξ,Pwie=∑iPi(Te)TΦT(ξi),Mθle=-∑lMl(Te)TΦT(ξl);Ke,1=GAkle∫01(Te)T(dΦTdξ)(dΦdξ)Tedξ,Ke,2=-GAk∫01(Te)T(dΦTdξ)ΦTedξ,Ke,3=(Ke,2)T,Ke,4=EIle∫01(Te)T(dΦTdξ)(dΦdξ)Tedξ+GAlek∫01(Te)TΦTΦTedξ.The element solving equation of vibration problem is(11)|K-λM|=0,
where M is the consistent mass matrix as follows:
(12)M=ρAle(Te)T(dΦTdξ)(dΦdξ)Te.K is the stiffness matrix as follows:
(13)K=[Ke,1Ke,2Ke,3Ke,4].ρ is density and λ is the vibration engine value.


### 2.4. The BSWI Mindlin Plate Element

According to the Mindlin plate theory, the displacements and rotations are interpolated independently, then the generalized potential energy function for the Mindlin plate bending problem is [32]
(14)Πp=∬Ω12κTDκdxdy-∬Ω12γTCγdxdy-∬Ω12q(x,y)wdxdy,
where q is the concentrated load. Consider
(15)κ={-∂θx∂x-∂θy∂y-(∂θx∂y+∂θy∂x)}T,D=Et312(1-μ2)[1μ0μ10001-μ2],γ={(∂w∂x-θx)(∂w∂y-θy)0}T,C=5Et12(1+μ)[100010000],
where θx and θy are rotations in directions x and y, t is the thickness of the Mindlin plate, and μ is the Poisson ratio. Other symbols are the same as in (7). Taking BSWI scaling functions as interpolating function is as follows:
(16)w(ξ,η)=ΦTewe,θ(ξ,η)=ΦTeθe.
Submitting (16) into (14), according to the generalized variational principle, the BSWI element formulation of the solving equation can be obtained as follows:
(17)[Ke,1Ke,2Ke,3Ke,4Ke,5Ke,6Ke,7Ke,8Ke,9][weθxeθye]=[Pe00],
where we denote the following:
(18)Γ10,0=(T1e)Tlex∫01Φ1TΦ1dξT1e,Γ10,1=(T1e)T∫01Φ1TdΦ1dξdξT1e,Γ11,0=(T1e)Tlex∫01dΦ1TdξΦ1dξT1e,Γ11,1=(T1e)T1lex∫01dΦ1TdξdΦ1dξdξT1e.Γ2i,j  (i,j=0,1) is similar to Γ1i,j if the lengths of element in x-axis lex, dξ, T1e, and Φ1 are replaced by ley, dη, T2e, and Φ2, respectively. Then the terms of (17) can be written as
(19)Ke,1=kEt2(1+μ)(Γ11,1⊗Γ20,0+Γ10,0⊗Γ21,1);Ke,2=-kEt2(1+μ)Γ11,0⊗Γ20,0;Ke,3=-kEt2(1+μ)Γ10,0⊗Γ21,0;Ke,4=(Ke,2)T;Ke,5=Et312(1-μ2){Γ11,1⊗Γ20,0+(1-μ)2Γ10,0⊗Γ21,1}+kEt2(1+μ)Γ10,0⊗Γ20,0;Ke,6=Et312(1-μ2)(μΓ11,0⊗Γ20,1+(1-μ)2Γ10,1⊗Γ21,0);Ke,7=(Ke,3)T;Ke,8=(Ke,6)T;Ke,9=Et312(1-μ2){Γ10,0⊗Γ21,1+(1-μ)2Γ11,1⊗Γ20,0}+kEt2(1+μ)Γ10,0⊗Γ20,0;Pe=(Te)Tlexley∬01f(ξ,η)(Φ1⊗Φ2)Tdξdη.The element solving equation of vibration problem is(20)|K-λM|=0,
where M is the following consistent mass matrix:
(21)M=ρlexleyt(Te)T(dΦTdξ)⊗(dΦdη)Te.K is the following stiffness matrix:
(22)K=[Ke,1Ke,2Ke,3Ke,4Ke,5Ke,6Ke,7Ke,8Ke,9].ρ is density and λ is the vibration engine value.


### 2.5. BSWI Plane Rigid Frame Element

The BSWI plane rigid frame element can be seen as the simple superposition of the axial rod element and the Euler beam element. The displacement of element nodes, solving domain, and corresponding coordinates are shown in Figure3.Figure 3
The solving domain and nodes displacement of plane rigid frame.It can be seen from Figure3 that there are 2j+1+4 DOFs of one element. The physical DOFs of BSWI plane rigid element in local coordinate (X,Y) are(23)δe=(u1w1M1u2w2M2⋯un-1wn-1Mn-1unwnMn)T,
while the corresponding DOFs in global coordinate (X-,Y-) are
(24)δ-e=(u-1w-1M-1u-2w-2M-2⋯u-n-1w-n-1M-n-1u-nw-nM-n)T.In order to construct the BSWI plane rigid frame element, the relationship between global and local coordinates should be obtained. Consider(25)[uiwiMi]=gi[u¯iw¯iM¯i],i=1,…,n,
where
(26)gi=[cos⁡αsinα0-sinαcos⁡α0001],i=1,…,n.
Therefore,
(27)δe=Geδ¯e.Therefore, the BSWI plane rigid frame element for bending problem analysis in global coordinate is(28)K¯δ¯=P¯.The BSWI plane rigid frame element for vibration problem analysis in global coordinate is(29)|K--λM-|=0.Then, by introducing the BSWI axial rod element and the BSWI Euler beam element [40], the BSWI plane rigid frame can be obtained.


## 2.1. One-Dimensional BSWI on the Interval[0,1]

Goswami et al. constructed the scaling functions of BSWI and gave themth order wavelet formulations at scale j=0 in 1995 [39]. To any scale j, the mth order BSWI scaling functions must satisfy the following condition [39] in order to have at least one inner wavelet on the interval [0,1]:
(1)2j≥2m-1.Thenj scale mth order BSWI (BSWImj) scaling functions ϕm,kj(ξ) and the corresponding wavelets ψm,kj(ξ) can be evaluated by the following formulas:(2)ϕm,kj(ξ)={ϕm,kl(2j-lξ),k=-m+1,…,-1(0boundaryscalingfunctions)ϕm,2j-m-kl(1-2j-lξ),k=2j-m+1,…,2j-1(1boundaryscalingfunctions)ϕm,0l(2j-lξ-2-lk),k=0,…,2j-m(innerscalingfunctions)ψm,kj(ξ)={ψm,kl(2j-lξ),k=-m+1,…,-1(0boundarywaveletfunctions)ψm,2j-m-kl(1-2j-lξ),k=2j-m+1,…,2j-1(1boundarywaveletfunctions)ψm,0l(2j-lξ-2-lk),k=0,…,2j-m(innerwaveletfunctions).Letj0 be the scale which satisfies (1); then let l=0 in (2), for any j>j0; according to the 0 scale mth order scaling functions and wavelets given in [39], we can get the j scale mth order scaling functions and wavelets. There are m-1 boundary scaling functions and wavelets, 2j-m+1 inner scaling functions, and 2j-2m+2 inner wavelets. The scaling functions and the wavelets on the interval [0,1] can be written in vector form as follows:
(3)Φ={ϕm,-m+1j(ξ)ϕm,-m+2j(ξ)⋯ϕm,2j-1j(ξ)},(4)Ψ={ψm,-m+1j(ξ)ψm,-m+2j(ξ)⋯ψm,2j-1j(ξ)}.
The corresponding eleven scaling functions and eleven wavelets of BSWI43 are given in Figure 1.BSWI43 on the interval [0,1]: (a) scaling functions and (b) wavelets.
(a)
(b)


## 2.2. Two-Dimensional Tensor Product BSWI on the Interval[0,1]

By means of tensor product, 2D BSWI can be easily extended from 1D BSWI. A separable wavelet semiorthonormal basis at scalej of L2(R2) is constructed with tensor products of the multiresolution approximation spaces Vj1 and Vj2, so the tensor product subspace can be expressed as Fj=Vj1⊗Vj1, and the scaling functions can be expressed as
(5)Φ=Φ1⊗Φ2,
where Φ1={ϕm,-m+1j(ξ)ϕm,-m+2j(ξ)⋯ϕm,2j-1j(ξ)}, indicating one row vector combined by the mth order scaling functions at scale j. Φ2={ϕm,-m+1j(η)ϕm,-m+2j(η)⋯ϕm,2j-1j(η)}, indicating another row vector combined by the mth order scaling functions at scale j. The wavelets functions are
(6)Ψ1=Φ1⊗ψ2Ψ2=ψ1⊗Φ2Ψ3=ψ1⊗ψ2,
where ψ1 and ψ2 are wavelets of BSWImj in (4). Figure 2(a) shows all the scaling functions of 2D tensor products BSWI, and the wavelets of tensor products BSWI are shown in Figures 2(b), 2(c), and 2(d).Two-dimensional tensor product BSWI: (a) scaling functionsΦ=Φ1⊗Φ2, (b) wavelets Ψ1=Φ1⊗ψ2, (c) wavelets Ψ2=ψ1⊗Φ2, and (d) wavelets Ψ3=ψ1⊗ψ2.
(a)
(b)
(c)
(d)


## 2.3. The BSWI Timoshenko Beam Element

Considering the influence of shear deformation, the generalized function of potential energy of the Timoshenko beam element can be given as [32]
(7)Πp(w,θ)=∫0LEI2(-dθdx)2dx+∫0LGA2k(dwdx-θ)2dx-∫0Lf(x)wdx-∑iPiw(xi)+∑iMiθ(xi),
where EI is the bending rigidity, θ is cross section rotation, G is shear modulus, A is the area of cross section, Pi is the lump force, Mi is the lump bending moment, and k is the shear deformation coefficient. According to the Timoshenko beam theory, w(ξ) and θ(ξ) can be independently interpolated by BSWI scaling functions as
(8)w(ξ)=ΦTewe,θ(ξ)=ΦTeθe.Translating (7) into interval [0,1] and taking (8) into (7), let δΠp=0, and the element solving equation of bending problem can be obtained as follows:
(9)[Ke,1Ke,2Ke,3Ke,4][weθe]=[Pwe0]+[PwieMθle],
where
(10)Pwe=le∫01(Te)TΦTf(ξ)dξ,Pwie=∑iPi(Te)TΦT(ξi),Mθle=-∑lMl(Te)TΦT(ξl);Ke,1=GAkle∫01(Te)T(dΦTdξ)(dΦdξ)Tedξ,Ke,2=-GAk∫01(Te)T(dΦTdξ)ΦTedξ,Ke,3=(Ke,2)T,Ke,4=EIle∫01(Te)T(dΦTdξ)(dΦdξ)Tedξ+GAlek∫01(Te)TΦTΦTedξ.The element solving equation of vibration problem is(11)|K-λM|=0,
where M is the consistent mass matrix as follows:
(12)M=ρAle(Te)T(dΦTdξ)(dΦdξ)Te.K is the stiffness matrix as follows:
(13)K=[Ke,1Ke,2Ke,3Ke,4].ρ is density and λ is the vibration engine value.


## 2.4. The BSWI Mindlin Plate Element

According to the Mindlin plate theory, the displacements and rotations are interpolated independently, then the generalized potential energy function for the Mindlin plate bending problem is [32]
(14)Πp=∬Ω12κTDκdxdy-∬Ω12γTCγdxdy-∬Ω12q(x,y)wdxdy,
where q is the concentrated load. Consider
(15)κ={-∂θx∂x-∂θy∂y-(∂θx∂y+∂θy∂x)}T,D=Et312(1-μ2)[1μ0μ10001-μ2],γ={(∂w∂x-θx)(∂w∂y-θy)0}T,C=5Et12(1+μ)[100010000],
where θx and θy are rotations in directions x and y, t is the thickness of the Mindlin plate, and μ is the Poisson ratio. Other symbols are the same as in (7). Taking BSWI scaling functions as interpolating function is as follows:
(16)w(ξ,η)=ΦTewe,θ(ξ,η)=ΦTeθe.
Submitting (16) into (14), according to the generalized variational principle, the BSWI element formulation of the solving equation can be obtained as follows:
(17)[Ke,1Ke,2Ke,3Ke,4Ke,5Ke,6Ke,7Ke,8Ke,9][weθxeθye]=[Pe00],
where we denote the following:
(18)Γ10,0=(T1e)Tlex∫01Φ1TΦ1dξT1e,Γ10,1=(T1e)T∫01Φ1TdΦ1dξdξT1e,Γ11,0=(T1e)Tlex∫01dΦ1TdξΦ1dξT1e,Γ11,1=(T1e)T1lex∫01dΦ1TdξdΦ1dξdξT1e.Γ2i,j  (i,j=0,1) is similar to Γ1i,j if the lengths of element in x-axis lex, dξ, T1e, and Φ1 are replaced by ley, dη, T2e, and Φ2, respectively. Then the terms of (17) can be written as
(19)Ke,1=kEt2(1+μ)(Γ11,1⊗Γ20,0+Γ10,0⊗Γ21,1);Ke,2=-kEt2(1+μ)Γ11,0⊗Γ20,0;Ke,3=-kEt2(1+μ)Γ10,0⊗Γ21,0;Ke,4=(Ke,2)T;Ke,5=Et312(1-μ2){Γ11,1⊗Γ20,0+(1-μ)2Γ10,0⊗Γ21,1}+kEt2(1+μ)Γ10,0⊗Γ20,0;Ke,6=Et312(1-μ2)(μΓ11,0⊗Γ20,1+(1-μ)2Γ10,1⊗Γ21,0);Ke,7=(Ke,3)T;Ke,8=(Ke,6)T;Ke,9=Et312(1-μ2){Γ10,0⊗Γ21,1+(1-μ)2Γ11,1⊗Γ20,0}+kEt2(1+μ)Γ10,0⊗Γ20,0;Pe=(Te)Tlexley∬01f(ξ,η)(Φ1⊗Φ2)Tdξdη.The element solving equation of vibration problem is(20)|K-λM|=0,
where M is the following consistent mass matrix:
(21)M=ρlexleyt(Te)T(dΦTdξ)⊗(dΦdη)Te.K is the following stiffness matrix:
(22)K=[Ke,1Ke,2Ke,3Ke,4Ke,5Ke,6Ke,7Ke,8Ke,9].ρ is density and λ is the vibration engine value.


## 2.5. BSWI Plane Rigid Frame Element

The BSWI plane rigid frame element can be seen as the simple superposition of the axial rod element and the Euler beam element. The displacement of element nodes, solving domain, and corresponding coordinates are shown in Figure3.Figure 3
The solving domain and nodes displacement of plane rigid frame.It can be seen from Figure3 that there are 2j+1+4 DOFs of one element. The physical DOFs of BSWI plane rigid element in local coordinate (X,Y) are(23)δe=(u1w1M1u2w2M2⋯un-1wn-1Mn-1unwnMn)T,
while the corresponding DOFs in global coordinate (X-,Y-) are
(24)δ-e=(u-1w-1M-1u-2w-2M-2⋯u-n-1w-n-1M-n-1u-nw-nM-n)T.In order to construct the BSWI plane rigid frame element, the relationship between global and local coordinates should be obtained. Consider(25)[uiwiMi]=gi[u¯iw¯iM¯i],i=1,…,n,
where
(26)gi=[cos⁡αsinα0-sinαcos⁡α0001],i=1,…,n.
Therefore,
(27)δe=Geδ¯e.Therefore, the BSWI plane rigid frame element for bending problem analysis in global coordinate is(28)K¯δ¯=P¯.The BSWI plane rigid frame element for vibration problem analysis in global coordinate is(29)|K--λM-|=0.Then, by introducing the BSWI axial rod element and the BSWI Euler beam element [40], the BSWI plane rigid frame can be obtained.


## 3. Stochastic Finite Element Based on B-Spline Wavelet on the Interval

The application of the MCS-FEM to large-scale structures was impossible due to its excessive computation cost a few years ago due to the huge calculation amount. Let us consider the event that the random variablex can take the discrete values x=1,…,N with probabilities p1,p2,…,pN, respectively. The corresponding probability distribution functions can be expressed as
(30)p(x)=∑i=1Npiδ(x-i),
where δ(x) is the Dirac function, and the corresponding cumulative distribution function can be expressed as
(31)P(x)={0,ifj<0,∑i=1jpj,if1≤j≤N,1,ifj>1.Since a discrete distribution is used to approximate the continuous distribution in the real events, the sampling counterN should be big enough.While BSWI-FEM has the property of high accuracy and efficiency [26], through inosculating this technology with MCS, the limitation of applying MCS-FEM to large-scale structures can be alleviated further.Therefore, the stochastic wavelet finite element based on B-spline wavelet on the interval is proposed by the aid of combining Monte Carlo method with BSWI-FEM. The first step is to select the initialization parameters, such as element type and element number. Then the sampling counterN should be decided, and the probability model of variables also needs to be selected. With these inputs, the calculation of BSWI-SFEM can be started as shown in Figure 4. Then the stiffness matrix K and the force vector P are calculated by BSWI-SFEM. Finally, through BSWI-SFEM, the results can be obtained. This process will be repeated again and again until converged. Considering the robustness of the procedure, the sampling counter N is set as a variable in course; this can guarantee the convergence.Figure 4
The flowchart of BSWI-SFEM.


## 4. Numerical Examples

To stochastic analysis problem, it is hard to define an exact solution. But in one sense, the solution obtained by a large number of samples could be a good approximation of the real exact solution. In this paper, we use this kind of sample exact solution as an index to evaluate the property of the constructed BSWI-FEM. Furthermore, we will give a traditional result of finite element analysis before the BSWI-FEM analysis to validate the accuracy of the proposed method.All the programs are implemented in the MATLAB environment. In all stochastic simulations, a computer with a 2.5 GHz Intel CPU and 2 gigabytes of memory is used, and the operating system is 32-bit Windows XP. We assume that all the parameters obey the normal distribution: the mean value of length is 1000 mm (1000 mm × 1000 mm when it is plate), the mean value of the cross section for beam is 20 mm × 12 mm, the mean value of plate thickness is 10 mm, the mean value amplitude of load is 100 N, the mean value of density is 7917 kg/m3, and the mean value of Young’s modulus is 206 MPa. The standard deviations of these parameters are 5% of their mean values.


### 4.1. Bending of the Timoshenko Beam

A clamp supported Timoshenko beam under a uniform load (as shown in Figure5) is considered in this example. The parameters are mentioned above. It is demonstrated in Table 1 that the proposed BSWI element has a good agreement with the exact solution [38].Table 1
Comparison of center deflection of clamp supported beam under uniform load.
Mesh method (DOFs)
Center deflection [mm]
80 traditional FEM (162)
0.14114
160 traditional FEM (322)
0.15399
320 traditional FEM (642)
0.15758
1 BSWI FEM (11)
0.15876
Exact solution [38]
0.15802Figure 5
The Clamp supported Timoshenko beam under a uniform load.


### 4.2. Stochastic Analysis of the Timoshenko Beam Bending

A stochastic analysis for a clamp supported Timoshenko beam under a uniform load as shown in Figure5 is considered in this numerical example. The results of BSWI-SFEM are shown in Table 2 compared with conventional FEM. It is evident that 1 BSWI-SFEM element can achieve a similar accuracy as 320 conventional Timoshenko elements. However, the proposed method only needs 0.078 s to establish the computation, while the conventional method needs 154.875 s as shown in Figure 6.Table 2
Comparison of center deflection of clamp supported beam under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [mm]
Time consumption [s]
80 traditional FEM (162)
0.14660
0.781000
160 traditional FEM (322)
0.16066
5.765000
320 traditional FEM (642)
0.16462
154.875000
1 BSWI FEM (11)
0.16593
0.078000
Sample exact solution
0.16368
—Figure 6
Time consumption comparison of center deflection for clamp supported beam under uniform load stochastic analysis.Figure7 presents the mean value of deflection at the center of beam which varies with the counter of sample. It can be seen that 1 BSWI-SFEM element obtains a good approximation of the center deflection compared with conventional method. Although the conventional method could approximate the sample exact solution with mesh refinement, the time consumption increases sharply. This kind of increase would result in unacceptable time consumption for the conventional method when it is used in large-scale structures. However, owing to the efficiency and accuracy of BSWI-SFEM, this problem can be greatly eased.Relationship between the mean value of beam center deflection and sampling number: (a) sampling number 1–20 and (b) sampling number 20–100.
(a)
(b)


### 4.3. Stochastic Analysis of the Timoshenko Beam Vibration

Figure8 shows the vibration analysis results of the clamp supported Timoshenko beam. BSWI results are compared with traditional finite element method. It can be seen that 1 BSWI-SFEM element can obtain a good approximation of natural frequencies compared with traditional method. Although the conventional method could approximate the sample exact solution with mesh refinement, the time consumption increases sharply. In this numerical example, the proposed BSWI method only needs 6.005 s to establish the computation with 1 element, while the conventional method needs 240.758 s to obtain the results in a similar accuracy with 320 elements.Relationship between the first three natural frequencies and sampling numbers of the Timoshenko beam. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


### 4.4. Bending of the Mindlin Plate

In this part, the Mindlin plate bending problem as presented in Figure9 is used to test the accuracy of the BSWI Mindlin plate element. The corresponding parameters are chosen as the mean value of these stochastic parameters mentioned above. The numerical results given in Table 3 indicate that the proposed method has a good agreement with the exact solution [38] once more.Table 3
Comparison of center deflection of the clamp and simply supported Mindlin plate under uniform load.
Mesh method (DOFs)
Center deflection [μm]
The clamp supported Mindlin plate
The simply supported Mindlin plate
100 traditional FEM (363)
7.9426
22.598
400 traditional FEM (1323)
7.9676
22.637
900 traditional FEM (2883)
7.9722
22.644
1 BSWI FEM (363)
7.9761
22.651
Exact solution [35]
7.9760
22.649The Clamp and simply supported Mindlin plate under a uniform load: (a) clamp supported and (b) simply supported.
(a)
(b)


### 4.5. Stochastic Analysis of the Mindlin Plate Bending

Stochastic analysis of a square Mindlin plate with clamp and simply supported boundary conditions subjected to a uniform load as shown in Figure9 is considered in this numerical example. The stochastic bending analysis of this plate is carried out by BSWI-SFEM incorporated with the Monte Carlo simulation. The center deflection of the two kinds of boundary conditions is presented in Tables 4 and 5. The results show that the BSWI-SFEM can give a better approximation for the sample exact solution compared with the conventional method. Furthermore, the time consumption as shown in Figures 10 and 11 indicates that BSWI-SFEM is more efficient than the conventional method, which allows the proposed approach to be used for the large-scale structures. The further comparisons are presented in Figures 12 and 13. Using the mean value of output variables, the reliability calculations can be thereafter performed.Table 4
Comparison of center deflection of the clamp supported Mindlin plate under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [μm]
Time consumption [s]
100 traditional FEM (363)
8.0250
6.766000
400 traditional FEM (1323)
8.0506
90.562000
900 traditional FEM (2883)
8.0553
824.875000
1 BSWI FEM (363)
8.0592
4.015000
Sample exact solution
8.0590
—Table 5
Comparison of center deflection of the simply supported Mindlin plate under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [μm]
Time consumption [s]
100 traditional FEM (363)
22.883
6.821000
400 traditional FEM (1323)
22.922
93.610000
900 traditional FEM (2883)
22.930
816.320000
1 BSWI FEM (363)
22.936
3.986000
Sample exact solution
22.934
—Figure 10
Time consumption comparison of center deflection for the clamp supported Mindlin plate under uniform load stochastic analysis.Figure 11
Time consumption comparison of center deflection for the simply supported Mindlin plate under uniform load stochastic analysis.Relationship between the mean value of plate center deflection and sampling number: (a) sampling number 10–35, (b) sampling number 35–100.
(a)
(b)Relationship between the mean value of plate center deflection and sampling number: (a) sampling number 15–50, (b) sampling number 50–100.
(a)
(b)


### 4.6. Stochastic Analysis of the Mindlin Plate Vibration

The vibration analysis results of the clamp supported Mindlin plate are shown in Figure14. The first three natural frequency results of BSWI element are compared with conventional element. It can be seen from the results that 1 BSWI element can obtain the results in a similar accuracy with 320 traditional elements. However, 1 BSWI method only needs 10.397 s to establish the computation, while the conventional method needs 4320.385 s to obtain the results in a similar accuracy with 320 elements.Relationship between the first three natural frequencies and sampling numbers of the Mindlin plate. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


### 4.7. Stochastic Analysis of Plane Rigid Frame Bending

As shown in Figure15, the plane rigid frame combined by 15 components is considered in this example. We assume that all the parameters obey the normal distribution. The corresponding material parameters are elastic modulus E=2.06×1011 N/m2, width b=0.012 m, height h=0.012 m, and density ρ=7917 kg/m3, respectively. The standard deviations of these parameters are 5% of their mean values.Figure 15
Coordinates, material parameters, and load parameters of plane rigid frame.As shown in Figure16, the bending analysis results of BSWI element on points A and B are compared with BEAM3 element. Although BEAM3 element can obtain the results with high accuracy by grid refinement, it will spend much more time on computation. However, the same solving domain of plane rigid frame is discreted by only 15 BSWI elements; it can achieve the results in a similar accuracy with 480 BEAM3 elements. So the efficiency and solving precision of BSWI-SFEM plane rigid element are verified.Relationship between the bending analysis results and sampling number of the plane rigid frame. (a) Point A and (b) point B.
(a)
(b)


### 4.8. Stochastic Analysis of Plane Rigid Frame Vibration

The stochastic vibration analysis of plane rigid frame is carried out by BSWI-SFEM incorporated with the Monte Carlo simulation in this example. The corresponding analysis results are shown in Figure17. BSWI element results are compared with BEAM3 element. It can be seen from the comparison that, although BEAM3 element can acquire the results with high accuracy by grid refinement, much more time is spent. BSWI-SFEM is superior, with only 15 elements, and the results with high efficiency and accuracy can be achieved.Relationship between the first three natural frequencies and sampling numbers of plane rigid frame. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


## 4.1. Bending of the Timoshenko Beam

A clamp supported Timoshenko beam under a uniform load (as shown in Figure5) is considered in this example. The parameters are mentioned above. It is demonstrated in Table 1 that the proposed BSWI element has a good agreement with the exact solution [38].Table 1
Comparison of center deflection of clamp supported beam under uniform load.
Mesh method (DOFs)
Center deflection [mm]
80 traditional FEM (162)
0.14114
160 traditional FEM (322)
0.15399
320 traditional FEM (642)
0.15758
1 BSWI FEM (11)
0.15876
Exact solution [38]
0.15802Figure 5
The Clamp supported Timoshenko beam under a uniform load.


## 4.2. Stochastic Analysis of the Timoshenko Beam Bending

A stochastic analysis for a clamp supported Timoshenko beam under a uniform load as shown in Figure5 is considered in this numerical example. The results of BSWI-SFEM are shown in Table 2 compared with conventional FEM. It is evident that 1 BSWI-SFEM element can achieve a similar accuracy as 320 conventional Timoshenko elements. However, the proposed method only needs 0.078 s to establish the computation, while the conventional method needs 154.875 s as shown in Figure 6.Table 2
Comparison of center deflection of clamp supported beam under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [mm]
Time consumption [s]
80 traditional FEM (162)
0.14660
0.781000
160 traditional FEM (322)
0.16066
5.765000
320 traditional FEM (642)
0.16462
154.875000
1 BSWI FEM (11)
0.16593
0.078000
Sample exact solution
0.16368
—Figure 6
Time consumption comparison of center deflection for clamp supported beam under uniform load stochastic analysis.Figure7 presents the mean value of deflection at the center of beam which varies with the counter of sample. It can be seen that 1 BSWI-SFEM element obtains a good approximation of the center deflection compared with conventional method. Although the conventional method could approximate the sample exact solution with mesh refinement, the time consumption increases sharply. This kind of increase would result in unacceptable time consumption for the conventional method when it is used in large-scale structures. However, owing to the efficiency and accuracy of BSWI-SFEM, this problem can be greatly eased.Relationship between the mean value of beam center deflection and sampling number: (a) sampling number 1–20 and (b) sampling number 20–100.
(a)
(b)


## 4.3. Stochastic Analysis of the Timoshenko Beam Vibration

Figure8 shows the vibration analysis results of the clamp supported Timoshenko beam. BSWI results are compared with traditional finite element method. It can be seen that 1 BSWI-SFEM element can obtain a good approximation of natural frequencies compared with traditional method. Although the conventional method could approximate the sample exact solution with mesh refinement, the time consumption increases sharply. In this numerical example, the proposed BSWI method only needs 6.005 s to establish the computation with 1 element, while the conventional method needs 240.758 s to obtain the results in a similar accuracy with 320 elements.Relationship between the first three natural frequencies and sampling numbers of the Timoshenko beam. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


## 4.4. Bending of the Mindlin Plate

In this part, the Mindlin plate bending problem as presented in Figure9 is used to test the accuracy of the BSWI Mindlin plate element. The corresponding parameters are chosen as the mean value of these stochastic parameters mentioned above. The numerical results given in Table 3 indicate that the proposed method has a good agreement with the exact solution [38] once more.Table 3
Comparison of center deflection of the clamp and simply supported Mindlin plate under uniform load.
Mesh method (DOFs)
Center deflection [μm]
The clamp supported Mindlin plate
The simply supported Mindlin plate
100 traditional FEM (363)
7.9426
22.598
400 traditional FEM (1323)
7.9676
22.637
900 traditional FEM (2883)
7.9722
22.644
1 BSWI FEM (363)
7.9761
22.651
Exact solution [35]
7.9760
22.649The Clamp and simply supported Mindlin plate under a uniform load: (a) clamp supported and (b) simply supported.
(a)
(b)


## 4.5. Stochastic Analysis of the Mindlin Plate Bending

Stochastic analysis of a square Mindlin plate with clamp and simply supported boundary conditions subjected to a uniform load as shown in Figure9 is considered in this numerical example. The stochastic bending analysis of this plate is carried out by BSWI-SFEM incorporated with the Monte Carlo simulation. The center deflection of the two kinds of boundary conditions is presented in Tables 4 and 5. The results show that the BSWI-SFEM can give a better approximation for the sample exact solution compared with the conventional method. Furthermore, the time consumption as shown in Figures 10 and 11 indicates that BSWI-SFEM is more efficient than the conventional method, which allows the proposed approach to be used for the large-scale structures. The further comparisons are presented in Figures 12 and 13. Using the mean value of output variables, the reliability calculations can be thereafter performed.Table 4
Comparison of center deflection of the clamp supported Mindlin plate under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [μm]
Time consumption [s]
100 traditional FEM (363)
8.0250
6.766000
400 traditional FEM (1323)
8.0506
90.562000
900 traditional FEM (2883)
8.0553
824.875000
1 BSWI FEM (363)
8.0592
4.015000
Sample exact solution
8.0590
—Table 5
Comparison of center deflection of the simply supported Mindlin plate under uniform load (stochastic analysis).
Mesh method (DOFs)
Center deflection [μm]
Time consumption [s]
100 traditional FEM (363)
22.883
6.821000
400 traditional FEM (1323)
22.922
93.610000
900 traditional FEM (2883)
22.930
816.320000
1 BSWI FEM (363)
22.936
3.986000
Sample exact solution
22.934
—Figure 10
Time consumption comparison of center deflection for the clamp supported Mindlin plate under uniform load stochastic analysis.Figure 11
Time consumption comparison of center deflection for the simply supported Mindlin plate under uniform load stochastic analysis.Relationship between the mean value of plate center deflection and sampling number: (a) sampling number 10–35, (b) sampling number 35–100.
(a)
(b)Relationship between the mean value of plate center deflection and sampling number: (a) sampling number 15–50, (b) sampling number 50–100.
(a)
(b)


## 4.6. Stochastic Analysis of the Mindlin Plate Vibration

The vibration analysis results of the clamp supported Mindlin plate are shown in Figure14. The first three natural frequency results of BSWI element are compared with conventional element. It can be seen from the results that 1 BSWI element can obtain the results in a similar accuracy with 320 traditional elements. However, 1 BSWI method only needs 10.397 s to establish the computation, while the conventional method needs 4320.385 s to obtain the results in a similar accuracy with 320 elements.Relationship between the first three natural frequencies and sampling numbers of the Mindlin plate. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


## 4.7. Stochastic Analysis of Plane Rigid Frame Bending

As shown in Figure15, the plane rigid frame combined by 15 components is considered in this example. We assume that all the parameters obey the normal distribution. The corresponding material parameters are elastic modulus E=2.06×1011 N/m2, width b=0.012 m, height h=0.012 m, and density ρ=7917 kg/m3, respectively. The standard deviations of these parameters are 5% of their mean values.Figure 15
Coordinates, material parameters, and load parameters of plane rigid frame.As shown in Figure16, the bending analysis results of BSWI element on points A and B are compared with BEAM3 element. Although BEAM3 element can obtain the results with high accuracy by grid refinement, it will spend much more time on computation. However, the same solving domain of plane rigid frame is discreted by only 15 BSWI elements; it can achieve the results in a similar accuracy with 480 BEAM3 elements. So the efficiency and solving precision of BSWI-SFEM plane rigid element are verified.Relationship between the bending analysis results and sampling number of the plane rigid frame. (a) Point A and (b) point B.
(a)
(b)


## 4.8. Stochastic Analysis of Plane Rigid Frame Vibration

The stochastic vibration analysis of plane rigid frame is carried out by BSWI-SFEM incorporated with the Monte Carlo simulation in this example. The corresponding analysis results are shown in Figure17. BSWI element results are compared with BEAM3 element. It can be seen from the comparison that, although BEAM3 element can acquire the results with high accuracy by grid refinement, much more time is spent. BSWI-SFEM is superior, with only 15 elements, and the results with high efficiency and accuracy can be achieved.Relationship between the first three natural frequencies and sampling numbers of plane rigid frame. (a) First-order circular frequency, (b) second-order circular frequency, and (c) third-order circular frequency.
(a)
(b)
(c)


## 5. Conclusions

Based on B-spline wavelet on the interval and the Monte Carlo simulation, the stochastic finite element method is developed in the paper. The proposed method can deal with the response variability resulting from the spatial variability of the structural material properties. Taking the Timoshenko beam and the Mindlin plate as examples, the corresponding model are constructed and solved. The numerical examples illustrate that the proposed method can achieve a high analytical accuracy and efficient convergence in solving the stochastic problems. Although only three kinds of stochastic BSWI elements are constructed in the paper, further extensive work is easy to be established for solving the buckling and vibration problems of structures based on the same interpolation functions and approaches. The BSWI-SFEM is feasible to carry out the reliability analysis of large-scale engineering structures.



---

*Source: 104347-2014-05-08.xml*