Site amplification effects of a radially multi-layered semi-cylindrical canyon on seismic response of an earth and rockfill dam

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

Contents lists available at ScienceDirect

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

Site amplification effects of a radially multi-layered semi-cylindrical canyon on seismic response of an earth and rockfill dam

T

⁎

Ning Zhanga, Yu Zhanga, Yufeng Gaoa, , Ronald Y.S. Pakb, Jun Yangc a

Key Laboratory of Ministry of Education for Geomechanics and Embankment Engineering, Hohai University, No. 1, Xikang Road, Nanjing 210098, China Department of Civil, Environmental and Architectural Engineering, University of Colorado at Boulder, CO 80309-0428, USA c Department of Civil Engineering, The University of Hong Kong, Hong Kong, China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Wave propagation and scattering Ground motion amplification Inhomogeneity, Canyon Earth and rockfill dam Varying shear wave velocity Wave concentration Site effects Wave function expansion

An analytical treatment of an earth and rockfill dam in a radially inhomogeneous multi-layered semi-cylindrical canyon in a half-space under obliquely incident plane SH waves is presented. In terms of a radial wave function expansion and a transfer matrix approach, a rigorous approach is formulated to consider inhomogeneity of the deformable canyon and associated phenomena of wave redirection, reflection, transmission, scattering and concentration. Upon confirmation of its accuracy with a past exact solution for the degenerated case of a dam in a homogeneous canyon, the multi-layer solution is extended to the case of a dam in a canyon covered by a surficial zone with a power-law shear-wave velocity model as a generalized class of smooth in-situ variations. A comprehensive set of numerical examples are presented to illustrate the sensitivity of the dam response to the inhomogeneity profile of the top canyon zone, the frequency content and the angle of the seismic wave incidence. As illustrated in both frequency and time domains, wave concentration is found to be the main reason for the high amplifications of the displacement within the dam.

1. Introduction The seismic response of earth and rockfill dams is a subject of considerable interest in earthquake engineering. A comprehensive review on the research progress in an early stage can be found in Gazetas [1] and Gazetas and Dakoulas [2]. Due to the complexity of the problem, a large number of theoretical and experimental studies have recently been conducted with emphasis on its different aspects, for example, the advanced nonlinear constitutive models of dam materials [3,4], the liquefaction of dam foundation [5], the seismically induced failure of the face plate of rockfill dams [6–8], the performance-based seismic analysis procedure [9] and the validation of dynamic analysis approaches by well-documented field performance data [10,11]. A common assumption in most of these studies is that the excitations at all point along the dam-canyon interface are identical (i.e., spatially uniform input ground motion). However, it has been understood for many years that the seismic motion in a canyon is spatially non-uniform [2,12]. Some numerical research has been directed toward the response of earth and rockfill dams to spatially non-uniform earthquake ground motion by using either deterministic methods by considering wave propagation along the dam foundation/canyon [13–17] or random

⁎

vibration analysis with aid of stochastic ground motion models [18]. Although numerical methods are more flexible and useful tools, analytical closed-form solutions are valuable in revealing the mathematical and physical nature of the problem as well as testing the numerical methods. To explore the effect of the spatial variation of the ground motion, Dakoulas and his co-workers were able to propose analytical solutions for the lateral response of earth and rockfill dams in rectangular, semi-cylindrical and semi-elliptical canyons under obliquely incident harmonic SH waves [19–21]. The developed analytical solutions were based on a generalized shear-beam model of dam and have been used to examine the effects of the canyon geometry, the dam-to-canyon impedance ratio, and the obliquity and frequency of the incident waves. More importantly, the results demonstrated the necessity of accounting for the wave scattering phenomenon associated with the dam-filled homogeneous canyon in an elastic half-space. Earth and rockfill dams are often constructed on weathered rock or even alluvial soil. A case in point is the Lechago dam in Spain whose foundation contains a 5–10 m thick weathered layer as well as alluvial soils [22]. For engineers in China, it is also common to encounter weathered rocks in a canyon site and to use them as dam foundation [23]. In engineering seismology, modification of ground motion due to wave scattering by local irregular topographies such as canyons is a

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

https://doi.org/10.1016/j.soildyn.2018.09.014 Received 6 June 2018; Received in revised form 12 August 2018; Accepted 12 September 2018 0267-7261/ © 2018 Elsevier Ltd. All rights reserved.

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

α in the x-z plane. Both a Cartesian and a cylindrical coordinate system are employed as indicated in Fig. 1(b).

subject of fundamental interest [24–29]. An analytical solution in terms of Bessel-Hankel series for SH wave scattering by a semi-cylindrical canyon with a surficial layer with a continuous power-law variation of its shear modulus on its curved boundary were recently realized by Zhang et al. [30]. It was found that the inhomogeneity profile of the finite surficial layer due to surficial weathering or soil sedimentation to a finite depth from water flows or other geological mechanisms can be very influential to the level and pattern of the ground motion amplifications. Since a dam canyon often contains multiple layers of different material properties due to weathering or soil sedimentation [22], a question that arises is whether and to what extent the site amplification from an inhomogeneous canyon affects the dynamic response of a dam built in it. To this end, this study develops an analytical solution to the lateral response of dams in a radially multi-layered semi-cylindrical canyon in a homogeneous half-space and addresses the engineering problem from the perspective of seismology with an emphasis on the wave propagation and scattering aspect. To illustrate the importance of incorporating not only the canyon's topology but also its specific depth-wise material variation in its surficial region in dam response assessments, the rigorous formulation is employed to examine the case where the shear wave velocity of a surficial zone is described by a canonical power-law profile. In what follows, it is first illustrated that the formulation is amenable to a rigorous analytical treatment by means of separation of variables and a new developed transfer matrix approach for radially layered media. Then, spatial distributions of displacement amplitudes across the crest as well as the maximum longitudinal section of the dam in an inhomogeneous canyon are studied in the frequency domain, relative to the case of a dam-homogeneous canyon system. Finally, key features of wave scattering and concentration phenomena are demonstrated in terms of synthetic seismograms and snapshots in the time domain.

3. Analytical formulation 3.1. Wave motion in the dam Under anti-plane shear wave excitation, the response of the dam idealized by a generalized shear beam model can be characterized fully by the displacement wd (r, θ, t) in the y-direction (e.g., Gazetas and Dakoulas [2]). For the time-harmonic problem, one may write in complex notation the displacement response wd(r, θ, t) as wd(r, θ)eiωt so that the governing equation of motion can be written as (e.g., Dakoulas [20])

∂2wd ∂w ∂2wd ∂w + 2r d + − tan θ d + r 2kd2 wd = 0, ∂r 2 ∂r ∂θ 2 ∂θ π π − ≤θ≤ , 2 2

r2

r ≤ H, (1)

and the traction-free boundary conditions along the dam crest as d τθy (r , θ) =

Gd* ∂wd (r , θ) π = 0, θ = ± , r ≤ H , r ∂θ 2

(2)

where kd = ω/ Vd* is the wave-number of the medium in the dam with Vd* = Gd*/ ρd and Gd* = Gd (1 + 2iβd ) being the complex shear wave velocity and the complex shear modulus in the dam, respectively. The solution of Eq. (1) satisfying Eq. (2) has been derived by Dakoulas [20] and can be written in the following form ∞

wd (r , θ) =

∑ An jn (kd r ) Φn (θ),

(3)

n=0

where An is the constant to be determined, jn (⋅) is the spherical Bessel functions of the first kind of order n, and

2. Simplified model for a dam-inhomogeneous canyon system under plane SH waves

Φn (θ) = inPn [ cos(θ + α )] + (−i)nPn [ cos(θ − α )],

(4)

with Pn being the Legendre polynomial.

A simplified model of an earth dam built in a homogeneous semicylindrical canyon has been developed by Dakoulas [20]. To seek a balance between mathematical convenience and physical reality, the dam was idealized as a “shear beam” which assumes either uniform or average response values for the upstream-downstream direction. To reveal the influences of input ground motion amplification from an inhomogeneous canyon on an earth dam, the dam-homogeneous canyon model by Dakoulas [20] is extended to a dam-inhomogeneous canyon model in this paper. Consider the SH wave propagation and scattering problem associated with an earth dam-filled multi-layered semi-cylindrical canyon. A perspective view of the geometry of the dam-layered canyon system in a half-space is given in Fig. 1(a). As shown in Fig. 1(b), both the dam and the canyon are assumed to have a semi-circular cross section in the x-z plane. The height of the dam or the depth of the canyon is denoted as H, and the length of the dam or the width of the canyon is L = 2H. The inhomogeneous surficial zone of the canyon is allowed to be composed of M semi-circular bonded layers which share the same geometric center at O. The inner and outer radii of the mth layer (m = 1, 2, …, M) are denoted as Hm-1 and Hm, respectively. The total thickness of the layered zone is S = HM – H0, with H0 being equal to H. The material properties of the media in each layer and in the half-space are all assumed to be homogeneous, isotropic and linearly elastic. The mass density, shear modulus, and shear wave velocity of the media in the mth layer are denoted as ρm, Gm and Vm, and those of the underlying elastic half-space as ρh, Gh and Vh, respectively. Following Dakoulas [20], the dam has a triangular cross section in the y-z plane (Fig. 1c) and is idealized as a generalized shear-beam model with homogeneous linearly hysteretic soil of mass density ρd, shear modulus Gd, shear wave velocity Vd, and material hysteretic damping ratio βd. The seismic excitation is taken to be a steady train of timeharmonic plane SH waves with a circular frequency ω, an incident angle

3.2. Wave motion in the radially layered canyon The steady state dynamic response in the mth layer of the canyon can also be characterized fully by the anti-plane horizontal displacement wm (r, θ) in the y-direction, with the governing equation of motion (Helmholtz equation) as

1 ∂wm 1 ∂2wm ∂2wm + + 2 + km2 wm = 0, 2 ∂r r ∂r r ∂θ 2 π π − ≤ θ ≤ , m = 1, 2, ... ,M 2 2

Hm − 1 ≤ r ≤ Hm, (5)

and the traction-free boundary conditions along the horizontal portion of the ground surface of the surficial layers as (m) τθy (r , θ) =

Gm ∂wm (r , θ) π = 0, θ = ± , Hm − 1 ≤ r ≤ Hm, m = 1, 2, ... r ∂θ 2 (6) ,M ,

where km = ω/ Vm is the wave-number of the medium in the mth layer. Appropriate for semi-annular domains with the conditions in Eq. (6), a general wave function expansion of wm(r, θ) in terms of Bessel series can be reduced from Eq. (21) in Zhang et al. [30] by setting β = 0. To facilitate the total solution for a dam-filled canyon, according to Dakoulas [20], it's necessary to reform the wave field wm(r, θ) by using the spherical Bessel series as ∞

wm (r , θ) =

∞

∑ Bn(m) jn (km r ) Φn (θ) + ∑ Cn(m) yn (km r ) Φn (θ), n=0

= 1, 2, …, M , where 146

Bn(m)

and

m

n=0

Cn(m)

(7)

are constants pertaining to the mth layer, and jn (⋅)

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 1. 3D model of a dam in a radially layered semi-cylindrical canyon subject to plane SH waves.

where kh = ω/ Vh is the shear wave-number of the medium. The solution is subject to the Sommerfeld radiation condition for r →∞ as well as the traction-free surface condition of (h) τθy (r , θ) =

Gh ∂wh (r , θ) π = 0, r ≥ HM , θ = ± . r ∂θ 2

(9)

With an incident plane SH wave excitation of constant amplitude with an arbitrary angle of incidence α in the r-θ plane, it is convenient to decompose the resultant free wave field in a cylindrically-gorged half-space into

wh (r , θ) = w f (r , θ) + w s (r , θ)

(10)

(r , θ) is the free-field motion corresponding to a flat halfwhere space without the layered canyon, while the second part is the scattered wave field w s (r , θ) representing the perturbation of w f (r , θ) due to the presence of the canyon and its surficial inhomogeneity. Consistent with the wave expansion format in Eqs. (33) and (34) in Dakoulas [20], the free-field and scattered wave fields can be listed as wf

Fig. 2. Shear wave velocity profiles for a dam in a canyon with a M-layer surficial zone with material parameters Vh/Vd = 4, V0/Vd = 4/3, ρ/ρd = 1.5, βd = 10% and geometric parameters S/H = 1.

and yn (⋅) are spherical Bessel functions of the first and second kinds of order n, respectively.

∞

w f (r , θ) =

∑ (2n + 1) jn (kh r ) Φn (θ),

(11)

n=0

∞

3.3. Wave motion in underlying gorged half-space

w s (r , θ) =

r ≥ Hm, −

π π ≤θ≤ , 2 2

(12)

n=0

Analogous to the treatment of the multi-layered canyon zone, the governing equation for the anti-plane horizontal displacement wh(r,θ) in the underlying gorged homogeneous half-space satisfies the Helmholtz equation

∂ 2 wh 1 ∂wh 1 ∂ 2 wh + + 2 + kh 2wh = 0, 2 ∂r r ∂r r ∂θ 2

∑ Dn hn(2) (kh r ) Φn (θ), hn(2)

where {Dn} are constant coefficients and (⋅) is the spherical Hankel function of the second kind of order n. When combined with the time factor eiωt, Eq. (12) represents outgoing waves and satisfies the required radiation condition as well as the traction-free surface condition in Eq. (9) identically on the horizontal ground surface of the gorged halfspace.

(8) 147

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 3. Displacement amplitude amplification factors along dam crest (−1 ≤ x/H ≤ 1) and along surface of half-space near dam (x/H < −1 & x/H > 1) for a damM-layer canyon system whose shear wave velocity profiles are shown in Fig. 2 under incident angles α of 0°, 30°, 60° and 90° at dimensionless frequency Ω = 12.

Fig. 4. Comparison of the present solution with previous exact solution (Dakoulas [20]) in terms of displacement amplitude amplification factors along dam crest (−1 ≤ x/H ≤ 1) and along surface of half-space near dam (x/H < −1 & x/H > 1) for ξ = 0, S/H = 1, and βd = 10% at Ω = 4.5.

conditions of both traction and displacement fields along (a) the semicircular interface at r = H between the dam and the 1st layer of the multi-layered zone, (b) all (M-1) semi-circular internal interfaces in the inhomogeneous zone, and (c) the semi-circular interface at r = HM between the Mth layer and the underlying rutted half-space. Expressed

3.4. Total solution for seismic problem of the dam-layered canyon By virtue of the representations of the wave fields wd in Eq. (3), {wm }, m= 1 to M in Eq. (7) and wh in Eq. (10), the remaining requirements for the realization of the full dynamic solution are the continuity 148

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 5. Displacement amplitude amplification factors along dam crest (−1 ≤ x/H ≤ 1) and along surface of half-space near dam (x/H < −1 & x/H > 1) for a canyon with a homogeneous (ξ = 0) or inhomogeneous (ξ = 0.5, 1) surficial zone with dimensionless thickness S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 3.

substituting Eqs. (3) and (7) into the continuity conditions of Eqs. (13) and (14) yields the following system of equations in matrix form

explicitly, they are

τryd (r ,

θ) =

τry(1)

π π (r , θ), r = H , − ≤ θ ≤ , 2 2

π π wd (r , θ) = w1 (r , θ), r = H , − ≤ θ ≤ , 2 2

(13)

yn (k1 H ) ⎤ ⎛ Bn(1) ⎞ ⎛ jn (kd H ) ⎞ ⎡ jn (k1 H ) ⎢ G1 j ′ (k1 H ) G1 y′ (k1 H ) ⎥ ⎜ (1) ⎟ = ⎜G *j ′ (kd H ) ⎟ An , n = 0, 1, 2, ... n ⎠ ⎣ n ⎦ ⎝Cn ⎠ ⎝ d n

(14)

(19)

π π τry(m) (r , θ) = τry(m + 1) (r , θ), r = Hm, − ≤ θ ≤ , m = 1, 2, …, M − 1, 2 2 (15) wm (r , θ) = wm + 1 (r , θ), r = Hm, −

π π ≤ θ ≤ , m = 1, 2, …, M − 1, 2 2 (16)

π π τry(M ) (r , θ) = τryh (r , θ), r = HM , − ≤ θ ≤ , 2 2 wM (r , θ) = wh (r , θ), r = HM , − where τryd = Gd*

∂wd ∂r

where jn′ (⋅) and yn′ (⋅) denote the derivatives of jn (⋅) and yn (⋅) , respectively. By using the Wronskian relationship for spherical Bessel functions [31], the aforementioned matrix can be simplified as

π π ≤θ≤ , 2 2

(1) ⎛ Bn ⎞ ⎛t1 ⎞ ⎜C (1) ⎟ = ⎜t2 ⎟ An = tAn , n = 0, 1, 2, ... ⎝ n ⎠ ⎝ ⎠

where

(17)

t1 = k1 H 2 [jn (kd H ) yn′ (k1 H ) − (18)

∂w ∂w ∂wd in the dam, τry(m) = Gm ∂rm = ρm (Vm)2 ∂rm ∂r ∂wh ∂wh 2 Gh ∂r = ρh (Vh) ∂r in the half-space.

= ρd (Vd*)2 τryh

(20)

Gd* j ′ (kd H ) yn (k1 H )], G1 n

t2 = k1 H 2 [−jn (kd H ) jn′ (k1 H ) +

= in the mth layer and By virtue of the orthogonality of the Legendre polynomial,

Gd* j ′ (kd H ) jn (k1 H )], G1 n

(21)

(22)

By means of Eq. (7) and the Wronskian relationship for spherical 149

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 6. Displacement amplitude amplification factors along dam crest (−1 ≤ x/H ≤ 1) and along surface of half-space near dam (x/H < −1 & x/H > 1) for a canyon with different dimensionless thickness of the surficial zone (S/H = 0, 0.1, 1, 2) with inhomogeneity exponents ξ = 0.5 under incident angles α of 0°, 30°, 60° and 90° at Ω = 4.5.

Bessel functions [31], Eqs. (15) and (16) translate to the matrix recurrence relationship

(m) T21 = km + 1 Hm2 [−jn (km Hm) jn′ (km + 1 Hm) +

Gm j ′ (km Hm) jn (km + 1 Hm)], Gm + 1 n (26)

(m) (m) (m + 1) (m) (m) ⎛B ⎞ ⎛ Bn ⎞ ⎡ T11 T12 ⎤ ⎛ Bn ⎞ = T (m) ⎜ n(m) ⎟, m = 1, 2, …, M − 1, ⎜C (m + 1) ⎟ = ⎢ (m) (m) ⎥ ⎜ (m) ⎟ ⎝Cn ⎠ ⎝ n ⎠ ⎣ T21 T22 ⎦ ⎝Cn ⎠ (23)

(m) T22 = km + 1 Hm2 [−yn (km Hm) jn′ (km + 1 Hm) +

(27) By Eqs. (7) and (10), the interface conditions of Eqs. (17) and (18) between the Mth layer and the half-space similarly yield

between the wavefield coefficents in the (m+1)th and mth layers, with the elements of the transfer matrix T (m) being (m) T11 = km + 1 Hm2 [jn (km Hm) yn′ (km + 1 Hm) −

(M ) (M ) (M ) (M ) ⎛⎜ 2n + 1⎟⎞ = ⎡ T11 T12 ⎤ ⎛ Bn ⎞ = T (M ) ⎛ Bn ⎞, ⎢ (M ) ⎜ (M ) ⎥ ⎜ (M ) ⎟ (M ) ⎟ Dn ⎠ ⎣ T21 T22 ⎦ ⎝Cn ⎠ ⎝ ⎝Cn ⎠

Gm j ′ (km Hm) yn (km + 1 Hm)], Gm + 1 n (24)

(m) T12 = km + 1 Hm2 [yn (km Hm) yn′ (km + 1 Hm) −

Gm y′ (km Hm) jn (km + 1 Hm)]. Gm + 1 n

where the elements of the transfer matrix (M ) T11

Gm y′ (km Hm) yn (km + 1 Hm)], Gm + 1 n (25)

(28) are

G 2 = ikh HM [jn (kM HM ) hn(2) ′ (kh HM ) − M jn′ (kM HM ) hn(2) (kh HM )], Gh

(29)

GM y′ (kM HM ) hn(2) (kh HM )], Gh n

(30)

(M ) 2 T12 = ikh HM [yn (kM HM ) hn(2) ′ (kh HM ) −

150

T (M )

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 7. Displacement amplitude amplification factors along dam crest (−1 ≤ x/H ≤ 1) and along surface of half-space near dam (x/H < −1 & x/H > 1) for homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 2/3, 1/3) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 6. (M ) 2 T21 = ikh HM [−jn (kM HM ) jn′ (kh HM ) +

GM j ′ (kM HM ) jn (kh HM )], Gh n

(31)

(M ) 2 T22 = ikh HM [−yn (kM HM ) jn′ (kh HM ) +

GM y′ (kM HM ) jn (kh HM )]. Gh n

(32)

By successive application of the recurrence formula in Eq. (23) and the relationship in Eq. (20), Eq. (28) can be re-arranged to relate the coefficients for the half-space directly to those for the dam via

⎛⎜ 2n + 1⎟⎞ = T (M ) T (M − 1)⋅⋅⋅T (1) tA n D ⎠ ⎝ n

tˆ = tˆAn = ⎜⎛ 1 ⎟⎞ An ˆ ⎝t2 ⎠

Bn(1) = (2n + 1) t1/ tˆ1,

(36)

Cn(1) = (2n + 1) t2/ tˆ1.

(37)

Dn = (2n + 1) tˆ2/ tˆ1,

(38)

By means of Eqs. (35)–(38), all the other unknown coefficients Bn(m + 1) and Cn(m + 1) for m = 1 to M-1 can be obtained via Eq. (23) as

(33)

(m) (m) (m) (m) Bn(m + 1) = T11 Bn + T12 Cn , m = 1, …, M − 1,

(39)

(m) (m) (m) (m) Cn(m + 1) = T21 Bn + T22 Cn , m = 1, …, M − 1.

(40)

where tˆ denotes the global transfer vector (M ) (M − 1) (M ) (M − 1) (1) T12 T (1) ⎤ ⎛t1 ⎞ ⎤ ... ... ⎡ T11 ⎛⎜tˆ1 ⎞⎟ = ⎡ T11 T12 ⎤ ⎡ T11 [ ... ... ] ⎢ (1) 12 ⎜ ⎟ ⎢ (M ) (M ) ⎥ ⎢ (M − 1) (M − 1) ⎥ (1) ⎥ t2 ˆt2 T22 ⎝ ⎠ ⎣ T21 T22 ⎦ ⎣ T21 ⎣ T21 T22 ⎦ ⎝ ⎠ ⎦

With all the necessary coefficients found in closed form, the wave fields for the dam, the multi-layered canyon and the underlying halfspace can be calculated efficiently with suitable truncation of the series solution in Eqs. (3), (7), (11) and (12) at the Nth term. By checking the continuity conditions of both traction and displacement fields in Eqs. (13)–(18) as well as convergence tests, N ≥ 35 is proved to be more than sufficient for all the cases studied herein.

(34)

relating the response of the gorged half-space to that of the 1st canyon layer. Once the vector tˆ is assembled for a specific multi-layer configuration, the unknown coefficients An , Bn(1) , Cn(1) and Dn can be solved in closed form as

An = (2n + 1)/ tˆ1,

(35) 151

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

2 ⎧ ρd Vd* , ⎪ G (r ) = ρV02 r H ⎨ ⎪ ρV 2, ⎩ h

( )

r
, H ≤ r ≤ HM . r > HM

(42)

For the power-law wave velocity profile, representation of the smoothly inhomogeneous model by a M-layer piecewise homogeneous one is achieved at varying frequencies and directions of the incoming waves. Defining the dimensionless excitation frequency [15,20] as

Ω = ωH / Vd = 2πH / λ

with λ being the wavelength of the incident waves, a M-layer discretization of the smoothly inhomogeneous zone with M = 100 are proved to be sufficient for the frequency range Ω ≤ 12π studied herein. Figs. 2 and 3 show an example of convergence tests in terms of M corresponding to the case for Vh/Vd = 4, V0/Vd = 4/3, ρ/ρd = 1.5, βd = 10% and S/H = 1 at Ω = 12. As can be seen from Fig. 2, the shear wave velocities in the dam (r/H ≤ 1) and the half-space (r/H > 2) are kept constant while the inhomogeneous zone of the canyon (1 ≤ r/H ≤ 2) is divided into M homogeneous layers with M varying from 2, 20–100. Fig. 3 shows that the surface displacement amplitude amplification factors denoted as |AF| (i.e. displacement amplitudes normalized by the free-field motion amplitude) for M = 20 match exactly those for M = 100. This implies that even a 20-layer discretization of the inhomogeneous zone is enough for this case. For higher frequencies of incident waves, finer discretization is needed. To be conservative, M = 100 is adopted in the following calculations. Furthermore, the approximation of the power-law inhomogeneous canyon to the layered homogeneous model is very efficient. For example, the maximum computing time of different cases in Fig. 3 is only several minutes by using a personal computer (Intel®Core™i7–[email protected] GHz&8.00 GB with Windows 7 Professional 64 bit). The high efficiency of the proposed solution facilitates a systematical parametric study and a detailed wave concentration analysis in Section 6.

Fig. 8. Displacement amplitude amplification factors at midcrest of dam (x/H = 0) versus dimensionless frequency Ω for homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 2/3, 1/3) canyons with S/H = 1.

4. Representation of a power-law inhomogeneous canyon profile by the layered homogeneous model In soil dynamics, a soil layer with a continuously varying shear wave velocity or modulus is of fundamental relevance [32–34]. A power-law variation of the shear wave velocity or modulus in the inhomogeneous zone is one particular case of interest because such material variation of the soil or rock on a canyon surface can represent surficial weathering to a finite depth from water flows or other geological mechanisms [35]. The ground motion amplification by a semicylindrical canyon covered by a power-law inhomogeneous layer has recently been investigated by some of the authors [30]. To go a further step, this study aims at revealing the influence of the ground motion amplification by an inhomogeneous canyon on an earth and rockfill dam. Due to the general usefulness of the transmission matrix approach to the cylindrical multi-layered canyon, excellent representations of the smoothly inhomogeneous canyon models can be achieved by means of a sufficiently fine multi-layer discretization of the surficial zone and assigning the material properties to each layer according to the target continuous material profile. The shear wave velocity of the semi-annular inhomogeneous zone is thus taken to be a power function of the radial distance from the axis of cylindrical geometry of a surficial canyon as described in Eq. (1a) in Zhang et al. [30]. Therefore, the shear wave velocity profile for the maximum cross section of the dam-inhomogeneous canyon system is in the form of

r HM ⎩

(43)

5. Validation of solution If the inhomogeneity exponent ξ is set to be zero and the materials of the inhomogeneous zone and the half-space are set to be identical and uniform (i.e., V0 = Vm = Vh, ρm = ρh, and Gm = Gh), the presented solution can be degenerated to the canonical model of an earth dam in a homogenous semi-cylindrical canyon of Dakoulas [20]. To check the present solution, the results given in Fig. 7 of Dakoulas [20] are taken as validation examples. Fig. 4 shows the comparison of the present results with previous ones by Dakoulas [20] in terms of displacement amplitude amplification factors along dam crest and along surface of half-space near dam for the inhomogeneity exponent ξ = 0 (i.e. V0/Vh = 1), dimensionless thickness of the inhomogeneity zone S/H = 1, and material damping of dam βd = 10% at dimensionless frequency Ω = 4.5 under vertical (α = 0°), oblique (α = 30°, 60°) and horizontal (α = 90°) SH wave incidences. One can see from the results that the current solution agree very well with the referenced solution. Results for cases with other parameters agree equally well but are not presented here for brevity.

( )

(41)

where Vd* is the complex shear wave velocity in the dam, V0 = V(H) and Vh = V(HM) are the shear wave velocity at the upper and lower surface of the inhomogeneous zone, respectively. In this study, ξ = ln (V0/Vh)/ln(H/HM) parameterizes the rate of variation of the material inhomogeneity with respect to r with V0/Vh serving as an indicator of the overall variation of the material wave speed going from the top to the bottom of the inhomogeneous zone. Taking the mass densities ρm and ρh to be identical and constant (i.e., ρm = ρh = ρ, m = 1 to M) henceforth, the shear modulus of the dam-inhomogeneous canyon system is a function of r via

6. Illustrative results and discussions 6.1. Parametric analysis of surface displacement amplitude amplification factor With the preceding confirmation of the correctness of the damlayered canyon solution, a comprehensive parametric study of the effects of the inhomogeneous surficial zone of the canyon via the powerlaw model with respect to the parameters ξ, S/H and V0/Vh of the canyon surficial zone, the angle α of the incident wave, and the dimensionless frequency Ω is presented in this section. The material 152

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 9. Contour plots of displacement amplitudes across the maximum longitudinal section of dams in homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 0.5) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 6.

example, one can note that the maximum amplification for ξ = 1 and 0.5 can reach as high as 650% and 600%, respectively, while it is around 550% for ξ = 0 (Fig. 5c). Another inhomogeneity parameter S/H also has an important role in the dynamic response of dams. Fig. 6 shows that no additional amplification occurs for a dam in a canyon with a very thin inhomogeneity zone (S/H = 0.1) relative to the reference case of a dam in a homogeneous canyon (S/H = 0) at a dimensionless frequency Ω = 4.5. With the definition of Ω = 2πH/λ in Eq. (42), one can see that S/H = 0.1 in conjunction with Ω = 4.5 represents the case of a long wavelength λ that is equal to about 14 times the thickness of the inhomogeneity zone. Therefore, the inhomogeneity zone is too thin to have an influence on the wave amplification in this case. According to Zhang et al. [30], in fact, it also makes sense because the dynamic input of the dam will not be amplified due to the ignorable site amplification effect of ground motion by a very thin surficial inhomogeneity zone covering an empty canyon. As the dimensionless thickness of the surfacial inhomogeneity zone of the canyon increases, the site amplification of ground motion along the canyon surface becomes stronger and therefore the dam response becomes more intense for the cases of S/H = 1 and 2, relative to the cases of S/H = 0 and 0.1 (Fig. 6).

configuration can be described by ρh/ρd = ρm/ρd = 1.5 and Vh/Vd = 3 (i.e., ρhVh/ρdVd = 4.5). Note that the effect of the dam to half-space material impedance ratio, defined as IR = ρhVh/ρdVd, has been studied by Dakoulas [20] and therefore kept constant in the following results. To illustrate how the presence of the surficial inhomogeneous canyon zone can influence the anti-plane response of the dam, the surface displacement amplitudes amplification factors for a canyon with varying inhomogeneity material parameters ξ and V0/Vh or varying geometric parameter S/H under incident SH waves at four incident angles (α = 0°, 30°, 60° and 90°) are shown in Figs. 5–7. The hysteretic damping of the dam material is set to be βd = 10% for all three figures. Based on the results, it is evident that additional amplification to the response is caused by the inhomogeneity parameters ξ, S/H and V0/Vh, especially within the dam crest for all angles of the incoming waves. For a dimensionless frequency Ω of 3, Fig. 5 shows that the level of the amplification of dam response is improved by an increase of the inhomogeneity exponent ξ, relative to the case of ξ = 0 (a dam in a homogeneous canyon). There are noticeable differences in the maximum amplification factor of the dam crest for different degrees of material inhomogeneity of the canyon. Taking the case of α = 60° as an 153

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 10. Contour plots of displacement amplitudes across the maximum longitudinal section of dams in homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 0.5) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 12.

especially for the first resonances. Taking the case of βd = 10% as an example, the first resonances for the cases of homogeneous (V0/Vh = 1), slightly inhomogeneous (V0/Vh = 2/3) and heavily inhomogeneous (V0/Vh = 1/3) canyons are characterized by |AF| = 5.3 at Ω = 3, |AF| = 5.7 at Ω = 2.9 and |AF| = 7.2 at Ω = 2.5, respectively. Comparing Fig. 8(b) to Fig. 8(a), one can see a similar trend with sharper resonances and stronger amplifications of midcrest motions. From Figs. 5–7, one can find that the influence of the angle of wave incidence on the dam motion is generally strengthened by the presence of the inhomogeneous zone. As the angle α increases, both the amplitude increase and the location shift (to the right direction) of the maximum response on the dam crest become more substantially for the case of an inhomogeneous canyon. The additional amplification of the dam crest caused by the increase of the inhomogeneity degree of the canyon is the largest for the case of horizontal incidence of waves (α = 90°) and is relatively weak for the cases of vertical (α = 0°) and oblique (α = 30°, 60°) incidence. For example, the maximum displacement amplitude amplification factor can reach as high as 1500% for α = 90° in Fig. 7(d) while it is not larger than 370%, 450% and 580% for α = 0°, 30° and 60° in Fig. 7(a)-(c), respectively. According to Zhang et al. [30], such a huge difference in additional amplification by

In the case of a higher excitation frequency (e.g., Ω = 6 in Fig. 7), the difference caused by V0/Vh is also apparent. V0/Vh = 1, 2/3 and 1/3 represent the cases of a homogeneous canyon, a slightly inhomogeneous and heavily inhomogeneous one, respectively. It is clear that the maximum amplification factors along the dam crest tend to increase as the value of V0/Vh decreases. For example, the maximum amplification factor can reach as high as 1500% for the case of a dam in a heavily inhomogeneous canyon under a horizontal incidence of waves, while the reference value of maximum amplification factor is around 600% for the case of a homogeneous canyon (Fig. 7d). To illustrate how the inhomogeneous canyon can affect the resonant characteristics of the dam, Fig. 8 shows the midcrest amplification factors as a function of the dimensionless frequency Ω for homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 2/3, 1/3) canyons with S/ H = 1. Note that the frequency variation of the midcrest amplification factor (r = 0) is identical for all angles of incidence of waves which is in consistent with Figs. 3 and 4 in Dakoulas [20]. Fig. 8(a) and (b) correspond to the cases of dam material hysteretic damping of βd = 10% and 5%, respectively. For both cases, the resonance frequencies become lower and the corresponding amplifications of midcrest motion become more substantial as the inhomogeneity degree of the canyon improves, 154

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 11. Contour plots of displacement amplitudes across the maximum longitudinal section of dams in homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 0.5) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 18.

waves. From Figs. 9–11, one can see clearly the signs of wave concentration phenomenon, i.e., convergence of wave energy in localized regions and therefore high amplification of wave fields (e.g., Davis et al. [36], Semblat et al. [37], Tsaur and Chang [38]). Fig. 9(a)-(d) correspond to the homogeneous case of V0/Vh = 1. For this case, the maximum displacement amplitudes for α = 0°, 30°, 60° and 90° are 8.3 at (x/H, z/H) = (0, 0.61), 6.4 at (x/H, z/H) = (0.2, 0.6), 7.8 at (x/H, z/H )= (0.45, 0) and 12.2 at (x/H, z/H) = (0.55, 0), respectively. In contrast, for the inhomogeneous case in Fig. 9(e)-(h), the maximum displacements for α = 0°, 30°, 60° and 90° are 12.9 at (x/H, z/H) = (0, 0.64), 15.5 at (x/H, z/H) = (0.33, 0.63), 12.5 at (x/H, z/ H) = (0.56, 0) and 25.8 at (x/H, z/H) = (0.65, 0), respectively. This demonstrates significant influences of the angles of wave incidence on the dam response for both homogeneous and inhomogeneous configurations. On the one hand, the maximum displacement amplitudes at the focuses (i.e., locations of wave concentrating) for α = 90° are about twice those for α = 0°, 30° and 60°. On the other hand, the focuses tend to shift upward and/or rightward from the central axis of the dam as α increases from 0° to 90°. Comparing Fig. 9(e)-(h) to Fig. 9(a)-(d), one can find an apparent enhancement of the wave concentrating intensity due to the inhomogeneous zone of the canyon. Furthermore, the

the inhomogeneity zone between the horizontally and non-horizontally incidence cases does not exist for a surficially inhomogeneous canyon filled with no dam. This implies that the dynamic amplification of waves in a dam-canyon system must involve some mechanism different from that for a pure canyon. As is illustrated later, this is because the horizontally incident wave tends to concentrate on the dam crest while the vertically and obliquely incident waves concentrate below the dam crest. That is to say, the maximum response of the dam under vertically and obliquely incoming waves may not occur on the dam crest but within the dam body. 6.2. Wave concentration in the frequency domain Figs. 9–11 show the spatial distributions of displacement amplitudes across the maximum longitudinal section of dams (see Fig. 1(b) for reference) in homogeneous (V0/Vh = 1) and moderately inhomogeneous (V0/Vh = 0.5) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90° at Ω = 6, 12 and 18, respectively. Different from the surface motion amplification factor which is normalized by the free-field motion amplitude, the displacement amplitudes in Figs. 9–11 are normalized by the unit amplitude of incident SH 155

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 12. Contour plots of displacement amplitudes across the maximum longitudinal section of dams or alluvial valleys in a homogeneous canyon (V0/Vh = 1, S/H = 1) under incident angles α of 0°, 30°, 60° and 90° at Ω = 18.

damping reduces from βd = 10% to βd = 5%. Furthermore, similar patterns of spatial distributions of displacement amplitudes can be seen across the maximum longitudinal section of dams for the two cases. When it comes to the elastic case of βd = 0 in Fig. 12(e)-(h), the dam vibrates more intense than the previous two viscoelastic cases. It is interesting to see multiple positions of wave concentrations across the longitudinal section of the elastic dam in contrast to the single wave concentration in the viscoelastic dams. Since wave phenomena between a dam and an alluvial valley are connected [19], a comparison of wave concentration between them is valuable. Trifunac [39] proposed an analytical solution to the scattering of SH waves by a semi-circular alluvial valley. Since then, many other authors had extended his work to more complicated models such as a circular-arc alluvial valley in a elastic [40] or poroelastic [41] halfspace. Liang and his coworkers did a systematical study on scattering of SH, SV and P waves by circular-arc layered alluvial valleys [42–44]. All of these studies are relevant to the dam-inhomogeneous canyon model studied herein. Because the longitudinal section of the dam studied herein is semi-circular, it is appropriate to select a semi-circular alluvial valley for comparison. The results in Fig. 12(i)-(l) are calculated on the basis of the canonical solution for a semi-circular alluvial valley by Trifunac [39]. Similar to results of the semi-circular elastic dam (Fig. 12e-h), the wave also tends to concentrate at multiple locations in an elastic alluvial valley (Fig. 12i-l). The patterns of the spatial distribution of the displacement amplitudes for the elastic dam and the alluvial valley are comparable as well. A major difference is that the maximum responses of the elastic dam is much larger than those of the

locations of the focuses for the inhomogeneous case (V0/Vh = 0.5) move downward and/or rightward relative to the homogeneous case (V0/Vh = 1) for all angles of wave incidence. In the cases of higher excitation frequencies (Ω = 12 in Fig. 10; Ω = 18 in Fig. 11), the oscillation of the dam response with locations becomes more frequently but the change in maximum displacement amplitude is not obvious. As one can notice in Figs. 10 and 11, the maximum response within the dam generally takes place below the crest for the cases of α = 0°, 30° and 60° and on the dam crest only for the limiting case of α = 90°. This reflects the necessity to calculate the earthquake responses for the entire longitudinal section of the dam rather than only those on the dam crest. In fact, huge differences may exist between the dam responses on the crest and those at the underground focuses. For example, there is no obvious drop in maximum displacement amplitudes as the frequency Ω increases from 6 to 12 and even 18 as mentioned earlier. However, the midcrest responses may drop quickly within this frequency band based on the results in Fig. 8(a). In addition, one can notice that the influences of α and V0/Vh on dam responses for Ω = 12 in Fig. 10 and Ω = 18 in Fig. 11 are in consistence with those for Ω = 6 in Fig. 9. The results in Figs. 9–11 correspond to the case of a dam with a fixed hysteretic damping βd = 10%. To gain more insights into the wave concentration phenomenon, other cases for different βds are shown in Fig. 12. Fig. 12(a)-(d) correspond to a viscoelastic dam with βd = 5% in a homogeneous canyon, while Fig. 12(e)-(h) are for an elastic dam with βd = 0. Comparing Fig. 12(a)-(d) to Fig. 11(a)-(d), one can find, as expected, an obvious increase in the dam response as the

156

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 13. Synthetic seismograms along dam crest (−1 km ≤ x ≤ 1 km) and along surface of half-space near dam for homogeneous (V0/Vh = 1) and inhomogeneous (V0/Vh = 0.5) canyons with S/H = 1 under incident angles α of 0°, 30°, 60° and 90°.

with the characteristic frequency fc = 1.5 Hz and the calculated frequencies ranging from 0 Hz to 6 Hz with 1/16 Hz intervals. For illustration, the characteristic parameters of the dam are set to be H = 1 km, Vd = 1 km/s and βd = 10%. The geometric and material properties for the surficial zone of the canyon and the half-space are S/H = 1, Vh/Vd = 3 with V0/Vh = 0 (homogeneous canyon) or 0.5 (inhomogeneous canyon), and ρm/ρd = ρh/ρd = 1.5. Fig. 13 shows the synthetic seismograms of a dam-inhomogeneous canyon system (V0/Vh = 0.5) under different incident angles of α = 0°, 30°, 60° and 90° versus those of a dam-homogeneous canyon system (V0/Vh = 1). Each synthetic seismogram contains 81 time series received from equally spaced positions located along the dam crest (|x| ≤ 1 km), the surface of the canyon (1 km ≤ |x| ≤ 2 km) and the surface of the half-space (2 km ≤ |x| ≤ 4 km). Fig. 13(a) corresponds to the case

alluvial valley which may be due to the 3D nature of the dam model.

6.3. Wave concentration in the time domain To show more clearly the wave propagation and concentration patterns in and around the dam-inhomogeneous canyon system, timedomain results in the form of synthetic seismograms and snapshots are given in Figs. 13–17. They are obtained from the frequency-domain displacement transfer functions and a Fast Fourier transform algorithm (see, e.g., Zhang et al. [27]). The incident time signal considered is a Ricker wavelet

R (t ) = (2π 2fc2 t 2 − 1) e−π

2 2 2 fc t

(44) 157

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 14. Snapshots at 12 specified times for a dam-inhomogeneous canyon system with βd = 10%, S/H = 1 and V0/Vh = 0.5 under incident angles α = 0°.

For the case of α = 90°, the amplitude of the direct wave at about t = 3.4 s for the inhomogeneous case is about twice as much as the homogeneous case, which is consistent with the frequency domain results in Figs. 9–11. Figs. 14–17 show the wave fields at 12 specified times in terms of underground snapshots corresponding to the synthetic seismograms in Fig. 13(e)-(h), respectively. Each snapshot is for a mesh with 161 × 81 equally spaced receivers located within the rectangle of −4 km ≤ x ≤ 4 km by 0 km ≤ z ≤ 4 km. For the case of vertical incidence in Fig. 14, the Ricker-type direct wave marked by D1 begins to propagate upward with a speed of Vh = 3 km/s from the reference line at the depth of z = 5.75 km at the time of t = 0 s and arrives at z = 2 km when t = 1.25 s (Fig. 14a). The direct wave D1 shows up on the ground surface of the half-space (2 km ≤ |x| ≤ 4 km, z = 0 km) at around t = 1.92 s (Fig. 14b). In contrast, it takes more time for the direct wave D1 to reach the ground surface within the inhomogeneous zone of the canyon (1 km < r < 2 km) because of the decreasing wave velocity of the medium in the upward path, i.e., the wave velocity decreases

of a dam in a homogeneous canyon under a vertical incidence of waves at α = 0°. Due to the symmetry of the dam-canyon model, the results for this case are symmetrical. The wave front of the first arrival direct wave forms a straight line on the ground surface outside the dam (|x| > 1 km) at the time of about t = 2 s, while it becomes a curve along the dam crest (|x| ≤ 1 km) from t = 2 s to t = 2.5 s. The surface motions within the dam crest are stronger than those in the half-space near the dam. From Fig. 13(b)-(d) for α = 30°, 60° and 90°, respectively, one can see that the surface motions of the dam tend to increase as the incident angle α increases which confirms the conclusion of the frequency domain results in Dakoulas [20]. Fig. 13(e)-(h) correspond to the case of a dam in an inhomogeneous canyon. Comparing them to Fig. 13(a)-(d), respectively, one can find that both the direct waves and the later arrival waves are amplified by the presence of the inhomogeneous surficial zone of the canyon for all angles of wave incidence. For the cases of α = 0°, 30° and 60°, the later arrival waves on the dam crest are amplified more seriously by the inhomogeneous canyon and their amplitudes can be comparable with the direct waves. 158

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 15. Snapshots at 12 specified times for a dam-inhomogeneous canyon system with βd = 10%, S/H = 1 and V0/Vh = 0.5 under incident angles α = 30°.

the other part moves up toward the dam crest. Interestingly, the upgoing wave is reflected and converges again near the focus of R2 as shown in Fig. 14(l). Comparing Fig. 14(l) to Fig. 14(e), one can notice the similarity in the pattern of the twice wave concentration phenomena. The difference is the lower intensity of the second concentration. In addition to the wave concentrating process, snapshots in Fig. 14 can further explain the synthetic seismogram in Fig. 13(e). Outside the dam, the scattered waves S2 and S3 appear one after another on the surface of the half-space as shown in Fig. 13(e). Based on the wave propagation process outside the dam in Fig. 14(d)-(i), it is clear that the reflected waves R2 on the left and right sides of the half-space cross each other at about t = 3.00 s and pass the canyon shoulders just after t = 3.75 s. The scattered wave S2 is therefore generated at the canyon shoulders and reaches the ground surface at t = 4.00 s. From the wave process in Fig. 14(g)-(j), one can see the scattered wave S3 is induced when the dissipated waves expand from the location at around x = 0 km and z = 1.4 km at t = 3.50 s and meet the canyon shoulders at about t = 4.25 s. Another important feature in the synthetic

gradually with r in the form of V(r) = V0r = Vhr/2. Similarly, the direct wave D1 will not reach the dam midcrest until t = 2.50 s with a slower speed of Vd = 1 km/s (Fig. 14c). Besides the time delay, the wave velocity change in the upward path shortens the wave length of the Ricker wavelet as shown in Fig. 14(c). Note that the reflected wave R1 is generated when the direct wave D1 enters the dam through the damcanyon interface at r = 1 km (Fig. 14b and c). After the wave R2 reflected by the ground surface moves downward, scattered wave S1 is generated from both shoulders of the canyon (x = ± 1 km, z = 0 km) as annotated in Fig. 14(c). Fig. 14(d) and (e) demonstrate the process of the concentrating of the reflected wave R2 within the dam due to the redirection effect of the cylindrical boundary of the dam. The reflected wave R2 keeps converging toward the focus located at the lower half of the central axis of the dam (Fig. 14f). It is worth noting that the location of the wave concentration agrees well with the frequency domain results in Figs. 9(e)-11(e). After the convergence, the focus acts as secondary wave sources and generates new waves expanding outward. Fig. 14(g)-(k) show that one part of secondary waves dissipates into the half-space after passing through the inhomogeneous canyon zone and 159

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 16. Snapshots at 12 specified times for a dam-inhomogeneous canyon system with βd = 10%, S/H = 1 and V0/Vh = 0.5 under incident angles α = 60°.

that the location of the first wave concentration moves upward and rightward further, i.e., closer to the dam crest. This may explain the strengthening of the surface motions on the right side dam crest as α increases (Fig. 13). A major difference between the limiting case of α = 90° and the other three cases of α = 0°, 30° and 60° is that the horizontally travelling direct wave concentrates right on the dam crest. This is the reason why the direct wave D1 on the right side dam crest is so intense during the period from t = 3–4 s (see Fig. 13(h)). From Figs. 17(b) and 13(h), one can see the reflected wave R1 occurs when the direct wave D1 enters the dam horizontally at about t = 2.00 s. The direct wave D1 is redirected gradually and converges on the right side of the dam crest as shown in Fig. 17(b)-(d). The wave R2 reflected by the ground surface propagates along the dam bottom boundary (Fig. 17(e)) and induces the scattered wave on the left side of the horizontal ground surface marked by SL1 (Fig. 17(h)). The other two noticeable scattered waves SL2 and SL3 are generated similarly when the dissipated waves expand and pass the left shoulder of the canyon (Fig. 17(e)-(k)). These waves are also marked in Fig. 13(h). From Fig. 13(h), one can see several scattered

seismogram is that the amplitude of later arrival waves is comparable with that of the direct wave. From the wave process within the dam in Fig. 14(e)-(j), it is evident that the later arrival waves are associated with the wave concentration phenomenon. In fact, these waves are surface motions induced by the up-going waves generated from the focus. The wave concentration also happens in the case of α = 30° as shown in Fig. 15. Similar to the case of α = 0°, the first convergence of the energy occurs when the waves marked by R2 reflected by the horizontal ground surface are redirected by the velocity contrast on the dam-canyon interface. Comparing Fig. 15(e) to Fig. 14(e), it looks like the wave fields for α = 30° are rotated to the right side of the dam in an angle of 30° with respect to those for α = 0°. The location of the first focus for α = 30° thus moves toward the right direction from the central axis of the dam, which agrees well with the frequency domain results in Figs. 9–11. Accordingly, the second wave concentration happens at the left side of the dam as can be seen in Fig. 15(l). Fig. 16 shows the process of wave propagation and concentration for the case of α = 60° which is very similar to that for α = 30°. Note 160

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

Fig. 17. Snapshots at 12 specified times for a dam-inhomogeneous canyon system with βd = 10%, S/H = 1 and V0/Vh = 0.5 under incident angles α = 90°.

cylindrical canyon under obliquely incident plane SH waves is presented. Founded on a wave function expansion method and a transfer matrix approach, the multi-layered solution for the dam-inhomogeneous canyon problem is verified by its comparison with a past exact solution for a dam-homogeneous canyon problem. To demonstrate the amplified sensitivity of seismic ground motion in or near the dam to the local inhomogeneity of the canyon in a half-space, the case of a dam in a canyon with a surficial zone that is described by a smoothly inhomogeneous power-law shear-wave velocity model is developed. Through a systematic parametric study in the frequency domain, the critical importance to consider the effect of the degree of inhomogeneity with the incidence obliquity and frequency of the seismic wave in dynamic response modeling at or near a dam is illustrated. Wave concentration is found to be the major mechanism of high amplification of displacement within the dam. Time domain results demonstrate the process of wave propagating, scattering and concentrating in the dam-inhomogeneous canyon system. The high amplification of displacements at the wave focuses within the dam should warrant careful engineering attention in its anti-earthquake design and

waves show up on the right side of the dam crest after SL1-SL3. In fact, as shown in Fig. 17(k) and (l), this is due to a similar process of wave propagation in an opposite direction as that from t = 4–5 s (see Fig. 17(f)-(h)). One should also notice that the dam response is getting weaker and weaker as time goes by because of the penetration of waves into the half-space. Based on the foregoing time domain analysis in Figs. 13–17, it should be evident that the amplification and deamplification of dynamic responses as shown in Figs. 5–11 are the result of the complex interference and concentration of direct/reflected waves and scattered waves. The velocity variation in the inhomogeneous zone of the canyon and the velocity contrast on the dam-canyon interface can redirect the direct/reflected waves and converge energy at a localized area within the dam. 7. Conclusions In this paper, a rigorous elastodynamic formulation for an earth and rockfill dam in a radially piecewise-homogeneous multi-layered semi161

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

the Fundamental Research Funds for the Central Universities (Grant No. 2018B14014), the National Key Basic Research Program of China (“973” Program) (Grant No. 2015CB057901), the Public Service Sector R&D Project of Ministry of Water Resource of China (Grant No. 201501035-03), and the 111 Projects (Grant No. B13024).

reinforcement. Acknowledgements This work was supported by the National Natural Science Foundation of China (Grant No. 41630638，51608172 and 51479050), Appendix A. Notations used in this paper

H = Height of the dam, depth of the semi-circular canyon L = Length of the dam, width of the semi-circular canyon S = Thickness of the surface zone of the canyon Hm = The outer radii of the mth layer ρm, Gm, Vm = The mass density, shear modulus, and shear wave velocity of the mth layer ρh, Gh, Vh = The mass density, shear modulus, and shear wave velocity of the half-space ρd, Gd, Vd = The mass density, shear modulus, and shear wave velocity of the dam βd, Gd* , Vd* = The hysteretic damping ratio, the complex shear modulus, and the complex shear wave velocity of the dam V0 = The shear wave velocity on the surface of canyon V(r) = The shear wave velocity profile for the dam-inhomogeneous canyon system ξ = Material inhomogeneity exponent for the surface zone of the canyon km, kh, kd = The wave-numbers of the mth layer, the half-space and the dam wm, wh, wd = The displacement response of the mth layer, the half-space and the dam α = Incident angles of the plane SH waves ω = Circular frequency Ω = Dimensionless frequency jn (⋅) , yn (⋅) = The spherical Bessel functions of the first and second kinds of order n

(A.1) (A.2) (A.3) (A.4) (A.5) (A.6) (A.7) (A.8) (A.9) (A.10) (A.11) (A.12) (A.13) (A.14) (A.15) (A.16) (A.17)

hn(2) (⋅) = The spherical Hankel function of the second kind of order n T(m) = Transfer matrix M = Total number of layers in the surface zone of the canyon N = Truncation number of the series solution

(A.18) (A.19) (A.20) (A.21)

earthquake ground motion. J Eng Mech-ASCE 2001;127:932–9. [19] Dakoulas P, Hashmi H. Wave passage effects on the response of earth dams. Soils Found 1992;32:97–110. [20] Dakoulas P. Response of earth dams in semicylindrical canyons to oblique SH waves. J Eng Mech-ASCE 1993;119:74–90. [21] Dakoulas P, Hsu C-H. Response of dams in semielliptical canyons to oblique SH waves. J Eng Mech-ASCE 1995;121:379–91. [22] Alonso EE, Olivella S, Soriano A, Pinyol NM, Esteban F. Modelling the response of Lechago earth and rockfill dam. Géotechnique 2011;61:387–407. [23] Wang S. Treatment of highly weathered rock mass at the abutment of earth-rock dam in Qiaoqi hydropower project. Eng Sci 2005;7:78–83. [24] Trifunac MD. Scattering of plane SH waves by a semi-cylindrical canyon. Earthq Eng Struct Dyn 1973;1:267–81. [25] Wong HL, Trifunac MD. Surface motion of a semi-elliptical alluvial valley for incident plane SH waves. Bull Seismol Soc Am 1974;64:1389–408. [26] Gao Y, Zhang N, Li D, Liu H, Cai Y, Wu Y. Effects of topographic amplification induced by a U-shaped canyon on seismic waves. Bull Seismol Soc Am 2012;102:1748–63. [27] Zhang N, Gao Y, Cai Y, Li D, Wu Y. Scattering of SH waves induced by a nonsymmetrical V-shaped canyon. Geophys J Int 2012;191:243–56. [28] Lee VW, Liu W-Y. Two-dimensional scattering and diffraction of P- and SV-waves around a semi-circular canyon in an elastic half-space: an analytic solution via a stress-free wave function. Soil Dyn Earthq Eng 2014;63:110–9. [29] Sanchez-Sesma FJ, Crouse CB. Effects of site geology on seismic ground motion: early history. Earthq Eng Struct Dyn 2015;44:1099–113. [30] Zhang N, Gao Y, Pak RYS. Soil and topographic effects on ground motion of a surficially inhomogeneous semi-cylindrical canyon under oblique incident SH waves. Soil Dyn Earthq Eng 2017;95:17–28. [31] Abramowitz M, Stegun IA. Handbook of mathematical functions, with formulas, graphs, and mathematical tables. New York: Dover; 1972. [32] Towhata I. Seismic wave propagation in elastic soil with continuous variation of shear modulus in the vertical direction. Soils Found 1996;36:61–72. [33] Vrettos C. Simple inversion procedure for shallow seismic refraction in continuously nonhomogeneous soils. Soil Dyn Earthq Eng 1996;15:381–6. [34] Travasarou T, Gazetas G. On the linear seismic response of soils with modulus varying as a power of depth-the Maliakos marine clay. Soils Found 2004;44:85–93. [35] Alvarez-Rubio S, Sanchez-Sesma FJ, Benito JJ, Alarcon E. The direct boundary element method: 2D site effects assessment on laterally varying layered media (methodology). Soil Dyn Earthq Eng 2004;24:167–80. [36] Davis PM, Rubinstein JL, Liu KH, Gao SS, Knopoff L. Northridge earthquake damage caused by geologic focusing of seismic waves. Science 2000;289:1746–50.

References [1] Gazetas G. Seismic response of earth dams: some recent developments. Soil Dyn Earthq Eng 1987;6:2–47. [2] Gazetas G, Dakoulas P. Seismic analysis and design of rockfill dams: state-of-the-art. Soil Dyn Earthq Eng 1992;11:27–61. [3] Rampello S, Cascone E, Grosso N. Evaluation of the seismic response of a homogeneous earth dam. Soil Dyn Earthq Eng 2009;29:782–98. [4] Elia G, Amorosi A, Chan AHC, Kavvadas MJ. Fully coupled dynamic analysis of an earth dam. Géotechnique 2011;61:549–63. [5] Sharp MK, Adalier K. Seismic response of earth dam with varying depth of liquefiable foundation layer. Soil Dyn Earthq Eng 2006;26:1028–37. [6] Dakoulas P. Nonlinear seismic response of tall concrete-faced rockfill dams in narrow canyons. Soil Dyn Earthq Eng 2012;34:11–24. [7] Liu J, Liu F, Kong X, Yu L. Large-scale shaking table model tests on seismically induced failure of concrete-faced rockfill dams. Soil Dyn Earthq Eng 2016;82:11–23. [8] Arici Y, Ozel HF. Comparison of 2D versus 3D modeling approaches for the analysis of the concrete faced rock-fill Cokal Dam. Earthq Eng Struct Dyn 2013;42:2277–95. [9] Sica S, Pagano L. Perfomance-based analysis of earth dams: procedures and application to a sample case. Soils Found 2009;49:921–39. [10] Zou D, Xu B, Kong X, Liu H, Zhou Y. Numerical simulation of the seismic response of the Zipingpu concrete face rockfill dam during the Wenchuan earthquake based on a generalized plasticity model. Comput Geotech 2013;49:111–22. [11] Han B, Zdravkovic L, Kontoe S, Taborda DMG. Numerical investigation of the response of the Yele rockfill dam during the 2008 Wenchuan earthquake. Soil Dyn Earthq Eng 2016;88:124–42. [12] Chopra AK, Wang J-T. Earthquake response of arch dams to spatially varying ground motion. Earthq Eng Struct Dyn 2010;39:887–906. [13] Abouseeda H, Dakoulas P. Response of earth dams to P and SV waves using coupled FE-BE formulation. Earthq Eng Struct Dyn 1996;25:1177–94. [14] Abouseeda H, Dakoulas P. Non-linear dynamic earth dam-foundation interaction using a BE-FE method. Earthq Eng Struct Dyn 1998;27:917–36. [15] Papalou A, Bielak J. Seismic elastic response of earth dams with canyon interaction. J Geotech Geoenviron Eng-ASCE 2001;127:446–53. [16] Papalou A, Bielak J. Nonlinear seismic response of earth dams with canyon interaction. J Geotech Geoenviron Eng-ASCE 2004;130:103–10. [17] Seiphoori A, Haeri SM, Karimi M. Three-dimensional nonlinear seismic analysis of concrete faced rockfill dams subjected to scattered P, SV, and SH waves considering the dam-foundation interaction effects. Soil Dyn Earthq Eng 2011;31:792–804. [18] Chen M-T, Harichandran RS. Response of an earth dam to spatially varying

162

Soil Dynamics and Earthquake Engineering 116 (2019) 145–163

N. Zhang et al.

[41] Zhou X-L, Jiang L-F, Wang J-H. Scattering of plane wave by circular-arc alluvial valley in a poroelastic half-space. J Sound Vib 2008;318:1024–49. [42] Liang JW, Zhang Y, Gu X, Lee VW. Surface motion of circular-arc layered alluvial valley for incident plane SH waves. Chin J Geotech Eng 2000;22:396–401. [43] Liang JW, Yan L, Lee VW. Scattering of plane P waves by circular-arc layered alluvial valleys: an analytical solution. Acta Seismol Sin 2001;23:167–84. [44] Liang J, Yan L, Qin D, Lee VW. Dynamic response of circular-arc sedimentary valley site under incident plane SV waves. China Civ Eng J 2003;36:74–82.

[37] Semblat JF, Dangla P, Kham M, Duval AM. Seismic site effects for shallow and deep alluvial basins: in-depth motion and focusing effect. Soil Dyn Earthq Eng 2002;22:849–54. [38] Tsaur DH, Chang KH. Scattering and focusing of SH waves by a convex circular-arc topography. Geophys J Int 2009;177:222–34. [39] Trifunac MD. Surface motion of a semi-cylindrical alluvial valley for incident plane SH waves. Bull Seismol Soc Am 1971;61:1755–70. [40] Yuan X, Liao Z-P. Scattering of plane SH waves by a cylindrical alluvial valley of circular-arc cross-section. Earthq Eng Struct Dyn 1995;24:1303–13.

163