Two-dimensional numerical approach for the vibration isolation analysis of thin walled wave barriers in poroelastic soils

Two-dimensional numerical approach for the vibration isolation analysis of thin walled wave barriers in poroelastic soils∗ J.D.R. Bordón, J.J. Aznárez, O. Maeso Instituto Universitario de Sistemas Inteligentes y Aplicaciones Numéricas en Ingeniería, Universidad de Las Palmas de Gran Canaria, Edificio Central del Parque Científico y Tecnológico del Campus Universitario de Tafira, 35017 Las Palmas de Gran Canaria, Spain {jdrodriguez,jaznarez,omaeso}@iusiani.ulpgc.es March 9, 2016 Abstract This paper is concerned with the vibration isolation efficiency analysis of total or partially buried thin walled wave barriers in poroelastic soils. A two-dimensional time harmonic model that treats soils and structures in a direct way by combining appropriately the conventional Boundary Element Method (BEM), the Dual BEM (DBEM) and the Finite Element Method (FEM) is developed to this aim. The wave barriers are impinged by Rayleigh waves obtained from Biot’s poroelasticity equations assuming a permeable free-surface. The suitability of the proposed model is justified by comparison with available previous results. The vibration isolation efficiency of three kinds of wave barriers (open trench, simple wall, open trench-wall) in poroelastic soils is studied by varying their geometry, the soil properties and the frequency. It is found that the efficiency of these wave barriers behaves similarly to these in elastic soils, except for high porosities and small dissipation coefficients. The efficiency of open trench-wall barriers can be evaluated neglecting their walls if they are typical sheet piles. This does not happen with walls of bigger cross-sections, leading in general to efficiency losses. Likewise, increasing the burial depth to trench depth ratio has a negative impact on the efficiency. Keywords: vibration isolation, poroelasticity, Rayleigh waves, thin structures, Dual Boundary Element Method, BEM-FEM coupling 1 Introduction The vibrations induced by machinery or vehicles can travel through the soil to nearby constructions, which can annoy people or cause the malfunction of devices located inside of these. In order to reduce the vibrations, a wave barrier can be installed at a point of the transmission path. The design of each vibration isolation system depends on the source of vibrations, the properties of the transmission path, and the isolation requirements. An open trench is a very efficient system because its stress-free boundaries act as perfect reflectors of elastic waves. Its efficiency greatly depends on the ratio between the Rayleigh wavelength and the trench depth. However, for soil stability reasons, especially in water saturated soils, a pure open trench can not be excavated to any desired depth. Thus, other systems such as in-filled trenches, or the installation of sheet piles or rows of piles, are often used. Another option is reinforcing the open trench by installing retaining sheet piles or concrete walls on both sides of the trench. This type of wave barrier is called an open trench-wall. ∗Draft of the paper originally published in Computers and Geotechnics 71 (2016) 168-179 http://dx.doi.org/10.1016/j.compgeo.2015.08.007 1 There exist a vast literature about the design and analysis of wave barriers. Before the numerical computing era, only experimental studies were performed in order to assess these problems, where the works by Barkan [7] and Woods [35] blaze a trail. Nowadays, analytical, semi-analytical and numerical methods, mainly the Boundary Element Method (BEM), are being used, although experimental methods are still being used to confirm and/or parametrize mathematical models, e.g. [23]. Three kinds of wave barriers in elastic soils have been extensively studied: open and in-filled trenches, and rows of piles. The open and in-filled trenches have been studied through two-dimensional BEM models by Emad et al. [20], Beskos et al. [9, 27], and even formulas for a simplified design have been given by Ahmad et al. [2]. They were studied in three-dimensional problems using BEM models by Banerjee et al. [6] and Dasgupta et al. [16]. The vibration isolation produced by rows of piles have been studied by Avilés et al. [4] analytically, and by Kattis et al. [25] using a three-dimensional BEM model. The open trench-wall systems have been rarely studied, to the authors’ knowledge only Tsai et al. [34] using a two-dimensional multidomain BEM model. When compared with elastic soils, much less works dealing with the efficiency of wave barriers in poroelastic soils exist. Cai et al. [13, 12] and Xu et al. [36] studied the isolation efficiency of rows of piles in poroelastic soils using semi-analytical methods, and Cao et al. [14] did the same for open trenches under a moving load. As it is seen, the BEM has been widely applied to study these types of problems because of its own capability to deal with unbounded regions. The Finite Element Method (FEM) has been used also, but mainly in combination with the BEM, being the FEM used for structural parts of the problem. Among other coupled BEM-FEM models used in this field, those developed for the study of the isolation of vibrations produced by moving loads (trains) are of great interest nowadays. To this end, the models developed by Andersen et al. [3] and François et al. [21] are great exponents. The aim of this paper is twofold. Firstly, to present a two-dimensional BEM-FEM dynamic model for soil-structure interaction analyses, where the structures are thin, and one or both of their faces can interact with the surrounding media. The initial idea of this model for fluid-structure analyses has already been presented [11]. Here, the model is expanded by considering a Biot’s poroelastic surrounding media. The second aim is to apply the proposed model to study the efficiency of wave barriers buried in this medium. Three kinds of wave barriers are studied: open trench, simple barrier (thin in-filled trench), and open trench-wall. They are impinged by surface waves generated by a Rayleigh incident wave field assuming a permeable free-surface. The rest of the paper is organized as follows. The Biot’s poroelasticity model is briefly described in Section 2.1. In Section 2.2, the Rayleigh waves on a permeable free-surface are discussed for this model. In section 2.3, the conventional BEM and the Dual BEM for the Biot’s poroelasticity are presented. The soil-structure coupling conditions are described in Section 2.4. In Section 3.1, results obtained from the proposed model are compared with published results. In Section 3.2, a study of the previously mentioned wave barriers under incident Rayleigh waves is presented. 2 Methodology 2.1 Biot’s poroelasticity A very general representation of soils can be done by the Biot’s poroelasticity model [10]. This model is able to represent a two-phase medium consisting of a solid frame saturated by a fluid. Let uiand τi j be the displacements and stresses of the solid phase, Uiand τthe displacements and equivalent stress of the fluid phase, and i,j∈[1,2]. The governing equations in the time domain can be written as: µ∇2u+∇£¡λ+µ+Q2/R¢(∇·u)+Q(∇·U)¤+X=ρ11 ¨ u+ρ12 ¨ U+b¡˙ u−˙ U¢(1) ∇[Q(∇·u)+R(∇·U)]+X′=ρ12 ¨ u+ρ22 ¨ U−b¡˙ u−˙ U¢(2) 2 and the stress-strain relationships as: τi j =δi j £¡λ+Q2/R¢(∇·u)+Q(∇·U)¤+µ¡ui,j+uj,i¢(3) τ=Q(∇·u)+R(∇·U)(4) where Xand X′are the body forces of the solid and fluid phases, respectively, λand µare the Lamé’s parameters of the solid phase, Qand Rare the Biot’s coupling parameters, bis the dissipation constant, and ρ11 =(1−φ)ρs+ρa,ρ12 =−ρa,ρ22 =φρf+ρa, being φthe porosity, ρsthe solid phase density, ρfthe fluid phase density, and ρathe additional aparent density. In the following, in order to avoid confusion, the subscripts 1 and 2 are used to denote solid phase and fluid phase variables, respectively, while x,yare used to denote coordinates. Using the Helmholtz decomposition: ux=∂ϕ1 ∂x+∂ψ1 ∂y,uy=∂ϕ1 ∂y−∂ψ1 ∂x,Ux=∂ϕ2 ∂x+∂ψ2 ∂y,Uy=∂ϕ2 ∂y−∂ψ2 ∂x(5) and considering null body forces, two decoupled sets of two equations are obtained from Eqs. (1-2): ϕ1,ϕ2(¡N=λ+2µ+Q2/R¢∇2ϕ1+Q∇2ϕ2=ρ11 ¨ ϕ1+ρ12 ¨ ϕ2+b¡˙ ϕ1−˙ ϕ2¢(6) Q∇2ϕ1+R∇2ϕ2=ρ12 ¨ ϕ1+ρ22 ¨ ϕ2−b¡˙ ϕ1−˙ ϕ2¢(7) ψ1,ψ2(µ∇2ψ1=ρ11 ¨ ψ1+ρ12 ¨ ψ2+b¡˙ ψ1−˙ ψ2¢(8) 0=ρ12 ¨ ψ1+ρ22 ¨ ψ2−b¡˙ ψ1−˙ ψ2¢(9) The first set is related with a rotational-free (P) displacement field due to scalar potentials ϕ1and ϕ2, and the second set with a divergence-free (S) displacement field due to scalar potentials ψ1and ψ2. In the time harmonic regime, these equations lead to the three well known bulk modes of wave propagation in Biot’s poroelasticity. Onwards, the circular frequency is denoted as ω, and the assumed time harmonic term is exp(iωt), which is omitted for brevity. If only the time harmonic potentials ϕi=Piexp(−ikPx) are considered, then the bulk P mode is obtained from: P1,P2(£ω2ˆ ρ11 −k2 PN¤P1+£ω2ˆ ρ12 −k2 PQ¤P2=0 (10) £ω2ˆ ρ12 −k2 PQ¤P1+£ω2ˆ ρ22 −k2 PR¤P2=0 (11) where ˆ ρ11 =ρ11 −ib/ω,ˆ ρ22 =ρ22 −ib/ωand ˆ ρ12 =ρ12 +ib/ω. The wavenumbers kPare obtained from its characteristic equation: kP=± 1 p2³a1±¡a2 1−4a0¢1/2´1/2 ,a1=ω2µˆ ρ22 R+ˆ ρ11 +ˆ ρ22Q2/R2−ˆ ρ122Q/R λ+2µ¶,a0=ω4ˆ ρ11 ˆ ρ22 −ˆ ρ2 12 R¡λ+2µ¢ (12) where two of the solutions are relevant incoming waves (Re(kP)>0). Hence, two P modes exist: the wavenumber associated with the fastest wave speed is kP1, while the wavenumber associated with the slowest wave speed is kP2. If only the time harmonic potentials ψi=Siexp(−ikSx) are considered, then the bulk S mode is obtained from: S1,S2(£ω2ˆ ρ11 −k2 Sµ¤S1+ω2ˆ ρ12S2=0 (13) ω2ˆ ρ12S1+ω2ˆ ρ22S2=0 (14) and the wavenumber kSis obtained from its characteristic equation: kS=±ωÈ ρ11 −ˆ ρ2 12/ˆ ρ22 µ!1/2 (15) where only one solution is a relevant incoming wave (Re(kS)>0). 3 2.2 Rayleigh waves on a permeable free-surface The Rayleigh waves are surface waves that exist when a half-space is in contact with the vacuum through its free-surface. For a half-space y≤0, three different cases can be considered at the free-surface y=0: permeable (τi j nj=0, τ=0), impermeable (τi j +τδi j )nj=0, (Uj−uj)nj=0) or partially permeable. In this paper, only the permeable case is considered. The potentials for the surface mode are composed by unknown functions Ri=Ri(y) and a wave propagating in the positive xdirection: ϕi=R(P) ie−ikRx,ψi=R(S) ie−ikRx(16) Once substituted into the governing equations, a set of four ordinary differential equations are obtained. It can be converted into a fourth order equation in terms of R(P) 1, and a second order equation in terms of R(S) 1. The solution of the fourth order equation leads to: R(P) 1=P11ekRP1y+P12ekRP2y(17) R(P) 2=P21ekRP1y+P22ekRP2y=D1P11ekRP1y+D2P12ekRP2y Dj=¡λ+2µ¢k2 P j −ω2¡ˆ ρ11 −Q/Rˆ ρ12¢ ω2¡ˆ ρ12 −Q/Rˆ ρ22¢ (18) where the wavenumbers kRPjare obtained from: kRPj =±qk2 R−k2 Pj (19) being physically meaningful only those with Re(kRPj)>0, i.e. those producing evanescent displacements when y→−∞. The solution of the second order equation leads to: R(S) 1=S1ekRSy(20) R(S) 2=S2ekRSy=−ˆ ρ12/ˆ ρ22ekRSy(21) where the wavenumber kSis: kRS =±qk2 R−k2 S(22) being meaningful only that with Re(kRS)>0. Therefore, the potentials are: ϕi=³Pi1ekRP1y+Pi2ekRP2y´e−ikRx,ψi=SiekRSye−ikRx(23) At this point, displacements and stresses can be written as functions of three amplitudes (P11,P12 and S1) and the Rayleigh wavenumber kR. Applying the permeable boundary conditions τxy =0, τyy =0 and τ=0 to the free-surface at y=0, one obtains the following set of three equations: h−2µikRkRP1iP11 +h−2µikRkRP2iP12 +hµ¡2k2 R−k2 S¢iS1=0 h2µk2 R−(N+QD1)k2 P1iP11 +h2µk2 R−(N+QD2)k2 P2iP12 +h2µikRkRSiS1=0 h−(Q+RD1)k2 P1iP11 +h−(Q+RD2)k2 P2iP12 =0 (24) 4 where N=λ+2µ+Q2/R. After some algebraic manipulations using the relationships between wavenumbers given by Eqs. (12), (15), (19) and (22), the characteristic equation associated with this homogeneous set of equations can be written in a similar fashion than that of the elastic case [1, Eq. (5.95)]: ¡2−r2¢2−4p1−r2µH2q1−G1r2−H1q1−G2r2¶=0 (25) where: r=kS kR ,Hj=£µ/(λ+2µ)¤k2 S−k2 Pj k2 P1 −k2 P2 ,Gj= k2 Pj k2 S (26) Eq. (25) is arranged in a new way which is more tractable than others previously obtained, e.g. [17, 37]. In fact, all terms are dimensionless, well behaved, and depend only on the bulk wavenumbers and Lamé’s parameters. It is direct to verify that this equation collapse into the elastic equation; if φ→0, then kP1 →0, kP2 →kelastic P,kS→kelastic Sand (λ+2µ)/µ→(kelastic S/kelastic P)2. 2.3 Conventional and Dual Boundary Element Method The poroelastic soil region is treated numerically using the BEM. Two classes of region boundaries are considered: ordinary and crack-like. A crack-like boundary is an oriented boundary composed by two boundaries sharing the same space but with opposite orientations. It represents the idealization of a null thickness discontinuity within the region: a crack (void) or an inclusion. The BEM relies on the discretization of the Boundary Integral Equations (BIE) used to build a solvable linear system of equations. The conventional BEM uses the Singular BIE (SBIE), and it is able to deal with ordinary boundaries but not with crack-like boundaries. This occurs because identical SBIE are obtained when collocating them on the crack-like boundary, leading to a singular linear system of equations. It can be solved by using the Dual BEM [24, 32] when collocating on a crack-like boundary. The DBEM is a proper combination of the SBIE and the Hypersingular BIE (HBIE), which is commonly used for crack analysis. However, in this work, it is used to couple a crack-like boundary with an elastic inclusion modeled as a structural member. Let Ωbe a Biot poroelastic region, and Γ=∂Ωits boundary with outward unit normal n. Using the weighted residual formulation proposed by Domínguez [18, 5], the SBIE for a collocation point xinot located at a crack-like boundary can be written as: ·Jci 00 0 0ci lk ¸½τi ui k¾+ZΓ"−(U∗ n00 +J X ′∗ jnj)t∗ 0k −U∗ nl0 t∗ lk #½τ uk¾dΓ=ZΓ·−τ∗ 00 u∗ 0k −τ∗ l0 u∗ lk ¸½Un tk¾dΓ Ciui+ZΓ T∗udΓ=ZΓ U∗tdΓ (27) where indicial notation l,k∈[1,2]with summation convention is used, and body forces has been neglected. The secondary variables at the boundary are the fluid normal displacement Un=Uknkand the solid traction tl=τlk nk. The vector ucontains all the primary variables, while tcontains all the secondary variables. The fundamental solution matrix U∗was obtained by Dominguez [18] using the Kupradze procedure [26], and was written in a compact form. However, in order to ease the developments done in the present work, we follow the idea of Maeso et al. [29] of writing the fundamental solution separately in a way that resembles the fundamental solutions of acoustics and elastodynamics: τ∗ 00 =1 2πη,η=1 k2 1−k2 2·α1K0(ik1r)−α2K0(ik2r)¸,αj=k2 j−µ λ+2µk2 3(28) u∗ 0k =− 1 2πΘr,k,Θ=µQ R−Z¶1 λ+2µ 1 k2 1−k2 2·ik1K1(ik1r)−ik2K1(ik2r)¸(29) 5 τ∗ l0 =1 2πJΘr,l(30) u∗ lk =1 2πµµψδlk −χr,lr,k¶ ψ=K0(ik3r)+1 ik3rK1(ik3r)−1 k2 1−k2 2·β1 1 ik1rK1(ik1r)−β2 1 ik2rK1(ik2r)¸ χ=K2(ik3r)−1 k2 1−k2 2·β1K2(ik1r)−β2K2(ik2r)¸,βj=µ λ+2µk2 j−k2 1k2 2 k2 3 (31) where r=|x−xi|is the distance between collocation and observation points, k1=kP1,k2=kP2,k3=kS, J=1/( ˆ ρ22ω2), Z=ˆ ρ12/ˆ ρ22, and Kn(z)is the modified Bessel function of the second kind, order nand argument z. By doing so, the fundamental solution matrix U∗is composed by four submatrices: 00,0k, l0 and lk; where the first index is associated with the load and the second index with the observation, being 0associated with the fluid phase and l,kwith the solid phase. Also, the diagonal submatrices 00 and lk have the same kind of singularities as the well known acoustics and elastodynamic problems, respectively. In fact, the free-terms ci 00 and ci lk are completely similar to that problems. The off-diagonal submatrices 0k and l0 associated with the coupling between phases have one lower order of singularity than the diagonal submatrices. Using Eqs. (2), (3) and (4), the fundamental solution matrix T∗can be obtained from: U∗ n00 +J X ′∗ jnj=−Jτ∗ 00,jnj−Zu∗ 0jnj(32) t∗ 0k =·λu∗ 0m,mδkj +µ³u∗ 0k,j+u∗ 0j,k´¸nj+Q Rτ∗ 00nk(33) U∗ nl0 =−Jτ∗ l0,jnj−Zu∗ ljnj(34) t∗ lk =·λu∗ lm,mδkj +µ³u∗ lk,j+u∗ lj,k´¸nj+Q Rτ∗ l0nk(35) which again can be written in a way that resembles the corresponding fundamental solutions of acoustics and elastodynamics, see A. The HBIE is obtained by establishing the secondary variables at the collocation point, which requires the SBIE and its derivatives with respect to the collocation point (ä,k= ∂ä/∂xi k): Ui n=Ui jni j=−Jτi ,jni j−Zui jni j(36) ti l=τi l j ni j=hλui m,mδlj +µ³ui l,j+ui j,l´ini j+Q Rτini l(37) where niis the unit normal at the collocation point. When considering a collocation point located at a boundary, corner points are excluded in order to avoid multivalued ni, thus Γ(xi)∈C1holds for the HBIE. Therefore, the HBIE can be written as: ci·1 0 0δlk ¸½Ui n ti k¾+ZΓ·−s∗ 00 s∗ 0k −s∗ l0 s∗ lk ¸½τ uk¾dΓ=ZΓ·−d∗ 00 d∗ 0k −d∗ l0 d∗ lk ¸½Un tk¾dΓ Citi+ZΓ S∗udΓ=ZΓ D∗tdΓ (38) where again the diagonal submatrices of D∗and S∗and the free-term resemble those of acoustics and elastodynamic problems. The matrices D∗and S∗are written in A. When xi∈Γ, both the SBIE and the 6 HBIE contain singular integrals, being at most weakly singular those associated with U∗, at most strongly singular (Cauchy Principal Value integrals) those associated with T∗and D∗, and at most hypersingular (Hadamard Finite Part integrals) those associated with S∗. The treatment of those integrals for this problem is analogous to that of acoustics [11] and elastodynamics [33, 15] problems. The treatment of the HBIE is based on a regularization process that requires the integrands (excluding the term r−2) belong to the Hölder function space C1,α[30]. To do so, the collocation point must be in a boundary point where the primary variables are differentiable, i.e. ui∈C1. The SBIE and HBIE shown in Eqs. (27) and (38), respectively, are valid for a collocation point located inside or outside the domain, and at an ordinary boundary, with the condition that Γ(xi)∈C1and ui∈ C1for the HBIE. However, for a collocation point located at a crack-like boundary, the BIEs have to be modified. A crack-like boundary is composed by two sub-boundaries, denoted as positive +and negative −faces. Hence, the integration domain associated with a crack-like boundary can be divided into these two faces, which are geometrically coincident but have opposite orientations. Taking into account this, the BIEs for a collocation point located at a crack-like boundary can be built using a limit to the boundary approach, see e.g. [11]. In this case, the SBIE and HBIE can be written as: 1 2·J0 0δlk ¸³ui++ui−´+ZΓ T∗udΓ=ZΓ U∗tdΓ(39) 1 2·1 0 0δlk ¸³ti+−ti−´+ZΓ S∗udΓ=ZΓ D∗tdΓ(40) where it has been assumed that Γ(xi)∈C1for both equations. Since each face has its own set of variables, (ui+,ti+) for the positive face and (ui−,ti−) for the negative face, and its own boundary conditions, two linearly independent BIEs are needed. By examining Eqs. (39) and (40), it is clear that neither the SBIE nor the HBIE are able to give independently enough conditions. Using both BIEs is the only way to directly get linearly independent equations for a crack-like boundary. Consequently, Eqs. (39) and (40) are called Dual BIEs [24], and their application to the BEM is called the Dual BEM [32]. In this work, the discretization of the BIEs is performed using quadratic elements. For ordinary boundaries, the SBIE is used in the conventional way using nodal collocation for all nodes, except at corner points where double nodes with non-nodal collocation are applied. For crack-like boundaries, the Dual BIEs are used applying the Multiple Collocation Approach [11, 22]. The Rayleigh incident field defined in Section 2.2 is introduced into the model by formulating the BIEs in terms of the scattered field variables [19]. 2.4 Soil-structure coupling The conventional BEM and the Dual BEM used to treat the poroelastic soil region, together with the structural FEM model of thin structures, and the appropriate coupling conditions, make possible building a general BEM-FEM model of soils interacting with plate-like structures. In this two-dimensional approach, thin structures of infinite depth are assumed, hence an equivalent two-dimensional beam model with a modified Young’s modulus Em/(1−ν2) can be used. The beam finite element has been already described [11, 31], but a brief summary is given here for the sake of completeness. Let’s consider an elastic material with Young’s modulus Em, Poisson’s ratio ν, and hysteretic damping coefficient ξ. Then, the effective beam Young’s modulus including the hysteretic damping is E=Em(1+i2ξ)/(1−ν2). The beam is straight, with a cross-section defined by its area Aand its inertia I, and follows the Euler-Bernoulli beam theory with added rotational inertia. It has three nodes and eight degrees of freedom: vertex nodes i=1,2 have translation u(i) 1,u(i) 2and rotation θ(i), while the central node i=3 has only translation u(3) 1,u(3) 2. Axial and lateral quadratic distributed loads over each 7 Ωp Γ+ Ωp Γ− 1 1 1 3 3 3 2 2 2 Ωs τ(2+),t(2+) k s′ 2 (2s) Υ+ j Υ− j Υs j s′ 1 (2s) θ(2s) u(2s) k U(2+) n,u(2+) k τ(2−),t(2−) k U(2−) n,u(2−) k τ(3+),t(3+) k s′ 2 (3s) s′ 1 (3s) u(3s) k U(3+) n,u(3+) k τ(3−),t(3−) k U(3−) n,u(3−) k τ(1+),t(1+) k s′ 2 (1s) s′ 1 (1s) θ(1s) u(1s) k U(1+) n,u(1+) k τ(1−),t(1−) k U(1−) n,u(1−) k n −n x′ 2 x′ 1 Ωp Γ 1 1 3 3 2 2 Ωs τ(2),t(2) k s′ 2 (2s) Υj Υs j s′ 1 (2s) θ(2s) u(2s) k U(2) n,u(2) k τ(3),t(3) k s′ 2 (3s) s′ 1 (3s) u(3s) k U(3) n,u(3) k τ(1),t(1) k s′ 2 (1s) s′ 1 (1s) θ(1s) u(1s) k U(1) n,u(1) k n x′ 2 x′ 1 Figure 1: Types of coupling. Left: poroelastic ordinary boundary element - beam finite element (BEMFEM coupling). Right: poroelastic crack-like boundary element - beam finite element (DBEM-FEM coupling). element are considered, being their nodal values s′ 1 (i)and s′ 2 (i), respectively. By doing so, a node-bynode correspondence between the beam finite element and the quadratic boundary element is achieved. Therefore, the element-wise FEM equation is: ¡K−ω2M¢·u=Q·s′(41) where uis the vector of kinematic variables in global coordinates, s′is the vector of distributed loads in local coordinates, Kand Mare the stiffness and mass matrices in global coordinates, respectively, and Qis the load matrix that converts distributed loads in local coordinates into equivalent nodal forces and moments in global coordinates. The beam finite element can be coupled with an ordinary boundary element (BEM-FEM coupling) or with a crack-like boundary element (DBEM-FEM coupling), see Fig. 1. Let Υs jbe a beam finite element belonging to a elastic region Ωs, and x′ 1,x′ 2its local axis. Let Γbe a boundary (ordinary or crack-like) of a poroelastic region Ωp, being the orientation of Γdefined by its outward unit normal vector n, and Υja boundary element of that boundary. Assuming an impermeable interface with perfect bonding between Γand ∂Ωs, the coupling conditions when Γis an ordinary boundary consist of the following compatibility and equilibrium equations: U(i) n=u(i) knk(42) u(i) l=u(is) l(43) τ(i)nl+t(i) l+s′ 1 (is)x′ 1l+s′ 2 (is)x′ 2l=0 (44) where indicial notation l,k∈[1,2]with summation convention is used, and iis the local index of a node. The coupling conditions when Γis a crack-like boundary consist of the following compatibility and equi8 librium equations: U(i+) n=u(i+) knk(45) U(i−) n=−u(i−) knk(46) u(i+) l=u(is) l(47) u(i−) l=u(is) l(48) τ(i+)nl+t(i+) l−τ(i−)nl+t(i−) l+s′ 1 (is)x′ 1l+s′ 2 (is)x′ 2l=0 (49) If, instead of a poroelastic soil, an elastic soil or an inviscid fluid is considered, these coupling equations lead directly to the appropriate ones by simply removing the equations and variables that do not exist in that kind of medium. 3 Results and discussion 3.1 Comparison with published results To the authors’ knowledge, there are not any published results where directly validate the proposed numerical approach. However, is possible to compare it with published results obtained by other models sharing some aspects. Thereby, in this section it is compared with the classical vibration isolation paper by Beskos et al. [9]. Beskos et al. [9] studied the vibration isolation of open and filled trenches using a two-dimensional conventional BEM elastodynamic model. The problem under consideration is a trench, open or filled with concrete, with a depth to width ratio d/w=10, impinged by waves coming from a footing 5dbehind the trench, vibrating with a frequency corresponding to a Rayleigh wavelength λRequal to the depth d. The elastic soil has a density ρ=1785 kg/m3, shear modulus µ=132 MPa, Poisson’s ratio ν=0.25 and hysteretic damping ξ=0.03. The equivalent poroelastic soil used in our model has a porosity φ=0.001, fluid density ρf=0.001 kg/m3, solid density ρs=1785 kg/m3, null additional aparent density, solid Lamé’s parameters µ=λ=132 MPa, solid phase hysteretic damping ξs=0.03, Biot’s parameters R=Q=0.1 MPa, and null dissipation coefficient. The concrete for the filled trench barrier has a density ρ=2449 kg/m3, shear modulus µ=4.52628 GPa, Poisson’s ratio ν=0.25 and hysteretic damping ξ=0.15. Fig. 2 shows a comparison between their results and our results using the vertical displacement amplitude reduction ratio Ay: Ay(x)=¯¯uy¡x,y=0¢¯¯ ¯¯uno barrier y¡x,y=0¢¯¯ (50) For the open trench, we have used an open trench with d/w=10 using a conventional BEM model, but also an open trench with the null width assumption (d/w→∞) using the Dual BEM. For the filled trench, we have used a filled trench with d/w=10 using a conventional multidomain BEM model, and also a filled trench using our DBEM-FEM model, i.e. from the soil point of view the trench has null thickness but preserves its structural behaviour. There exist differences between the results of Beskos et al. and our results, although the main tendencies are similar. It is probably due to the fact that they used constant boundary elements and an important truncation of the free-surface mesh, which was also noticed by Ahmad et al. [2]. In both problems, the differences between the results using the real geometry (conventional BEM) and the results using the null width assumption (DBEM) are very small. Therefore, it is justified using the proposed DBEM-FEM model for thin structures (d/w≤10) in these kinds of problems. 9 d/h=100, l/d=0.50 ω∗ 1.41.21.00.80.6 d/h=20, l/d=0.50 ω∗ 1.41.21.00.80.6 d/h=10, l/d=0.50 ω∗ ¯ Ay 1.41.21.00.80.6 0.3 0.2 0.1 0.0 d/h=100, l/d=0.25d/h=20, l/d=0.25d/h=10, l/d=0.25 ¯ Ay 0.3 0.2 0.1 0.0 d/h=100, l/d=0d/h=20, l/d=0 h/t=20 h/t=6 h/t=1 Open trench d/h=10, l/d=0 ¯ Ay 0.3 0.2 0.1 0.0 Figure 10: ¯ Aycomparison between open trench and open trench-wall for different d/h,l/dand h/t ratios (φ=0.20, b∗=0.2, νs=0.30, d/w=2) Each column corresponds to a different d/hratio, and each row correspond to a different l/dratio. Four curves are drawn on each graph, one corresponding with the open trench, and three corresponding with h/t={1,6,20}. It is seen that the open trench-wall converges to the open trench as d/hand h/tincrease, as it should be. The cross-sections corresponding with a plate with uniform thickness (h/t=1) have a considerable impact on ¯ Ay, increasing the efficiency for ω∗<1 and l/d=0, but decreasing it in the rest of the cases. The cross-sections associated with the sheet pile idealization have a small influence on the efficiency when compared with the open trench. 4 Conclusions It has been developed a two-dimensional time harmonic model combining the BEM and the FEM for the isolation efficiency analysis of total or partially buried thin walled wave barriers in poroelastic soils. The SBIE and the HBIE needed for the conventional BEM and the Dual BEM for poroelastic regions have been obtained. Also, we describe the coupling conditions between a beam finite element and an ordinary boundary (BEM-FEM), and between a beam finite element and a crack-like boundary (DBEM-FEM), being the latter a new type of coupling element. The considered ground vibrations are Rayleigh waves propagating on a permeable free-surface, which have been obtained, and whose characteristic equation is written is simple new form. The open trench, simple wall and open trench-wall are studied varying their geometry, soil properties and frequency. The soil is assumed to be a sandstone following a linear relationship between porosity and solid dry bulk modulus. In the study, several values of porosity φ, Poisson’s ratio νsand dimensionless dissipation coefficient b∗are considered. From the point of view of isolation efficiency of all wave barri16 ers, it is found that the porosity φis relevant when is near the critical porosity φcr and the dimensionless dissipation coefficient is b∗<5. Also, results do not vary significantly beyond b∗>5, and Poisson’s ratio νsbecomes relevant only for dimensionless frequency ω∗<0.8. Qualitatively, the open trench and the simple wall (thin in-filled trench) behave similarly to those in elastic soils, except for high porosities and small dimensionless dissipation coefficients. For the evaluation of the isolation efficiency of an open trench-wall, it is found that the influence of the walls can be ignored if they are typical sheet piles, and if the dimensionless frequency ω∗lies between 0.5 and 1.5. This is not the case when walls with bigger cross-sections are used, leading in general to an efficiency loss. Wall burial depths l/d>0 lead to efficiency losses, especially for high porosities and low dimensionless frequencies ω∗<1. Acknowledgements This work was supported by the Subdirección General de Proyectos de Investigación of the Ministerio de Economía y Competitividad (MINECO) of Spain and FEDER through research Project BIA2014-57640R and also by the Agencia Canaria de Investigación, Innovación y Sociedad de la Información (ACIISI) of the Government of the Canary Islands and FEDER through research Project ProID20100224. J.D.R. Bordón was recipient of the fellowship TESIS20120051 from the Program of predoctoral fellowships of the ACIISI until September 2014, and currently is recipient of the research fellowship FPU13/01224 from the Ministry of Education, Culture and Sports of Spain. The authors are grateful for this support. A Fundamental solution matrices T∗, D∗and S∗ Elements of the fundamental solution matrix T∗: U∗ n00 +J X ′∗ jnj=1 2πW0 ∂r ∂n(54) W0=µZΘ−J∂η ∂r¶(55) t∗ 0k =1 2π·T01r,k ∂r ∂n+T02nk¸(56) T01 =−2µµ∂Θ ∂r−1 rΘ¶(57) T02 =−λµ∂Θ ∂r+1 rΘ¶−2µ1 rΘ+Q Rη(58) U∗ nl0 =1 2πµ ·W1r,l ∂r ∂n+W2nl¸(59) W1=Zχ−µµ∂Θ ∂r−1 rΘ¶(60) W2=−Zψ−µ1 rΘ(61) t∗ lk =1 2π·T1r,lr,k ∂r ∂n+T2µ∂r ∂nδlk +r,knl¶+T3r,lnk¸(62) 17 T1=−2µ∂χ ∂r−2 rχ¶(63) T2=∂ψ ∂r−1 rχ(64) T3=λ µµ∂ψ ∂r−∂χ ∂r−1 rχ¶−2 rχ+Q R J Θ(65) Elements of the fundamental solution matrix D∗: d∗ 00 =1 2πJW0 ∂r ∂ni(66) d∗ 0k =1 2πµ µ−W1r,k ∂r ∂ni+W2ni k¶(67) d∗ l0 =1 2πJµ−T01r,l ∂r ∂ni+T02ni l¶(68) d∗ lk =1 2π·T1r,lr,k ∂r ∂ni−T2µ−∂r ∂niδlk +r,lni k¶−T3r,kni l¸(69) Elements of the fundamental solution matrix S∗: s∗ 00 =1 2π·Q1 ∂r ∂n ∂r ∂ni+Q2njni j¸(70) Q1=Z2 µχ−2Zµ∂Θ ∂r−1 rΘ¶+Jµ∂2η ∂r2−1 r ∂η ∂r¶(71) Q2=Z2 µψ+2Z1 rΘ−J1 r ∂η ∂r(72) s∗ 0k =1 2π½S01r,k ∂r ∂n ∂r ∂ni+S02nk ∂r ∂ni+S03 ·ni k ∂r ∂n+r,knjni j¸¾ (73) S01 =−2Zµ∂χ ∂r−2 rχ¶−2µ·−∂2Θ ∂r2+3 rµ∂Θ ∂r−1 rΘ¶¸ (74) S02 =Z·λ µµ∂ψ ∂r−∂χ ∂r−1 rχ¶−2 rχ¸+λ·∂2Θ ∂r2+1 rµ∂Θ ∂r−1 rΘ¶¸+2µ1 rµ∂Θ ∂r−1 rΘ¶+Q R J µZΘ−J∂η ∂r¶(75) S03 =−Zµ∂ψ ∂r−1 rχ¶−2µ1 rµ∂Θ ∂r−1 rΘ¶(76) s∗ l0 =1 2π½−S01r,l ∂r ∂n ∂r ∂ni+S02ni l ∂r ∂n−S03 ·−nl ∂r ∂ni+r,lnjni j¸¾ (77) s∗ lk =µ 2π½S1·r,lni k ∂r ∂n−r,knl ∂r ∂ni−∂r ∂n ∂r ∂niδlk +r,kr,lnjni j¸+S2µr,kni l ∂r ∂n−r,lnk ∂r ∂ni¶+S3r,lr,k ∂r ∂n ∂r ∂ni +S4hnjni jδlk +nlni ki+S5nkni lo 18 (78) S1=−∂2ψ ∂r2+1 rµ∂ψ ∂r+3∂χ ∂r−6 rχ¶(79) S2=2λ µ·−∂2ψ ∂r2+∂2χ ∂r2+1 rµ∂ψ ∂r−2 rχ¶¸+4 rµ∂χ ∂r−2 rχ¶−2Q R J µ∂Θ ∂r−1 rΘ¶(80) S3=4·−∂2χ ∂r2+1 rµ5∂χ ∂r−8 rχ¶¸ (81) S4=−2 rµ∂ψ ∂r−1 rχ¶(82) S5=λ2 µ2·−∂2ψ ∂r2+∂2χ ∂r2+1 rµ−∂ψ ∂r+2∂χ ∂r¶¸−λ µ 4 rµ∂ψ ∂r−∂χ ∂r−1 rχ¶+4 r2χ+Q R 1 J·−2λ µµ∂Θ ∂r+1 rΘ¶−4 rΘ+Q Rµη¸ (83) References [1] J. D. Achenbach. Wave propagation in elastic solids, volume 16 of Applied mathematics and mechanics. North-Holland, 1973. [2] S. Ahmad and T. M. Al-Hussaini. Simplified design for vibration screening by open and in-filled trenches. Journal of Geotechnical Engineering, 117:67–88, 1991. [3] L. Andersen and S. R. K. Nielsen. Reduction of ground vibration by means of barriers or soil improvement along a railway track. Soil Dynamics and Earthquake Engineering, 25:701–716, 2005. [4] J. Avilés and F. J. Sánchez-Sesma. Foundation isolation from vibrations using piles as barriers. Journal of Engineering Mechanics, 114:1854–1870, 1988. [5] J. J. Aznárez, O. Maeso, and J. Domínguez. BE analysis of bottom sediments in dynamic fluidstructure interaction problems. Engineering Analysis with Boundary Elements, 30:124–136, 2006. [6] P. K. Banerjee, S. Ahmad, and K. Chen. Advanced application of BEM to wave barriers in multilayered three-dimensional soil media. Earthquake Engineering and Structural Dynamics, 16:1041– 1060, 1988. [7] D. D. Barkan. Dynamics of bases and foundations. McGraw-Hill Book Company, Inc., 1962. [8] J. G. Berryman. Confirmation of Biot’s theory. Applied Physics Letters, 37(4):382–384, 1980. [9] D. E. Beskos, B. Dasgupta, and I. G. Vardoulakis. Vibration isolation using open or filled trenches. Part 1: 2-D homogeneous soil. Computational Mechanics, 1:43–63, 1986. [10] M. A. Biot. Theory of propagation of elastic waves in a fluid-saturated porous solid. I. Low-frequency range. Journal of Acoustical Society of America, 28(2):168–178, 1956. 19 [11] J. D. R. Bordón, J. J. Aznárez, and O. Maeso. A 2D BEM-FEM approach for time harmonic fluidstructure interaction analysis of thin elastic bodies. Engineering Analysis with Boundary Elements, 43:19–29, 2014. [12] Y. Q. Cai, G. Y. Ding, and C. J. Xu. Amplitude reduction of elastic waves by a row of piles in poroelastic soil. Computers and Geotechnics, 36:463–473, 2009. [13] Y. Q. Cai, G. Y. Ding, C. J. Xu, and J. Wang. Vertical amplitude reduction of Rayleigh waves by a row of piles in a poroelastic half-space. International Journal for Numerical and Analytical Methods in Geomechanics, 33:1799–1821, 2009. [14] Z. Cao and Y. Q. Cai. Isolation of train-induced ground-borne vibration by trenches on a poroelastic half-space. Journal of Engineering Mechanics, 139:580–593, 2013. [15] F. Chirino and R. Abascal. Dynamic and static analysis of cracks using the hypersingular formulation of the boundary element method. International Journal for Numerical Methods in Engineering, 43:365–388, 1998. [16] B. Dasgupta, D. E. Beskos, and I. G. Vardoulakis. Vibration isolation using open or filled trenches. Part 2: 3-D homogeneous soil. Computational Mechanics, 6:129–142, 1990. [17] H. Deresiewicz. The effect of boundaries on wave propagation in a liquid-filled porous solid: IV. Surface waves in a half-space. Bulletin of the Seismological Society of America, 52(3):627–638, 1962. [18] J. Domínguez. Boundary element approach for dynamic poroelastic problems. International Journal for Numerical Methods in Engineering, 35:307–324, 1992. [19] J. Domínguez. Boundary Elements in Dynamics. International Series on Computational Engineering. Computational Mechanics Publications, 1993. [20] K. Emad and G. D. Manolis. Shallow trenches and propagation of surface waves. Journal of Engineering Mechanics, 111:279–282, 1985. [21] S. François, M. Schevenels, P. Galvín, G. Lombaert, and G. Degrande. A 2.5D coupled FE-BE methodology for the dynamic interaction between longitudinally invariant structures and a layered halfspace. Computer Methods in Applied Mechanics and Engineering, 199:1536–1548, 2010. [22] R. Gallego and J. Domínguez. Hypersingular BEM for transient elastodynamics. International Journal for Numerical Methods in Engineering, 39:1681–1705, 1996. [23] P. Galvín and J. Domínguez. Experimental and numerical analyses of vibrations induced by highspeed trains on the Córdoba–Málaga line. Soil Dynamics and Earthquake Engineering, 29:641–657, 2009. [24] H.-K Hong and J.T. Chen. Derivations of integral equations of elasticity. Journal of Engineering Mechanics, 114(6):1028–1044, 1988. [25] S. E. Kattis, D. Polyzos, and D. E. Beskos. Vibration isolation by a row of piles using a 3-D frequency domain BEM. International Journal for Numerical Methods in Engineering, 46:713–728, 1999. [26] V. D. Kupradze. Three-dimensional problems of the mathematical theory of elasticity and thermoelasticity, volume 25 of Applied mathematics and mechanics. North-Holland, 1979. 20 [27] K. L. Leung, D. E. Beskos, and I. G. Vardoulakis. Vibration isolation using open of filled trenches: Part 3: 2-D non-homogenous soil. Computational Mechanics, 7:137–148, 1990. [28] C. H. Lin, V. W. Lee, and M. D. Trifunac. The reflection of plane waves in a poroelastic half-space saturated with inviscid fluid. Soil Dynamics and Earthquake Engineering, 25:205–223, 2005. [29] O. Maeso, J. J. Aznárez, and F. García. Dynamic impedances of piles and groups of piles in saturated soils. Computers and Structures, 83:769–782, 2005. [30] P. A. Martin, F. J. Rizzo, and T. A. Cruse. Smoothness-relaxation strategies for singular and hypersingular integral equations. International Journal for Numerical Methods in Engineering, 42:885–906, 1998. [31] L. A. Padrón, J. J. Aznárez, and O. Maeso. BEM-FEM coupling model for the dynamic analysis of piles and pile groups. Engineering Analysis with Boundary Elements, 31:473–484, 2007. [32] A. Portela, M. H. Aliabadi, and D. P. Rooke. The dual boundary element method: effective implementation for crack problems. International Journal for Numerical Methods in Engineering, 33:1269– 1287, 1992. [33] A. Sáez, R. Gallego, and J. Domínguez. Hypersingular quarter-point boundary elements for crack problems. International Journal for Numerical Methods in Engineering, 38:1681–1701, 1995. [34] P. H. Tsai and T. S. Chang. Effects of open trench siding on vibration-screening effectiveness using the two-dimensional boundary element method. Soil Dynamics and Earthquake Engineering, 29:865–873, 2009. [35] R. D. Woods. Screening of elastic surface waves by trenches. PhD thesis, University of Michigan, 1967. [36] B. Xu, J. F. Lu, and J. H. Wang. Numerical analysis of the isolation of the vibration due to Rayleigh waves by using pile rows in the poroelastic medium. Archive of Applied Mechanics, 80:123–142, 2010. [37] J. Yang. A note on Rayleigh wave velocity in saturated soils with compressible constituents. Canadian Geotechnical Journal, 38(6):1360–1365, 2001. 21