A RIM-based Time-domain Boundary Element Method for Three-Dimensional Non-homogeneous Wave Propagations

Liu Liqiand Wang Haitao

1 Introduction

The boundary element method(BEM)offers an efficient way to solve wave propagation problems.When establishing the integral equations,the dynamic fundamental solutions,transformed-domain fundamental solutions or static fundamental solutions can be used.The first approach,generally represented as TD-BEM(TD for time-domain),can preserve the causality condition.Since the radiation condition is satisfied automatically,it is suitable for in finite and semi-in finite problems,see Beskos(1987,1997),Manolis and Beskos(1988).However,the storage of matrices generated at each time-step should be mentioned as one of the main bottlenecks against implementing this method to solve large-scale problems.In order to improve its efficiency,several fast algorithms have been implemented in this field.Of particular interest are the fast multipole method(FMM)[Rokhlin(1985);Greengard and Rokhlin(1997a,1997b)]and the adaptive cross approximation method(ACA)[Bebendorf and Rjasanow(2003);Kolk,Weber and Kuhn(2005);Kolk and Kuhn(2006);Benedetti,Aliabadi and Dav(2008);Benedetti,Milazzo and Aliabadi(2009)].Furthermore,the conventional dynamic fundamental solutions employed in TD-BEM are limited on linear elastic homogeneous problems[Sohrabi-Bidar,Kamalian and Jafari(2010);Wang and Yao(2013);Galvin and Romero(2014a)]or layered homogeneous problems[Birgisson and Crouch(1998)],while specific fundamental solutions for other types of problems,for example,half-space problems[Johnson(1974);Galvin and Romero(2014b)],poroelastic problems[Chen(1994a,1994b);Schanz(2001a,2001b)],non-homogeneous problems[Sanchez-Sesma,Madariaga and Irikura(2001);Luzon,Sanchez-Sesma,Perez-Ruiz,Ramirez-Guman and Pech(2009)],etc.,are still unavailable in closed form.

The second approach can be performed either in the Fourier-domain[Dominguez and Roesset(1978);Mansur,Abreu and Carrer(2004)]or in the Laplace-domain[Cruse and Rizzo(1968)].Tadeu and his co-workers proposed an iterative procedure coupling the normal derivative of the integral equation(TBEM)and the method of fundamental solutions(MFS)to solve the transient acoustic[Tadeu,António and Ferreira(2013)]and elastic[António,Tadeu and Ferreira(2013)]wave propagation problems in the presence of multiple inclusions in frequency domain.In order to obtain the time-domain results,the inverse transformation is required.A new approach combining time and Laplace domains is to use the convolution quadrature method(CQM)developed by Lubich(1988a,1988b),for example,Schanz and Antes(1997),Schanz(2001b),Abreu,Carrer and Mansur(2003),Antes,Schanz and Alvermann(2004)and Abreu,Mansur and Carrer(2006).It utilizes Laplace-domain fundamental solutions and provides solutions directly in time-domain.Li,Zhang,Xie,Zheng and Guo(2014)developed a general method to tackle arbitrary,non-null initial conditions in the analysis of three-dimensional elastodynamic problems using CQM.The main merits of CQM are that only the Laplace-domain fundamental solutions are needed so that it can be applied to problems where time-domain fund a mental solutions are not available[Carrer and Mansur(2010)]and that it shows better stability behavior[Banjai and Schanz(2012)].

The third approach is based on the static fundamental solutions(Kelvin fundamental solutions).In this case,domain integrals occur in the boundary-domain integration equation(BDIE).Generally,domain discretization is required to evaluate these domain integrals(D-BEM)[Carrer and Mansur(2010);Dong,Zhang,Xie,Lu,Han and Wang(2015);Dong,Zhang,Xie,Lu,Li,Han and Li(2015)].Such treatment may offset the main advantage of BEM that only boundary discretization is required.In order to avoid domain discretization,various techniques have been developed to transform these domain integrals into boundary integrals,for example,dual reciprocity method(DRM)[Albuquerque,Sollero and Aliabadi(2002);Vanegas,Patino and Vargas(2014);Useche and Albuquerque(2015);Vanegas and Patino(2015)]introduced in Ref.[Nardini and Brebbia(1983)]and explained in detail by Partridge,Brebbia and Wrobel(1992),and Radial Integration Method(RIM)proposed by Gao(2002a,2002b).In DRM,the integrands are approximated with a series of prescribed basis functions and particular solutions are derived based on these basis functions and the differential operator of the problem.Then the domain integrals are transformed to the boundary employing the particular solutions.The main challenge to implement DRM is that it may be difficult to obtain the particular solutions for some 3D problems.The RIM can transform any domain integrals into a boundary integral and a radial integral which is independent of geometry.It has been employed in nonlinear and nonhomogeneous elastic problems[Gao,Zhang and Guo(2007)],fracture analysis[Gao,Zhang,Sladek and Sladek(2008)]and transient heat conduction problems[Yao,Yu,Gao and Gao(2014)].

After space discretization,a set of second-order ordinary differential equations with respect to time variable are derived.These equations can be solved using time-stepping algorithms such as the Wilson-θmethod,the Newmark method and the Houbolt method.For a long time the Houbolt method was probably the most popular and reliable time integration method used in the D-BEM[Hatzigeorgiou and Beskos(2002);Soares,Carrer and Mansur(2005);Pereira,Karam,Carrer and Mansur(2012)]or DR-BEM formulations[Albuquerque,Sollero and Aliabadi(2002,2004);Useche and Albuquerque(2015)];mainly,due to the fact that unwanted complex higher frequencies that cause numerical instability are damped out.However,it provides numerical results with considerable numerical damping and one has to select a small time step to obtain accurate results.The Newmark method[Newmark(1959)]has been widely used in FEM formulations,presenting over the Houbolt method a better control of the stability and accuracy,according to the values of the parametersβandγ,see Bathe(1996)and Cook(2002).Alternative time-marching schemes are also developed,for example,Carrer and Mansur(2004)and Souza,Carrer and Martins(2004).

To the best knowledge of the authors,few BEM work is found for analyzing tran-sient response of3-D non-homogeneous body in time domain.In this paper,a RIM-based time-domain BEM using static fundamental solutions is developed to treat this problem.The entire domain is assumed to be continuously non-homogeneous and the Kelvin fundamental solutions are employed to form the boundary-domain integration equation.The domain integrals resulting from the inertia term and the non-homogeneity term are transformed to boundary integrals by RIM.Consequently,only the boundary is discretized and8-node discontinuous elements from Mi and Aliabadi(1992)are employed.The Wilson-θmethod is implemented to perform the marching in time and its performance is tested in terms of accuracy and stability.The paper is organized as follows:Firstly,the BDIE for 3D non-homogeneous domain utilizing Kelvin fundamental solutions is derived.Secondly,the domain integrals in the BDIE are converted to boundary integrals with RIM,resulting in a pure boundary integration equation(BIE).Then,a set of second-order ordinary differential equations are obtained by numerical implementation of the aforementioned BIE and solved by the Wilson-θmethod.Finally,numerical examples are given in order to demonstrate accuracy,capacity to treat non-homogeneity and large-scale potential of the proposed method.

2 The boundary-domain integration equation for3-D non-homogeneous problems

We consider a 3D domainVbounded by boundaryS.In the absence of body force,the governing equation for linear elastodynamic problems can be written as

whereujis displacement at pointxand timet;σijis stress andρis density.It is assumed that the shear modulusμ(x)and densityρ(x)vary with Cartesian coordinates while Poisson’s ratioνis constant.In this case,the elasticity tensorCijkl(x)can be written as

with

According to the generalized Hooke’s law,the relationship betweenσijand displacement gradientsuk，lis

A weak form of Eq.(1)can be written as

whereUmjis the weight or test function.Herein,the elastic displacement fundamental solution is chosen as the weight function,which satisfies the following partial differential equation:

where Δ(x，y)is the Dirac delta function;μ0is the shear modulus of a homogeneous domain.Since the value ofμ0will not affect the result,it is set to be 1.Uijis written as

whereris the distance between the source pointyand field pointx.Substitution of Eq.(4)into Eq.(5)and application of Gauss’s divergence theorem yield

wheretj=σijnjis boundary traction,njis the outward unit normal to the boundaryS,and

Substituting Eq.(6)into Eq.(8),we have,

The boundary integration equation can be obtained by taking the source pointyto the boundary as,

wherec(y)=0.5 ifyis on the smooth boundary.

3 Transforming domain integrals to the boundary with RIM

As can be seen from Eq.(12),two domain integrals caused by inertia term and non-homogeneity term respectively are included.A conventional way to treat such domain integrals is to discretize the domain into cells,which to some extent offsets the advantage of BEM.In order to keep only boundary discretization,the Radial Integration Method(RIM)proposed by Gao(2002a,2002b)is adopted herein to transform the domain integrals to equivalent boundary integrals.

The basic idea of RIM is that,the domain integral of a general known functionf(x)withxrepresenting the Cartesian coordinates,

can be transformed into a boundary integral and a radial integral which is independent of the geometry:

wherer(z)is the distance between the origin and boundary pointz.

As for the unknown displacementsujand accelerationjat fi eld pointx,they are approximated by a series of prescribed basis functions in the first place so that RIM can be implemented.It is recommended that a combination of the radial basis functions(RBFs)and the polynomials in terms of global coordinates can give promising results[Golberg,Chen and Bowman(1999)]:

In the above,φA(R)is the RBF;is the distance from theA-th Application PointxAto the field pointx.NAis the total number of Application Points.(A=1，···，NA，i=1，2，3)are undetermined coefficients which are independent of coordinates.In this paper,the compactly supported 4th-order spline RBF is adopted,i.e.,

in whichdAis the radius of the support region at theA-th Application Point.

Substituting Eq.(15)into the domain integrals in Eq.(12)and transforming them into boundary integrals,we have

wherer(y，z)is the distance from source pointyto boundary pointz;and

Eqs.(17)and(18)can be expressed in a simple form as:

where

The following relationships are used in Eq.(17)and Eq.(18):

Substituting Eq.(17)and Eq.(18)into Eq.(12),a pure boundary integral equation can be obtained:

4 Discretization of boundary integration equation

The boundarySis discretized intoNeboundary elements andNbboundary nodes,

Also,NIinternal nodes are casted into the domainV,making the total number of nodes to beN.

In order to replace the undetermined coefficient vectors{α}and{β}in Eq.(35)with the Application Points in Eq.(15)are collocated at all nodes,which meansNA=N.Then the displacement and acceleration ofk-th node can be expressed as

Rewrite Eq.(37)and Eq.(38)in matrix form,

where

whereIis the 3×3 unit matrix.

If no two nodes share the same coordinates,the matrixis invertible and thereby

According to Eq.(40),the matrixcan be expressed in the block form as

Then Eq.(42)can be rewritten as

Substitution of Eq.(44)into Eq.(35)leads to

where

After collocating the source pointyat every boundary and internal node,a set of second-order partial differential equations are formed,which can be expressed in matrix form:

The subscriptsb1andb2mean that the variables are on known displacement boundary and known traction boundary respectively.The subscriptImeans the variables belong to internal nodes.

5 Time-marching scheme

Eq.(47)can be expressed in a more general form:

In order to solve Eq.(48),the Wilson-θmethod is implemented.A brief review of this method is given in the following and readers are referred to Ref.[Bathe(1996)]for detailed descriptions.In Wilson-θmethod,a linear variation of acceleration is assumed during the time interval[t，t+θΔt],whereθ≥1.0.Letτ，0≤τ≤θΔtdenote the increase in time;then from timettot+θΔt,

For unconditional stability,θ≥1.37 is required in FEM and usuallyθ=1.4 is employed.In the following section,the sensitivity to this parameter will be examined.

6 Numerical Examples

A C++code has been developed based on the proposed method.In the code,each entry is calculated and stored as aneight-bytevalue.The program runs on a laptop with Intel Core(TM)Duo T9550 processor and 4GB physical memory.In order to demonstrate accuracy,capacity to treat non-homogeneity and large-scale potential of the proposed method,a number of numerical tests are carried out,including a stability analysis on the time-step size and parameterθin Wilson-θmethod,an accuracy analysis on the number of internal nodes and a comparison of storage between the presented method and conventional TD-BEM.

6.1 Test 1:Stability analysis on the time-step size and parameter θ.

In this example,the classical 1D numerical model is considered.This model consists of a one-dimensional bar of lengthLunder Heaviside-type loadpz=p0H(t-0),wherep0=5.Since the scheme presented herein is derived for 3D case,a cuboid domain in dimension 1×1×2 is analyzed.Fig.1 shows the definition of the model and boundary conditions.

The boundary mesh is illustrated in Fig.2,consisting of 40 boundary elements and 63(3×3×7)internal nodes.Eight-node discontinuous elements from Mi and Aliabadi(1992)as shown in Fig.3 are adopted to treat edge and corner nodes,where the positioning parameter 0＜λ＜1 of collocation nodes stands for the degree of continuity.It is taken to be 2/3 in the analyses.Poisson’s ratio is set to be 0.In order to represent non-homogeneity,the shear modulus and mass density areassumed to vary quadratically and linearly alongzdirection respectively:

Figure 1: One dimensional bar:geometry and loading definition.

In order to measure the time-step size Δt,a dimensionless variable,sayβΔt,is usually adopted[Mansur(1983);Junior(2007);Oyarzun,Loureiro,Carrer and Mansur(2011)],

wherelis the smallest boundary element length andcdis the wave velocity.S-ince the wave velocity is a variable in this case,the maximum ofcdis chosen to determine the time-step.

In order to study the influence of time-step size and parameterθon stability,various combinations listed in Tab.1 are chosen in the analyses.The BEM results are compared with graded finite element method proposed by Santare,Thamburaj and Gazonas(2003)in 1D.

Figs.4-6 show the displacementuzin the middle of the bar and the tractiontzat the fixed end under different parameters.As can be seen from Fig.4(case 1,2 and 3)that results become unstable whenθ=1.3whileθ=1.4and1.8give good stability within the time range considered.This is in accordance with the conclusion in FEM thatθ≥1.37 is recommended.Fig.5 shows that only the smallest time step length

(βΔt=0.2 in case 4)leads to instability whenθ=1.4.Fig.6 shows that all three values ofβΔtgive satisfactory stabilities whenθ=1.8.In general values of 0.3 and 0.5 is recommended forβΔt.

Figure 2: One dimensional bar:boundary discretization and internal nodes.

Figure 3: Eight-node discontinuous element.

Table 1: Parameters for Wilson-θ method.

Figure 4: Displacement component uzat middle point and traction component tzat fi xed end:Wilson-θ analysis with βΔt=0.3.

Figure 5: Displacement component uzat middle point and traction component tzat fi xed end:Wilson-θ analysis with θ=1.4.

Figure 6: Displacement component uzat middle point and traction component tzat fi xed end:Wilson-θ analysis with θ=1.8.

6.2 Test 2:Sensitivity of accuracy to the number of internal nodes.

In this example,the influence of number of internal nodes on accuracy is studied.The model and boundary mesh as given in test 1 are used,while 0,3 and 63 internal nodes are scattered into the domain respectively for comparison.For 0-internalnodes case,no internal nodes are used.For 3-internal-nodes case,three internal nodes are equally spaced along the longitude axis of the rod.For 63-internal-nodes case,the pattern is shown in Fig.2.The parameters for time marching scheme are chosen to beβΔt=0.3，θ=1.4,which have been verified in the first test.

The BEM results of displacement and traction in terms of various numbers of internal nodes are shown in Fig.7 and compared with FEM.It is observed that the number of internal nodes has little impact on the accuracy and that the BEM results with all the three cases are in good agreement with the FEM results.This indicates that only a few or even no internal nodes are required in the proposed method for the 3-D non-homogeneous wave propagations.

Figure 7: Displacement component uzat middle point and traction component tzat fixed end:Wilson-θ analysis.

6.3 Test 3:Storage analysis:a comparison with conventional TD-BEM.

In this example,the storage requirement versus the number of DOFs is studied.The model in Test 1 is meshed with different element sizes in order to generate different scales.The parameters for time marching scheme are chosen to beβΔt=0.3，θ=1.4.Total time length considered is set to be 15.The storage versus the number of DOFs by using the proposed RIM-based BEM is outputted as shown inFig.8 and compared with that from the conventional TD-BEM(based on dynamic fundamental solutions and direct time stepping scheme,see Manolis and Beskos(1988)).As this is a comparison only on storage,for the conventional TD-BEM,homogeneous material is assumed withμ=1000 andρ=1000.

It is shown that the physical memory consumed by conventional TD-BEM grows dramatically with the increase of number of DOFs.When the number of DOF-s reaches 3840,the physical memory is projected to be 10.4 GB,which already exceeds the maximum physical memory of the computer.The physical memory needed for 6000 number of DOFs projects to 36.9 GB.This indicates that the conventional TD-BEM meets serious challenges when treating large-scale problems.The reason behind this is that the storage of conventional TD-BEM depends not only on the number of DOFs,but also on the number of coefficient matrices,which depends strongly on the total time steps required for a wave travelling through the entire domain.

By contrast,the physical memory used by the proposed method(herein presented as RIM-based BEM)just rises slightly with the increase of number of DOFs.Only about 1.5 GB is required for 6000 number of DOFs with 354 time steps.This is because the number of matrices arising from this method is time-independent and can be used for all time steps.This advantage makes this method potential for the wave propagation analysis of large-scale non-homogeneous problems.

Figure 8: Storage requirement:a comparison.

7 Conclusions

In this paper,a RIM-based time domain BEM was proposed to solve transient dynamic problems in 3D non-homogeneous domains.The Kelvin fundamental solution was adopted to form the boundary-domain integration equations.The use of the Radial Integration Method(RIM)kept the advantage of BEM that the discretization and integration were restricted to the boundary.The Wilson-θmethod was adopted as time-marching scheme.The main advantages of the proposed method are that it can treat 3-D non-homogeneous wave propagations with only boundary discretization and that the associated storage is independent of the number of time steps,making it suitable for large-scale problems.Numerical results demonstrated 1)good stability,2)satisfactory accuracy and 3)large-scale potential of the proposed method.

Acknowledgement:Financial supports from the National S&T Major Project(Grant No.2011ZX06002-10)and the National Natural Science Foundation of China(Grant No.11072128)are gratefully acknowledged.

Abreu,A.I.;Carrer,J.A.M.;Mansur,W.J.(2003):Scalar wave propagation in 2D:a BEM formulation based on the operational quadrature method.Engineering Analysis with Boundary Elements,vol.27,pp.101-105.

Abreu,A.I.;Mansur,W.J.;Carrer,J.A.M.(2006):Initial conditions contribution in a BEM formulation based on the convolution quadrature method.International Journal for Numerical Methods in Engineering,vol.67,pp.417-434.

Albuquerque,E.L.;Sollero,P.;Aliabadi,M.H.(2002):The boundary element method applied to time dependent problems in anisotropic materials.International Journal of Solids and Structures,vol.39,pp.1405-1422.

Albuquerque,E.L.;Sollero,P.;Aliabadi,M.H.(2004):Dual boundary element method for anisotropic dynamic fracture mechanics.International Journal for Numerical Methods in Engineering,vol.59,pp.1187-1205.

Antes,H.;Schanz,M.;Alvermann,S.(2004):Dynamic analyses of plane frames by integral equations for bars and Timoshenko beams.Journal of Sound and Vibration,vol.276,pp.807-836.

António,J.;Tadeu,A.;Ferreira,P.(2013):Iterative coupling between the TBEM and the MFS Part II-Elastic wave propagation.CMES:Computer Modeling in Engineering&Sciences,vol.91,no.5,pp.337-354.

Banjai,L.;Schanz,M.(2012):Wave propagation problems treated with convolu-tion quadrature and BEM.Lecture Notes in Applied&Computational Mechanics,pp.145-184.

Bathe,K.(1996):Finite Element Procedures.New Jersey,Prentice Hall Inc.

Bebendorf,M.;Rjasanow,S.(2003):Adaptive low-rank approximation of collocation matrices.Computing,vol.70,pp.1-24.

Benedetti,I.;Aliabadi,M.H.;Dav,G.(2008):A fast 3D dual boundary element method based on hierarchical matrices.International Journal of Solids&Structures,vol.45,pp.2355–2376.

Benedetti,I.;Milazzo,A.;Aliabadi,M.H.(2009):A fast dual boundary element method for 3D anisotropic crack problems.International Journal for Numerical Methods in Engineering,vol.80,pp.1356-1378.

Beskos,D.E.(1987):Boundary Element Methods in Dynamic Analysis.Applied Mechanics Reviews,vol.40,pp.1-23.

Beskos,D.E.(1997):Boundary Element Methods in Dynamic Analysis:Part II(1986-1996).Applied Mechanics Reviews,vol.50,pp.149-197.

Birgisson,B.;Crouch,S.L.(1998):Elastodynamic boundary element method for piecewise homogeneous media.International Journal for Numerical Methods in Engineering,vol.42,pp.1045-1069.

Carrer,J.A.M.;Mansur,W.J.(2004):Alternative time-marching schemes for elastodynamic analysis with the domain boundary element method formulation.Computational Mechanics,vol.34,pp.387-399.

Carrer,J.A.M.;Mansur,W.J.(2010):Scalar wave equation by the boundary element method:A D-BEM approach with constant time-weighting functions.International Journal for Numerical Methods in Engineering,vol.81,pp.1281-1297.

Chen,J.(1994a):Time-domain fundamental solution to biot complete equations of dynamic poroelasticity 1:2-dimensional solution.International Journal of Solids and Structures,vol.31,pp.1447-1490.

Chen,J.(1994b):Time-domain fundamental solution to biot complete equations of dynamic poroelasticity 2:3-dimensional solution.International Journal of Solids and Structures,vol.31,pp.169-202.

Cook,R.(2002):Concepts and Applications of Finite Element Analysis.New York,Wiley.

Cruse,I.T.A.;Rizzo,F.J.(1968):A direct formulation and numerical solution of the general transient elastodynamic problem.Journal of Mathematical Analysis and Applications,vol.22,pp.244-259.

Dominguez,J.;Roesset,J.M.(1978):Dynamic stiffness of rectangular foundations.Nasa Sti/recon Technical Report N,vol.79.

Dong,Y.Q.;Zhang,J.M.;Xie,G.Z.;Lu,C.J.;Han,L.;Wang,P.(2015):A general algorithm for the numerical evaluation of domain integrals in 3D boundary element method for transient heat conduction.Engineering Analysis with Boundary Elements,vol.51,pp.30-36.

Dong,Y.Q.;Zhang,J.M.;Xie,G.Z.;Lu,C.J.;Li,Y.;Han,X.;Li,G.Y.(2015).A General Algorithm for Evaluating Domain Integrals in 2D Boundary Element Method for Transient Heat Conduction.International Journal of Computational Methods,vol.12.

Galvin,P.;Romero,A.(2014a):AMATLAB toolbox for soil-structure interaction analysis with finite and boundary elements.Soil Dynamics and Earthquake Engineering,vol.57,pp.10-14.

Galvin,P.;Romero,A.(2014b):A 3D time domain numerical model based on half-space Green’s function for soil-structure interaction analysis.Computational Mechanics,vol.53,pp.1073-1085.

Gao,X.W.(2002a):The radial integration method for evaluation of domain integrals with boundary-only discretization.Engineering Analysis with Boundary Elements,vol.26,pp.905-916.

Gao,X.W.(2002b):A boundary element method without internal cells for two dimensional and three-dimensional elastoplastic problems.Journal of Applied Mechanics-Transactions of the Asme,vol.69,pp.154-160.

Gao,X.W.;Zhang,C.;Guo,L.(2007):Boundary-only element solutions of 2D and 3D nonlinear and nonhomogeneous elastic problems.Engineering Analysis with Boundary Elements,vol.31,pp.974-982.

Gao,X.W.;Zhang,C.;Sladek,J.;Sladek,V.(2008):Fracture analysis of functionally graded materials by a BEM.Composites Science and Technology,vol.68,pp.1209-1215.

Golberg,M.A.;Chen,C.S.;Bowman,H.(1999):Some recent results and proposals for the use of radial basis functions in the BEM.Engineering Analysis with Boundary Elements,vol.23,pp.285-296.

Greengard,L.;Rokhlin,V.(1997a):A new version of the fast multipole method for the Laplace equation in three dimensions.Acta.Numer.,vol.6,pp.229-269.

Greengard,L.;Rokhlin,V.(1997b):A fast algorithm for particle simulations.J.Comput.Phys.,vol.135,pp.280-292.

Hatzigeorgiou,G.D.;Beskos,D.E.(2002):Dynamic elastoplastic analysis of 3-D structures by the domain/boundary element method.Computers&Structures,vol.80,pp.339-347.

Johnson,L.R.(1974):Green’s Function for Lamb’s Problem.Geophysical Journal of the Royal Astronomical Society,vol.37,pp.99–131.

Junior,D.S.(2007):A time-marching scheme based on implicit Green’s functions for elastodynamic analysis with the domain boundary element method.Computational Mechanics,vol.40,pp.827-835.

Kolk,K.;Weber,W.;Kuhn,G.(2005):Investigation of 3D crack propagation problems via fast BEM formulations.Computational Mechanics,vol.37,pp.32-40.

Kolk,K.;Kuhn,G.(2006):The advanced simulation of fatigue crack growth in complex 3D structures.Archive of Applied Mechanics,vol.76,pp.699-709.

Li,Y.,Zhang,J.,Xie,G.,Zheng,X.,Guo,S.(2014).Time-Domain BEM Analysis for Three-dimensional Elastodynamic Problems with Initial Conditions.CMES:Computer Modeling in Engineering&Sciences,101(3),187-206.

Lubich,C.(1988a):Convolution qadrature and discretized operational calculus.1.Numerische Mathematik,vol.52,pp.129-145.

Lubich,C.(1988b):Convolution quadrature and discretized operational calculus 2.Numerische Mathematik,vol.52,pp.413-425.

Luzon,F.;Sanchez-Sesma,F.J.;Perez-Ruiz,J.A.;Ramirez-Guman,L.;Pech,A.(2009):In-plane seismic response of inhomogeneous alluvial valleys with vertical gradients of velocities and constant Poisson ratio.Soil Dynamics and Earthquake Engineering,vol.29,pp.994-1004.

Manolis,G.D.;Beskos,D.E.(1988):Boundary element methods in elastodynamics.London,Unwin Hyman.

Mansur,W.J.(1983):A time-stepping technique to solve wave propagation problems using the boundary element method.PhD,University of Southampton.

Mansur,W.J.;Abreu,A.I.;Carrer,J.A.M.(2004):Initial conditions contribution in frequency-domain BEM analysis.Cmes-Computer Modeling in Engineering&Sciences,vol.6,pp.31-42.

Mi,Y.;M.Aliabadi,H.(1992):Dual boundary element method for3-Dimensional fracture-mechanics analysis.Engineering Analysis with Boundary Elements,vol.10,pp.161-171.

Nardini,D.;Brebbia,C.(1983):A new approach to free vibration analysis using boundary elements.Applied mathematical modelling,vol.7,pp.157-162.

Newmark,N.(1959):A method of computation for structural dynamics.ASCE J.Eng.Mech.Div,vol.85,pp.67-94.

Oyarzun,P.;Loureiro,F.S.;Carrer,J.A.M.;Mansur,W.J.(2011):A timestepping scheme based on numerical Green’s functions for the domain boundary element method:The ExGA-DBEM Newmark approach.Engineering Analysis with Boundary Elements,vol.35,pp.533-542.

Partridge,P.W.;Brebbia,C.A.;Wrobel,L.C.(1992):The dual reciprocity boundary element method,Computational Mechanics Publications Southampton.

Pereira,W.L.A.;Karam,V.J.;Carrer,J.A.M.;Mansur,W.J.(2012):A dynamic formulation for the analysis of thick elastic plates by the boundary element method.Engineering Analysis with Boundary Elements,vol.36,pp.1138-1150.

Rokhlin,V.(1985):Rapid solution of integral equations of classical potential theory.J.Comput.Phys.,vol.60,pp.187-207.

Sanchez-Sesma,F.J.;Madariaga,R.;Irikura,K.(2001):An approximate elastic two-dimensional Green’s function for a constant-gradient medium.Geophysical Journal International,vol.146,pp.237-248.

Santare,M.H.;Thamburaj,P.;Gazonas,G.A.(2003):The use of graded finite elements in the study of elastic wave propagation in continuously nonhomogeneous materials.International Journal of Solids and Structures,vol.40,pp.5621-5634.

Schanz,M.;Antes,H.(1997):Application of’Operational Quadrature Methods’in time domain Boundary Element Methods.Meccanica,vol.32,pp.179-186.

Schanz,M.(2001a):Wave propagation in viscoelastic and poroelastic continua:a boundary element approach.Berlin,Springer.

Schanz,M.(2001b):Application of 3D time domain boundary element formulation to wave propagation in poroelastic solids.Engineering Analysis with Boundary Elements,vol.25,pp.363-376.

Soares,D.;Carrer,J.A.M.;Mansur,W.J.(2005):Non-linear el as to dynamic analysis by the BEM:an approach based on the iterative coupling of the D-BEM and TD-BEM formulations.Engineering Analysis with Boundary Elements,vol.29,pp.761-774.

Sohrabi-Bidar,A.;Kamalian,M.;Jafari,M.K.(2010):Seismic response of 3-D Gaussian-shaped valleys to vertically propagating incident waves.Geophysical Journal International,vol.183,pp.1429-1442.

Souza,L.A.;Carrer,J.A.M.;Martins,C.J.(2004):A fourth order finite difference method applied to elastodynamics:Finite element and boundary element formulations.Structural Engineering and Mechanics,vol.17,pp.735-749.

Tadeu,A.;António,J.;Ferreira,P.(2013):Iterative coupling between the TBEM and the MFS Part I-Acoustic wave propagation.CMES:Computer Modeling in Engineering&Sciences,vol.91,no.3,pp 153-176.

Useche,J.;Albuquerque,E.L.(2015):Transient dynamic analysis of shear deformable shallow shells using the boundary element method.Engineering Structures,vol.87,pp.1-7.

Vanegas,J.D.;Patino,I.D.;Vargas,C.A.(2014):Boundary element approaches for filling simulations of anisotropic reinforced preforms used in the resin transfer molding process.Journal of Composite Materials,vol.48,pp.2603-2629.

Vanegas,J.D.;Patino,I.D.(2015):Filling simulation of the RTM process in isotropic homogeneous/non-homogeneous media using the boundary element method.Advanced Composite Materials,vol.24,pp.113-139.

Wang,H.;Yao,Z.(2013):ACA-accelerated Time Domain BEM for Dynamic Analysis of HTR-PM Nuclear Island Foundation.Cmes-Computer Modeling in Engineering&Sciences,vol.94,pp.507-527.

Yao,W.A.;Yu,B.;Gao,X.W.;Gao,Q.(2014):A precise integration boundary element method for solving transient heat conduction problems.International Journal of Heat and Mass Transfer,vol.78,pp.883-891.