3D scattering of obliquely incident plane SV waves by an alluvial valley embedded in a fluid-saturated, poroelastic layered half-space

1.   Introduction

Earthquake investigations in the field show that local valley sites significantly influence the amplification of ground motion. For example, during the 1985 earthquake in Mexico, structures built on the alluvial basin suffered heavy damage, whereas the damage to buildings situated on bedrock on the outskirts of the Mexico City was significantly less.

Since the pioneering studies of Aki and Larner (1970) and Trifunac (1971), a great deal of work has been reported on wave scattering by valleys. The resulting problem can be solved analytically or numerically. Analytical methods (Aki and Larner 1970; Trifunac 1971; Wong and Trifunac 1974; Yuan and Liao 1995; Liang et al. 2003) have the advantage in that they explore the physical nature of the particular problem; however, numerical methods (Kawase and Aki 1989; Dravinski and Mossessian 1987; Sainchez-Sesma et al. 1993; Liang and Liu 2010; Ba and Liang 2011) are more suitable for dealing with problems that have complex geometries and varying material properties.

Although the majority of studies are still limited to a 2D case in a uniform half-space, few articles consider obliquely incident waves or a layered site. In spite of this, some studies have found that soil layers have important effects on both the amplitude and the spectrum of ground motion (Ba and Liang 2011a, b), and that the 3D wave scattering by an alluvial valley is fundamentally different from that of a 2D case (Liang et al. 2009, 2010). Articles (Liang et al. 2009, 2010; De Barros and Luco 1995; Khair et al. 1989, 1991; Liu et al. 1991) studied the 3D wave scattering problem by an infinitely long valley or valleys, but the solutions reported are still restricted to the elastic case.

In reality, soil is not homogeneously elastic but poroelastic and often saturated with fluid, especially in coastal areas. In this paper, we extend the method presented in article (Ba and Liang 2011) to apply it to a fluid-saturated, poroelastic medium to study the 3D scattering by a valley in a saturated layered site. By using Green's functions for moving loads as the fundamental solutions, the presented method yields very accurate results with low computational cost.

2.   Calculation model and method

The model used in this paper makes use of an infinitely long alluvial valley with an arbitrary cross section embedded in a site composed of saturated soil layers and underlying bedrock. The excitation is represented by plane SV waves, which incident obliquely from the bedrock. The incoming wave arrives with an azimuth θ_h relative to the axis of the valley (Fig. 1a) and an incidence θ_v to the vertical axis (Fig. 1b). The profile of the valley cross section is shown in Fig. 1c.

Figure  1.  a Horizontal plane, definition of angle θ_h, b vertical plane, definition of angle θ_v, and c cross section of the saturated valley

DownLoad: Full-Size Img PowerPoint

The total wave filed near the valley can be expressed as the summation of the free wave field and the scattered wave filed. The free wave field is the motion produced by the obliquely incident SV waves in the absence of the valley and the scattered wave field is induced by the alluvial valley. At first, the free wave field is calculated using the direct stiffness method. Then, the scattered wave fields inside and outside of the valley are simulated by the dynamic responses produced by moving fictitious distributed loads on the interface of the valley. The densities of the moving distributed loads can be determined through the continuity conditions associated with interface S. Finally, the total motion is obtained by adding the scattered wave filed to the free wave filed.

2.1   Biot's equation for a poroelastic medium

The dynamic equilibrium equations in a fluid-saturated, poroelastic medium (Biot 1956a, b) can be written as:

where U and w are, respectively, the solid frame's displacement vector and the displacement vector of the pore fluid with respect to the solid frame; w = n(W - U), W is the displacement vector of the pore fluid; ρ = (1 - n)ρ_s + nρ_f is the total density, ρ_s and ρ_f are, respectively, the solid grain density and the pore fluid density, n is the porosity; b is a parameter that accounts for internal fractionation; m = ρ₂₂/n² = (nρ_f + ρ_a)/n², ρ_a is the coupled density between the solid grain and the pore fluid; G and λ are the two Lame's constants corresponding to the solid frame, λ_c = λ + α²M; α and M are the Biot's parameters, which account for compressibility of a two phase material.

2.2   Free-field response

We use the direct stiffness method to calculate the free wave field. The key step within the method is to build an exact 3D dynamic stiffness matrix of the saturated layered site. As each saturated soil layer contains the upward traveling and downward traveling P1, P2, SV, and SH waves, the displacement and stress amplitudes at the top and bottom interfaces of each soil layer can be written as:

where U_x1, U_y1, U_z1, w_z1 and U_x2, U_y2, U_z2, w_z2 are the displacement amplitudes of the solid frame and of the pore fluid at the top and bottom interfaces of the saturated soil layer. τ_zx1, τ_zy1, σ_z1, p_f1 and τ_zx2, τ_zy2, σ_z2, p_f2 are the stress and pore pressure amplitudes at the top and bottom interfaces. A_P1, B_P1, A_P2, B_P2, A_SV, B_SV, A_SH, and B_SH are the amplitudes of the upward and downward traveling waves associated with P1, P2, SV, and SH waves, respectively. Eliminating the amplitudes of the upward and downward traveling waves in Eq. (2) results in a dynamic stiffness matrix of the 3D saturated soil layer, which can be written as [S_P1-P2-SV-SH^L]. For the saturated half-space, only the downward traveling waves exist, and the matrix of the saturated half-space [S_P1-P2-SV-SH^R] can be derived by letting A_P1 = A_P2 = A_SV = A_SH = 0. By introducing the continuity condition at the soil layer interfaces, we can obtain the global dynamic stiffness matrix [S_P1-P2-SV-SH] by assembling [S_P1-P2-SV-SH] for each soil layer and [S_P1-P2-SV-SH^R] for the half-space.

Once the 3D global dynamic stiffness matrix is obtained, the dynamic equation of the fluid-saturated poroelastic layered site can be expressed as

where {U} is the vector of the amplitude of displacement at the solid layer interfaces and {Q} is the external load vector.

For plane waves obliquely incident to the bedrock, the last four elements of vector {Q} can be determined by Eq. (4) with other elements equal to zero. Solving Eq. (3), the displacement at the interface of each soil layer is obtained. Then, the amplitudes of the upward and downward traveling waves in each soil layer can be determined by Eq. (2a). Finally, the displacement and stress amplitudes for any location can be obtained by Eqs. (2a) and (2b).

where U_x0, U_y0, U_z0, and w_z0 are the out-cropping motion. More details about the free wave field can be found in article (Ba 2008).

2.3   Green's functions for moving loads in a fluid-saturated, poroelastic layered half-space

The problem of 3D scattering of obliquely incident plane waves by a 2D valley is somewhat simpler than the full 3D problem. For any two cross sections, the wave field will be identical but phase shifted in the frequency domain. Therefore, dynamic responses of distributed loads (pore pressure) acting on inclined lines that move parallel to the y axis can be used to simulate the scattered wave field. The dynamic responses of moving distributed loads (pore pressure) are usually called moving Green's functions for moving loads.

The moving Green's functions are calculated in the following way. First, Green's function for moving distributed loads (pore pressure) acting on inclined lines in a layered site is derived. Then a Fourier integral along the axis of the valley (y axis) is used to calculate moving Green's functions. The moving velocity is c_s^R/(cosθ_hcosθ_v) from Fig. 1, where c_s^R is the S-wave velocity of the bedrock.

To derive Green's function for distributed loads (pore pressure) acting on inclined lines, the distributed loads are transformed into the wave-number domain by a Fourier transform. Then, calculations are performed to obtain Green's function for the wave-number domain. Finally, an inverse Fourier transform is used to transform Green's functions back into the space domain. In the wave-number domain, two fictitious interfaces are introduced to fix the distributed loads. Corresponding reaction forces are calculated, and then reaction forces with reversed directions are applied as loads on the saturated layered site to calculate the global response. The total response can be obtained by adding the results of the fixed layer and the results of the reversed reaction forces.

If [g_u(s)] and [g_t(s)] are the moving Green's functions for displacement and stress in a drained free surface, the resulting displacement (include the displacement of pore fluid with respect to the soil frame) and stress (include pore pressure) at the valley interface S can be expressed as

where the superscript 'p' in the displacement terms {U_x^p, U_y^p, U_z^p, w_n^p}^T and in the stress terms {t_x^p, t_y^p, t_z^p, t_f^p}^T indicates the results attributable to the moving fictitious distributed loads or pore pressure. The elements in the vector {p_x0, p_y0, p_z0, p_f0}^T are the fictitious distributed loads in the x-, y-, and z-directions and the fictitious pore pressure.

Similarly, if [g_uu(s)] and [g_ut(s)] are the Green's functions corresponding to an undrained free surface, the resulting displacement (include the displacement of pore fluid with respect to the solid frame) and stress (include pore pressure) at the valley interface S can be expressed as

2.4   Boundary conditions

If the free surface of the saturated layered site and the valley interface are both permeable (i.e., a drained boundary), the continuity conditions at the displacement interface (including the displacement of pore fluid with respect to the solid frame) and stress (including pore pressure) can be expressed as

with [W(s)] is the weighting function. If [W(s)] takes the value 1 in a single element and 0 in all the others, the integral can be evaluated over each element separately. Substituting Eqs. (5) and (6) into Eq. (9) gives:

where

where g_u^L(s) and g_u^V(s) are the Green's functions pertaining to the displacement of the layered saturated site and of the valley, g_t^L(s) and g_t^V(s) are the Green's functions pertaining to the stress of the layered saturated site and of the valley. Elements in the vector {p_1x, p_1y, p_1z, p_1f}^T and {p_2x, p_2y, p_2z, p_2f}^T are the distributed load and pore pressure applied on the valley interface to calculate the moving Green's functions for the layered half site and for the valley, respectively. Combining the solutions of Eqs. (10) and (3), the surface displacement for the valley can be written as:

If the free surface of the layered site and the interface of the valley are both impermeable (i.e., an undrained boundary), the boundary condition dictates that the displacement of the soil frame and the stress are continuous, but the displacement of pore fluid with respect to the solid frame is equal to zero. The displacement in the valley can then be obtained using the same method that was used to solve the drained boundary case.

3.   Verification of the method

When θ_h = 90° the present method will be reduced to the 2D case. We first compared our results with those of Liang and Liu (2010) for the 2D case of a semi-circular valley in a uniform saturated half-space. The material parameters of the saturated half-space and of the valley are shown in Table 1 for a porosity of 0.3; the superscript 'H' represents the half-space and the superscript 'V' represents the valley. The dimensionless frequency is defined as η = ωa/πc_β = 0.5, where a is the radius of the valley and c_β is the shear wave velocity of the saturated half-space. It is evident from Fig. 2 that our results are consistent with those of Liang and Liu (2010).

Table  1.  Parameters for the saturated layer, saturated half-space, and saturated valley (Lin et al. 2005)

 | Show Table

DownLoad: CSV

Figure  2.  Comparison our results with those of article (Liang and Liu 2010)

DownLoad: Full-Size Img PowerPoint

When ρ_f = 0 the presented method will reduce to the dry soil case. We also compared our results with those of De Barros and Luco (1995) for 3D scattering caused by a 2D semi-circular valley in an elastic layered site. The parameters are as follows: Poisson's ratio ν = 1/3, the dimensionless incident frequency η = ωa/π c_s^L = 0.5, the damping ratio ζ = 0.005, the ratio of the shear wave velocity of the valley to that of the half-space c_s^V/c_s^L = 0.5, the mass density ratio ρ_s^V/ρ_s^L = 2/3, the horizontal incident angle θ_h = 45°, and the vertical incident angle θ_v = 30°, 60°, and 90°. Figure 3 also shows that our results are relatively consistent with those of De Barros and Luco (1995).

Figure  3.  Comparison of our results with those of article (De Barros and Luco 1995)

DownLoad: Full-Size Img PowerPoint

4.   Numerical results

Figure 4 shows the surface displacement amplitudes of a uniform saturated half-space under drained boundary conditions for different obliquely incident angles. The parameters of the half-space and of the valley are shown in Table 1 for a porosity of 0.3. We define the dimensionless incident frequency as where ω is the incident circular frequency and a is the radius of the valley. Other parameters are as follows: η = 0.5 and 1.0, θ_v = 45°, θ_h = 0°, 30°, 60°, and 90°. |U_x|, |U_y|, and |U_z| in Fig. 4 are the solid frame displacement amplitudes and |A_SV| is the amplitude of the incident plane SV wave.

Figure  4.  Displacement amplitudes corresponding to different obliquely incident angels

DownLoad: Full-Size Img PowerPoint

Figure 4 shows that the incidence angle has a significant effect on the displacement amplitude. As the horizontal incident angle increases, the amplitudes in the x-direction increase gradually, while the amplitudes in the y-direction decrease gradually. This is probably due to an increase in the vibrational component in the x-direction as the incidence direction tends to the x axis.

Analysis of the results in Fig. 4 also indicates that the 3D dynamic response (θ_h ≠ 90°) is fundamentally different from the 2D case (θ_h = 90°). As a result, the obliquely incident SV waves should not simply be decomposed into the cross section of the valley and into the axis of the valley, and then calculated as the summation of the in-plane results and the anti-plane results.

Figure 5 shows the amplitudes of the surface displacement around a valley embedded in a uniform saturated half-space under undrained boundary conditions with porosities n = 0.1, 0.3, and 0.34. The material parameters used are given in Table 1, where superscript 'H' represents the half-space and superscript 'V' represents the valley. The dimensionless incident frequency η = 0.5 is defined as in Fig. 4, the horizontal incident angle is θ_h = 45°, and the vertical incident angle is θ_v = 5°, 30°, 60°, and 90°.

Figure  5.  Displacement amplitudes for different porosities

DownLoad: Full-Size Img PowerPoint

Figure 5 shows that porosity has a significant effect on wave scattering throughout the valley. Large differences are observed related to changes in porosity for the surface displacement amplitudes and large phase shifts. For low porosity (n = 0.1), the amplitudes corresponding to a saturated case are almost identical to those of a dry poroelastic case. The amplitudes of the horizontal displacement decrease gradually with an increase in porosity, which may be due to the smaller horizontal amplitudes of the free wave field corresponding to larger porosities. It can also be seen from Fig. 5 that the displacement amplitudes in the x- and z-directions are more complex than the amplitudes in the y-direction. This is probably because the valley is infinite in the y-direction.

To study the effect of a saturated soil layer, Fig. 6 depicts the simplest model, which consists of a single layer overlying elastic bedrock. The applied material parameters are given in Table 1 for a porosity of 0.3, with superscript "L" representing the saturated layer and superscript "V" representing the valley. The parameters of the elastic bedrock are as follows: mass density ρ^R = 1, 855 kN/m³, Poisson's ratio v^R = 0.25, damping ratio ζ^R = 0.02, and the ratio of the shear module of the bedrock to the layer is = 2.0 and 5.0. The thickness of the layer is H/a = 2.0. The dimensionless frequency is defined as η= with chosen values of η = 0.25 and 0.75. The horizontal incident angle is θ_h = 45°, the vertical incident angles are θ_v = 5°, 30°, 60°, and 90°.

Figure  6.  Displacement amplitudes for a saturated layered site with different soil layer rigidities

DownLoad: Full-Size Img PowerPoint

This study shows that surface displacement amplitudes can differ greatly between a layered site and a uniform half-space. This is because the dynamic responses of the layered case are determined by the dynamic characteristics of the valley and the soil layer, whereas the dynamic responses of the uniform case are determined solely by the valley. Furthermore, the differences between the uniform case and the layered case cannot be solely attributed to the differences in the free wave field; there is an interaction between the valley and the soil layer, which can be called the soil layer-valley interaction. When the incident frequency is equal to the fundamental frequency of the layered site (η = 0.25), the displacement amplitudes are significantly larger than those of other incident frequencies. It can also be seen from Fig. 6 that, although the shear module ratio between the bedrock and the saturated layer has changed, the distribution of the amplitudes still looks similar. This is most likely due to the unchanging self-dynamic characteristic of the layered site as the thickness of the soil layer is constant.

To study the influence of soil layer thickness, Fig. 7 shows the surface displacement amplitudes for additional soil layer thicknesses H/a = 1.0 and 4.0. The same model is used as in Fig. 6 with = 5.0 and η = 0.5. All other parameters are the same as those in Fig. 6.

Figure  7.  Displacement amplitudes for a saturated layered site with different soil layer thicknesses

DownLoad: Full-Size Img PowerPoint

A comparison of the results in Figs. 6 and 7 shows that soil layer thickness has a significant effect not only on amplitude but also on the distribution of surface displacement. For the same dimensionless incident frequency η = 0.5, the distribution of surface displacement is entirely different for different soil layer thicknesses H/a = 1.0, 2.0, and 4.0. This may be due to changes in the soil layer thickness, which result in a change in the layered site's self-vibrational characteristics. The displacement amplitudes decrease gradually with the increase of the soil layer thickness, and this is because material damping is considered in the presented method and the control point selected at the outcrop of the bedrock.

5.   Conclusions

This paper examines the 3D scattering of obliquely incident plane SV waves by an alluvial valley embedded in a fluid-saturated, poroelastic layered half-space. The direct stiffness method is used to calculate the free wave field, and Green's functions corresponding to the movement of distributed loads and pore pressure acting on inclined lines are used to simulate the scattered wave filed. Numerical results and analyses are performed by taking examples of a valley embedded in a uniform saturated half-space and in a single saturated soil layer overlying elastic bedrock. We conclude the following:

(1) The oblique incident angel has a significant effect on surface displacement amplitudes. In general, the surface displacement amplitudes in the x-direction increase gradually with an increase in the horizontal incident angle (the incident direction gradually tends to the x-direction).

(2) Porosity has an important effect on surface displacement amplitudes. Large differences in the displacement amplitudes and large phase shifts are observed for different porosities.

(3) Total dynamic responses are determined by the soil layer-valley interaction. The stiffness ratio of the bedrock to the saturated layer affects the amplitudes, whereas the soil layer thickness affects not only the amplitudes but also the spatial distribution of the surface displacement.