Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
Contents lists available at ScienceDirect
Soil Dynamics and Earthquake Engineering journal homepage: www.elsevier.com/locate/soildyn
On effective characteristic of Rayleigh surface wave propagation in porous fluid-saturated media at low frequencies Yu Zhang a,b,n, Yixian Xu c,d,e, Jianghai Xia c,e, Shuangxi Zhang a,b, Ping Ping c a
School of Geodesy and Geomatics, Wuhan University, Wuhan 430079, China Key Laboratory of Geospace Environment and Geodesy, Ministry of Education, Wuhan 430079, China c Institute of Geophysics and Geomatics, China University of Geosciences, Wuhan 430074, China d State Key Laboratory of Geological Processes and Mineral Resources, China University of Geosciences, Wuhan 430074, China e Subsurface Imaging and Sensing Laboratory, China University of Geosciences, Wuhan 430074, China b
art ic l e i nf o
a b s t r a c t
Article history: Received 29 April 2013 Received in revised form 8 November 2013 Accepted 18 November 2013 Available online 6 December 2013
The analytical dispersion and waveform solutions of Rayleigh surface wave for the Biot fluid-saturated model are obtained at low frequencies for a homogeneous half space. The equivalent solution is also obtained by the equivalent-viscoelastic representation of the fluid-saturated model based on a single viscoelastic element for each wave modulus. The effective characteristics of the validations and limitations for the equivalent-viscoelastic model are analyzed by comparison of the numerical solutions of the fluid-saturated model and the equivalent model for the surface wave propagations. Our calculations show that the free boundary effects on the frequency dependent dispersion and time dependent dynamical waveforms of the surface wave in the Biot model are well fitted in a relative narrow low frequency band by the Zener elements in case of the frequency is much lower than the critical frequency f c of the porous material. The effective characteristics for air filled cases with a higher f c show a better result. Furthermore, if the critical frequency f c is low, always with high permeability κ under near surface condition, at low frequencies (e.g. the seismic frequency band o200 Hz) the surface fluid drainage conditions influence Rayleigh-wave propagations obviously. The frequency range must hence be carefully checked for the viscoelastic representations. When the validated frequency range is defined, the viscoelastic elements can solve the transient surface wave propagation in porous media effectively. The convolution integral in wave modeling can be replaced by memory variables, which makes the field quantities calculated at every time step need not be stored. The effective representation saves the consumptions of computer time and storages, and supplies a more convenient approach to apply the surface wave considering poroelasticity. & 2013 Elsevier Ltd. All rights reserved.
Keywords: Surface wave Biot's theory Poroelastic Viscoelastic Numerical approximation
1. Introduction Essential is understanding of surface-wave propagation along a free surface for near-surface surface wave applications, e.g., Spectral Analysis of Surface Waves (SASW, [41]) and Multi-channel Analysis of Surface Waves (MASW, [42,58]), etc. The characteristics of surface-wave propagation bring useful information of the subsurface earth media. In general geological conditions, subsurface earth media possess pores, which make it necessary to investigate the surface wave propagation in porous media. The well-known Biot theory [10,11] provided a theoretical framework to investigate the wave propagations in porous media. The global solid and fluid viscous coupling interactions make the energy
n Corresponding author at: School of Geodesy and Geomatics, Wuhan University, Wuhan 430079, China. Tel.: þ86 27 68778371. E-mail addresses:
[email protected],
[email protected] (Y. Zhang).
0267-7261/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.soildyn.2013.11.007
dissipations of wave. Most studies about the surface wave propagation along different interfaces in porous media have been progressing on the basis of Biot's theory, e.g. [27,28,30,1,2,3,38]. For near surface air/solid configurations, Rayleigh waves, formed by interfering of compressional and shear waves on the free surface, are always adopted to express the surface wave characteristics along the earth surface. The dispersion and dynamic solutions of Rayleigh wave in porous media have been being studied in the Biot fluid-saturated model for decades. Jones [33] first discovered that only one type of non-dispersive Rayleigh wave exists in a poroelastic half space if fluid viscosity is ignored. Deresiewicz [22] investigated the velocity dispersion of Rayleigh wave, and found that a local minimum of the velocity is present in the dispersion curves considering pore fluid viscosity. Paul [43,44] firstly analyzed the displacement disturbance under a surface load in a poroelastic half space by using Helmholtz decomposition. Tajuddin [56] examined the effects of different fluid drainage on the free surface for the non-dispersive surface-wave propagation
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
in the same model of Jones [33]. [31] obtained the harmonic response in a poroelastic half space under a steady-state surface traction. Philippacopoulos [45] studied classical Lamb problem in poroelastic media. He used the displacement potentials to decouple the Biot's equations and obtained the dynamic Green's functions in a half poroelastic space. Senjuntichai and Rajapakse [55] studied the dynamic response in a two-dimensional half plane. Two types of dispersive surface waves were verified in a widefrequency range at an impermeable half space in a simple mixture poroelastic model [25,26,4]. Zhang et al. [62,63] studied the sensitivity of the different physical properties on Rayleigh surface-wave propagations under different free surface fluid drainage conditions. Some researches concerned modeling of an unsaturated or a doublepore medium based on Biot's theory (e.g., [7,53,54,8,9,49,50,36]). In these models, surface-wave propagations were also studied (e.g., [21,20,37]). It is well known that the presence of the Biot slow compressional wave makes Biot's differential equations stiff at low frequencies. Therefore, the numerical approaches require a very large amount of storage and computational time [15,16], which always face unstable problems on wave simulations. In order to circumvent this difficulty, to seek an effective solution for wave propagation in porous media is a necessary approach to compress these large consumptions and unstable problems. Global solid and fluid interaction effects are modeled by a viscoelastic equivalent element for the body waves in porous media, as first proposed by [29]. They showed the attenuation quality factor (Q) of the fast compressional wave can be analytically approximated by a single linear standard or a Zener viscoelastic relaxation mechanism [5], if the quality factor satisfies Q4 5. Carcione [17] successfully applied the Zener equivalent representations to replace Biot's and squirt flow energy dissipation mechanism [24], and modeling the elastic wave fields containing fast compressional wave and shear wave in a equivalent-viscoelastic model effectively, which allows the convolution integral to be replaced by memory equations. Since the Biot attenuation mechanism is not related to bulk deformation, the Zener model, which generalizes compressibility to a relaxation function, is directly match the wave dispersion and attenuation by using relaxation functions associated with each wave mode (see [17]). This approximate method is also applied to describe the wave scattering problems [40]. Recent studies also show that the Zener viscoelastic model can effectively represent the attenuation of body waves in unsaturated media (e.g., [52,46,47,19]), doublepore media (e.g., [34,35]), and multiphase solid media (e.g., [48]). These procedures all show that the convolution operators in complex media can be replaced by the memory variables in equivalent representations. Although the attenuation in porous media is not caused by the bulk deformation, the direct fits on velocity and attenuation dispersion by using relaxation functions associated with each wave mode also shows that a single-phase equivalent viscoelastic representation can actually represent the macroscopic attenuation of wave propagation in porous medium reliable with appropriate viscoelastic relaxation mechanism (e.g., Zener relaxation function) for each wave modulus (also see [18]). Because of the free surface, the numerical approaches for Rayleigh surface-wave propagations faced more harsh unstable problems, even in elastic media. We must solve the differential equation in very small grid [59,60]. The appearance the pore fluid phase makes the free surface effects on the surface-wave propagation more complicate [62,63], and undoubtedly aggravates this difficulty. These difficulties prevent us to introduce poroelasticity to the surface wave application in an economic way. Thus, in parallel with the above issue for seeking the effective body wave representations in porous media, after we explicitly analyzed the characteristics of Rayleigh surface waves in porous media [62,63], the questions as whether Rayleigh surface wave can be effectively
95
represented by an equivalent viscoelastic model, and how to solve the effective surface wave propagation with the free surface still need to be suggested. To the best of our knowledge, there is no paper providing explicit analysis about the effective viscoelastic characteristics of Rayleigh surface waves in fluid-saturated porous media. Especially, how to solve Rayleigh wave propagation effectively represented in viscoelastic mechanisms for fluid-saturated media, which is the main concern of our present study, is still lack at present. Therefore, we apply an equivalent-viscoelastic model to investigate the effective nature of Rayleigh surface-wave propagations in a fluid-saturated model and propose several numerical examples for the simplified approximate method to solve the dynamic problems in poroelastic media and assess its validity and limitations. In this paper, we first derive the phase velocity and the quality factor expressions in a fluid-saturated Biot's model in Sec II, and develop the analytical dynamic transient waveform solutions of Rayleigh waves subjected to a vertical stress source in Sec III, respectively. Second, we obtain the equivalent-viscoelastic model for Rayleigh waves by replacing the matrix complex modulus by a Zener relaxation function for each wave modulus, respectively. Meanwhile, we derive the analytical transient solutions in this effective model in Sec IV. The effective characteristics of Rayleigh wave for several typical rock samples, i.e., velocity and attenuation dispersion, dynamical waveform as well, are compared and analyzed between the poroelastic and viscoelastic models in Sec V. Finally, discussions and conclusions are briefly summarized.
2. A porous fluid-saturated model A framework of solid grains and a connected pore space filled entirely by fluid comprise the fluid-saturated porous medium. The volume ratio of the pore space is defined as porosity, ϕ. Two displacement vectors, u and U, describe the motions of the multiphase of the solid and the fluid, respectively. Usually, w ¼ ϕðU uÞ is introduced to characterize the relative fluid-solid displacement. The stress tensor for porous media is a volume average effect of the solid frame and pore fluid. For a general isotropic case, the relations of stress and strain can be drawn from Hooke's law [13]. τij ¼ λδij eþ 2μeij αδij p;
ð1Þ
p ¼ αe þ Mς;
ð2Þ
where τij is the stress tensor of the porous medium, and p is the pore fluid pressure. eij is the solid strain tensor, and e ¼ ∇ u, ς ¼ ∇ w. λ ¼ K m ð2=3Þμ, K m is the bulk modulus of the drained porous media, and μ is the shear modulus of the elastic frame. The coefficients α and M can be written by [12] α ¼ 1
Km 1 αϕ ϕ ; and ¼ þ M Ks Kf Ks
ð3Þ
where K s and K f are the solid grain bulk modulus and the fluid modulus, respectively. By combining the momentum conservation and the stress– strain relation, the equations governing the motion of porous media are written by μ∇2 u þ ðλC þ μÞ∇ð∇ uÞ þ αM∇ð∇ wÞ ¼ ∂tt ðρu þ ρf wÞ; and αM∇ð∇ uÞ þ M∇ð∇ wÞ ¼ ∂tt ðρf u þmwÞ þ b∂t w;
ð4Þ
where ρ ¼ ð1 ϕÞρs þ ϕρf is the mass density of the fluid-saturated porous media, and ρs and ρf are the densities of the solid grains and the pore fluid. λC ¼ λ þ α2 M, and m ¼ Cρf =ϕ, where C denotes the tortuosity of the pore space. It has a theoretical estimated value of C ¼ ð1=2Þð1 þ 1=ϕÞ for spherical solid grain constituent [6].
96
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
The viscous drag-force coupling is represented by the damping coefficient b0 ¼ ðη=κÞ, and b ¼ b0 F r , where η is the fluid dynamic viscosity, and κ denotes the intrinsic dynamic permeability of the porous medium. Johnson et al. [32] introduced the viscous correction factor F r in the frequency domain can be written by F r ðf Þ ¼ ð1 þ ði=2Þ M s f =f C Þ1=2 , where M s is the shape factor (always equal to 1), and the reference critical frequency is defined as f c ¼ ðηϕ=2πρf κÞ. To simplify these equations, let us express the displacements in terms of scalar and vector potentials [13] u ¼ ∇ φ1 þ ∇ ψ 1 ; and w ¼ ∇ φ2 þ ∇ ψ 2 ;
ð5Þ
where φ1 and φ2 are the scalar potential, and ψ 1 and ψ 2 are the vector potential of the solid- and fluid-displacements, respectively. By applying Fourier transform for the time derivation of Eq. (4), and substituting Eq. (5), we obtain ε1 ∇2 φ^ 1 þ ε2 φ^ 1 φ^ 2 ¼ 0; ∇4 φ^ 1 þ ε3 ∇2 φ^ 1 þ ε4 φ^ 1 ¼ 0 ϑ1 ψ^ 1 ψ^ 2 ¼ 0; and 2
ð6Þ
where the mark ^ represents the parameters in the frequency domain. With the following notations, λ þ 2μ αðmω2 ibωÞ ρf ω2 ðρ αρf Þω2 αðmω2 ibωÞ ρf ω2
ε3 ¼
αMðρ αρf Þ þ ρf ðλ þ 2μÞ þ ðλC þ 2μÞ½αðm ib=ωÞ ρf 2 ω αMðλ þ 2μÞ
ε4 ¼
ρf ðρ αρf Þ þ ρ½αðm ib=ωÞ ρf 4 ω ; and αMðλ þ 2μÞ ρf ω2 ; mω2 ibω
where ω is the angular frequency, and ω ¼ 2πf . The complex wavenumbers for the body waves in the fluid-saturated medium can be written from the Helmholtz decomposition of Eq. (6) by kpi ¼
ks ¼ ω
ε3 7ðε23 4ε4 Þ1=2 2
;
i ¼ 1; 2;
ð7Þ
ð8Þ
where the frequency-dependent kpi , i¼1, 2 corresponds two types of compressional waves (P1 and P2 waves), and ks corresponds the shear wave in the porous medium, respectively. Therefore, there exists a slow compressional wave, named Biot's slow wave or P2 wave, which is highly dissipative and dispersive. This P2 wave was predicted by Biot's poroelastic theory [10,11], and then confirmed in laboratory experiments [51]. In the general geological conditions, the reference critical frequency f C is very high. Therefore, the low frequency assumption does not affect the most practical applications in seismology. The phase velocity and the quality factor can be drawn from the relations ki ðωÞ ¼ ω=vi ðωÞ, i ¼ p1; p2; s. We use the concept of complex velocity to obtain the phase velocity and the loss angle as a function of angular frequency. They are given by i ¼ p1; p2; s;
ð9Þ
ð10Þ
The relation between the loss angle and the standard definition of quality factor for an elasticity is Q 1 ¼ tan θ. 3. Rayleigh surface waves in a porous fluid-saturated medium To investigate the propagation of Rayleigh surface waves, it is necessary to solve the boundary value problem related to free surface conditions. We consider a three-dimensional fluid-saturated half-space in a cylindrical coordinate system ðr; θ; zÞ, and z 4 0 represents the medium. With a vertical external force on the solid frame at (0, 0, 0), the free surface boundary conditions at the interface z ¼0 for the solid frame can be written by FðωÞδðrÞ ; and 2πr
τ^ rz ðr; z; ωÞjz ¼ 0 ¼ 0
ð11Þ
where FðωÞ means the Fourier transform of the source time function FðtÞ, δ denotes Dirac delta function. The boundary for the pore fluid on the free surface can be written by [23,14,62,63] ð12Þ
where χ denotes the surface stiffness, which represents the free surface drainage conditions for pore fluids. When χ-1, the pores on the free surface are drained, referred as an “open-pore” condition; when χ ¼ 0, the pores on the free surface are undrained, referred as a “closed-pore” condition. Besides, Eq. (12) represents an intermediate state. Although, under some condition, different drainage conditions can affect the surface-wave propagations [62,63], in the present discussion, the “open-pore” condition, as ^ z; ωÞjz ¼ 0 ¼ 0, is assumed for typical real world conditions. pðr; Following Zhang et al. [63], the general transient solutions for Eq. (6) in cylinder coordinate are shown to be Z 1 u^ z ¼ ð A0 γ 1 e γ 1 z B0 γ 2 e γ 2 z þ C 0 ξe ςz ÞξJ 0 ðξxÞdξ; ð13Þ u^ r ¼
!1=2
ρ þ ϑ1 ρf 1=2 ; μ
1 k ðωÞ vi ¼ ℜ i ; ω
i ¼ p1; p2; s;
^ z; ωÞ χ pðr; ^ z; ωÞjz ¼ 0 ¼ 0; ½iω wðr;
ε2 ¼
ϑ1 ¼
ℑðk2 ðωÞÞ 1 1 i θ ¼ tan 1 ; ¼ tan ℜðk2 ðωÞÞ Q i
τ^ zz ðr; z; ωÞjz ¼ 0 ¼
∇2 ψ^ 1 þ ks ψ^ 1 ¼ 0;
ε1 ¼
and
^z¼ w ^r¼ w
Z
0
1 0
Z
ð A0 ξe γ 1 z B0 ξe γ 2 z þC 0 ςe ςz ÞξJ 1 ðξxÞdξ;
1 0
Z
1 0
ð14Þ
ð A0 δ1 γ 1 e γ 1 z B0 δ2 γ 2 e γ 2 z þ C 0 δ3 ξe ςz ÞξJ 0 ðξxÞdξ;
ð15Þ
ð A0 δ1 ξe γ 1 z B0 δ2 ξe γ 2 z þ C 0 δ3 ςe ςz ÞξJ 1 ðξxÞdξ;
ð16Þ
where ξ denotes the horizontal wavenumber, which can describe 2 the propagation of Rayleigh surface waves. γ 2i ¼ ξ2 kpi ; i ¼ 1; 2 2 and ς2 ¼ ξ2 ks . J 0 and J 1 are the first kind Bessel functions of the zero and the first orders. The amplitude ratios for the fluid phase in pores can be drawn from the eigenvalue problems of Eq. (6) as ðλ þ 2μÞðγ 2i ξ2 Þ þ ρω2 αρf ω2 ; αðmω2 ibωÞ ρf ω2 ρf ω δ3 ¼ ib þ mω δi ¼
i ¼ 1; 2
The stress components can be obtained by the solutions of displacements by Z 1 2 ½ð λkp1 þ2μγ 21 ÞA0 e γ 1 z τ^ zz ¼ 0
^ 0 ðξxÞdξ; þð λkp2 þ 2μγ 22 ÞB0 e γ 2 z 2μξςC 0 e ςz αpξJ 2
ð17Þ
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
Z τ^ zr ¼
1 0
μ½2ξγ 1 A0 e γ 1 z þ2ξγ 2 B0 e γ 2 z þ ðξ2 þ ς2 ÞC 0 e ςz ξJ 1 ðξxÞdξ; ð18Þ
p^ ¼
Z
1 0
½ðα þ δ1 ÞMkp1 A0 e γ 1 z þ ðα þ δ2 ÞMkp2 B0 e γ 2 z ξJ 0 ðξxÞdξ; 2
2
ð19Þ
The coefficients A0 C 0 can be determined by the free surface boundary conditions of Eqs. (11) and (12) as M 2 2 FðωÞ ; ðξ þ ς2 Þðα þ δ2 Þkp2 R1 2π M 2 FðωÞ ; and B01 ¼ ðξ2 þ ς2 Þðα þ δ1 Þkp1 R1 2π 2ξM 2 2 FðωÞ C0 ¼ ½ðα þ δ2 Þγ 1 kp2 ðα þ δ1 Þγ 2 kp1 R1 2π A0 ¼
ð20Þ
where R1 ¼ ðξ2 þ ς2 Þðs′1 g′2 s′2 g′1 Þ 4ξ2 ςðγ 1 g′2 γ 2 g′1 Þ; 2
s′i ¼ λkpi þ 2μγ 2i ; g′i ¼
i ¼ 1; 2; and
2 Mðα þδi Þkpi ;
i ¼ 1; 2
ð21Þ
The equation R1 also represents Rayleigh function in porous fluid-saturated media. By solving the equation of R1 ¼ 0, the reasonable complex root of ξ represents the wavenumber of Rayleigh surface waves in porous fluid-saturated media. 3.1. An equivalent-viscoelastic model for Rayleigh surface waves Based on the porous fluid-saturated model, the equivalentviscoelastic model is obtained by replacing the complex modulus in the constitutive equations by the viscoelastic model [18]. In a single-phase viscoelastic model, two types of body waves can be effectively represented in the relaxation mechanism associated with the velocity and the attenuation of each mode. The effective characteristics of Rayleigh surface waves can be straightforwardly represented by the combining of these two body-wave relaxation mechanisms. We use the standard viscoelastic model (Zener model) to represent the relaxed mechanism of the fast waves (the P1 wave and the shear wave) in a porous medium, and investigate the effective characteristics and the discrepancy of the representations for the fluid-saturated model. The complex modulus of a Zener viscoelastic element can be written by [18] 1 þ iωτε MðωÞ ¼ M 0 ; ð22Þ 1 þ iωτs where M 0 is the relaxed modulus, τε and τs are the strain relaxation time and the stress relaxation time of the Zener element, respectively. Fitting the relaxed modulus of the element in Eqs. ((7) and (8)) for the fast waves, the effective wavenumbers can be written by 1=2 ρ n kp ¼ kp1 ¼ ω ; ð23Þ M p ðωÞ n
ks ¼ ks ¼ ω
ρ M s ðωÞ
97
described by these complex modulus. Moreover, only one viscoelastic relaxation element is enough to represent the each body wave modulus in the porous medium. Therefore, one set of relaxation time τε and τs need to be determined for each wave modulus. For Zener model, they can be written by qffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 ð Q 20 þ 1 þ 1Þ; and τε ¼ ω0 Q 0 qffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 τs ¼ ð Q 20 þ 1 1Þ; ð25Þ ω0 Q 0 where ω0 ¼ 2πf 0 corresponds the peak attenuation of the wave loss mechanism, Q 0 is the quality factor at the frequency f 0 . Hence, the multiphase fluid-saturated medium becomes to a single-phase viscoelastic effective medium. The boundary problems for Rayleigh surface waves can be similarly derived by the free surface boundary conditions as Eq. (11). Nevertheless, the different pore fluid drainage conditions on the surface is ignored and effectively represented by the complex relaxations in the single-phase solid. The equivalent representations can be easily obtained from Eq. (25) by replacement the two wavenumbers with f 0i and Q 0i ; i ¼ p; s for the fast P1 wave and the shear wave without considering the slow P2 wave generated in the pore fluid phase. The effective displacement representations can be written by omitting the slow wave elements in Eqs. (13)–(16) as Z 1 n n n u^ z ¼ ð An0 γ np e γ p z þ C n0 ξe e ς z Þξe J 0 ðξe xÞdξe ; ð26Þ 0
n u^ r ¼
^ nz ¼ w ^ nr ¼ w
Z
1 0
Z
1 0
Z
1 0
n
n
ð An0 ξe e γ p z þ C n0 ςn e ς z Þξe J 1 ðξe xÞdξe ;
ð27Þ
n
n
ð28Þ
n
n
ð29Þ
ð δ1 An0 γ np e γ p z þ δ3 C n0 ξe e ς z Þξe J 0 ðξe xÞdξe ; ð δ1 An0 ξe e γ p z þ δ3 C n0 ςn e ς z Þξe J 1 ðξe xÞdξe ; n2
n2
where γ pn2 ¼ ξ2e kp and ςn2 ¼ ξ2e ks . ξe denotes the effective horizontal wavenumber. The coefficients An0 and C n0 can be determined by the free surface conditions as An0 ¼
ðξ2e þ ςn2 Þ FðωÞ ; M s R′ 2π
ð30Þ
C n0 ¼
2ξe γ np FðωÞ ; M s R′ 2π
ð31Þ
R′ ¼ ðξ2e þ ςn2 Þ2 4ξ2e γ np ςn ;
ð32Þ
where
To represent the effective characteristics of Rayleigh surface waves, Zener elements should be chosen to make Eqs. (23) and (24) over a full frequency range. However, the wave pulse used in seismic is usually through several hundred hertz, a relatively narrow band. Therefore, the effective frequency-dependent material properties in a limit frequency range determine the effective characteristics of the surface-wave propagations which we concern here.
1=2 ;
ð24Þ
where kp1 and ks are the complex wavenumbers of the fast P1 wave and the shear wave, respectively. M p ðωÞ and M s ðωÞ are the complex modulus for the fast P1 wave and the shear wave, respectively. For the P wave M 0p ¼ K G þ ð4=3Þ μ, and for the shear wave M 0s ¼ μ, where K G is Gassmann bulk modulus, and K G ¼ K m þ α2 M. The Biot dissipation mechanism can be approximately
3.2. Effective characteristics of Rayleigh surface waves We consider specific examples to validate the frequency range of applicability and the consistency of the effective characteristics of the equivalent-viscoelastic models for the fluid-saturated models. We choose representative samples of fluid-saturated materials (unconsolidated sand, Massilon sandstone, and clay) containing either water or air. Unconsolidated sand of small bulk modulus and high intrinsic
98
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
permeability represents typical geological material in very shallow aquifers. Malison sandstone and clay typify important consolidated porous media near surface. The elastic and hydraulic properties of the fluids (water and air) and the sedimentary materials are listed in Table 1 [39]. The phase velocity V R ðf Þ and the attenuation loss angle θR ðf Þ derived from the complex wavenumbers of Eqs. (9) and (10), the frequency-dependent dispersion curves of the velocity and attenuation of the surface wave can be computed from Rayleigh equations in the fluid-saturated model and the equivalentviscoelastic model at different frequencies. In addition, the time dependent transient dynamic waveforms are also computed, using Ricker-type source wavelets. The Ricker wavelet time function is given by FðtÞ ¼ ½1 2π 2 f s ðt t 0 Þ2 eπ 2
2 2 f s ðt t 0 Þ2
;
ð33Þ
We choose the peak frequencies f s of 5, 50 and 200 Hz, and t 0 ¼ 0:1 s for general seismic applications. The oscillatory integrals in Eqs. (13)–(16) and (26)–(29) can be implemented by the fast coverage algorithms, e.g., PTAM [61] and the residue integration [57]. We take the vertical solid displacements as examples for the typical discussion. In the following figures, the solid lines are the results for the different fluid-saturated models, and the dashed lines are the results for the effective equivalent-viscoelastic materials with a single Zener relaxation function for each body wave modulus. At the peak frequencies of each body wave modulus, the quality factors are set to be equal to those of the fluid-saturated models. The differences between the solid lines and dashed lines show the effective characteristics of the representations by Zener elements. Fig. 1 shows the dispersion curves of Rayleigh waves for the unconsolidated sand saturated by water (Fig. 1a) and air (Fig. 1b) over the whole frequency range. For the water-saturated case, the
Table 1 Material properties for three different pore fluids and porous media [39]. Parameter Pore fluids Bulk modulus of fluid (GPa) Material density of fluid (kg/m3) Viscosity (10 3N s/m2)
Kf ρf η
Porous media Bulk modulus of porous media (GPa) Shear modulus of porous media (GPa) Bulk modulus of solid (GPa) Material density of solid (kg/m3) Porosity Intrinsic permeability (10 13 m2)
Km μ Ks ρs ϕ κ
Water
Air
2.234 997 1
1.45 10 4 1.2 1.8 10 3
Unconsolidated sand
Massilon sandstone
Clay
0.125 0.06 35 2650 0.44 1100
1.021 1.441 35 2650 0.23 9
2.55 1.62 20 2650 0.4 3
Fig. 1. The effective dispersion curves of the phase velocity (left) and the loss angle (right) of Rayleigh waves for the unconsolidated sand saturated by water (a) and air (b). The critical frequency f c equals to 638.5 Hz for the water case (a), and 9549.3 Hz for the air case (b). The two equivalent-viscoelastic models have the peak frequencies of 436.9 Hz, 421.0 Hz and the quality factors of 17.5, 13.3 of the fast body waves for the water case (a), and the peak frequencies of 1.2489 104 Hz, 5840.9 Hz and the quality factors of 4488.5, 9200.6 for the air case (b), respectively. The equivalent representation underestimates the attenuation for the water case at low frequencies.
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
critical frequency f c equals to 638.5 Hz for the Biot dissipation mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼436.9 Hz and the quality factor Q 0p ¼ 17:5 for the fast P1 wave, and f 0s ¼ 421:0 Hz and Q 0s ¼ 13:3 for the shear wave. Significant discrepancy can be found in the velocity and the loss angle for frequencies around and greater than f c . Although the equivalent-viscoelastic model can approximate the body-wave propagations in fluid-saturated model in great agreements [17,18], these equivalent representations cannot represent the local minimum of the velocity dispersion of Rayleigh surface waves at the relative low frequencies, which is perturbed by the slow P2 wave. Meanwhile, the equivalent representations underestimate the attenuation of Rayleigh waves. It is because that the unconsolidated sand bears high permeability, Rayleigh surface waves present great attenuation. Furthermore, the frequency 1=2 range of f dependent diffusion at low frequencies for the slow P2 waves is very narrow in the fluid-saturated model, which makes the high velocity propagatory P2 waves in the seismic frequency band (o200 Hz) [62]. This phenomena strongly amplify the attenuation of Rayleigh surface waves in the fluid-saturated model (Fig. 1a, right), which the equivalent-viscoelastic model cannot represent. For the air-saturated case, the critical frequency f c equals to 9549.3 Hz for the Biot dissipation mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼1.2489 104 Hz and the quality factor Q 0p ¼4488.5 for the fast P1 wave, and f 0s ¼5840.9 Hz and Q 0s ¼ 9200.6 for the shear wave. In this case, the velocity dispersions are unconspicious, and the attenuation coefficients are very low. The dispersion curves of the velocity and the loss angle in the equivalent-viscoelastic model show effective consistency to the fluid-saturated model. Fig. 2 shows the computed waveforms of the vertical solid displacements of the surface wave propagations at offset 100 m on the free surface for the unconsolidated sand with the source peak frequencies f s of 5 Hz (left diagrams), 50 Hz (center diagrams), and 200 Hz (right diagrams). The wave modes are labeled on the
99
seismograms by the arrival time. In the water-saturated case (Fig. 2a), the waveforms for f s ¼50 Hz and f s ¼ 200 Hz show poor matches between the equivalent-viscoelastic model and the water-saturated model. There is a fairly well match between the equivalent representations and the fluid-saturated models for f s ¼5 Hz. In the air-saturated case (Fig. 2b), all the results show acceptable agreements. Fig. 3 shows the dispersion curves of Rayleigh waves for Malison sandstone saturated by water (Fig. 3a) and air (Fig. 3b) over the whole frequency range. For the water-saturated case, the critical frequency f c equals to 4.0795 104 Hz for the Biot dissipation mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼1.6403 104 Hz and the quality factor Q 0p ¼103.3 for the fast P1 wave, and f 0s ¼1.5550 104 Hz and Q 0s ¼51.9 for the shear wave. The equivalent-viscoelastic model also cannot effectively represent the local minimum for the velocity dispersion, and the viscoelastic mechanism also underestimates the attenuation for the water-saturation model around f c . For the air-saturated case, the critical frequency f c equals to 6.1009 105 Hz for the Biot dissipation mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼2.3542 105 Hz and the quality factor Q 0p ¼9.169 104 for the fast P1 wave, and f 0s ¼2.2816 105 Hz and Q 0s ¼3.954 104 for the shear wave. The Rayleigh surface-wave dispersion curves of the velocity and the loss angle show very good consistency. The Zener model effectively represents the surfacewave propagation in the fluid saturated model. Fig. 4 shows the computed waveforms of the vertical solid displacements of the surface wave propagations at offset 300 m on the free surface for Malison sandstone with the source peak frequencies f s of 5 Hz (left diagrams), 50 Hz (center diagrams), and 200 Hz (right diagrams). For the water-saturated case (Fig. 4a), the waveforms for f s ¼ 200 Hz show a poor match between the Zener model and the water-saturated model. The others are in acceptable agreements between the equivalent representations and the fluid-saturated models.
Fig. 2. The effective analytical transient waveforms of the surface wave propagations of the unconsolidated sand for the vertical solid displacements for the water (a) and the air (b) saturated models. The source pulses have the peak frequencies of 5 Hz (left), 50 Hz (center), and 200 Hz (right). The offset is 100 m. Except for 50 Hz and 200 Hz for the water case (a), the waveforms show a fairly well agreement with the fluid-saturated model.
100
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
Fig. 3. The effective dispersion curves of the phase velocity (left) and the loss angle (right) of Rayleigh waves for Malison sandstone saturated by water (a) and air (b). The critical frequency f c equals to 4.0795 104 Hz for the water case (a), and 6.1009 105 Hz for the air case (b). The two equivalent-viscoelastic models have the peak frequencies of 1.6403 104 Hz, 1.5550 104 Hz and the quality factors of 103.3, 51.9 of the fast body waves for the water case (a), and the peak frequencies of 2.3542 105 Hz, 2.2816 105 Hz and the quality factors of 9.169 104, 3.954 104 for the air case (b), respectively. The equivalent representation underestimates the attenuation for the water case at low frequencies.
Fig. 4. The effective analytical transient waveforms of the surface wave propagations of the unconsolidated sand for the vertical solid displacements for the water (a) and the air (b) saturated models. The source pulses have the peak frequencies of 5 Hz (left), 50 Hz (center), and 200 Hz (right). The offset is 300 m. Except for 200 Hz for the water case (a), the waveforms show a fairly well agreement with the fluid-saturated model.
Fig. 5 shows the dispersion curves of Rayleigh waves for the clay saturated by water (Fig. 5a) and air (Fig. 5b) over the whole frequency range. For the water-saturated case, the critical frequency f c equals to 2.1285 105 Hz for the Biot dissipation
mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼1.6437 105 Hz and the quality factor Q 0p ¼ 5.383 104 for the fast P1 wave, and f 0s ¼1.2889 105 Hz and Q 0s ¼16.4 for the shear wave. For the air-saturated case, the critical
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
101
Fig. 5. The effective dispersion curves of the phase velocity (left) and the loss angle (right) of Rayleigh waves for the clay saturated by water (a) and air (b). The critical frequency f c equals to 2.1285 105 Hz for the water case (a), and 3.1831 106 Hz for the air case (b). The two equivalent-viscoelastic models have the peak frequencies of 1.6437 105 Hz, 1.2889 105 Hz and the quality factors of 5.383 104, 16.4 of the fast body waves for the water case (a), and the peak frequencies of 1.8635 106 Hz, 1.8194 106 Hz and the quality factors of 1.365 104, 1.160 104 for the air case (b), respectively.
Fig. 6. The effective analytical transient waveforms of the surface wave propagations of the unconsolidated sand for the vertical solid displacements for the water (a) and the air (b) saturated models. The source pulses have the peak frequencies of 5 Hz (left), 50 Hz (center), and 200 Hz (right). The offset is 300 m. Except for 200 Hz for the water case (a), the waveforms show a fairly well agreement with the fluid-saturated model.
frequency f c equals to 3.1831 106 Hz for the Biot dissipation mechanism. The Zener fits are obtained with the peak attenuation frequency f 0p ¼1.8635 106 Hz and the quality factor Q 0p ¼ 1.3654 104 for the fast P1 wave, and f 0s ¼ 1.8194 106 Hz and Q 0s ¼ 1.160 104 for the shear wave. Except for the local minimum
of the velocity dispersion in the water-saturation models, the Rayleigh surface-wave dispersion curves of the velocity and the loss angle show very good consistency. Fig. 6 shows the computed waveforms of the vertical solid displacements of the surface wave propagations at offset 300 m on
102
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
the free surface for the clay with the source peak frequencies f s of 5 Hz (left diagrams), 50 Hz (center diagrams), and 200 Hz (right diagrams). For the water-saturated case (Fig. 6a), the waveform for f s ¼200 Hz shows a deficient representation of the Zener model for the water-saturated model. The others are in acceptable consistency between the equivalent representations and the fluid-saturated models.
4. Discussions and conclusions The comparisons of the effective characteristics for the fluidsaturated model provide a means for checking the validations of the equivalent Zener representations. Although the slow wave effect in the dispersion is ignored, Rayleigh surface waves in the fluid-saturated models can be effectively approximated by the Zener equivalent-viscoelastic model in case of the frequency range is much lower than the critical frequency f c . Therefore, the fluidsaturated materials with high critical frequency f c can be more effectively represented by the equivalent-viscoelastic models. Moreover, it is noted that the surface drainage effects on the surface waves can be ignored in case of high f c (high coupling damping b) in the seismic frequency band ( o200 Hz) [62,63]. The porous materials of high f c bear the effective characteristics can be represented more universally. The numerical methods for modeling the seismic waves propagation in porous media based on the equivalent representations by the Zener elements have progressed, e.g., Carcione [17,18], Picotti et al. [46,47], and Liu et al. [34,35]. By considering the free surface boundary conditions, these methods can extent to the surface wave propagations. In summary, we have analyzed the effective characteristics of Rayleigh surface wave propagations in fluid-saturated media for viscoelastic representations. The velocity and attenuation dispersions, as well as the transient displacement waveforms in a homogeneous half-space of the fluid-saturated models and the equivalent-viscoelastic models expressed by the Zener elements are analytically investigated. By fitting the attenuation coefficients of the fluid-saturated model at the peak frequencies for the fast body waves, the equivalent-viscoelastic model with a single Zener element for each fast wave modulus can effectively represent the surface wave propagations for the velocity and attenuation dispersions, and the transient waveforms in the fluid-saturated model. It should be noted that the restriction of the effective characteristics is that the frequency is much lower than the critical frequency f c of the fluid saturated material. Although, the equivalent-viscoelastic representation do not completely reflect the local minimum in the velocity dispersion perturbed by the slow P2 wave in the fluid-saturated model, which may also be shown by the underestimated attenuation, the transient waveform in the fluid-saturated model at low frequencies, especially for Rayleigh wave can be effectively represented. The porous materials with high critical frequency f c show a more universal effective characteristics and better representation by the equivalentviscoelastic models at low frequencies. For the porous materials with low critical frequency f c (usually with high permeability κ), we must be careful to check the validated frequency range for the effective characteristics. Notably, our discussion is based on the sample materials in this paper. It is important to check the validity for the actual materials before using the equivalent representations in practical applications. The significant advantage of applying the equivalent-viscoelastic model is that they allow memory equations to replace the convolution integral in the numerical modeling for saving the field quantities storing consumptions, which also presents a simple way to introduce a constrained model considering poroelasticity for surface wave inversion.
Acknowledgments Gratitude is expressed for the financial support by the Natural Basic Research Program of China (the “973 Project”, Grant No. 2013CB733303), the National Natural Science Foundation of China (NSFC, Grant No. 41304077, 40974079), Postdoctoral Science Foundation of China (Grant No. 2013M531744), and Key Laboratory of Geospace Environment and Geodesy (Grant No. 12-02-03). References [1] Allard JF, Jansen G, Vermeir G, Lauriks W. Frame-borne surface waves in airsaturated porous media. J Acoust Soc Am 2002;111:690–6. [2] Allard JF, Henry M, Glorieux C, Petillon S, Lauriks W. Laser induced surface modes at an air-porous medium interface. J Appl Phys 2003;93:1298–304. [3] Allard JF, Henry M, Glorieux C, Lauriks W, Petillon S. Laser induced surface modes at water-elastic and poroelastic interfaces. J Appl Phys 2004;95: 528–35. [4] Albers B, Wilmanski K. Monochromatic surface waves on impermeable boundaries in two-component poroelastic media. Continuum Mech Thermodyn 2005;17:269–85. [5] Ben-Menahem A, Singh SG. Seismic waves and sources. New York: SpringerVerlag New-York Inc; 1981. [6] Berryman JG. Confirmation of Biot's theory. Appl Phys Lett 1980;37(4):382–4. [7] Berryman JG, Thigpen L, Chin RCY. Bulk elastic wave propagation in partial saturated porous solids. J Acoust Soc Am 1988;84(1):360–73. [8] Berryman JG, Pride SR. Volume averaging, effective stress rules, and inversion for microstructural response of multicomponent porous media. Int J Solids Struct 1998;35:4811–43. [9] Berryman JG. Extension of poroelastic analysis to double-porosity materials: new technique in microgeomechanics. J Eng Mech 2002;128(8):840–7. [10] Biot MA. Theory of propagation of elastic waves in fluid filled porous solid I. Low frequency range. J Acoust Soc Am 1956;28:168–78. [11] Biot MA. Theory of propagation of elastic waves in fluid filled porous solid II. Higher frequency range. J Acoust Soc Am 1956;28:179–91. [12] Biot MA, Willis DG. The elastic coefficients of the theory of consolidation. J Appl Mech 1957;33:1482–98. [13] Biot MA. Mechanics of deformation and acoustic propagation in porous media. J Appl Phys 1962;33:1482–98. [14] Bourbié T, Coussy O, Zinszner B. Acoustics of porous media. Paris: Gulf Publishing Company; 1987. [15] Carcione JM, Quiroga-Goode G. Some aspects of the physics and numerical modelling of Biot compressional waves. J Comput Acoust 1995;3(4):261–80. [16] Carcione JM, Quiroga-Goode G. Full frequency-range transient solution for compressional waves in a fluid-saturated viscoacoustic porous medium. Geophys Prospect 1996;44:99–129. [17] Carcione JM. Viscoelastic effective rheologies for modelling wave propagation in porous media. Geophys Prospect 1998;46:249–70. [18] Carcione J. Wave fields in real media: wave propagation in anisotropic, anelastic, porous, and electromagnetic media. 2nd ed.. Amsterdam: Elsevier; 2007. [19] Carcione JM, Gei D, Picotti S, Michelini A. Cross-hole electromagnetic and seismic modeling for CO2 detection and monitoring in a saline aquifer. J Petrol Sci Eng 2012;100:162–72. [20] Chao G, Smeulders DM, van Dongen MEH. Dispersive surface waves along partially saturated porous media. J Acoust Soc Am 2006;119(3):1347–55. [21] Dai Z, Kuang Z, Zhao S. Rayleigh waves in a double porosity half-space. J Sound Vib 2006;298:319–32. [22] Deresiewicz H. The effect of boundaries on wave propagation in a liquid filled porous media: IV. Surface waves in a half-space. Bull Seismol Soc Am 1962;52 (3):627–38. [23] Deresiewicz H, Skalak R. On uniqueness in dynamic poroelasticity. Bull Seismol Soc Am 1963;53(4):783–8. [24] Dvorkin J, Nolen-Hoeksema R, Nur A. The squirt-flow mechanism: macroscopic description. Geophysics 1994;59(3):428–38. [25] Edelman I, Wilmanski K. Asymptotic analysis of surface waves at vacuum/ porous medium and liquid/porous medium interfaces. Continuum Mech Thermodyn 2002;14:25–44. [26] Edelman I. Surface waves at vacuum/porous medium interface: low frequency range. Wave Motion 2004;39:111–27. [27] Feng S, Johnson DL. High-frequency acoustic properties of a fluid/porous solid interface. I. New surface mode. J Acoust Soc Am 1983;74:906–14. [28] Feng S, Johnson DL. High-frequency acoustic properties of a fluid/porous solid interface. II. The 2D reflection Green's functions. J Acoust Soc Am 1983;74: 915–24. [29] Geertma J, Smit DC. Some aspects of elastic wave propagation in fluidsaturated porous solids. Geophysics 1961;26:169–81. [30] Gubaidullin AA, Kuchugurina OY, Smeulders DM, Wisse CJ. Frequencydependent acoustic properties of a fluid/porous solid interface. J Acoust Soc Am 2004;116:1474–80. [31] Halpern MR, Christiano P. Response of poroelastic halfspce to steady-state harmonic surface traction. Int J Nemer Anal Mech Geomech 1986;10:609–32.
Y. Zhang et al. / Soil Dynamics and Earthquake Engineering 57 (2014) 94–103
[32] Johnson DL, Koplik J, Dashen R. Theory of dynamic permeability and tortuosity in fluid-saturated porous-media. J Fluid Mech 1987;176:379–402. [33] Jones J. Rayleigh wave in a porous elastic saturated solid. J Acoust Soc Am 1961;33:959–62. [34] Liu X, Greenhalgh S, Zhou B. Transient solution for poro-viscoacoustic wave propagation in double porosity media, and its limitations. Geophys J Int 2009;178:375–93. [35] Liu X, Greenhalgh S, Zhou B. Approximating the wave moduli of double porosity media at low frequencies by a single Zener or Kelvin-Voigt element. Geophys J Int 2010;181:391–8 (1117). [36] Lo W, Sposito G, Majer E. Wave propagation through elastic porous media containing two immiscible fluids. Water Resour Res 2005;12(3):1105–17. [37] Lo W. Propagation and attenuation of Rayleigh waves in a semi-infinite unsaturated poroelastic medium. Adv Water Resour 2008;31(10):1399–410. [38] Markov MG. Low-frequency Stoneley wave propagation at the interface of two porous half-spaces. Geophys J Int 2009;177:603–8. [39] Mavko G, Mukerji T, Dvorkin J. The rock physics handbook: tools for seismic analysis of porous media. Cambridge: Cambridge University Press; 1998. [40] Morochnik V, Bardet JP. Viscoelastic approximation of poroelastic media for wave scattering problems. Soil Dyn Earthq Eng 1996;15:337–46. [41] Nazarian S, Stokoe II KH, Hudson WR. Use of spectral analysis of surface waves method for determination of moduli and thicknesses of pavement systems. Transp Res Record 1983;930:38–45. [42] Park CB, Miller RD, Xia J. Multi-channel analysis of surface waves (MASW). Geophysics 1999;64(3):800–8. [43] Paul S. On the displacement produced in a porous elastic half-space by an impulsive line load (non-dissipative case). Pure Appl Geophys 1976;114: 605–14. [44] Paul S. On the disturbance produced in a semi-infinite poroelastic medium by a surface load. Pure Appl Geophys 1976;114:615–27. [45] Philippacopoulos AJ. Lamb's problem for fluid saturated porous media. Bull Seism Soc Am 1988;78:908–23. [46] Picotti S, Carcione JM, Rubini G, Santos JE, Cavallini F. A viscoelstic representation of wave attenuation in porous media. Comp Geosci 2010;36:44–53. [47] Picotti S, Carcione JM, Gei D, Rossi G, Santos JE. Seismic modeling to monitor CO2 geological storage—the Atzbach-Schwanenstadt gas field. J Geophys Res 2012;117(6):B06103. [48] Picotti S, Carcione JM, Santos JE. Oscillatory numerical experiments in finely layered anisotropic viscoelastic media. Comp Geosci 2012;43:83–9.
103
[49] Pride SR, Berryman JG. Linear dynamics of double-porosity and dualpermeability materials I. Governing equations and acoustic attenuation. Phys Rev E 2003;68:036603. [50] Pride SR, Berryman JG. Linear dynamics of double-porosity and dualpermeability materials II. Fluid transport equations. Phys Rev E 2003;68: 036604. [51] Plona TJ. Observation of a second bulk compressional wave in a porous medium at ultrasonic frequencies. Appl Phys Lett 1980;36:259–61. [52] Rubino JG, Ravazzoli CL, Santos JE. Equivalent viscoelastic solids for heterogeneous fluid-saturated porous rock. Geophysics 2009;74:N1–13. [53] Santos JE, Corbero JM, Douglas J. Static and dynamic behavior of a porous solid saturated by a two-phase fluid. J Acoust Soc Am 1990;87(4):1428–38. [54] Santos JE, Corbero JM, Douglas J. A model for wave propagation in a porous medium saturated by a two-phase fluid. J Acoust Soc Am 1990;87(4):1439–48. [55] Senjuntichai T, Rajapakse RKND. Dynamic Green's functions of homogeneous poroelastic half-plane. J Eng Mech ASCE 1994;120:2381–404. [56] Tajuddin M. Rayleigh waves in a poroelastic half-space. J Acoust Soc Am 1984;75:682–4. [57] van Dalen KN, Drijkoningen GG, Smeulder DM. On wavemodes at the interface of a fluid and a fluid-saturated poroelastic solid. J Acoust Soc Am 2010;127: 2240–51. [58] Xia J, Miller RD, Park CB. Estimation of near-surface shear-wave velocity by inversion of Rayleigh wave. Geophysics 1999;64:691–700. [59] Xu Y, Xia J, Miller RD. Numerical investigation of implementation of air–earth boundary by acoustic–elastic boundary approach. Geophysics 2007;72(5): SM147–SM153. [60] Zeng C, Xia J, Miller RD, Tsoflias GP. Application of the multiaxial perfectly matched layer (M-PML) to near-surface seismic modeling with Rayleigh waves. Geophysics 2011;76(3):T43–52. [61] Zhang H, Chen X, Chang S. An efficient numerical method for computing synthetic seismogram for a layered half-space with sources and receivers at close or same depths. Pure Appl Geophys 2003;160:467–86. [62] Zhang Y, Xu Y, Xia J. Analysis of dispersion and attenuation of surface waves in poroelastic media in the exploration seismic frequency band. Geophys J Int 2011;187:871–88. [63] Zhang Y, Xu Y, Xia J. Wave fields and spectra of Rayleigh waves in poroelastic media in the exploration seismic frequency band. Adv Water Resour 2012;49:62–71.