地球物理学报  2021, Vol. 64 Issue (8): 2766-2779   PDF    
流体饱和半空间二维地形三分量弹性波散射间接边界元模拟
梁建文1,2, 吴孟桃2, 巴振宁1,2     
1. 天津大学水利工程仿真与安全国家重点实验室, 天津 300350;
2. 天津大学土木工程系, 天津 300350
摘要:无限长局部地形地震波斜入射响应问题称为二维三分量问题,在计算量远小于三维的情况下,一定程度上反映了近地表场地的三维动力响应特征.基于天然土体的成层性及固液两相耦合特性,以层状多孔介质内部移动线荷载(孔隙水压)动力格林函数作为基本解,开展流体饱和半空间二维地形三分量弹性波散射的2.5维间接边界元模拟研究.总场响应由自由波场和散射波场叠加构成,前者可由直接刚度法求得,后者则通过施加移动虚拟均布荷载和移动虚拟孔隙水压所产生的动力响应来模拟.该方法优势在于离散仅限于地形底边界(无须离散自由地表),格林函数计算不存在奇异性(荷载可直接加在边界上),容易控制计算精度,对复杂边界条件具有很强的适应性.在退化验证和精度比较的基础上,以梯形凹陷和半椭圆沉积地形为例,模拟了时域和频域的流体饱和半空间三维弹性波散射响应.研究表明:局部地形的地震动响应依赖于入射频率、入射角度、边界透水条件、土层刚度和土层厚度等,入射波、反射波和散射波相互干涉,极大延长了位移的振动持续时间.
关键词: 流体饱和半空间      地震波传播      二维三分量      边界元法      移动线源     
IBEM simulation of three-component scattering of elastic waves in a fluid-saturated half-space with 2D topography
LIANG JianWen1,2, WU MengTao2, BA ZhenNing1,2     
1. State Key Laboratory of Hydraulic Engineering Simulation and Safety, Tianjin University, Tianjin 300350, China;
2. Department of Civil Engineering, Tianjin University, Tianjin 300350, China
Abstract: Seismic response of infinite local geometry subjected to obliquely incident waves is called a 2-D three-component problem, which reflects the 3D dynamic response characteristics of the near-surface site to some extent under the premise that calculation amount is much less than a 3D case. Considering the stratification and solid-liquid coupling properties of natural soil, a 2.5D IBEM (indirect boundary element method) is used to solve three-component scattering in a fluid-saturated, poroelastic layered half-space with 2D topography. In this IBEM dynamic the Green functions are used for moving forces (pore pressure) acting on an internal inclined line as its fundamental solutions. In the calculation, the total wave fields are decomposed into a free field and a scattered field, the former can be determined by the direct stiffness method, and the latter is simulated by the dynamic response of a moving fictitious distributed load as well as fictitious moving pore pressure applied to each discrete inclined line element. The advantage of this IBEM is that the discretization is only confined to the geometry boundary (no need to discretize the free surface), and no singularity in the calculation of Green's function (the load can be directly applied to the boundary), so it is easy to meet the calculation accuracy and complex boundary conditions. On the basis of degradation verification and precision comparison, taking a trapezoidal canyon and semi-elliptical valley as examples, the scattering effect of 3-D elastic waves in fluid-saturated half-space both in the frequency and time domains is simulated. The results show that the earthquake ground motion of the local topography depends on the incident frequency, incident angle, boundary drainage condition, soil layer stiffness and soil layer thickness. Due to the existence of topography, incident, reflected and scattered waves interfere with each other, which greatly prolongs the vibration duration of displacement.
Keywords: Fluid-saturated half-space    Seismic wave propagation    2D three-component    Boundary element method    Moving line source    
0 引言

震害调查和强震观测表明:局部地形(如凹陷、沉积和凸起)显著影响地震动的幅值及空间分布(张建毅等, 2012).而近年来,众多超高复杂建筑、重大交通工程和大型水利设施等不断出现,这些位于场地条件复杂地区的防震减灾工作亟需精确可靠的地震动参数.因此开展场地地震动局部地形效应研究具有重要的理论价值和巨大的社会效应.

局部场地地震效应,实质上是场地对弹性波的散射问题,近40年来,引起了国内外学者的重点关注.研究手段主要包括波函数展开解析解(Lee, 1984袁晓铭和廖振鹏, 1993)和各种数值求解方法,例如有限元法(景立平等, 2005;宋贞霞和丁海平, 2013)、有限差分法(廖振鹏等, 1981Zhang and Chen, 2006)、谱单元法(Komatitsch et al., 2004周红等, 2010)、边界元法(Fu and Wu, 2001Ge et al., 2005梁建文和巴振宁, 2007)、广义R/T递推传播矩阵法(Chen, 1990, 1995Ge and Chen, 2007)、离散波数法(Bouchon et al., 1989Fu and Bouchon, 2004)和各种逼近方法(Fu, 2005Hu et al., 2009Yu et al., 2010管西竹等, 2011Yu and Fu, 2012, 2014陈高祥等, 2017).

从地形散射机理来看,依据入射波方向与地形轴线的关系可分为二维散射(入射波场与地形垂直)和三维散射(入射波场与地形不垂直).前者自Trifunac(1971)开创性地给出半圆凹陷地形对均匀半空间SH波的散射后,Wong和Jennings(1975)Sánchez-Sesma和Campillo(1991)Yuan和Men(1992)杜修力和熊建国(1993)刘晶波(1996)Fu和Wu(2001)Fu等(2002)李伟华和赵成刚(2004)梁建文等(2006)Tsaur(2011)Gao等(2012)陈少林等(2014)Lee和Liu(2014)巴振宁等(2017)Shyu等(2018)高玉峰(2019)等相继对常规局部地形地震散射影响进行了大量深入研究.相对于二维散射,地形对三维地震波散射的研究成果还较少,主要有Sánchez-Sesma等(1983)Mossessian和Dravinski(1989)Sánchez-Sesma和Luzón(1995)Chávez-García(2003)赵成刚和韩铮(2007)梁建文等(2014)王铭锋等(2017)Ba和An(2019)等采用各类方法求解了均匀或层状半空间中局部地形的三维响应问题,获得了二维模型无法反映的复杂场地三维动力响应特征.三维模型与二维相比有着本质差别,能够全面反映地震动的空间分布特征,但对计算内存需求大、求解耗时长,尤其对于时域响应解.

在不失一般性的情况下,建立二维地形的三分量散射模型(即2.5维模型)是一种计算量远小于三维,并能在一定程度上反映场地三维动力响应特征的有效方法.Pedersen等(1994)采用基本解为全空间格林函数的间接边界元方法(IBEM)研究了二维沉积谷地的三维散射问题,但该方法需要离散自由地表.Stamos和Beskos(1996)采用直接边界元方法(DBEM)求解了二维凹陷地形的三分量弹性波散射,其解也仅适用于均匀半空间情况.Ba和Liang(2010a, b)采用基本解为半空间格林函数的IBEM求解了层状半空间凹陷地形的2.5维散射问题.以上研究仅限于单相介质,而实际地表土层很多是饱和多孔介质,尤其是滨海地区的沉积谷地中沉积层多为饱和土.根据Biot理论(Biot, 1956),波在流体饱和介质中的传播与弹性介质中的传播特征截然不同.目前关于流体饱和半空间二维地形的三维地震散射问题研究,仍鲜有报道.

本文将Ba和Liang(2010a, b)给出的方法进一步拓展到层状流体饱和半空间情况,推导建立了饱和土层及下卧基岩的三维精确动力刚度矩阵和移动线源下的动力格林函数,进而求解了流体饱和半空间中二维地形三分量弹性波散射问题.文中通过与已发表结果的比较验证了方法的可靠性和求解精度,并以梯形凹陷和半椭圆沉积地形为算例,在空间域内给出了一些时频数值结果,以期为揭示层状介质中局部地形地震波的三维散射机理提供一定参考.

1 理论模型及求解

图 1所示,流体饱和层状场地中存在一无限长任意截面形状局部地形,场地模型由水平成层的饱和土及其下的基岩半空间组成,局部地形与土层的交界面为S.地震体波(P、SV或SH波)由下卧面入射,入射波波阵面的法向量在水平面上的投影与y轴的夹角为θh(水平入射角),在竖直面内与z轴的夹角为θv(竖向入射角).

图 1 流体饱和半空间二维地形三分量弹性波散射模型示意图 Fig. 1 Sketch of three-component scattering of elastic waves by 2-D topography in a fluid-saturated half-space

由于弹性波斜入射下,半空间中局部地形沿其y轴的两个相邻截面的响应只存在相位差(因截面位置不同产生),那么可仅选取任一截面进行离散和求解,然后通过沿地形轴向的波数拓展求得其他截面位置的动力响应.本文采用以移动斜线荷载(孔隙水压)动力格林函数为基本解的2.5维IBEM进行问题求解.具体求解时,通过波场分离将整体模型分解为开口的层状半空间(波场1)和闭合的局部地形(波场2),并假定波场1中存在着自由波场和散射波场,波场2中仅存在散射波场.首先提出层状饱和场地自由波场(无地形存在)求解的精确动力刚度矩阵方法,然后通过在各离散斜线单元上施加一组沿轴线方向(y轴)移动的三坐标方向虚拟均布荷载和虚拟孔隙水压,求解位移和应力的格林影响函数来模拟散射波场.虚拟移动均布荷载和孔隙水压的密度,可通自由地表和地形的边界条件(位移、应力连续条件和透水条件)求得.对于时域响应则可通过对频域结果进行快速Fourier逆变换得到.

1.1 Biot理论及其一般解

根据文献(Biot, 1962),流体饱和多孔介质体系的土体本构方程和渗流连续方程可表示为

(1a)

(1b)

其中,i, j=x, y, zδij为克罗内克函数,εij为土骨架应变,σijpf分别为饱和介质总应力和孔隙水压,eξ分别为土骨架和流体的体积应变,μ*λ*为骨架材料的两个复Lamé常数.对两相多孔介质而言,阻尼作用不仅源于材料本身应力-应变关系,流-固黏性耦合效应还会导致额外的阻尼作用,进一步加剧能量的耗散.根据对应原理(Liang et al., 2006),土骨架的材料阻尼可以通过Lamé常数的复数形式引入:μ*=μ(1+2ζi),(λ+2μ)*=(λ+2μ)(1+2ζi),这里ζ为土体阻尼比,为虚数单位.

以土骨架位移U和流体相对位移w表示的饱和土动力控制方程为

(2a)

(2b)

其中,ρ=(1-n)ρs+fn为孔隙率,ρρsρf分别为两相土、土骨架和流体的质量密度;bmαM均为饱和相关参数.

函数f(x, y, z, t)对时间t以及水平坐标xy的Fourier变换及其逆变换可表示为

(3a)

(3b)

其中,上标“=”表示对空间分量(x, y)的Fourier变换,上标“~”表示对时间t的Fourier变换.对式(2)进行式(3a)的Fourier变换并进行求解,可得到土骨架位移以及流体相对位移在频率-波数域中解为

(4a)

(4b)

(4c)

(4d)

其中,, , , , AP1BP1AP2BP2ASVBSVASHBSH分别为上行和下行P1、P2、SV和SH波的幅值.kP1kP2kS分别为P1、P2和S波的波数,可由下式确定:

(5)

其中,

(6a)

(6b)

1.2 饱和介质三维精确动力刚度矩阵及自由波场求解

考虑上行和下行的波幅值,每一饱和土层上下表面处的土骨架位移和流体相对于土骨架位移可由式(7a)表示.结合式(1)和式(4),同时令饱和土层上表面处的外力

,令饱和土层下表面处的外力,则饱和土层上下表面的外力可由式(7b)表示:

(7a)

(7b)

其中,下标“1”、“2”分别表示饱和土层上下界面的相关变量.由式(7)可求得饱和土层上下表面处外力和位移之间的关系,即频率-波数域内饱和土层的精确动力刚度矩阵,记为.考虑基岩半空间内只存在下行波,同样可求得刚度矩阵.组集后即为频率-波数域内流体饱和半空间整体的三维刚度矩阵,这样可得到如下的层状饱和半空间的动力平衡方程:

(8)

(9)

对于由下卧基岩入射的弹性波,的最后四个元素可由式(9)求得,而其他元素为零,其中为地震波在基岩面处的输入位移.求解式(8)和式(9),由Fourier逆变换可得到频率-空间域内地形边界面各离散斜线单元上的土骨架位移和流体相对位移向量,以及各离散斜线单元上的应力和孔隙水压向量.

1.3 移动线荷载(孔隙水压)动力格林函数及散射波场模拟

流体饱和半空间层内部作用的移动斜线荷载(孔隙水压)如图 2所示,荷载作用层厚度为dpx0py0pz0分别为沿坐标xyz方向的外荷载幅值,pf0为孔隙水压幅值.该移动线源在时间-空间域中展开形式为:

(10)

图 2 流体饱和层中移动斜线均布荷载(孔隙水压) Fig. 2 Moving distributed forces (pore pressure) acting on an inclined line in fluid-saturated layers

其中,p0为线源幅值的一般表示,δ为狄拉克函数,θ(0°<θ<180°)为线源与x轴的夹角,线源的移动速度c应取为入射波沿y轴的视速度,(入射波为P1波)或(入射波为SV或SH波),cP1RcSR分别为下卧饱和半空间的P1波和S波波速.对式(10)采用式(3b)表示的Fourier变换,得到移动线源在频率-波数域中的表达式为

(11)

图 2,在节点1、2处引入附加交界面以使荷载作用介质层固定,总响应涉及的特解、齐解和反力解可由直接刚度法方便求得.流体饱和半空间动力格林函数的详细推导过程可参考层状弹性介质的思路(Ba and Liang, 2010a),不再赘述.

为频率-波数域中动力响应,并令ky=ω/c,则由Fourier逆变换可求得频率-空间域内的动力响应(位移或应力),即移动斜线荷载(孔隙水压)格林函数可写作:

(12)

利用该基本解,我们得到散射波场在地形边界S各单元上产生的土骨架位移以及流体相对位移、应力和孔隙水压计算公式为

(13a)

(13b)

其中,分别为位移和应力的移动荷载动力格林函数;下标“g”表示由三坐标方向移动斜线荷载和移动孔隙水压产生的动力响应;为在凹陷地形表面S各斜线单元上,为模拟散射波场而施加的三坐标方向的虚拟荷载密度.

1.4 边界条件及求解

求解中涉及到的边界条件,按地形类型可分为:(1)零应力边界条件;(2)位移和应力连续边界条件.按是否透水则有:(a)排水边界(透水情况);(b)不排水边界(不透水情况).注意到,对于透水情况,孔隙流体自由流动,孔隙水压为零;对于不透水情况,孔隙流体无法自由流动,流体较土骨架的相对位移为零.而无论是否透水,求解步骤均相同,以下给出一般性边界条件:

(1) 零应力边界条件(如凹陷地形,以“Sc”表示)

(14)

其中[W(s)]为权函数,一般取为单位矩阵.将式(13)代入式(14)得:

(15)

其中,

(16a)

(16b)

(2) 位移和应力连续边界条件(如沉积地形,以“Sv”表示)

(17a)

(17b)

同样,将式(13)代入式(17)得:

(18a)

(18b)

其中,

(19a)

(19b)

(19c)

(19d)

求解式(15)和(18),可得层状饱和半空间中局部地形附近地表任意位置的位移和应力(孔隙水压)幅值,对于时域动力响应可通过快速Fourier逆变换求得.

2 退化验证和精度比较

凹陷地形的二维散射是本文方法在入射波场与凹陷垂直时一种特殊情况.图 3呈现了本文方法与Liang等(2006)给出的流体饱和半空间中半圆凹陷地形对入射SV波散射结果的比较.地形和土层参数包括:n=0.375,μ=183 MPa,ρs=2650 kg·m-3ρf=1000 kg·m-3v=0.25,α=0.99,b=1.0 MPa·m-4m=1000 kg·m-3M=5600 MPa,ζ=0.01.计算中无量纲频率取为,其中a为凹陷半径,SV波入射角度ψ=60°.图中标注“透水/不透水”为自由地表和凹陷表面均为透水/不透水情况.由图可见,本文结果与Liang等(2006)给出结果吻合良好,从一方面验证了方法的正确性.

图 3 本文结果与Liang等(2006)给出结果的比较 Fig. 3 Comparison of the solutions given by the present study and Liang et al.(2006)

沉积地形的二维散射计算同样是本文方法的退化解,本文通过与Chen等(2008)给出的弹性半空间中半圆沉积谷地对入射SH波散射结果的比较来验证方法的正确性.计算中,令所有饱和相关参数取值为零即可退化为单相弹性介质情况.沉积层与土层的剪切刚度比为μL/μv=104,密度比为ρL/ρv=2/3,泊松比为vL=vv=1/3.无量纲频率定义为η=ωa/(csLπ),其中a为谷地半径,csL为土层波速.图 4给出了η=2.0时,不同入射角度(ψ=0°、60°)下沉积附近地表的幅值变化.由图可见,本文结果与Chen等(2008)的结果也十分一致,从另一方面验证了方法的正确性.

图 4 本文结果与Chen等(2008)给出结果的比较 Fig. 4 Comparison of the solutions given by the present study and Chen et al.(2008)

de Barros和Luco(1995)采用边界积分方程法首次给出弹性半空间中二维半圆沉积谷地三分量弹性波散射结果的数值解.在本文方法退化为单相介质情况的基础上,分别以均匀地形和层状地形对P波散射的求解模型为例,图 5比较了本文方法和文献方法求得的三个坐标方向地表位移幅值在不同斜入射角度下的计算结果.相关参数取值如下:αv=2βv=βLαL=2βLαR=2βR=4βLρv=(2/3)ρLρR=(4/3)ρLζR=ζL=ζv=0.005,其中αβ分别为压缩波波速和剪切波波速,上标“v”、“L”和“R”分别表示与沉积内、沉积外和下卧基岩有关的变量.无量纲频率取为η=ωa/πβL=0.5,地震波入射角度取为θh=0°和θv=30°.由图 5可见,无论在均匀介质还是在层状介质中,两种方法的计算结果均具有良好的一致性,较好地验证了本文方法的求解精度.

图 5 本文结果与De Barros和Luco(1995)给出结果的比较 Fig. 5 Comparison of the solutions given by the present study and Barros & Luco (1995)
3 数值算例 3.1 凹陷地形对弹性波的三维散射

为研究凹陷地形对弹性波的三维散射,以均匀饱和土层中梯形凹陷为算例模型,并假定斜入射波为SH波进行研究.材料参数如下:n=0.3、μ=37 MPa、ρs=2650 kg·m-3ρf=1000 kg·m-3M=60.72 MPa、m=7222 kg·m-3v=0.25、α=0.8287和b=0,不考虑材料阻尼.定义无量纲频率η=ωL/(πμ/ρs),其中L表示梯形凹陷的底边宽度.我们比较了地表透水和地表不透水时,水平方向±4L,地表下2L范围孔隙水压云图结果,如图 6图 7所示.图中竖向入射角度为θv=30°,水平入射角度为θh=0°和90°(0°代表垂直情况),入射频率为η=0.5、1.0和2.0.

图 6 地表透水时梯形凹陷附近的孔压变化 (a) 透水:η=0.5, n=0.3, θh=0°, θv=30°; (b) 透水:η=0.5, n=0.3, θh=90°, θv=30°; (c) 透水:η=1.0, n=0.3, θh=0°, θv=30°; (d) 透水:η=1.0, n=0.3, θh=90°, θv=30°; (e) 透水:η=2.0, n=0.3, θh=0°, θv=30°; (f) 透水:η=2.0, n=0.3, θh=90°, θv=30°. Fig. 6 Variation of pore pressure near the trapezoidal canyon with drained condition
图 7 地表不透水时梯形凹陷附近的孔压变化 (a) 不透水:η=0.5, n=0.3, θh=0°, θv=30°; (b) 不透水:η=0.5, n=0.3, θh=90°, θv=30°; (c) 不透水:η=1.0, n=0.3, θh=0°, θv=30°; (d) 不透水:η=1.0 n=0.3, θh=90°, θv=30°; (e) 不透水:η=2.0, n=0.3, θh=0°, θv=30°; (f) 不透水:η=2.0, n=0.3, θh=90°, θv=30°. Fig. 7 Variation of pore pressure near the trapezoidal canyon with undrained condition

结果显示,饱和排水情况(地表透水)和饱和不排水情况(地表不透水)的场地动力响应显著不同,且与入射波的频率、角度及观测点位置密切相关.总体上,饱和不排水时凹陷表面的孔隙水压幅值更大,原因在于孔隙流体的存在对波的传播有着耗能的作用,而随着观测点逐渐远离凹陷地形(x/a逐渐增大),两者的幅值差异逐渐减小,地表是否透水的影响减弱.随着入射频率的增大,半空间中孔隙水压幅值的变化及分布逐渐变得复杂.η=0.5(较低)时,幅值响应主要由其自由波场决定,η=2.0(较高)时,起主导作用的则为凹陷地形引起的附加散射波场.随着水平入射角度的变化,凹陷周围土层孔隙水压幅值的空间分布变化较为明显.θh=0°时,孔隙水压幅值关于凹陷轴线对称,θh=90°时,孔隙水压幅值在波入射侧较大,且分布较为复杂.

事实上,与经典的半圆和半椭圆凹陷地形相比(Wong and Jennings, 1975袁晓铭和廖振鹏, 1993Liang et al., 2006),梯形凹陷的土层和地形放大效应研究还较少,本文仅给出了一个均匀凹陷三维散射问题的简单结果,后续可针对成层场地中更多复杂截面形状凹陷地形开展研究.已有结果表明,对于给定的入射频率和入射角度,层状饱和半空间中凹陷附近动力响应由层状场地自身的动力特性(成层场地自由波场)和凹陷截面形状(凹陷地形形成的散射波场)共同决定(巴振宁和梁建文, 2013),而均匀饱和半空间情况仅取决于凹陷截面形状.因此为更为准确的确定凹陷地形附近的地震响应特征,需综合考虑凹陷附近的土层分布情况及其截面形状.

3.2 沉积地形对弹性波的三维散射

为研究沉积地形对弹性波的三维散射,以基岩上单一饱和土层中半椭圆沉积谷地为算例模型进行计算.材料参数如表 1所示,地表和沉积表面均假定为不透水边界.谷地深度与半径比为Rb/Ra=0.5,单一土层厚度为H/Ra=1.0.定义无量纲频率η=ωRa/(πcsL),其中csL表示饱和土层的波速.我们比较了斜入射P1波和SV波下,水平方向±2Ra,地表下1Ra范围孔隙水压云图结果,如图 8图 9所示.固定的入射角度θh(θv)=45°,变化的入射角度为θh(θv)=30°、60°和90°,入射频率为η=1.0.图中无量纲孔隙水压幅值定义为pf/kiAi,其中ki为入射P1或SV波的波数,Ai为入射P1或SV波的幅值.

表 1 基岩上单一饱和土层中半椭圆沉积谷地材料参数 Table 1 Material parameters of semi-elliptical alluvial valleys in a single saturated soil layer overlying the bedrock half-space
图 8 P1波入射下半椭圆沉积附近的孔压变化 (a) θv=45°, θh=30°;(b) θh=45°, θv=30°;(c) θv=45°, θh=60°;(d) θh=45°, θv=60°;(e) θv=45°, θh=90°;(f) θh=45°, θv=90°. Fig. 8 Pore pressure amplitudes around the semi-elliptical valley under P1 waves incidence
图 9 SV波入射下半椭圆沉积附近的孔压变化 (a) θv=45°, θh=30°;(b) θh=45°, θv=30°;(c) θv=45°, θh=60°;(d) θh=45°, θv=60°;(e) θv=45°, θh=90°;(f) θh=45°, θv=90°. Fig. 9 Pore pressure amplitudes around the semi-elliptical valley under SV waves incidence

结果显示,沉积地形的存在(图中虚线位置)显著改变了饱和层状半空间中孔隙水压振幅的分布.对于入射P1波,孔隙水压幅值的空间分布十分依赖于竖向入射角θv,但对水平入射角θh相对不敏感,其幅值大小似乎随着竖向入射角的增大而逐渐增大.孔隙水压幅值的最大值出现在靠近谷地与土层交界面的右侧,这意味着这些区域具有较大的应变响应.当θv=90°(P1波垂直入射)时,幅值的分布关于谷地轴线对称,响应与水平角度无关.对于入射SV波,水平和竖向角度都对孔隙水压的变化(尤其是幅值)有着显著影响.比较图 8图 9,可以明显看出SV波入射下谷地内部孔隙压力的空间分布更复杂,幅值也更大.其原因在于不同类型入射波幅值的差异主要是自由波场的差异,孔隙流体承担传递来的激励能量改变,孔隙水压自然不同.同样,当θv=90°(SV波垂直入射)时,幅值的分布具有对称性,这也为我们的结果提供了额外的信心.

为研究沉积地形在地震波作用下的时域响应,分别以均匀饱和半空间中半椭圆沉积谷地(图 10)和基岩上单一饱和土层中半椭圆沉积谷地(图 11abcd)为例,求解了Ricker时程(SV波形式)斜入射下沉积地形附近的地表位移幅值响应.入射Ricker时程U(τ)=ASV(2π2ηc2τ2-1)exp(-π2ηc2τ2),特征频率ηc定义为ηc=ωa/(πβ),无量纲时间定义为τ=/a,其中,a为沉积半径,β为均匀半空间或单一饱和土层的剪切波速,ASV为入射SV波的幅值.下文分析中,均取特征频率ηc=1.0,无量纲时间τ=0~8,无量纲计算频率范围η=0~5,共115个频率点(利用分段高斯积分提高积分精度).Ricker子波的斜入射角度取θh=45°和θv=45°,沉积地形及附近表面所有观测点均匀分布在-3ax≤3a范围内.

图 10 均匀饱和半空间中沉积附近地表位移幅值时域结果 Fig. 10 Time domain solutions for surface displacement amplitude around the valley in a uniform half-space
图 11 层状饱和半空间中沉积附近地表位移幅值时域结果 (a) , H/a=2.0; (b) , H/a=2.0; (c) , H/a=1.0; (d) , H/a=4.0. Fig. 11 Time domain solutions for surface displacement amplitude around the valley in a layered half-space

图 10展示了均匀半空间情况沉积附近地表位移时域响应,图中|Ux/ASV|、|Uy/ASV|和|Uz/ASV|分别表示三个坐标方向的无量纲位移.可以看到,沉积地形的存在使得谷地附近的位移幅值与自由波场(无沉积地形存在)存在明显差异,且由于局部不规则形成的反射波和散射波与入射波相互干涉,使得谷地附近位移幅值的持续时间明显长于自由波场情况.另外,SV波在整个半空间的传播过程在时域内可更为清楚的观察到,如谷地左侧(x/a≤-1.0)首先接收到直达波,随后接收到谷地左角点(x/a=-1.0)产生的散射波,最后接收到谷地右角点(x/a=1.0)产生的散射波.同时,入射波在较软的沉积介质中经过多次反射,导致谷地内部位移相较于外部的持续时间明显更长,且幅值更大.

图 11进一步示出了层状半空间情况对应的位移时程,其中,下卧基岩与土层的剪切刚度比取为和5.0,单一土层厚度考虑为H/a=1.0、2.0和4.0.从图中可以发现:1)总体上,与均匀半空间相比,单一饱和土层沉积地形附近位移时程的空间分布更为复杂,且位移幅值整体上大于均匀半空间情况,另外从图中还可以看出,由于下卧基岩层的存在,层状半空间中不同相速度的散射波相互作用,位移的振动持续时间也明显长于均匀半空间情况.2)比较图 11a11b发现,土层厚度相同时,随着剪切刚度比的增加,地表位移逐渐增大,两个水平方向波的到达延时缩短.这是由于波从下卧基岩折射到土层时,波的入射角度随剪切刚度比的增加而增大所致.3)比较图 11b11c11d发现,剪切刚度比相同时,随着土层厚度的增加,地表位移逐渐减小,总体到时逐渐延长.其原因在于,由于阻尼的作用(波在土层中的往复次数较无阻尼时显著减少),波的传播路径的增加,使得散射波逐渐衰减、到达延时变长.

4 结论

基于Biot理论,建立了层状多孔介质内部移动线荷载(孔隙水压)格林函数,进而提出了求解流体饱和半空间中二维局部地形三分量弹性波散射问题的间接边界元方法(IBEM).该方法优势在于离散仅限于地形边界(无须离散自由地表),且格林函数计算不存在奇异性(荷载直接加在边界上),因而计算精度容易控制,且对复杂边界条件具有很强的适应性.

文中以给出的方法具体研究了梯形凹陷地形和半椭圆沉积地形的三维地震效应,分别在频域和时域内进行了计算分析.频域结果表明,局部地形的存在显著改变了饱和层状半空间中孔隙水压的幅值及空间分布,且依赖于地震波类型、入射频率、入射角度和边界透水条件等影响因素.地形附近动力响应由半空间动力特性和地形截面形状共同决定.时域结果表明由于局部不规则形成的反射波和散射波与入射波相互干涉,使得地形附近位移幅值的持续时间明显长于自由波场情况.与均匀半空间相比,基岩的存在的位移时程的空间分布更为复杂,振动时长增加.层状半空间,土层自身特性(土层刚度、厚度等)对波的传播路径和持续时间有重要影响.

本文建立的以移动斜线荷载(孔隙水压)格林函数为基本解的IBEM在近地表复杂场地对地震波三维散射和土-结构相互作用方面有着重要应用前景.后续研究将进一步考虑有实际观测数据的场地,以期为区域地震危险性分析、烈度区划以及局部地形附近建筑物的抗震设防等工作提供一定理论基础.

References
Ba Z N, Liang J W. 2010a. 2.5 D scattering of incident plane SH waves by a canyon in layered half-space. Earthquake Science, 23(1): 25-33. DOI:10.1007/s11589-009-0085-3
Ba Z N, Liang J W. 2010b. 2.5 D scattering of incident plane SV waves by a canyon in layered half-space. Earthquake Engineering and Engineering Vibration, 9(4): 587-595. DOI:10.1007/s11803-010-0040-2
Ba Z N, Liang J W. 2013. Three-dimensional scattering by a two-dimensional canyon in a fluid-saturated, poroelastic layered half-space for obliquely incident plane P1 waves. China Civil Engineering Journal (in Chinese), 46(S2): 166-171.
Ba Z N, Huang D Y, Liang J W, et al. 2017. Scattering and diffraction of plane SH-waves by periodically distributed hill topographies. Chinese Journal of Geophysics (in Chinese), 60(3): 1039-1052.
Ba Z N, An D H. 2019. Seismic response of a 3-D canyon in a multilayered TI half-space modelled by an indirect boundary integral equation method. Geophysical Journal International, 217(3): 1949-1973. DOI:10.1093/gji/ggz122
Biot M A. 1956. Theory of propagation of elastic waves in a fluid-saturated porous solid. I. Low-frequency range. The Journal of the Acoustical Society of America, 28(2): 168-178. DOI:10.1121/1.1908239
Biot M A. 1962. Mechanics of deformation and acoustic propagation in porous media. Journal of Applied Physics, 33(4): 1482-1498. DOI:10.1063/1.1728759
Bouchon M, Campillo M, Gaffet S. 1989. A boundary integral equation-discrete wavenumber representation method to study wave propagation in multilayered media having irregular interfaces. Geophysics, 54(9): 1134-1140. DOI:10.1190/1.1442748
Chávez-García F J. 2003. Site effects in Parkway Basin: comparison between observations and 3-D modelling. Geophysical Journal International, 154(3): 633-646. DOI:10.1046/j.1365-246X.2003.02055.x
Chen G X, Fu L Y, Yu G X, et al. 2017. A quantitative analysis method for the seismic geological complexity of near surface. Chinese Journal of Geophysics (in Chinese), 60(3): 1062-1072. DOI:10.6038/cjg20170319
Chen J T, Chen P Y, Chen C T. 2008. Surface motion of multiple alluvial valleys for incident plane SH-waves by using a semi-analytical approach. Soil Dynamics and Earthquake Engineering, 28(1): 58-72. DOI:10.1016/j.soildyn.2007.04.001
Chen S L, Zhang L L, Li S Y. 2014. Numerical analysis of the plane SH waves scattering by semi-cylindrical alluvial valley. Engineering Mechanics (in Chinese), 31(4): 218-224.
Chen X F. 1990. Seismogram synthesis for multi-layered media with irregular interfaces by global generalized reflection/transmission matrices method. I. Theory of two-dimensional SH case. Bulletin of the Seismological Society of America, 80(6A): 1696-1724.
Chen X F. 1995. Seismogram synthesis for multi-layered media with irregular interfaces by global generalized reflection/transmission matrices method. Ⅱ. Applications for 2D SH case. Bulletin of the Seismological Society of America, 85(4): 1094-1106.
De Barros F C P, Luco J E. 1995. Amplification of obliquely incident waves by a cylindrical valley embedded in a layered half-space. Soil Dynamics and Earthquake Engineering, 14(3): 163-175. DOI:10.1016/0267-7261(94)00047-K
Du X L, Xiong J G, Guan H M. 1993. The boundary integration equation methods in scattering of plane SH waves. Acta Seismologica Sinica (in Chinese), 15(5): 331-338. DOI:10.1007%2FBF02650400
Fu L Y, Wu R S. 2001. A hybrid BE-GS method for modeling regional wave propagation. Pure and Applied Geophysics, 158(7): 1251-1277. DOI:10.1007/PL00001222
Fu L Y, Wu R S, Campillo M. 2002. Energy partition and attenuation of regional phases by random free surface. Bulletin of the Seismological Society of America, 92(5): 1992-2007. DOI:10.1785/0120000292
Fu L Y, Bouchon M. 2004. Discrete wavenumber solutions to numerical wave propagation in piecewise heterogeneous media-I. Theory of two-dimensional SH case. Geophysical Journal International, 157(2): 481-498. DOI:10.1111/j.1365-246X.2004.02135.x
Fu L Y. 2005. Rough surface scattering: comparison of various approximation theories for 2D SH waves. Bulletin of the Seismological Society of America, 95(2): 646-663. DOI:10.1785/0120040081
Gao Y, Zhang N, Li D, et al. 2012. Effects of topographic amplification induced by a U-shaped canyon on seismic waves. Bulletin of the Seismological Society of America, 102(4): 1748-1763. DOI:10.1785/0120110306
Gao Y F. 2019. Analytical models and amplification effects of seismic wave propagation in canyon sites. Chinese Journal of Geotechnical Engineering (in Chinese), 41(1): 1-25.
Ge Z X, Fu L Y, Wu R S. 2005. P-SV wave-field connection technique for regional wave propagation simulation. Bulletin of the Seismological Society of America, 95(4): 1375-1386. DOI:10.1785/0120040173
Ge Z X, Chen X F. 2007. Wave propagation in irregularly layered elastic models: a boundary element approach with a global reflection/transmission matrix propagator. Bulletin of the Seismological Society of America, 97(3): 1025-1031. DOI:10.1785/0120060216
Guan X Z, Fu L Y, Tao Y, et al. 2011. Boundary-volume integral equation numerical modeling for complex near surface. Chinese Journal of Geophysics (in Chinese), 54(9): 2357-2367. DOI:10.3969/j.issn.0001-5733.2011.09.019
Hu S Z, Fu L Y, Yao Z X. 2009. Comparison of various approximation theories for randomly rough surface scattering. Wave Motion, 46(5): 281-292. DOI:10.1016/j.wavemoti.2009.03.001
Jing L P, Zhuo X Y, Wang X J. 2005. Effect of complex site on seismic wave propagation. Earthquake Engineering and Engineering Vibration (in Chinese), 25(6): 16-23.
Komatitsch D, Liu Q, Tromp J, et al. 2004. Simulations of ground motion in the Los Angeles basin based upon the spectral-element method. Bulletin of the Seismological Society of America, 94(1): 187-206. DOI:10.1785/0120030077
Lee V W. 1984. Three-dimensional diffraction of plane P, SV & SH waves by a hemispherical alluvial valley. International Journal of Soil Dynamics and Earthquake Engineering, 3(3): 133-144. DOI:10.1016/0261-7277(84)90043-3
Lee V W, Liu W Y. 2014. 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 Dynamics and Earthquake Engineering, 63: 110-119. DOI:10.1016/j.soildyn.2014.02.005
Li W H, Zhao C G. 2004. Scattering of plane SV waves by circular-arc alluvial valleys with saturated soil deposits. Chinese Journal of Geophysics (in Chinese), 47(5): 911-919.
Liang J W, Zhang Y S, Lee V W. 2006. Surface motion of a semi-cylindrical hill for incident plane SV waves: analytical solution. Acta Seismologica Sinica (in Chinese), 28(3): 238-249.
Liang J W, You H B, Lee V W. 2006. Scattering of SV waves by a canyon in a fluid-saturated, poroelastic layered half-space, modeled using the indirect boundary element method. Soil Dynamics and Earthquake Engineering, 26(6-7): 611-625. DOI:10.1016/j.soildyn.2006.01.012
Liang J W, Ba Z N. 2007. Surface motion of an alluvial valley in layered half-space for incident plane SH waves. Journal of Earthquake Engineering and Engineering Vibration (in Chinese), 27(3): 1-9.
Liang J W, Qi X Y, Ba Z N. 2014. Seismic response analysis of 3D canyon based on the viscous-spring boundary. Earthquake Engineering and Engineering Dynamics (in Chinese), 34(4): 21-28.
Liao Z P, Yang B P, Yuan Y F. 1981. Effects of three-dimensional topography on earthquake ground motion. Earthquake Engineering and Engineering Vibration (in Chinese), 1(1): 56-77.
Liu J B. 1996. The effect of local irregular topography on seismic ground motion. Acta Seismologica Sinica (in Chinese), 18(2): 239-245.
Mossessian T K, Dravinski M. 1989. Scattering of elastic waves by three-dimensional surface topographies. Wave Motion, 11(6): 579-592. DOI:10.1016/0165-2125(89)90028-0
Pedersen H A, Sanchez-Sesma F J, Campillo M. 1994. Three-dimensional scattering by two-dimensional topographies. Bulletin of the Seismological Society of America, 84(4): 1169-1183.
Sánchez-Sesma F J, Campillo M. 1991. Diffraction of P, SV, and Rayleigh waves by topographic features: a boundary integral formulation. Bulletin of the Seismological Society of America, 81(6): 2234-2253.
Sánchez-Sesma F J, Luzón F. 1995. Seismic response of three-dimensional alluvial valleys for incident P, S, and Rayleigh waves. Bulletin of the Seismological Society of America, 85(1): 269-284.
Sánchez-Sesma F J, Pérez-Rocha L E, Chávez-Pérez S. 1983. Diffraction of elastic waves by three-dimensional surface irregularities. Part Ⅱ. Bulletin of the Seismological Society of America, 73(6A): 1621-1636.
Shyu W S, Teng T J, Chou C S. 2018. Effect of geometry on in-plane responses of a symmetric canyon subjected by P waves. Soil Dynamics and Earthquake Engineering, 113: 215-229. DOI:10.1016/j.soildyn.2018.06.003
Stamos A A, Beskos D E. 1996. 3-D seismic response analysis of long lined tunnels in half-space. Soil Dynamics and Earthquake Engineering, 15(2): 111-118. DOI:10.1016/0267-7261(95)00025-9
Trifunac M D. 1971. Surface motion of a semi-cylindrical alluvial valley for incident plane SH waves. Bulletin of the Seismological Society of America, 61(6): 1755-1770.
Tsaur D H. 2011. Scattering and focusing of SH waves by a lower semielliptic convex topography. Bulletin of the Seismological Society of America, 101(5): 2212-2219. DOI:10.1785/0120100324
Wang M F, Zheng A, Zhang W B. 2017. Effect of local mountain topography on strong ground motion. Chinese Journal of Geophysics (in Chinese), 60(12): 4655-4670. DOI:10.6038/cjg20171210
Wong H L, Jennings P C. 1975. Effects of canyon topography on strong ground motion. Bulletin of the Seismological Society of America, 65(5): 1239-1257.
Yu G X, Fu L Y, Yao Z X. 2010. Comparison of different BEM+Born series modeling schemes for wave propagation in complex geologic structures. Geophysics, 75(3): T71-T82. DOI:10.1190/1.3392175
Yu G X, Fu L Y. 2012. Rytov series approximation for rough surface scattering. Bulletin of the Seismological Society of America, 102(1): 42-52. DOI:10.1785/0120110083
Yu G X, Fu L Y. 2014. Convergence analyses of different modeling schemes for generalized Lippmann-Schwinger integral equation in piecewise heterogeneous media. Soil Dynamics and Earthquake Engineering, 63: 150-161. DOI:10.1016/j.soildyn.2014.03.004
Yuan X M, Liao Z P. 1993. Series solution for scattering of plane SH waves by a canyon of circular-arc cross section. Earthquake Engineering and Engineering Vibration (in Chinese), 13(2): 1-11.
Yuan X M, Men F L. 1992. Scattering of plane SH waves by a semi-cylindrical hill. Earthquake Engineering & Structural Dynamics, 21(12): 1091-1098.
Zhang J Y, Bo J S, Wang Z Y, et al. 2012. Influence of local topography on seismic ground motion in Wenchuan earthquake. Journal of Natural Disasters (in Chinese), 21(3): 164-169.
Zhang W, Chen X F. 2006. Traction image method for irregular free surface boundaries in finite difference seismic wave simulation. Geophysical Journal International, 167(1): 337-353. DOI:10.1111/j.1365-246X.2006.03113.x
Zhao C G, Han Z. 2007. Three-dimensional scattering and diffraction of plane Rayleigh-waves by a hemispherical alluvial valley with saturated soil deposits. Chinese Journal of Geophysics (in Chinese), 50(3): 255-264.
Zhou H, Gao M T, Yu Y X. 2010. A study of topographical effect on SH waves. Progress in Geophysics (in Chinese), 25(3): 775-782.
巴振宁, 梁建文. 2013. 平面P1波斜入射下层状饱和半空间中凹陷地形周围的三维散射. 土木工程学报, 46(S2): 166-171.
巴振宁, 黄棣旸, 梁建文, 等. 2017. 层状半空间中周期分布凸起地形对平面SH波的散射. 地球物理学报, 60(3): 1039-1052.
陈高祥, 符力耘, 于更新, 等. 2017. 近地表地震地质复杂性的一种定量分析方法. 地球物理学报, 60(3): 1062-1072. DOI:10.6038/cjg20170319
陈少林, 张莉莉, 李山有. 2014. 半圆柱型沉积盆地对SH波散射的数值分析. 工程力学, 31(4): 218-224.
杜修力, 熊建国, 关慧敏. 1993. 平面SH波散射问题的边界积分方程分析法. 地震学报, 15(5): 331-338.
高玉峰. 2019. 河谷场地地震波传播解析模型及放大效应. 岩土工程学报, 41(1): 1-25.
管西竹, 符力耘, 陶毅, 等. 2011. 复杂地表边界元-体积元波动方程数值模拟. 地球物理学报, 54(9): 2357-2367. DOI:10.3969/j.issn.0001-5733.2011.09.019
景立平, 卓旭炀, 王祥建. 2005. 复杂场地对地震波传播的影响. 地震工程与工程振动, 25(6): 16-23. DOI:10.3969/j.issn.1000-1301.2005.06.004
李伟华, 赵成刚. 2004. 饱和土沉积谷场地对平面SV波的散射问题的解析解. 地球物理学报, 47(5): 911-919. DOI:10.3321/j.issn:0001-5733.2004.05.025
梁建文, 张彦帅, Lee V W. 2006. 平面SV波入射下半圆凸起地形地表运动解析解. 地震学报, 28(3): 238-249. DOI:10.3321/j.issn:0253-3782.2006.03.003
梁建文, 巴振宁. 2007. 弹性层状半空间中沉积谷地对入射平面SH波的放大作用. 地震工程与工程振动, 27(3): 1-9. DOI:10.3969/j.issn.1000-1301.2007.03.001
梁建文, 齐晓原, 巴振宁. 2014. 基于黏弹性边界的三维凹陷地形地震响应分析. 地震工程与工程振动, 34(4): 21-28.
廖振鹏, 杨柏坡, 袁一凡. 1981. 三维地形对地震地面运动的影响. 地震工程与工程振动, 1(1): 56-77.
刘晶波. 1996. 局部不规则地形对地震地面运动的影响. 地震学报, 18(2): 239-245.
王铭锋, 郑傲, 章文波. 2017. 局部山体地形对强地面运动的影响研究. 地球物理学报, 60(12): 4655-4670. DOI:10.6038/cjg20171210
袁晓铭, 廖振鹏. 1993. 圆弧形凹陷地形对平面SH波散射问题的级数解答. 地震工程与工程振动, 13(2): 1-11.
张建毅, 薄景山, 王振宇, 等. 2012. 汶川地震局部地形对地震动的影响. 自然灾害学报, 21(3): 164-169.
赵成刚, 韩铮. 2007. 半球形饱和土沉积谷场地对入射平面Rayleigh波的三维散射问题的解析解. 地球物理学报, 50(3): 255-264.
周红, 高孟潭, 俞言祥. 2010. SH波地形效应特征的研究. 地球物理学进展, 25(3): 775-782. DOI:10.3969/j.issn.1004-2903.2010.03.005