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 80 (2017) 116–126

Contents lists available at ScienceDirect

Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound

A stochastic perturbation edge-based smoothed finite element method for the analysis of uncertain structural-acoustics problems with random variables

MARK

⁎

F. Wua, L.Y. Yaoa, , M. Hua, Z.C. Heb a b

College of Engineering and Technology, Southwest University, Chongqing 400715, PR China State Key Lab of Advanced Technology for Vehicle Body Design & Manufacture, Hunan University, Changsha 410082, PR China

A R T I C L E I N F O

A BS T RAC T

Keywords: Stochastic analysis Edge-based smoothed finite element (ES-FEM) SP-ES-FEM Uncertainties Structural-acoustics problems

Among the current methods in predicting the response of structural-acoustics problems in mid-frequency regime, some problems such as low accuracy and inability to deal with the uncertainties still need to be solved. To eliminate these issues, a novel stochastic perturbation edge-based smoothed FEM method (SP-ES-FEM) is proposed for the analysis of structural-acoustics problems in this work. The edge-based smoothing technique is applied in the standard FEM approach to soften the over-stiff behavior of structural-acoustics problems aiming to improve the accuracy of deterministic response predictions. Then, this approach, for the first time, intends to introduce the first-order perturbation technique into the edge-based smoothed FEM theory frame especially for the probabilistic analysis of structural-acoustics problems. The response of the coupled systems can be expressed simply as a linear function of all the pre-defined input variables by using the change of variable techniques. Due to the linear relationships of variables and response, the probability density function and cumulative probability density function of the response can be obtained based on the simple mathematical transformation of probability theory. The proposed approach not only improves the numerical accuracy of deterministic output quantities with respect to a given random variable, but also can handle the randomness well in the systems. Two numerical examples for frequency response analysis of random structural-acoustics are presented and verified by Monte Carlo simulation, to demonstrate the effectiveness of the present method.

1. Introduction The structural-acoustics problems are commonly encountered in numerous engineering systems whenever there is a dynamic loading, such as land vehicles, sea vessels and the aircrafts. The frequency interesting of the structural-acoustics systems range often from the “low-frequency” up to “high-frequency” regime, depends on the actual application. Different types of numerical simulation techniques have been developed for predictions, aiming to cover different frequency ranges. Most of these numerical prediction techniques can be categorized largely in two big groups: either deterministic or statistical methods. The former is suitable for lower frequency range, and the later for higher. Finite element method (FEM) [1–3], a typical deterministic method, is mostly successfully used to model low frequency vibrational and acoustical behavior of structural-acoustic systems. However, the prediction for the frequency response of a complex structural-acoustic system becomes exponentially difficult as frequency increases. In

⁎

practical level, two causes contribute to the difficulties in the higher frequency prediction. One main cause lies in the fact that a precise mathematical model can hardly be devised to capture the detailed deformation of the system, due to the decreasing wavelengths of deformation and the increasing model size. In the context of FEM, there is a “rule of the thumb” applied in controlling the elements sizes in order to obtain a relatively proper numerical solution [3]. Therefore, the computational effort of FEM will increase exponentially with the increasing of the frequency. The other cause is that the system response becomes more and more sensitive to small imperfections in the system with the increasing of the frequency, so that small manufacturing uncertainties may lead to huge variability in the frequency response. That means even a very detailed deterministic mathematical model cannot yield a reliable response prediction in a relatively higher frequency regime [4]. In recent years, the concern of developing a suitable estimation technique for structural-acoustic problems has shifted to the midfrequency regime, where structural-acoustic systems exhibit the lot of

Corresponding author. E-mail address: [email protected] (L.Y. Yao).

http://dx.doi.org/10.1016/j.enganabound.2017.03.008 Received 11 September 2016; Received in revised form 23 January 2017; Accepted 4 March 2017 0955-7997/ © 2017 Elsevier Ltd. All rights reserved.

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

instable in the calculation in the dynamic analysis [31]. Then, some researcher constructed the smoothing domain based on the edges of the elements and established a new method called ES-FEM [32–36]. In comparison with the traditional FEM or the above mentioned NS-FEM, the ES-FEM possesses suitable soften effects and help to build a reasonable system stiffness. And these effects enable the ES-FEM model to show neither “overly soft” nor “overly stiff” features. Therefore, the ES-FEM is able to reduce the desperation errors and thus enhance the accuracy of the prediction results. And the ES-FEM results are often found super convergence and ultra-accurate [32–36]. Given the superior performance of the edge-based smoothing technique (ES), it is natural to expect that the ES-FEM will greatly reduce the numerical dispersion error and obtain accurate results for acoustic problems. Thus, the embedding of the edge-based smoothing technique (ES) and change-of-variable technique into a similar hybrid stochastic perturbation FEM frame (SP-ES-FEM)presents a promising combination for solving structural-acoustic problem. The newly constructed ESSP-FEM will be expected to efficiently and effectively evaluate the response probability density function and cumulative density functions. This is also the essential motivation of this work. Based upon the above analysis, a stochastic perturbation ESFEM(SP-ES-FEM), which combines the first order stochastic perturbation technique with the edge-based smoothed finite element approach is proposed for structural-acoustic coupling problems. The paper is organized as follows: in the Section 2, the basic principles of probability theory and the ES-FEM for 3D acoustics problem are briefly described, respectively. In the Section 3, the hybrid SP-ES-FEM equations are constructed in order to achieve a first order perturbation stochastic approach for structural-acoustic coupling problems problem. In Section 4, numerical examples and applications are presented to demonstrate the performance of the SP-ES-FEM for the frequency response analysis of structural-acoustic coupling problems. Finally, a summary is given in Section 4 to conclude this work.

uncertainties. Considering the difficulties involved in developing a reliable deterministic method for the prediction of higher frequency domain, lots of research efforts are contributed to develop an alternative that takes into account model uncertainties and variability. For dealing with the uncertainties involved in the structural-acoustics problems, probabilistic methods are most commonly chosen with acceptable efficiency, where the uncertainty is expressed using different type random variables [5]. The most robust and simplest probabilistic method is the well-known Monte Carlo simulation method (MCM). It often can offer the most direct and reliable result at the expense of the large computational load, and therefore MCM results are chosen as the reference result in this work [6,7]. In the context of finite element method, there are another two classical probabilistic variants: the spectral stochastic finite element method (SSFEM) [8] and stochastic perturbation finite element method approach (SPFEM) [9,10]. In SSFEM, several of random Hermite polynomials is used to approximate the response quantity of the systems [8]. In SPFEM, all the response vectors or matrixes related to the pre-defined random input variables are needed to be expanded based the Taylor series expansion [9,10]. All the quantities can be expanded at the first order or higher order approximations according to the Taylor series. In this work, for a low variability level of the design parameters, only the first order expansion is applied in the process without losing much accuracy, due to the its feasible strategy [11,12]. If only the first order approximations are applied, the relationship between the response and the random variables can be approximated as the linear functions. Therefore, the computational process can be largely simplified, and the expectations and the standard deviations of the response can be easily obtained through the traditional SPFEM. However, the calculation of the probability distributions of the response is commonly ignored, unless the random responses is able be approximated as Gaussian random fields [13]. In this work, several different random input variables with different distribution are defined, where the random responses is apparent non-Gaussian. In order to observe the distributions of the random response, a change-of-variable technique is introduced [14,15]. There are two important steps for the calculation of the probability density functions (PDF) and cumulative density functions (CDF): firstly, the first order Taylor expansion technique helps to produce the linear functions between the inputs and responses. Secondly, the change-of variable can help to produce the PDF and CDF of the responses. Thus, the combination of classical SPFEM and the change-of-variable technique provide a promising approach for the comprehensive analysis of uncertainties. The deterministic method (FEM) is the key part for guaranteeing the accuracy in the probabilistic analysis. In the context of FEM, some more efficient techniques have been proposed to extend the frequency range of deterministic method(FEM) that are usual meant for low frequency problems. This category includes: stabilized methods [16,17],enriched methods [18,19], reduction techniques [20,21], higher order techniques [22,23],wave finite elements [24], and the ultraweak variational formulation [25]. However, the use of the deterministic method (FEM)always entails some inherent drawbacks, which are closely associated with the well-known “overly-stiff” feature of FEM and its sensitivity to numerical pollution errors. Some researches [26– 29] shows that the fundamental solution for eliminating the pollution errors is to soften the “stiffness” of the numerical structural-acoustic systems. Actually, a series of “soften” or “smoothing” gradient techniques have already been applied successfully into the traditional FEM frame [30–36], which are originally and extensively studied in meshfree methods [37,38] Smoothing gradient operations are performed on the newly constructed element integrate domains (smoothing domains), which can be built upon the nodes or edges of element. If the smoothing domain is constructed based on the nodes of elements, the node-based smoothed finite element method (NS-FEM) [30] can be built. However, research found that NS-FEM has overly soften the system stiffness, leading to an “overly-softy” feature. It is temporal

2. Basic principles for deterministic structural-acoustic systems Consider a structural-acoustic problem with thin elastic plate and an acoustic cavity, as shown in Fig. 1. The structural-acoustic system can be divided into a plate subsystem and an acoustic subsystem as shown in Fig. 1. Two subsystems are coupled exclusively through the interface ΓSd which is described by a set of degrees of freedom. The vibrating structure is modeled using the ES-FEM method [37,38], where the 3-node linear triangle plate elements are based on the Reissner-Mindlin plate theory. For overcoming the shear locking problem existing in the Reissner-Mindlin plate theory, the well-known discrete shear gap (DSG) technique is applied [39-40]. The ES-FEM formulation for structural domain then can be written as:

Ms u¨ + Ks u = Fs , in Ωs

(1)

in which u is defined as the unknown field variable displacements, Ms ,Ks is the lumped mass matrix and the edge-based smoothed

Fig. 1. Thestructural-acoustic coupling system.

117

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

stiffness matrix of plate, respectively Fs is the external force. More detailed derivation can be referred to Ref [37,38].

Cell i

2.1. ES-FEM formulations for 3D acoustic problems In this section, the ES-FEM formulations for 3D acoustic problems are briefly introduced.

n n

2.1.1. GS-Galerkin weak form and discretied system equations For acoustic domain, we first assume that the fluid is homogeneous, inviscid, and compressible and only undergoes small translational movement. Considering an enclosed cavity Ωf with Neumann boundary ΓN ,let p denote the acoustic pressure and k represent the wave number. Then the well-known governing equation for acoustic domain can be written as below:

Δp + k 2p = 0, in Ωf

Ff =

= Np (4)

∫Ω ∇N⋅∇NpdΩ + ∫Ω N⋅NpdΩ − jρω ∫Γ

N⋅vn dΓ = 0

N

v (x k ) =

∫Ω ∇N⋅∇NpdΩ + k 2 ∫Ω N⋅NpdΩ − jρω ∫Γ

N⋅vn dΓ = 0

N

v (x k ) = −

(5)

v (x k ) = −

∫Ω

v (x)dΩ

(12)

k

∫Ω dΩ denotes the volume of smoothing domain for edge k. k

1 jρωVk

∫Ω

∇pdΩ = −

k

1 ∑ Bi (xk) pi jρω I ∈ M k

1 jρωVk

∫Γ

k

p⋅n dΓ

(13)

(14)

where Mk represents the total number of nodes in the smoothing domain of edge k. Bi is defined as:

(6)

T

Bi (xk) = [ bi1 bi2 bi3 ]

bip =

(7)

1 Vk

∫Γ

(15)

Ni (x) np (x)dΓ

k

(16)

Finally, the smoothed stiffness matrix shown in Eq. (8) can be assembled based on the smoothed B as:

∫Ω (∇N)T ∇NdΩ

The smoothed acoustical stiffness matrix

K (k ) =

(8)

Mf =

(11)

Substituting the field variable (acoustic pressure) interpolation in form of Eq. (4) into Eq. (13), the smoothed velocity for edge k can be denoted as the following matrix form.

where:

Κf =

Nodal acoustic pressure in the domain

(10)

According to Green's theorem, the smoothed velocity can be expressed in terms of acoustic pressure:

Finally, the discretized system equations in Eq. (6)can be written in following matrix form:

M f {¨} p + K f {p} = −jρω {Ff }

1 Vk

where Vk =

In smoothed finite element formulations, the gradient component ∇N is replaced by a new smoothed item ∇N based on the edges of elements, for the purpose of “softening”the whole system stiffness. Then, the smoothed Garlerkin weak formulation for acoustic problem can be simply rewritten as:

−

The vector of nodal acoustic forces

2.1.2. Edge-based gradient smoothing operation for 3D acoustic domain In this section, the edge-based gradient smoothing technique for 3D acoustic problem is briefly introduced. First the acoustic domain is discretized using four node tetrahedral elements as that of standard FEM. Based on the edges of these tetrahedral elements, we can construct the edge-based gradient smoothing domains, which are also serving as integration domains. As shown in Fig. 2, the sub-smoothing domain of edge k in celli is created by connecting the centroid of celli to the two end-nodes of the edge k and the related surface triangles. Over each edge-based smoothing domain, the velocity v deduced by the gradient of acoustic pressure can be smoothed by gradient smoothing techniques as follows:

where pi is the unknown nodal pressure needed to be determined. Ni is the shape function in node i. N represents the generated shape function. p denotes the vector of generated pressure for each tetrahedron element. In the standard Galerkin weak formulations, the shape function is usually chosen as the weight function. Substituting the Eq. (4) into Eq. (2), then the standard Galerkin weak form for acoustic problem without acoustical damping can be expressed as:

−

T

N νn dΓ

Note that the lumped matrix is chosen during the analysis.

m

k2

∫Γ

{p}T = { p1, p2 , …, pn}

(3)

i =1

n

N

where j = −1 , ρ is known as the density of medium. vn is normal velocity on the boundary. In FEM, the continuous field variable pressure can be approximated using the shape function, and expressed as:

∑ Ni pi

n

n

Fig. 2. 3D edge-based smoothing domains constructed by connecting the centroid of cell i to end-nodes of the edge k and the related surface triangles.

where Δ denotes the Laplace operator, the wave number is calculated by:k = ω / c , in which ω is the angular frequency of the pressure oscillation, and c is the speed of sound traveling in the acoustic fluid. The boundary of the acoustic domain consists of three distinct boundary conditions: Dirichlet, Neumann and Robin conditions. In this work, the Neumann boundary condition is mainly used. The governing equation for the Neumann boundary of the interior acoustic domain can be expressed as below:

p=

n

Edge k

(2)

∇p⋅n = −jρωvn on ΓN

n

n

T

∫Ω N NdΩ

The acoustical mass matrix

∫Ω

k

T

B BdΩ =

Ns

∑ Vk BTB k =1

(17)

This assembled smoothed stiffness matrix in Eq. (17) is symmetric and banded, thanks to the compact supports of FEM shape functions.

(9) 118

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

n

Therefore the whole system equations can be solved efficiently. Without considering any structural damping, the discrete system of equations for a structural-acoustic system subjected to harmonic external excitations can be expressed as:

∑

D (b) = D (b) b = E (b) +

i =1

n

∑

F (b) = F (b) b = E (b) +

⎡ K f − ω 2 Mf ρω 2 HT ⎤ ⎡ p ⎤ ⎡ Ff ⎤ ⎢ ⎥⎢ ⎥ = ⎢ ⎥ ⎢⎣ Ks − ω 2 Ms ⎥⎦ ⎣ u ⎦ ⎣ Fs ⎦ −H

i =1

∂D (b) ∂bi ∂F (b) ∂bi

(bi − E (bi )) = D 0 + ΔD (23)

b = E (b)

(bi − E (bi )) = F 0 + ΔF (24)

b = E (b)

(18) Substituting Eqs. (25)–(26) into Eq. (20), we arrive at

where H is the spatial coupling matrix which states mutual transformation relationship between the interface DOFs u and the DOFs p of the model of the two different subsystems. Considering the velocity continuity at the interface, spatial coupling matrix can be written as:

H=

∫ NTs nf Nf dΓ

U = (D0 + ΔD)−1 (F 0 + ΔF)

(25)

According to Neumann expansions theory, in the case of spectral −1 −1 radius of ((D0) ΔD) is not big than 1, the expression of ((D0 + ΔD) ) can be expanded using Neumann series, thus we arrive at:

(19) −1

(D0 + ΔD)

where Ns ,Nf are the linear shape function of the interface and acoustical domain, respectively.nf is the outward normal vector along the interface of the acoustical domain. For easy reference, the Eq. (18) can be rewritten as:

−1

−1

−1

= (D0) −(D0) ΔD (D0)

(26)

Substituting Eq. (28) into Eq. (27), different perturbation term of U in Eq. (27) can be further written as: −1

U 0 = (D 0 ) F 0

(20)

DU = F

−1

ΔU1 = −(D0) (ΔDU 0−ΔF)

⎡ K f − ω 2 Mf ρω 2 HT ⎤ ⎥ is defined as the dynamic With D = ⎢ ⎢⎣ −H Ks − ω 2 Ms ⎥⎦ ⎡ p⎤ stiffness matrix of total system,U = ⎢ ⎥ is the vector of degrees of ⎣u⎦ ⎡ Ff ⎤ freedom of whole vibro-acoustic system,F = ⎢ ⎥ represents the total ⎣ Fs ⎦ external load vector.

−1

−1

ΔU 2 = −(D0) ΔDU1ΔUn = −(D0) ΔDUn −1

−1

n

∑ i =1

1 2!

n

n

i =1 j =1

F (b) = F (b) b = E (b) +

∑ i =1

1 + 2!

n

n

b = E (b)

⎛ n ⎜∑ ∂D (b) ⎜ ⎝ i =1 ∂bi

∑∑ i =1 j =1

⎞ (bi − E (bi )) ⎟ × ((D (b) b = E (b) )−1F (a) b = E (b) ⎟ ⎠ b = E (b)

(29)

Based on the Eq. (31), a simple mathematical transformation is carried on, then the unknown vector U can be written in a matrix and vector forms as follows:

U = (D (b) b = E (b) )−1F (a) b = E (b) −

(bi − E (bi ))

−1

∑ (D (b) b=E (b) ) i =1

b = E (b)

⎛ n ⎜∑ ∂F (b) ⎜ ⎝ i =1 ∂bi

(bi − E (bi ))(bj − E (bj )) + ⋯ b = E (b)

n

−

∑

b = E (bi)=1

∂D (b) ∂bi ⎛

n

+

−

b = E (b)

(bi − E (bi ))(bj − E (bj )) + ⋯

⎞ (D (b) b = E (b) )−1F (a) b = E (b) ⎟ E (bi ) ⎟ ⎠ b = E (b)

n

∑ (D (b) b=E (b) )−1 ⎜⎜∑ ∂F (b) ⎝ i =1

i =1

(bi − E (bi ))

∂F (b) ∂F (b) ∂bi ∂bj

(28)

U = ((D (b) b = E (b) )−1F (a) b = E (b) + (D (b) b = E (b) )−1 ⎛ n ⎞ ⎜∑ ∂F (b) (bi − E (bi )) ⎟ − ((D (b) b = E (b) )−1 ⎜ ⎟ ⎝ i =1 ∂bi ⎠

While random external force vector can be expressed as:

∂F (b) ∂bi

−1

Substituting Eqs. (25)–(26) into Eq. (30), we arrive at:

(21)

n

−1

n

∂D (b) ∂bi

∂D (b) ∂D (b) ∂bi ∂bj

∑∑

−1

U ≈ U0 + ΔU1 = (D0) F 0 + (D0) ΔF − (D0) ΔD (D0) F 0

The uncertainties of parameters are commonly encountered in numerous structural-acoustic system, and should be carefully handed in computational modeling and simulation. In this subsection, the SPES-FEM is constructed based on the traditional stochastic perturbation method, ES-FEM and the change-of variables techniques. Let us first introduce a given set of random field b(x) with its probability density function (PDF) and the expected value as pi (br ),Ei (br ), respectively, r=1,2…,R, i=1,2. Firstly, the Taylor expansions operations are performed on the dynamic stiffness matrix of whole structural-acoustic system and random external force vector at the expectation of the random vector. The whole dynamic stiffness matrix can be rewritten as follow:

D (b) = D (b) b = E (b) +

(27)

In order to simplify the computational complexity and eliminate computational cost, just as the as the low order Taylor expansions, only the first order perturbation of Neumann series is kept in the calculation while all the higher order perturbation terms are ignored artificially. Thus, the unknown vector U is can be written as:

3. Formulation for hybrid SP-ES-FEM

+

−1 −2

−1

− (D0) (ΔD (D0) ) ⋯

∂D (b) ∂bi

∂bi

b = E (b)

⎞ (D (b) b = E (b) )−1F (a) b = E (b) ⎟ bi = B0 + ⎟ b = E (b) ⎠

n

∑ Bi bi i =1

(30)

In which, the detailed expressions of B0 and Bi can be written as:

b = E (b)

n

(22)

B0 = (D (b) b = E (b) )−1F (a) b = E (b) −

−1

∑ (D (b) b=E (b) ) i =1

Taking into account the high order of the complexity of the calculation process and the calculation of the cost, the higher order of Taylor expansions terms are ignored for simplification. Therefore, the first order of Taylor expansions of matrix and vector can be expressed as:

⎛ n ⎜∑ ∂F (b) ⎜ ⎝ i =1 ∂bi

n

−

∑

b = E (b) i =1

∂D (b) ∂bi

⎞ (D (b) b = E (b) )−1F (a) b = E (b) ⎟ E (bi ) ⎟ ⎠ b = E (b) (31)

119

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

⎛ n ∂F (b) Bi = (D (b) b = E (b) )−1 ⎜∑ ⎜ ⎝ i =1 ∂bi

− b = E (b)

∂D (b) ∂bi

(D (b) b = E (b) )−1 b = E (b)

⎞ F (a) b = E (b) ⎟ ⎟ ⎠

(32)

The unknown vector U can be written in a more simple and apparent vector forms as follows:

⎧ B0,1 ⎫ ⎧ Bi,1 ⎫ ⎧ u1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⋮ ⋮ ⎪ ⎪ ⎪ ⋮ ⎪ ⎪ ⎪ U = ⎨ uk ⎬; B0 = ⎨ B0, k ⎬; Bi = ⎨ Bi, k ⎬ ⎪⋮⎪ ⎪ ⋮ ⎪ ⎪ ⋮ ⎪ ⎪u ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ m⎭ ⎩ B0, m ⎭ ⎩ Bi, m ⎭

n

U = B0 +

∑ Bi bi i =1

(33)

In which uk can be calculated as follows: Fig. 3. An flexible steel plate backed by a box of air.

n

∑ Bi,k bi

uk = B0, k +

(34)

i =1

As previous deviations, the relationship between the response of U and the input random b(x) has already been transformed to the linear functions as shown in Eq. (36). In order to obtain the response probability density function of the structural-acoustics, the changeof-variable technique is applied in this work. Because of the linear relationship between the response and inputs, the probability density function of response can be calculated steps by steps. If there are only two inputs, according to the change-of-variable technique, the probability density function of can be expressed as

fy1, r y2, r ( y1, r , y2, r ) = fy1, r ( y1, r ) fy2r ( y2, r ) = fy1, r ( y1, r ) fy2r (ur − B0, r − y1, r )

(35) Fig. 4. Three randomly chosen response points.

where y1, r = B1, k b1, y2, r = B2, k b 2 , fy1, r ( y1, r ), fy2r ( y2, r ) are the probability density of the inputs. If we assume that then using the change-of-variable technique, the probability density function of can be expressed as According to the change-of variable technique, the detailed derivations of the probability density function are expressed as follows: ∞

∞

fur (ur ) = ∫ y

1, r =−∞

∞

= ∫ y

1, r =−∞

=

∫y

2, r =−∞

Table 1 The pre-defined random variables of the structural-acoustic coupling model.

fy1, r y2r ( y1, r , y2, r ) dy1, r dy2, r

fy1, r ( y1, r ) fy2, r ( pr − B0, r − y1, r ) dy1, r

⎛ ur − B0, r − y1, r b1 ⎞ ∞ 1 f (b1) fa2 ⎜ ⎟ db1, r = 1, …, L ∫ B2, r b1=−∞ b1 B2, r ⎠ ⎝

1 B2, r

∞

∫x =−∞ fb (b1) fb 1

1

2

⎛ Y1, r − B1, r b1 ⎞ ⎟ db1 ⎜ B2, r ⎠ ⎝

Distribution type

Parameters of random variable

Modulus of elasticityE Poisson ratio ν

Gaussian distribution

μE = 2.1 × e11

σE = 5 × e9

Log-normal distribution Gaussian distribution Gumbel distribution

μν = −1.2

σν = 0.033

μt = 0.0015 μFs = 1

σt = 0.00015 σ Fs = 0.4

Thickness t The external force Fs

(36) 4. Numerical application of the hybrid SP-ES-FEM

If more than two random inputs are involved, we assume Y1, r = y1, r + y2, r , Thus, the probability density function of these two inputs are calculated just as the above case as follows:

fY1, r (Y1, r ) =

Random variables

In order to analysis the properties of the SP-ES-FEM, two numerical examples will be presented to study the frequency response of the systems. The first example is a simple rectangular plate coupled with closed acoustic cavity subjected to concentrated dynamic load traction at the center of the plate. For the second problem, the hybrid SP-ES-FEM is applied to automotive industry for higher accurate and more efficient mid-frequency solutions.

(37)

Then we assume Y2, r = Y1, r + y3, r ,y3, r = B3, k b3, based on the changeof variable technique, the PDF of Y2, r can be calculated similarly as: ∞

fY2, r (Y2, r ) =

∫b =−∞ fb 3

3

(b3) fY1, r (Y2, r − B3, r b3) db3

(38)

4.1. Validation of SP-ES-FEM

After several repeat calculations, the final PDF can be expressed as:

fYz, r (Yz, r ) =

∫b

The hybrid SP-ES-FEM approach is first applied to a standard example to demonstrate the validity and efficiency of this method. The structural-acoustics model consists of a flexible plate( ρ =7800 kg / m3, ν =0.3 and E=210 Gpa)and an acoustic cavity( ρ =1.225 kg /m3 and c = 340m / s ) as shown in Fig. 3. The closed acoustic cavity has dimensions of 500×245×200 mm and is divided by 1646 tetrahedron elements and 460 nodes, while the plate has dimensions of 500×245 mm, and is divided by 148 triangle elements and 184 nodes

∞

z +1=−∞

fbz +1 (bz +1) fYz −1, r (Yz, r − Bz +1, r bz +1) dbz +1 fur ∞

(u r ) =

∫b =−∞ fb n

n

(bn ) fYn −2, i (ur − B0, r − Bn, r bn ) dbz , r = 1, …, L (39)

Once the PDF is obtained, the CDF results can be easily obtained by simple mathematical transform. 120

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

0.8

F(X)

0.6

SP-ES-FEM MCM

0.4 0.2 0 -6

-4

-2

0 X

2

4

6

Fig. 5. The probability density functions of the frequency response (100 Hz) of Node 1calculated by MCM and SP-ES-FEM.

1

CFD

0.8 SP-ES-FEM MCM

0.6 0.4 0.2 0 -4

-3

-2

-1

0 X

1

2

3

4

Fig. 6. The cumulative distribution functions of the frequency response (100 Hz) of Node 1calculated by MCM and SP-ES-FEM.

3 2.5

SP-ES-FEM MCM

F(X)

2 1.5 1 0.5 0 -1

-0.8

-0.6

-0.4

-0.2

0 X

0.2

0.4

0.6

0.8

1

Fig. 7. The probability density functions of the frequency response (200 Hz) of Node 1calculated by MCM and SP-ES-FEM.

1

CFD

0.8

MCM GP-ES-FEM

0.6 0.4 0.2 0

-0.4

-0.2

0

X

0.2

0.4

0.6

0.8

Fig. 8. The cumulative density functions of the frequency response (200 Hz) of Node 1calculated by MCM and SP-ES-FEM.

these comparisons will be carried on the next real engineering example, and the conclusions are quite similar. In the random structuralacoustic systems, modulus of elasticityE, Poisson ratio ν , the thickness of the plate t and the external force Fs are considered as the random variables with some specified probability density distributions as shown in the Table 1. In order to further investigate the effectiveness of the proposed SPES-FEM, a Monte Carlo simulation method is used to calculate the frequency response of the above mentioned three random points. In

with a thickness of 1.5 mm. The walls (except the coupling wall) of cavity are assumed to be acoustically rigid with the surface velocity fixed at v=0. To illustrate the effectiveness of the SP-ES-FEM, three response points (Node 1, Node2, Node3) are randomly selected and observed as shown in the Fig. 4. Three response point coordinate position is as follow: Node 1 (46,443,111) mm, Node 2 (130,197,128) mm, Node 3 (105,239,200) mm. For the accuracy comparison, it is unnecessary to repeat the plots for accuracy comparison of frequency response, as

121

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

0.8 SP-ES-FEN MCM

F(X)

0.6 0.4 0.2 0 -5

-4

-3

-2

-1

0 X

1

2

3

4

5

Fig. 9. The probability density functions of the frequency response (100 Hz) of Node 2calculated by MCM and SP-ES-FEM.

1

CDF

0.8 SP-ES-FEM MCM

0.6 0.4 0.2 0 -3

-2

-1

0

1 X

2

3

4

5

Fig. 10. The cumulative density functions of the frequency response (100 Hz) of Node 2calculated by MCM and SP-ES-FEM.

3 2.5

SP-ES-FEM MCM

F(X)

2 1.5 1 0.5 0 -1.5

-1

-0.5

0

X

0.5

1

Fig. 11. The probability density functions of the frequency response (200 Hz) of Node 2calculated by MCM and SP-ES-FEM.

1

CFD

0.8

MCM SP-ES-FEM

0.6 0.4 0.2 0

-0.8

-0.6

-0.4

-0.2

X

0

0.2

0.4

0.6

Fig. 12. The cumulative density functions of the frequency response (200 Hz) of Node 2calculated by MCM and SP-ES-FEM.

It is found again, at Node 2, all the PDF or CDF calculated by SPES-FEM also agree well with the reference one (MCM). The follow Figs. 13–16 shows the frequency response results of the Node 2 at 100 Hz or 200 Hz. As shown in Figs. 5–16, it is found no matter in different points or at different interesting frequency, the proposed SP-ES-FEM always provide high accuracy results, which agree well with the MCM. Therefore it can be concluded that proposed SP-ES-FEM can deal with the randomness existing in structural-acoustic systems very well. To

MCM, 10,000 times calculation are carried on to obtain reliable results. The acoustic pressure in frequency 100 Hz, 200 Hz in three different points calculated by SP-ES-FEM. All these results are compared with the MCM results. The main numerical procedure can be referred to EqS. (27)–(39). The follow Figs. 5–8 shows the frequency response results of the Node 1 at 100 Hz or 200 Hz. As shows, at Node 1, all the PDF or CDF calculated by SP-ES-FEM agree well with the reference one (MCM). The follow Figs. 5–12 show the frequency response results of the Node 2 at 100 Hz or 200 Hz. 122

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

0.8

F(X)

0.6

SP-ES-FEM MCM

0.4 0.2 0 -5

-4

-3

-2

-1

0 X

1

2

3

4

5

Fig. 13. The probability density functions of the frequency response (100 Hz) of Node 3calculated by MCM and SP-ES-FEM.

1

CDF

0.8 SP-ES-FEM MCM

0.6 0.4 0.2 0 -3

-2

-1

0

1 X

2

3

4

5

Fig. 14. The cumulative density functions of the frequency response (100 Hz) of Node 3 calculated by MCM and SP-ES-FEM.

3 2.5

SP-ES-FEM MCM

F(X)

2 1.5 1 0.5 0 -1.5

-1

-0.5

0

X

0.5

1

Fig. 15. The probability density functions of the frequency response (200 Hz) of Node 3calculated by MCM and SP-ES-FEM.

1

CFD

0.8

MCM SP-ES-FEM

0.6 0.4 0.2 0

-0.8

-0.6

-0.4

-0.2

X

0

0.2

0.4

0.6

Fig. 16. The cumulative density functions of the frequency response (200 Hz) of Node 3 calculated by MCM and SP-ES-FEM.

4.2. Automobile passenger compartment with a flexible coping

further examine the efficiency of SP-ES-FEM, the CPU time of the SPES-FEM and the Monte Carlo simulation (MCM) codes calculated in the same computer is recorded. The CPU time of MCM (10,000 times) is about 18,235 s while SP-ES-FEM only need about 2.5 s. Through these direct comparison, without losing the computational accuracy, the efficiency of SP-ES-FEM is far higher than MCM does.

In automotive industry, continuous vehicles assembled from the same production line often exhibit huge variability in their measured responses. Meanwhile, a detailed FEM model with millions of degrees of freedom is usually required to predict the higher frequency response of vehicle, leading to severe computational inefficiency [39]. In this part, the hybrid SP-ES-FEM is applied to automotive industry for higher accurate and more efficient mid-frequency solutions. The 123

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

generated by the interaction of the floor panel and the passenger compartment cavity is a great concern, and it contributes strongly to the interior sound pressure level (SPL) in the automobile passenger compartment. As shown in Fig. 17, a simply FEM car model is established. The automobile passenger compartment is divided using 34,999 tetrahedron elements with 7942 nodes for ES-FEM. The results are compared against the reference result that is calculated using traditional coupled FEM model based on MCM. The vibration coming from the engine with the velocity of 5 m/s is applied to the floor panel of the passenger compartment. Let us now focus on the accuracy of the ESFEM. For the structural-acoustic problems, the frequency response of the middle points in the bottle face of the air box is computed from traditional FEM and the ES-FEM, where the same reasonable same discretization mesh is chosen. It is well known that if the FE mesh is enough small, the most correct result can be obtain through conventional FEM. Therefore, for the purpose of the comparison, the reference result are calculated using traditional FEM based on a very fine mesh. (34,999 tetrahedron elements with 7942 nodes). While the regular mesh is constructed based on 6122 tetrahedron elements with 1589 nodes. The frequency response between (100–290 Hz) of the plate is calculated using different methods as shown in Fig. 18. Obviously, as shown in the Fig. 18, ES-FEM is able to offer much accuracy results compared to FEM in frequency response calculation with the increasing of the interesting frequency. The performance of ES-FEM is excellent, where the ES-FEM's solutions are consistently more accurate than those methods mentioned above for structuralacoustic systems. From 200 Hz, the accuracy of FEM decreases with the increase of the frequency. FEM suffers from the dispersion errors mainly caused by the “overly-stiff” property. The result obtained from ES-FEM shows excellent agreements with the reference result, and it demonstrates that the pollution errors are great eliminated by the “soften effects” of ES-FEM. Even applied in the model with much less DOFs(1,589acoustic DOFs),ES-FEM is still able to provide the result in almost the same accuracy level as the result of FE-SEA using fine mesh (7942 acoustic DOFs). To sum up, the high accuracy and computational efficiency of ES-FEM for mid-frequency solution is illustrated in this example. In the further probabilistic analysis, the modulus of elasticityE, Poisson ratio ν , the thickness of the plate t and the external force Fs are considered as the random variables with some specified probability density distributions as shown in the Table 2. In order to further investigate the effectiveness of the proposed SPES-FEM, a Monte Carlo simulation method is used to calculate the frequency response of the driver's ear points at 100 Hz. In MCM, 2000 times calculation are carried on to obtain a reliable results. The acoustic pressure in frequency 100 Hz in driver ear left or right points calculated by SP-ES-FEM. All these results are compared with the MCM results. The main numerical procedure can be referred to Eqs. (27)–(39). The follow Figs. 19 and 20 shows the frequency response results of the left ear point at 100 Hz.

Fig. 17. An SP-ES-FEM model combined by the flexible coping and the passenger compartment.

180

FEM ES-FEM Reference

170

Response

160 150 140 130 120 110 100 100

150

200 Frequency

250

Fig. 18. The accuracy comparison of two methods. Table 2 The pre-defined random variables of the structural-acoustic coupling model. Random variables

Distribution type

Parameters of random variable

Modulus of elasticityE Poisson ratio ν

Gaussian distribution

μE = 2.1 × e11

σE = 5 × e9

Log-normal distribution Gaussian distribution Gumbel distribution

μν = −1.2

σν = 0.03

μt = 0.001 μFs = 5

σt = 0.00015 σ Fs = 2

Thickness t The external force Fs

vehicle body is made of panels, and is usually welded with numerous thin steel plates, among which, the floor structure of an automobile is one of the largest plates in the vehicle. The vibration and noise 8

x 10

7 6

SP-ES-FEM MCM

F(X)

5 4 3 2 1 0

-300

-200

-100

0

X

100

200

300

Fig. 19. The probability density functions of the frequency response (200 Hz) of left ear point calculated by MCM and SP-ES-FEM.

124

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al. 1 0.8 MCM SP-ES-FEM

CFD

0.6 0.4 0.2 0 -300

-200

-100

0 X

100

200

300

Fig. 20. The cumulative density functions of the frequency response (200 Hz) of left ear point calculated by MCM and SP-ES-FEM.

9

x 10

8 7

F(X)

6

SP-ES-FEM MCM

5 4 3 2 1 0 -200

-150

-100

-50

0

50 X

100

150

200

250

300

Fig. 21. The probability density functions of the frequency response (200 Hz) of right ear point calculated by MCM and SP-ES-FEM. 1 0.8 MCM SP-ES-FEM

CDF

0.6 0.4 0.2 0 -200

-150

-100

-50

0

50 X

100

150

200

250

300

Fig. 22. The cumulative density functions of the frequency response (200 Hz) of right ear point calculated by MCM and SP-ES-FEM.

of the present method.

The follow Figs. 21 and 22 shows the frequency response results of the right ear point at 100 Hz. The results for this complicated vehicle example re-enforces the finding from the previous simple example. The proposed SP-ES-FEM can offer good PDF and CDF results and it agrees well with the MCM results in the real car example. Also the efficiency of the SP-ES-FEM is very obvious, The CPU time of SP-ES-FEM is only about 150,088 s while MCN consume about 74 s.

(1) Due to the imbedded ES (Edge-based smoothing technique), the hybrid SP-ES-FEM significantly eliminates the dispersion error in comparison with conventional FEM in deterministic calculation. Also, in the probabilistic calculations, the higher accuracy in the deterministic computations apparently helps to enhance the accuracy of final results. (2) In the probabilistic analysis, the probability density functions and cumulative distribution functions of the structural-acoustic systems calculated by SP-ES-FEM agrees well the MCM does, moreover, it is found that SP-ES-FEM holds huge advantages in computation efficiency in comparison with MCM results without loss of much accuracy. (3) SP-ES-FEM is constructed based on the first order of Taylor expansions, if the further increase of the accuracy is needed, the higher order perturbation technique is recommended to be applied. However, higher order perturbation will significantly increase the computational burden.

5. Conclusions In this work, a hybrid SP-ES-FEM framework is proposed for predicting the average response of the structural-acoustic system in the mid-frequency regime. The proposed approach takes the best advantages of ES (Edge-based smoothing technique) to “soften” the whole stiffness of coupled structural-acoustic systems. Thus it significantly reduces the dispersion errors that are often observed in the standard FEM solution which is also become severer with the increase of frequency. What is more, the imbedded stochastic perturbation technique has dealt with the uncertainties well. Numerical examples of structural acoustic problems have demonstrated the following features 125

Engineering Analysis with Boundary Elements 80 (2017) 116–126

F. Wu et al.

Acknowledgements [19]

The author wishes to thank the support of the National Natural Science Foundation of China (NSFC) (Grant No 51605391), And also wants to thank the support of the Chongqing Science & Technology Commission (CSTC) (Grant No. 2016jcyjA0176) and Fundamental Research Funds for the Central Universities and Doctoral Fund(Grant No. XDJK2017B060 and No. SWU116017).

[20]

[21]

[22]

References

[23] [24]

[1] Everstine GC. Finite element formulatons of structural acoustics problems. Comput Struct 1997;65.3:307–21. [2] Liu GR, Quek SS. The finite element method: A practical course. ButterworthHeinemann; 2003. [3] Zienkiewicz C, Taylor RL. The finite element method: basic concepts and linear applications, 3rd ed. London: McGraw Hill; 1989. [4] Langley RS, Bremner P. A hybrid method for the vibration analysis of complex structural-acoustic systems. J Acoust Soc Am 1999;105:1657. [5] Stefanou George. The stochastic finite element method: past, present and future. Computer Methods in Applied Mechanics and Engineering. 198.9. 2009. p. 1031– 1051. [6] Boyaval Sébastien. A fast Monte–Carlo method with a reduced basis of control variates applied to uncertainty propagation and Bayesian estimation. Comput Methods Appl Mech Eng 2012;241:190–205. [7] Papadrakakis Manolis, Papadopoulos Vissarion. Robust and efficient methods for stochastic finite element analysis using Monte Carlo simulation. Comput Methods Appl Mech Eng 1996;134.3:325–40. [8] Ghanem Roger G, Pol D Spanos. Stochastic finite elements: a spectral approach. Cour Corp 2003. [9] Contreras Humberto. The stochastic finite-element method. Comput Struct 1980;12.3:341–8. [10] Liu Wing, Kam Ted, Belytschko , Mani A. Probabilistic finite elements for nonlinear structural dynamics. Comput Methods Appl Mech Eng 1986;56.1:61–81. [11] Papadimitriou Costas, Lambros S Katafygiotis, James L Beck. Approximate analysis of response variability of uncertain linear systems. Probab Eng Mech 1995;10.4:251–64. [12] Xia Baizhan, Yu Dejie, Liu Jian. Transformed perturbation stochastic finite element method for static response analysis of stochastic structures. Finite Elem Anal Des 2014;79:9–21. [13] Culla Antonio, Carcaterra Antonio. Statistical moments predictions for a moored floating body oscillating in random waves. J Sound Vib 2007;308.1:44–66. [14] Papoulis Athanasios, Pillai S Unnikrishna. Probability, random variables, and stochastic processes. Tata McGraw-Hill Education; 2002. [15] Xia Baizhan, Yu Dejie. Change-of-variable interval stochastic perturbation method for hybrid uncertain structural-acoustic systems with random and interval variables. J Fluids Struct 2014;50:461–78. [16] Harari Isaac, Magoulès Frédéric. Numerical investigations of stabilized finite element computations for acoustics. Wave Motion 2004;39.4:339–49. [17] Braack Malte, et al. Stabilized finite element methods for the generalized Oseen problem. Comput Methods Appl Mech Eng 2007;196.4:853–66. [18] Franca LEOPOLDOP, JAIRO VA Ramalho, Valentin F. Enriched finite element

[25] [26]

[27] [28]

[29]

[30] [31]

[32]

[33] [34]

[35]

[36] [37] [38]

[39] [40]

126

methods for unsteady reaction–diffusion problems. Commun Numer Methods Eng 2006;22.6:619–25. Farhat Charbel, Harari Isaac, Franca Leopoldo P. The discontinuous enrichment method. Comput Methods Appl Mech Eng 2001;190.48:6455–79. Moyroud Francois, Jacquet-Richardet Georges, Fransson Torsten. A comparison of two finite element reduction techniques for mistuned bladed disks. J Eng Gas Turbines Power 2002;124.4:942–52. Bennighof JK. et al. Meeting the NVH computational challenge: automated multilevel substructuring. In: Proceedings of the 18th international modal analysis conference. San Antonio, Texas; 2000. Craig Roy R. Structural dynamics: an introduction to computer methods, Vol. 40. New York: Wiley; 1981. Cohen G. Higher-order numerical methods for transient wave equations[M]. Springer Science & Business Media; 2013. Mencik J-M, Ichchou MN. Wave finite elements in guided elasto dynamics with internal fluid. Int J Solids Struct 2007;44.7:2148–67. Huttunen Tomi, Monk Peter, Kaipio Jari P. Computational aspects of the ultraweak variational formulation. J Comput Phys 2002;182.1:27–46. He ZC, et al. Dispersion error reduction for acoustic problems using the edge-based smoothed finite element method (ES-FEM). Int J Numer Methods Eng 2011;86.11:1322–38. Liu GR. On G space theory. Int J Comput Methods 2009;6.02:257–89. Liu GR. AG space theory and a weakened weak (W2) form for a unified formulation of compatible and incompatible methods: Part I theory. Int J Numer Methods Eng 2010;81.9:1093–126. Liu GR. AG space theory and a weakened weak (W2) form for a unified formulation of compatible and incompatible methods: Part II applications to solid mechanics problems. Int J Numer Methods Eng 2010;81.9:1127–56. Liu GR, et al. A node-based smoothed finite element method (NS-FEM) for upper bound solutions to solid mechanics problems. Comput Struct 2009;87.1:14–26. Zhang Zhi-Qian, Liu GR. Temporal stabilization of the node-based smoothed finite element method and solution bound of linear elastostatics and vibration problems. Comput Mech 2010;46.2:229–46. He ZC, et al. Coupled analysis of 3D structural–acoustic problems using the edgebased smoothed finite element method/finite element method. Finite Elem Anal Des 2010;46.12:1114–21. Cui XY, Li GY. Metal forming analysis using the edge-based smoothed finite element method. Finite Elem Anal Des 2013;63:33–41. Liu GR, Nguyen-Thoi T, Lam KY. An edge-based smoothed finite element method (ES-FEM) for static, free and forced vibration analyses of solids. J Sound Vib 2009;320.4:1100–30. He ZC, et al. An edge-based smoothed finite element method (ES-FEM) for analyzing three-dimensional acoustic problems. Comput Methods Appl Mech Eng 2009;199.1:20–33. Liu GR. Meshfree methods: moving beyond the finite element method[M]. Taylor & Francis; 2009. Chen Jiun-Shyan, et al. A stabilized conforming nodal integration for Galerkin mesh-free methods. Int J Numer Methods Eng 2001;50.2:435–66. Wang Dongdong, Wu Junchao. An efficient nesting sub-domain gradient smoothing integration algorithm with quadratic exactness for Galerkin meshfree methods. Comput Methods Appl Mech Eng 2016;298:485–519. Wu F, et al. A coupled ES-BEM and FM-BEM for structural acoustic problems. Noise Control Eng J 2014;62.4:196–209. Wu F, et al. A new hybrid smoothed FEM for static and free vibration analyses of Reissner–Mindlin Plates. Comput Mech 2014;54.3:865–90.