Plane strain dynamic response of a transversely isotropic multilayered half-plane

Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

Contents lists available at ScienceDirect

Soil Dynamics and Earthquake Engineering journal homepage: www.elsevier.com/locate/soildyn

Plane strain dynamic response of a transversely isotropic multilayered half-plane Zhi Yong Ai n, Yi Fan Zhang Department of Geotechnical Engineering, Key Laboratory of Geotechnical and Underground Engineering of Ministry of Education, College of Civil Engineering, Tongji University, Shanghai,China

art ic l e i nf o

a b s t r a c t

Article history: Received 19 November 2014 Received in revised form 2 April 2015 Accepted 13 April 2015

A semi-analytical method is developed to analyze the plane strain dynamic response of a transversely isotropic multilayered half-plane subjected to a time-harmonic surface or buried load. On the basis of the governing equations of motion in Cartesian coordinates, the analytical layer-elements of a single layer with a finite thickness and a half-plane are obtained through the Fourier transform and the corresponding algebraic operations. The analytical layer-element solution for the multilayered half-plane in the transformed domain can be derived in combination with the continuity conditions between two adjacent layers. After the boundary conditions are introduced, the corresponding solution in the frequency domain is recovered by the inverse Fourier transform. The comparison with an existing solution for an isotropic half-plane confirms the accuracy of the proposed method. Several examples are given to portray the influence of material anisotropy, the depth of external load, material stratification and the frequency of excitation on the vertical displacement and vertical normal stress. & 2015 Elsevier Ltd. All rights reserved.

Keywords: Analytical layer-element Plane strain Dynamic response Transversely isotropy Multilayered half-plane

1. Introduction The research of dynamic response problems is of great importance to the studies related to elastic wave propagation in soils caused by external loadings such as transportation, machine and pile driving. As the pioneer, Lamb [1] first studied the response of an isotropic elastic half-space subjected to a time-harmonic surface load. The follow-up studies for dynamic response problems dealing with an isotropic half-space or full-space can be seen in references (Achenbach [2], Miklowitz [3], Pak [4], Pak and Ji [5], Guzina and Pak [6], etc.). Scholars (Aspel and Luco [7], Pak and Guzina [8], Xu et al. [9], etc.) got analytical solutions for dynamic response of an isotropic multilayered half-space. However, soils in geotechnical engineering are generally transversely isotropic due to long-term sedimentation processes. Many researches, as in Pan and Chou [10,11], Yue et al. [12], Liao and Wang [13], and Wang and Liao [14] provided fundamental solutions for a transversely isotropic half-space subjected to different kinds of static loadings. Apart from transversely isotropy, soils also take on the phenomenon of layering. Solutions for a transversely isotropic multilayered medium under static loads can be found n Correspondence to: Department of Geotechnical Engineering, College of Civil Engineering, Tongji University, 1239 Siping Road, Shanghai 20092, China. Tel.: þ86 21 65982201; fax: þ 86 21 65985210. E-mail address: [email protected] (Z.Y. Ai).

http://dx.doi.org/10.1016/j.soildyn.2015.04.010 0267-7261/& 2015 Elsevier Ltd. All rights reserved.

in references (Small and Booker [15,16], Singh [17], Pan [18,19], Ai et al. [20], etc.). Therefore, it is more realistic to regard soils as a transversely isotropic multilayered medium and to study theirs dynamic response problems. Regarding the dynamic response of a transversely isotropic medium, Stoneley [21] was the earliest researcher focusing on wave propagation in a transversely isotropic medium. Synge [22] and Buchwald [23] studied the propagation of Rayleigh waves in a transversely isotropic medium. Payton [24] presented a time domain solution for displacements and stresses in a transversely isotropic full-space loaded by an instantaneously applied point force. Later, Payton [25] summarized the dynamic problems of a transversely isotropic elastic half-space under surface loads in his book published in 1983. Rajapakse and Wang [26,27] gave out the Green’s functions for a 2-D transversely isotropic half-plane and a non-axisymmetrical transversely isotropic half-space subjected to an interior time-harmonic load. The 3-D time harmonic Green’s function for a transversely isotropic medium was given by Zhu [28] and Yang et al. [29], respectively. Eskandri-Ghadi [30] introduced two potential functions as a general solution for a transversely isotropic medium. With the aid of the potential functions presented by Eskandri-Ghadi [30], Rahimian et al. [31] and Khojasteh et al. [32,33] achieved more subsequent studies in dynamic response problems of a transversely isotropic medium. As for a multilayered system, Khojasteh et al. [34] used the method of displacement potentials to obtain 3-D dynamic Green’s

212

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

functions in a two-layered transversely isotropic half-space. Later, they [35] extended the 3-D dynamic Green’s functions to a transversely isotropic multilayered medium. Ai et al. [36], Ai and Li [37] utilized the analytical layer-element method to derive the solutions for a transversely isotropic multilayered half-space under vertical and horizontal time-harmonic loads, respectively. From above, it can be found that most of the researches concern about a transversely isotropic half-space, while the researches of a multilayered system, especially the ones in Cartesian coordinates, are quite limited. In practical engineering, strip foundations, embankments and dams are generally considered as the plane strain problems, which is more suitable for Cartesian coordinates rather than cylindrical coordinates, therefore the dynamic response of a multilayered transversely isotropic medium in Cartesian coordinates should receive more attention. The purpose of this paper is to extend the analytical layerelement method to solve the plane strain dynamic response of a transversely isotropic multilayered half-plane subjected to a timeharmonic surface or buried load. Compared with the work of Refs. [34,35], the presented method only involves negative exponentials functions in the stiffness matrices of soils, which not only simplifies the calculation processes but also improves the numerical efficiency and stability. With the application of the Fourier transform and the corresponding algebraic operations, the analytical layer-elements which describe the relationship between stresses and displacements of a single layer with a finite thickness and a half-plane are obtained. According to the continuity conditions between adjacent layers, the global stiffness matrix equation is further achieved. After the boundary conditions are introduced, the solution in the frequency domain is achieved by taking the inversion of the Fourier transform. Selected numerical results are performed to demonstrate the accuracy of present method, and to discuss the influence of material anisotropy, material stratification, the depth of load and the frequency of excitation.

2. The analytical layer-elements for a single layer and a halfplane In a Cartesian coordinate system, defined that the z-axis is normal to the plane of isotropy, the governing equations of motion in the absence of body forces for an elastic body can be expressed as follows: ∂σ x ∂τxz ∂ 2 ux þ ¼ρ 2 ∂x ∂z ∂t

ð1aÞ

∂τxz ∂σ z ∂2 uz þ ¼ρ 2 ∂x ∂z ∂t

ð1bÞ

where σ x and σ z represent the normal stress components in the x and z directions, respectively; τxz stands for the shear stress component in the planes xz; ux and uz are the displacement components in the x and z directions, respectively; ρ denotes the density of the material; t is the time variable. The constitutive equations of a transversely isotropic body, which have five independent elastic parameters, can be written in terms of displacements as follows: ∂u ∂u σ x ¼ c11 x þ c13 z ∂x ∂z

ð2aÞ

∂u ∂u σ y ¼ c12 x þc13 z ∂x ∂z

ð2bÞ

∂ux ∂uz þ c33 ∂x ∂z

ð2cÞ

σ z ¼ c13

τxz ¼ c44




∂ux ∂uz þ ∂z ∂x

 ð2dÞ

where c11 ¼ λnð1  nμ2vh Þ, c12 ¼ λnðμh þ nμ2vh Þ, c13 ¼ λnμvh ð1 þ μh Þ, c33 ¼ λð1  μ2h Þ and c44 ¼ Gv are the five independent elastic parameters, in which n ¼ Eh =Ev , λ ¼ Ev =½ð1 þ μh Þð1  μh  2nμ2vh Þ. Here, Ev , Eh and Gv are the vertical Young’s modulus, horizontal Young’s modulus and shear modulus, respectively. In addition, μvh and μh are Poisson’s ratios characterizing horizontal strain due to parallel and normal stresses acting on the plane, respectively. We assume the load is time-harmonic with the circular frequency ω, so the displacement components may express in the form of ux ðx; z; tÞ ¼ ux ðx; zÞeiωt and uz ðx; z; tÞ ¼ uz ðx; zÞeiωt , and the harmonic time factor eiωt is suppressed. Substitution of Eqs. (2) into Eqs. (1) leads to the following equations: c11

∂ 2 ux ∂ 2 ux ∂ 2 uz ∂2 ux ¼ρ 2 þc44 2 þ ðc13 þ c44 Þ 2 ∂x∂z ∂x ∂z ∂t

ð3aÞ

∂ 2 ux ∂ 2 uz ∂ 2 uz ∂ 2 uz þ c33 2 þ c44 2 ¼ ρ 2 ∂x∂z ∂z ∂x ∂t

ð3bÞ

ðc13 þ c44 Þ

The integral transformation approaches are employed to reduce the partial differential equations mentioned above into ordinary differential equations. According to Sneddon [38], a Fourier integral transform is taken. The Fourier transform with respect to the variable x and its inversion are defined as Z þ1 1 ðux ; uz ; σ z ; τxz Þ ¼ ðiux ; uz ; σ z ; iτxz Þeiξx dx ð4aÞ 2π  1 ðux ; uz ; σ z ; τxz Þ ¼

Z

þ1 1

ð  iux ; uz ; σ z ;  iτxz Þeiξx dξ

ð4bÞ

where ξ is the Fourier pffiffiffiffiffiffiffiffi transform parameter with respect to the variable x, and i ¼  1. Eqs. (3) are treated by the Fourier transform Eq. (4a), then we have: ! 2 d du ð5aÞ ρω2  ξ2 c11 þ c44 2 ux  ðc13 þ c44 Þξ z ¼ 0 dz dz ! 2 dux d 2 þ ρω2  ξ c44 þ c33 2 uz ¼ 0 ðc44 þ c13 Þξ dz dz Eqs. (5) may be recast into:   2 ρω2  ξ2 c11 d ux ðc13 þ c44 Þ duz ¼ ξ  ux c44 c44 dz dz2   2 ρω2  ξ2 c44 d uz ðc44 þ c13 Þ dux ¼ ξ  uz c33 c33 dz dz2

ð5bÞ

ð6aÞ

ð6bÞ

In order to simplify the analysis, several variables are defined as follows:   h   0  iT W ξ; z ¼ U ξ; z ; U ξ; z ð7aÞ       T U ξ; z ¼ ux ξ; z ; uz ξ; z

ð7bÞ

    T   dux ξ; z duz ξ; z 0 ; U ξ; z ¼ dz dz

ð7cÞ

With the aid of Eqs. (7), Eqs. (6) take the following form:     dW ξ; z ¼ AðξÞW ξ; z dz

ð8Þ

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

where

"

2

0

6 0 6 6 2 AðξÞ ¼ 6 ξ c11  ρω2 c44 6 4 0

0

1

0

0

0

0

ξ2 c44  ρω2 c33



0

H ð1Þ 21

3

1

7 7

c44

7 5

ξðc13 þ c44 Þ

0

c33

H ð1Þ 22 ¼

det½AðξÞ  λI 44  ¼ 0

ð9Þ

where λ is the eigenvalue of AðξÞ; I 44 is the identity matrix of order 4  4. By solving Eq. (9), matrix AðξÞ has four eigenvalues 7 λ1 and 7 λ2 : qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi a1 þ a21 a2 2 a1  a21  a2 2 λ1 ¼  ; λ2 ¼  ð10Þ 2c33 c44 2c33 c44   where the intermediate variables a1 ¼ c213  c11 c33 þ 2c13c44 2 2 2 2 2 2 ξ þ ðc33 þ c44 Þρω and a2 ¼ 4c33 c44 c11 ξ  ρω c44 ξ  ρω . In this paper, the variables are as assumed to be in complex field. Therefore, there are only a21  a2 ¼ 0ðς ¼ 1Þ and a21  a2 a 0ðς ¼ 2Þ two different situations. When a21  a2 ¼ 0ðς ¼ 1 we have λ ¼ λ1 ¼  λ2 p ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  Þ, ¼ a1 =ð2c33 c44 Þ, so W ξ; z can be expressed as   W ξ; z ¼ ðc1 þ c2 zÞeλz þ ðc3 þ c4 zÞe  λz ð11Þ where c1 ; c2 ; c3 ; c4 are four coefficient vectors. Substituting Eq. (11) into Eq. (8), we can obtain: c1 ¼ m1 l1 þ m2 l2

ð12aÞ

c2 ¼ m1 l2

ð12bÞ

c3 ¼ m3 l3 þ m4 l4

ð12cÞ

c4 ¼ m3 l4

ð12dÞ

" m1 ¼ 1;

λ2 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 ; λ; ξðc13 þ c44 Þλ ξðc13 þ c44 Þ

"

#T ð13aÞ

#T

λ2 c44  ρω2 þ ξ2 c11 2λc44 m2 ¼ 0; ; 1; 2 ξ ð c 13 þc44 Þ ξðc13 þ c44 Þλ

ð13bÞ

"

λ2 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 m3 ¼ 1;  ;  λ; ξðc13 þ c44 Þλ ξðc13 þ c44 Þ "

λ2 c44  ρω2 þ ξ2 c11 2λc44 m4 ¼ 0; ; 1;  ξðc13 þ c44 Þ ξðc13 þ c44 Þλ2

#T ð13cÞ

#T

"

H ð1Þ 11 ¼ " H ð1Þ 12 ¼

ð1Þ H21

and

eλz

Hð1Þ 22 zeλz

d1 eλz

ðd2 þzd1 Þeλz

e  λz

ze  λz

d 3 e  λz

ðd4 þ zd3 Þe  λz

2

ð15aÞ # ð15bÞ

# ð15cÞ #

ð1  λzÞe  λz

ð15dÞ

½d3 ð1  λzÞ  d4 λe  λz

ρω  ξ c11 , d1 ¼ λ c44ξðþ c13 þ c44 Þλ

where

2

2

d2 ¼ λ

2

c44  ρω2 þ ξ c11 , ξðc13 þ c44 Þλ2 2

d3 ¼  d1 ,

d4 ¼ d2 .   When a21 a2 a 0ðς ¼ 2Þ, W ξ; z can be expressed as   W ξ; z ¼ c1 eλ1 z þ c2 e  λ1 z þ c3 eλ2 z þ c4 e  λ2 z

and

ð16Þ

where c1 ; c2 ; c3 ; c4 are also coefficient vectors. Similarly, substituting Eq. (16) into Eq. (8), we can obtain: c1 ¼ m1 l1

ð17aÞ

c2 ¼ m2 l2

ð17bÞ

c3 ¼ m3 l3

ð17cÞ

c4 ¼ m4 l4

ð17dÞ

" m1 ¼ 1;

λ21 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 ; λ1 ; 1 ξðc13 þc44 Þλ1 ξðc13 þc44 Þ

#T ð18aÞ

"

λ2 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 m2 ¼ 1;  1 ;  λ1 ; 1 ξðc13 þ c44 Þλ1 ξðc13 þc44 Þ "

λ2 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 m3 ¼ 1; 2 ; λ2 ; 2 ξðc13 þc44 Þλ2 ξðc13 þc44 Þ "

#T ð18bÞ

#T

λ2 c44 þ ρω2  ξ2 c11 λ2 c44 þ ρω2  ξ2 c11 m4 ¼ 1;  2 ;  λ2 ; 2 ξðc13 þ c44 Þλ2 ξðc13 þc44 Þ

ð18cÞ #T ð18dÞ

Applying Eqs. (17) to Eq. (16), we can also have:     W ξ; z ¼ Hð2Þ ξ; z L ð19Þ 2 3 ð2Þ Hð2Þ H11 12 5. Sub-matrices Hð2Þ , Hð2Þ , Hð2Þ and where Hð2Þ ðξ; zÞ ¼ 4 ð2Þ 11 12 21 H21 Hð2Þ 22 Hð2Þ 22 are of order 2  2. " # e  λ1 z eλ1 z ð2Þ H 11 ¼ d1 eλ1 z d2 e  λ1 z "

ð13dÞ

are of order 2  2. #

 λ e  λz

 λd3 e  λz

Hð2Þ 12 ¼

where l1 ; l2 ; l3 ; l4 are four arbitrary constants. Applying Eqs. (12) to Eq. (11), we have:     W ξ; z ¼ H ð1Þ ξ; z L ð14Þ 2 3 ð1Þ ð1Þ H11 H12

T 5 where L ¼ l1 ; l2 ; l3 ; l4 and Hð1Þ ðξ; zÞ ¼ 4 ð1Þ ð1Þ . Sub-matrices H21 H22 ð1Þ H12 ,

λeλz ðλz þ 1Þeλz ¼ λ z d1 λe ½d1 ðλz þ 1Þ þ d2 λeλz "

ξðc13 þ c44 Þ 7 7:

The eigenvalues of AðξÞ can be obtained by solving the characteristic equation as follow:

ð1Þ H11 ,

213

" Hð2Þ 21 ¼ " Hð2Þ 22 ¼ where

eλ2 z

e  λ2 z

d3 eλ2 z

d4 e  λ2 z

λ1 eλ1 z d1 λ1 eλ1 z λ2 eλ2 z d3 λ2 eλ2 z d1 ¼

ð20aÞ

#

 λ 1 e  λ1 z

ð20bÞ # ð20cÞ

 d 2 λ 1 e  λ1 z  λ 2 e  λ2 z

# ð20dÞ

 d 4 λ 2 e  λ2 z

λ21 c44 þ ρω2  ξ2 c11 ξðc13 þ c44 Þλ1 ,

d2 ¼  d1 ,

d3 ¼

λ22 c44 þ ρω2  ξ2 c11 ξðc13 þ c44 Þλ2

and

d4 ¼  d3 . Expressing the stress components in the form of σ z ðx; z; tÞ ¼ σ z ðx; zÞeiωt , τxz ðx; z; tÞ ¼ τxz ðx; zÞeiωt and applying the Fourier transform to Eq. (2c) and Eq. (2d), the relationship between displacements and stresses are established as 



σ z ξ; z ¼ c13 ξux þ c33

duz dz

ð21aÞ

214

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219 ðςÞ ðςÞ ðςÞ ðςÞ K ðhςÞ ðξ; zÞ ¼ ½  M 1 ðξ; zÞ  M 2 ðξ; zÞ ½ H 11 ðξ; zÞ H 12 ðξ; zÞ  þ 1 , the symbol þ 1 represents the generalized inverse matrix. Here, KðhςÞ ðξ; zÞ is a symmetric matrix of order 2  2, which is the analytical layer-element of a half-plane. The elements of the matrix are listed in Appendix B after sufficient simplifications.

Fig. 1. Stresses and displacements of a single layer with a finite thickness.





τxz ξ; z ¼ c44

3. The analytical layer-element solution for a multilayered half-plane

dux  c44 ξuz dz

ð21bÞ   It can be seen From Eqs. (21) that variables σ z ξ; z and τxz ξ; z     are constituted of variables ux ξ; z , uz ξ; z and their first derivatives, so the matrix form of Eqs. (21) is     V ξ; z ¼ RðξÞW ξ; z ð22Þ

RðξÞ ¼



      T V ξ; z ¼ τxz ξ; z ; σ z ξ; z

where

"



0

ξc13

 ξc44

c44

0

0

0

c33

and

#

.

With the aid of Eqs. (14) and (19), Eq. (22) can be rewritten as  V ξ; z ¼ MðςÞ ðξ; zÞL ð23Þ 

ðςÞ

where M ðςÞ ðξ; zÞ ¼ RðξÞH ðςÞ ðξ; zÞ ¼ ½ M 1

M 2ðςÞ ; sub-matrices Mð1ςÞ

and M2ðςÞ are of order 2  2. Combining Eqs. (14) and (19) with Eq. (23), the relationship between displacements and stresses of a single layer with a finite thickness is established in the following matrix form:  # # " "   V ξ; 0 U ξ; 0 ðςÞ     ¼ K ðξ; zÞ ð24Þ V ξ; z U ξ; z where K

ðςÞ

2

ðξ; zÞ ¼ 4

 M ð1ςÞ ðξ; 0Þ M 1ðςÞ ðξ; zÞ

 M ð2ςÞ ðξ; 0Þ M 2ðςÞ ðξ; zÞ

32 54

ςÞ H ð11 ðξ; 0Þ ðςÞ ðξ; zÞ H 11

ðςÞ H 12 ðξ; 0Þ

ςÞ H ð12 ðξ; zÞ

Under the assumption that surfaces at the top and bottom of each layer are horizontal, an n-layered transversely isotropic elastic system with an underlying half-plane is illustrated in Fig. 2. The thickness of the ith layer is hi ¼ H i H i  1 , where H i and H i  1 are the depths from the surface to the bottom and top of the ith layer, respectively. A strip time-harmonic uniform loading pðx; H i Þeiωt of intensity p, circular frequency ω and width 2a is applied at the depth of H i . Supposing the surface of the multilayered half-plane is free, we get:

σ z ðξ; 0Þ ¼ 0

ð26aÞ

τxz ðξ; 0Þ ¼ 0

ð26bÞ

Considering the continuity conditions between adjacent layers, we have: ux ðξ; H i Þ ¼ ux ðξ; H iþ Þ

ð27aÞ

uz ðξ; H i Þ ¼ uz ðξ; H iþ Þ

ð27bÞ 

σ z ðξ; Hi Þ ¼ σ z ðξ; Hiþ Þ þp ξ; Hi



ð27cÞ

31

τxz ðx; Hi Þ ¼ τxz ðx; Hiþ Þ

5

and ux ðξ are the transformed horizontal where ux ðξ displacements in the x direction at the depth z ¼ H i of the ith layer and the (i þ1)th layer, respectively. The meanings of other parameters can be understood in the same way except that   p ξ; H i is the applied load in the Fourier transformed domain iωt after  the  harmonic time factor  e is suppressed. It turns  out that p ξ; H i ¼ ðp sin ðξaÞ=πξÞ, p ξ; H i ¼ ðP=2πÞ and p ξ; H i ¼ 0, for a strip uniform loading p, a concentrated loading P and no force at the depth z ¼ H i , respectively. Applying Eq. (24) to each finite layer and Eq. (25) to the underlying half-plane, the global stiffness matrix of the multilayered half-plane is assembled in the form of

.

Here, KðςÞ ðξ; zÞ is a symmetric matrix of order 4  4, which is the so-called analytical layer-element for the plane strain dynamic response problem of a single layer. The specific elements of the matrix are listed in Appendix A after sufficient simplifications. The analytical layer-element associates the displacements and stresses of z ¼ 0 and arbitrary depth z in the Fourier transformed domain, as shown in Fig. 1. Besides, the relationship between displacements and stresses of a half-plane can also be established with a consideration of the regularity condition at infinity, which is in the following matrix

; H i Þ

ð27dÞ ; H iþ Þ

(28)

form:

 



 V ξ; z ¼ KðhςÞ ðξ; zÞ Uðξ; zÞ

where

ð25Þ

  where KiðςÞ ¼ K ξ; hi and KhðςÞ represent the analytical layerelement of the ith   layer and the underlying half-plane; and Fðξ; H i Þ ¼ ½0; p ξ; H i T is the load vector in the transformed domain.

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

215

Fig. 5. Comparison of the vertical stress along the z-axis with Ref. [26]. Fig. 2. A transversely isotropic multilayered half-plane.

Fig. 6. Influence of modulus ratio n on the vertical displacement along the z-axis.

Fig. 3. Comparison of the vertical displacement along the z-axis with Ref. [26].

4. Numerical examples and discussion

Fig. 4. Comparison of the vertical displacement along the x-axis with Ref. [26].

By virtue of the boundary conditions Eqs. (26), the unknown displacement vector ½Uðξ; H 0 Þ; Uðξ; H 1 Þ; ⋯; Uðξ; H i Þ; ⋯; Uðξ; H n  1 Þ; Uðξ; H n ÞT in the Fourier transformed domain can be derived by solving the system of linear equations Eq. (28). Then, the stress vector of each layer in the Fourier transformed domain can be obtained by substituting the corresponding displacement vector into Eq. (24). At last, the displacements and stresses in the frequency domain can be further obtained by taking the inversion of the Fourier transform.

Examples are performed to verify the accuracy of the proposed method and to discuss the influence of modulus ratio, the depth of load, stratified character and the frequency of excitation on the vertical displacement and vertical normal stress. It has been mentioned that the stress and displacement with the real and imagine parts in the following examples are solved in the frequency domain, so the physical meanings of the stress and displacement are not same as those in the time domain. When taking the inversion of the Fourier transform, numerical integrations are adopted. Since the results for displacements in the frequency domain are given in terms of infinite integrals with a complex-valued integrand, they cannot be carried out in closed forms due to the presence of singularities within the range of integration. The singularities may affect the accuracy of the integrations vitally, so it is significant to employ a suitable numerical quadrature scheme to evaluate the infinite integrals. Generally, there are two ways to deal with the singularities. One way is to add a small number of damping terms to the parameters of the medium (Aspel and Luco [7], Rajapakse and Wang [26,27], etc.). Since the artificially-added damping terms, such a medium isn’t strictly elasticity, while the influence of the small damping term are proved be negligible (Ai et al. [36]). The other way is to apply the method of residue to numerical integrations (Eskandri-Ghadi et al. [30], Rahimian et al. [31], etc.). In this paper, a small number of damping terms are added to the parameters of the medium to avoid the singularities   ~ ~ in numerical integrations.   We define:  c 11 ¼  c11 1 þ iξs , c 12 ¼ c12 1 þ iξs , c~ 13 ¼ c13 1 þ iξs , c~ 33 ¼ c33 1 þ iξs and c~ 44 ¼ c44 ð1 þ iξs Þ. The value of ξs ¼ 0:01 is the same as that of Rajapakse and

216

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

Wang [26,27]. For convenience of analysis, several variables are defined: the character of transversely isotropy n ¼ h =E v and pEffiffiffiffiffiffiffiffiffiffiffiffi m ¼ Gv =Ev , the non-dimensional frequency ω ¼ aω ρ=c44 , the vertical displacement u0z ¼ Gv uz =ðpaÞ and the vertical stress σ 0z ¼ σ z =p. In practical problems, if the displacement and stress within points in a layer are required, then it is convenient to define a set of fictitious planes through these points and to consider them as additional layers. On the other hand, if the external load is applied within a layer, then a fictitious interface is considered at the loading level. 4.1. Verification The following example is used to verify the feasibility of the proposed theory. The example given by Rajapakse and Wang [26] is about an isotropic half-plane subjected to a strip time-harmonic uniform buried loading. The time-harmonic uniform vertical strip loading of intensity p, circular frequency ω and width 2a is acting at a depth z0 ¼ a below the free surface. The parameters of the medium are μvh ¼ μh ¼ 0:25, Eh ¼ 25 MPa, Ev ¼ 25 MPa, Gv ¼ 10 MPa, and ω ¼ 1:0. The vertical displacement along the z-axis and x-axis, and the vertical stress along the z-axis are plotted in Figs. 3–5, respectively. It shows that the results from this paper are in good correspondence with those of Rajapakse and Wang [26]. Consequently, the accuracy of the method in this study can be confirmed.

applied at the surface. Poisson’s ratios is μh ¼ μvh ¼ 0:25, the density of the medium is ρ ¼ 2:0  103 kg=m3 and the non-dimensional frequency is ω ¼ 1:0. A series of n ¼ 1:0; 2:0; 3:0 and m ¼ 0:3 are selected for comparison. The real and imaginary parts of the vertical displacement along the z-axis and x-axis, and the vertical stress along the z-axis are plotted in Figs. 6–8, respectively. It shows that the vibration amplitudes of the vertical displacement and the vertical stress for both the real and imaginary parts decrease when modulus ratio n increases. What’s more, the vibration amplitudes of all curves also decrease with the increase of depth or horizontal distance. In particular, the vertical displacement and stress tend to zero when the depth or horizontal distance is large enough. However, the influence of modulus ratio n on the vertical displacement is limited.

4.2. The influence of character of transversely isotropy n In order to illustrate the influence of modulus ratio n, a degenerated example of a transversely isotropic half-plane is carried out here. A time-harmonic uniform loading as the previous example is

Fig. 9. Influence of modulus ratio m on the vertical displacement along the z-axis.

Fig. 10. Influence of modulus ratio m on the vertical displacement along the x-axis. Fig. 7. Influence of modulus ratio n on the vertical displacement along the x-axis.

Fig. 8. Influence of modulus ratio n on the vertical stress along the z-axis.

Fig. 11. Influence of modulus ratio m on the vertical stress along the z-axis.

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

4.3. The influence of character of transversely isotropy m Another instance is presented to study the influence of modulus ratio m. Cases of m ¼ 0:2; 0:3; 0:4 and n ¼ 2:0 are discussed below, when the dynamic load and other parameters of the medium still remain consistent with those in the previous instance. Figs. 9–11 depict the real and imaginary parts of the vertical displacement along the z-axis, the vertical displacement along x-axis and the vertical stress along the z-axis, respectively. Figs. 9–11 show that modulus ration m has a remarkable impact on the results. With the increase of modulus ration m, both the real and imaginary parts of the vibration amplitude increase. Besides, the curves of the stress and the displacement become more oscillatory. Both the real and imaginary parts also tend to zero when the depth or horizontal distance is large enough. Therefore, the character of transverse isotropy has a dramatically influence on the dynamic response of the medium.

217

at the depth z0 , which is consistent with the fact that the loading is applied at this level. From Fig. 13, we can find that the real and the imagine parts of the surface displacement deviate considerably when x=a r2, while the two parts are similar when x=a Z 2.

Table 1 The vertical Young’s modulus and shear modulus of three cases.

Case1 Case2 Case3

Gv

Ev1

Ev2

Ev3

2 2 2

6 12 6

3 3 6

2 2 4

4.4. The influence of the depth of load z0 In this example, a transversely isotropic half-plane subjected to a time-harmonic uniform loading is studied to verify the influence of the depth of load z0 . The character of transversely isotropy is n ¼ 2:0 and m ¼ 0:3, while Poisson’s ratios, the density of the medium and the non-dimensional frequency is unchanged. Three cases of different depth z0 =a ¼ 0; 1; 2 are studied below. The calculated results for vertical displacement are plotted in Figs. 12 and 13. As is illustrated in Fig. 12, both the real and the imagine parts of the displacement vibration amplitude decrease when the load moves down. Besides, the real part of the displacement has a kink Fig. 14. Influence of stratification on the vertical displacement along the z-axis.

Fig. 12. Influence of the depth of load z' on the vertical displacement along the z-axis.

Fig. 13. Influence of the depth of load z' on the vertical displacement along the x-axis.

Fig. 15. Influence of stratification on the vertical stress along the z-axis.

Fig. 16. Influence of frequency on the vertical displacement along the z-axis.

218

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

Fig. 17. Influence of frequency on the vertical stress along the z-axis.

4.5. The influence of stratification An example is considered to reveal regularities of a multilayered transversely isotropic half-plane subjected to a time-harmonic uniform surface load. A three-layered medium is chosen here and the ratio of the vertical Young’s modulus and shear modulus of three different cases are listed in Table 1. Other parameters are n1 ¼ n2 ¼ n3 ¼ 2, m1 ¼ m2 ¼ m3 ¼ 0:3, μvh1 ¼ μh1 ¼ μvh2 ¼ μh2 ¼ μvh3 ¼ μh3 ¼ 0:25, ρ1 ¼ ρ2 ¼ ρ3 ¼ 2:0  103 kg=m3 , ω ¼ 1:0 and h1 : h2 : a ¼ 4 : 3 : 1. The real and imagine parts of the vertical displacement and stress along the z-axis are illustrated in Figs. 14 and 15, respectively. It can be seen from Figs. 14 and 15, the curves of case 2 which has a stiffer first layer is rather different from those of case 1. As for case 3 which has a stiffer lower layer, we can observe that case 3 has small difference with case 1. It can be concluded that the displacement is more affected by the characteristic of upper layer.

medium are added to deal with the presence of singularities. The proposed method is verified to be precise by comparing with the existing solution for an isotropic half-plane. The parametric studies indicate that the degree of material transversely isotropy, stratified character, the depth of load and the frequency of excitation significantly influences the dynamic response. As to the transversely isotropic character, the vibration amplitude of the vertical displacement and stress decreases when n increases. Whereas the vibration amplitude increases with the increase of m. As to the depth of load, the displacement amplitude decreases when the load moves down and the real part of the displacement has a kink at the depth of the applied load. As to stratified character, the displacement and stress are more easily influenced by the first layer when the loading is applied on the surface. Furthermore, the decay of the displacement and stress with depth is smoother at low frequencies, but becomes increasingly oscillatory as frequency of excitation increases. The availability of an exact analytical solution and an accurate numerical procedure for evaluation enables the solution to serve as the basis of more complicated problems related to dynamic soils-structure interaction. Appendix A. The elements of KðςÞ (a symmetric matrix of order 4  4) When a21 a2 ¼ 0, the elements of the matrix Kð1Þ are: k11 ¼ 2λA1 c44 ðA1 A5 þ λ A6 c44 Þ=B; 2

k12 ¼ ξððA8  A2 ÞA9  λ A1 ðA2 þ A7 Þðc13  c44 ÞÞc44 =B 2

k13 ¼ 4e  zλ λA1 c44 ðA1 ðA4 A10 þ A11 Þ þ λ ðA10 þ A4 A11 Þc44 Þ=B; 2

  3 2 k14 ¼ 4e  zλ zλ ð1 A4 Þc33 c44 ðA1 þ λ c44 Þ2 = ξðc13 þ c44 ÞB ; k21 ¼ ξððA8  A2 ÞA9  λ A1 ðA2 þ A7 Þðc13  c44 ÞÞc44 =B; 2

4.6. The influence of the frequency of loading

k22 ¼ 2λ c33 c44 ðA1 ð4A4 zλ þ A3 Þ þ λ ð  A3 þ4zλA4 Þc44 Þ=B;

A three-layered transversely isotropic half-plane subjected to a time-harmonic uniform load is selected to investigate the influence of the frequency of excitation. The parameters except for the nondimensional frequency are the same as those of case 1 in the last example. Since most machine foundation vibrations are in the range of ω ¼ 0:5 3:0 (Gazetas [39]), the frequency range ω ¼ 0:5; 1:0; 3:0 is selected. The vertical displacement and stress along the z-axis are plotted in Figs. 16 and 17, respectively. Figs. 16 and 17 show that the influence of the frequency of the time-harmonic loading is clearly evident on these solutions for vertical displacement and stress. At low frequencies, the vertical displacement and stress of both real and imaginary parts decays smoothly with depth, while a gradually decaying oscillatory variation is noted at high frequencies. The curves become increasingly oscillatory and become close to the z-axis constantly with the increase of frequency.

  3 2 k23 ¼ 4e  zλ zλ ð1 A4 Þc33 c44 ðA1 þ λ c44 Þ2 = ξðc13 þ c44 ÞB ;

3

2

k24 ¼ 4e  zλ λ c33 c44 ðA1 ðA11 A4 þA10 Þ þ λ ðA10 A4 þA11 Þc44 Þ=B; 3

2

k31 ¼ k13 ; k32 ¼ k23 ; k33 ¼ k11 ; k34 ¼  k12 ; k41 ¼ k14 ; k42 ¼ k24 ; k43 ¼ k34 ; k44 ¼ k22 ; pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 where λ ¼  a1 =ð2c33 c44 Þ; A1 ¼ ρω2  ξ c11 ; A2 ¼ 1 þ e  4zλ ;  4zλ  2zλ A3 ¼ 1  e ; A4 ¼ e ; A5 ¼  A3 þ 4zλA4 ; A6 ¼ A3 þ 4zλA4 ; 2 2 4 A7 ¼ 2ð  1 þ 2z2 λ ÞA4 ; A8 ¼ 2ð1 þ 2z2 λ ÞA4 ; A9 ¼ A21  λ c13 c44 ; 2 A10 ¼ 1 þ zλ; A11 ¼ 1 þ zλ; B ¼ 2λ A1 ðA2 þA7 Þc44 þ ð  A2 þ 4 A8 ÞðA21 þ λ c244 Þ. When a21 a2 a 0, the elements of the matrix Kð2Þ are: k11 ¼ A1 c44 ðA3 A4 A8 λ1 þA2 A5 A7 λ2 Þðλ1  λ2 Þ=B; 2

2

k12 ¼ ξc44 ððA3 A5  A6 ÞA12 λ1 λ2  A2 A4 ðA11  2A10 λ1 λ2 ÞÞ=B; 2 2

5. Conclusions

k13 ¼ 2A1 c44 ðλ2  λ1 Þðe  zλ2 A2 A7 λ2 þ e  zλ1 A4 λ1 A8 Þ=B; 2

With the use of analytical element-layer method, the solution for a transversely isotropic multilayered half-plane subjected to a time-harmonic load is presented. The present method has the advantage that the stiffness matrices involve only numericallystable negative exponentials. Numerical integrations are adopted to evaluate the infinite integrals involved in the Fourier inverse transform. Negligible damping terms to the parameters of the

2

k14 ¼ 2ξA1 ðe  zλ1 A5  e  zλ2 A3 Þc44 ðc13 þ c44 Þλ1 λ2 ðλ1  λ2 Þ=B; 2

2

k21 ¼ ξc44 ððA3 A5  A6 ÞA12 λ1 λ2  A2 A4 ðA11  2A10 λ1 λ2 ÞÞ=B; 2 2

k22 ¼ c33 c44 λ1 λ2 ðA2 A5 A8 λ1 þ A3 A4 A7 λ2 Þðλ1  λ2 Þ=B; 2

2

Z.Y. Ai, Y.F. Zhang / Soil Dynamics and Earthquake Engineering 75 (2015) 211–219

k23 ¼  2ξA1 ðe  zλ1 A5  e  zλ2 A3 Þc44 ðc13 þc44 Þλ1 λ2 ðλ1  λ2 Þ=B; 2

2

k24 ¼  2c33 c44 λ1 λ2 ðe  zλ2 A2 A8 λ1 þe  zλ1 A4 A7 λ2 Þðλ1  λ2 Þ=B; 2

2

k31 ¼ k13 ; k32 ¼ k23 ; k33 ¼ k11 ; k34 ¼ k12 ; k41 ¼ k14 ; k42 ¼ k24 ; k43 ¼ k34 ; k44 ¼ k22 ; rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 

where λ1 ¼  a1 þ a21  a2 =ð2c33 c44 Þ; rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 

λ2 ¼

 a1 

a21  a2 =ð2c33 c44 Þ;

A1 ¼ ρω2  ξ c11 , A2 ¼ e  2zλ1  1; A3 ¼ e  2zλ1 þ 1; A4 ¼ 1  e  2zλ2 ; 2

A5 ¼ e  2zλ2 þ 1; A6 ¼ 4e  zðλ1 þ λ2 Þ ; A7 ¼ A1 þ c44 λ1 ; A8 ¼ A1 þ c44 λ2 ; 2

2

A9 ¼ ðA7 þc44 λ1 Þλ2 þ λ1 ðA8 þ c44 λ2 Þ;

2 2 2 2 A10 ¼ A1 ðc13  c44 Þ; 2 2 2 2 2 A11 ¼ ð 1 þ 2 ÞðA1  c13 c44 1 2 Þ; 2 2 2 2 A12 ¼  2A21 þ 2c13 c44 1 2 þ A10 ð 1 þ 2 Þ; 2 2 2 2 B ¼  2ðA3 A5 A6 ÞA7 A8 1 2  A2 A4 ðc244 1 2 ð 1 þ 2 Þ þ A1 A9 Þ.

λ

λ

λ λ λ λ λ λ

λ

λ

λ λ λ

λ

Appendix B. The elements of KðhςÞ (a symmetric matrix of order 2  2) When a21  a2 ¼ 0, the elements of the matrix Khð1Þ are: k11 ¼ 2Aλc44 =B; k12 ¼ k21 ¼ ξc44 ðA þ c13 λ Þ=B; 2

k22 ¼  2λ c33 c44 =B; pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 where λ ¼  a1 =ð2c33 c44 Þ; A ¼ ρω2  ξ c11 ; B ¼ A  c44 λ . ð2Þ 2 When a1  a2 a 0, the elements of the matrix Kh are: 3

k11 ¼ Ac44 ðλ1 þ λ2 Þ=B; k12 ¼ k21 ¼ ξc44 ðA þ c13 λ1 λ2 Þ=B; k22 ¼  λ1 λ2 c33 c44 ðλ1 þ λ2 Þ=B; where λ1 ¼ rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 

λ2 ¼

 a1 

a21  a2 =ð2c33 c44 Þ;

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi   a1 þ

a21  a2 =ð2c33 c44 Þ; A ¼ ρω2  ξ c11 , 2

B ¼ A  c44 λ1 λ2 . References [1] Lamb H. On the propagation of tremors over the surface of an elastic solid. Philos Trans R Soc Lond Ser A 1904;203:1–42. [2] Achenbach JD. Wave propagation in elastic solids. Amsterdam: North-Holland Publishing Company; 1973. [3] Miklowitz J. The theory of elastic waves and waveguides. Amsterdam: NorthHolland Publishing Company; 1984. [4] Pak RYS. Asymmetric wave propagation in an elastic half-space by a method of potentials. J Appl Mech 1987;54(1):121–6. [5] Pak RYS, Ji F. Axisymmetric stress-transfers from an embedded elastic cylindrical shell to a half-space. Proc R Soc Lond Ser A 1993;441(1912):237–59. [6] Guzina BB, Pak RYS. Static fundamental solutions for a bi-material full-space. Int J Solids Struct 1999;36(4):493–516. [7] Apsel RJ, Luco JE. On the Green’s functions for a layered half space, Part II. Bull Seismol Soc Am 1983;73(4):931–51. [8] Pak RYS, Guzina BB. Three-dimensional Green’s functions for a multilayered half-space in displacement potentials. J Eng Mech 2002;128(4):449–61.

219

[9] Xu J, Davies TG, Pan E. Efficient and accurate multi-layered elastostatic Green’s functions via the bi-material Green’s function. Eng Anal Bound Elem 2007;31 (8):683–91. [10] Pan YC, Chou TW. Point force solution for an infinite transversely isotropic solid. J Appl Mech 1976;43(4):608–12. [11] Pan YC, Chou TW. Green’s function solutions for semi-infinite transversely isotropic materials. Int J Eng Sci 1979;17(5):545–51. [12] Yue ZQ, Xiao HT, Tham LG, et al. Stresses and displacements of a transversely isotropic elastic half-space due to rectangular loadings. Eng Anal Bound Elem 2005;29(6):647–71. [13] Liao JJ, Wang CD, Jong JH. Elastic solutions for a transversely isotropic halfspace subjected to a point load. Int J Numer Anal Methods Geomech 1998;22 (6):425–47. [14] Wang CD, Liao JJ. Elastic solutions for a transversely isotropic half-space subjected to buried asymmetric-loads. Int J Numer Anal Methods Geomech 1999;23(2):115–39. [15] Small JC, Booker JR. Finite layer analysis of layered elastic materials using a flexibility approach. Part 1—strip loadings. Int J Numer Methods Eng 1984;20 (6):1025–37. [16] Small JC, Booker JR. Finite layer analysis of layered elastic materials using a flexibility approach. Part 2—circular and rectangular loadings. Int J Numer Methods Eng 1986;23(5):959–78. [17] Singh SJ. Static deformation of a transversely isotropic multilayered half-space by surface loads. Phys Earth Planet Inter 1986;42(4):263–73. [18] Pan E. Static response of a transversely isotropic and layered half-space to general surface loads. Phys Earth Planet Inter 1989;54(3):353–63. [19] Pan E. Static Green’s functions in multilayered half spaces. Appl Math Model 1997;21(8):509–21. [20] Ai ZY, Cang NR, Han J. Analytical layer-element solutions for a multi-layered transversely isotropic elastic medium subjected to axisymmetric loading. J Zhejiang Univ SCIENCE A 2012;13(1):9–17. [21] Stoneley R. The seismological implication of aeolotropy in continental structures. Geophys Suppl Mon Not R Astron Soc 1949;5(8):343–53. [22] Synge JL. Elastic waves in anisotropic media. J Math Phys 1957;35(35):323–34. [23] Buchwald VT. Rayleigh waves in transversely isotropic media. Q J Mech Appl Math 1961;14(4):293–317. [24] Payton RG. Green’s tensor for a constrained transversely isotropic elastic solid. Q J Mech Appl Math 1975;28(4):473–81. [25] Payton RG. Elastic wave propagation in transversely isotropic media. The Netherlands: Martinus Nijhoff; 1983. [26] Rajapakse R, Wang Y. Elastodynamic Green’s functions of orthotropic half plane. J Eng Mech 1991;117(3):588–604. [27] Rajapakse R, Wang Y. Green’s functions for transversely isotropic elastic half space. J Eng Mech 1993;119(9):1724–46. [28] Zhu H. A method to evaluate three-dimensional time-harmonic elastodynamic Green’s functions in transversely isotropic media. J Appl Mech 1992;59 (2S):S96–101. [29] Yang B, Pan E, Tewary VK. Three-dimensional Green’s functions of steady-state motion in anisotropic half-space and bimaterials. Eng Anal Bound Elem 2004;28(9):1069–82. [30] Eskandri-Ghadi M. A complete solutions of the wave equations for transversely isotropic media. J Elast 2005;81(1):1–19. [31] Rahimian M, Eskandari-Ghadi M, Pak RYS, et al. Elastodynamic potential method for a transversely isotropic solid. J Eng Mech 2007;133(10):1134–45. [32] Khojasteh A, Rahimain M, Eskandari M, et al. Asymmetric wave propagation in a transversely isotropic half-space in displacement potentials. Int J Eng Sci 2008;46(7):690–710. [33] Khojasteh A, Rahimain M, Pak RYS. Three-dimensional dynamic Green’s functions in transversely isotropic bi-materials. Int J Solids Struct 2008;45 (18  19):4952–72. [34] Khojasteh A, Rahimain M, Pak RYS, et al. Asymmetric dynamic Green’s functions in a two-layered transversely isotropic half-space. J Eng Mech 2008;134(9):777–87. [35] Khojasteh A, Rahimain M, Eskandari M, et al. Three-dimensional dynamic Green’s functions for a multilayered transversely isotropic half-space. Int J Solids Struct 2011;48(9):1349–61. [36] Ai ZY, Li ZX, Cang NR. Analytical layer-element solution to axisymmetric dynamic response of transversely isotropic multilayered half-space. Soil Dyn Earthq Eng 2014;60:22–30. [37] Ai ZY, Li ZX. Time-harmonic response of transversely isotropic multilayered half-space in a cylindrical coordinate system. Soil Dyn Earthq Eng 2014;66:69–77. [38] Sneddon IN. The use of integral transforms. New York: McGraw-Hill; 1972. [39] Gazetas G. Analysis of machine foundation vibrations: state of the art. Int J Soil Dyn Earthq Eng 1983;2(1):2–42.