VTI介质中多波多分量数值模拟及PML吸收边界条件

VTI介质中多波多分量数值模拟及PML吸收边界条件

张为彪;洪宇;郑洁;夏弋峻;郭飞;邹清文

【摘要】为了充分认识多波多分量地震勘探技术中各种波及其分量的传播特点,以便深入挖掘来自地下复杂介质信号中的有效信息,解决实际勘探中遇到的复杂问题,选择了一种特殊的各向异性介质——具有垂直对称轴的横向各向同性介质(VTI)中传播的波及其分量作为研究对象,对其波场进行了数值模拟,分析总结了其传播特点.为了解决边界反射的干扰问题,构建并应用了相应的PML吸收边界条件,模拟结果表明,吸收效果较好.

【期刊名称】《石家庄经济学院学报》

【年(卷),期】2017(040)005

【总页数】6页(P7-12)

【关键词】多波;多分量;各向异性;吸收边界条件;数值模拟

【作者】张为彪;洪宇;郑洁;夏弋峻;郭飞;邹清文

【作者单位】中海石油(中国)有限公司深圳分公司,广东深圳 518054;中海石油(中国)有限公司深圳分公司,广东深圳 518054;中海石油(中国)有限公司深圳分公司,广东深圳 518054;中海油能源发展股份有限公司工程技术分公司,广东深圳518054;中海石油(中国)有限公司深圳分公司,广东深圳 518054;中海石油(中国)有限公司深圳分公司,广东深圳 518054

【正文语种】中文

【中图分类】TE132

随着油气勘探程度的增加，勘探难度不断加大，出现了很多具有挑战性的地质勘探问题，目前最为常用、理论和技术发展也最为成熟的地震勘探方法—纵波勘探已无法完全解决这些问题，而多波多分量地震勘探技术在解决这些问题上有其独特的优势，如：消除模糊带、识别真假亮点、探测裂缝、获得与横波速度有关的各种岩性参数等。

波场数值模拟是人们理解波在不同介质中传播规律，认识和正确解释采集波信号携带的重要信息的重要手段[1-2]。为了理解多波多分量在地下介质中的传播特性，需对其波场进行数值模拟，模拟时，为了避免人工边界产生的反射干扰介质内的波场，必须引入吸收边界条件来处理边界。Cerjan等人[3]通过在人工边界附近引入损耗介质来衰减向外传播的波，该方法的缺点是增加了较大的计算量。另一种比较著名的吸收边界条件是Clayton吸收边界条件[4]，该吸收边界条件在旁轴近似理论的基础上导出，在特定的入射角和频率范围内，具有较好的吸收效果。Berenger[5]针对电磁波传播情况，给出了一种高效的完全匹配层吸收边界条件，并在理论上证明该方法可以完全吸收来自各个方向、各种频率的电磁波，而不产生任何反射。这种方法很快被引入到地震波场的模拟中，Frank（1996）把PML应用到一阶弹性波方程波场的模拟中[6]， Liu Q H（1999）把PML应用到球面坐标系和圆柱坐标系中的弹性波方程的波场的模拟中[7]，王守东（2003）针对声波方程，给出了完全匹配层吸收边界方法[8]。

为了充分认识多波多分量地震勘探技术中各种波及其分量的传播特点，以便深入挖掘来自地下复杂介质信号中的有效信息，解决实际勘探中遇到的复杂问题，需要对地下复杂介质中的各种波及其分量进行数值模拟，本文选择了一种特殊的各向异性介质—具有垂直对称轴的横向各向同性介质（VTI）中传播的波及其分量作为研究对象，对其波场进行了数值模拟和分析。为了解决边界反射的干扰问题,构建并应用了相应的PML吸收边界条件。

均匀各向异性介质中的弹性波方程[9]为：

其中，x1，x2，x3是质点在空间坐标系中的三个坐标分量，u1，u2，u3是质点

沿空间坐标系三个坐标轴方向的位移，Cijkl是弹性系数。

观察（1）式可知，共有81个弹性参数。但是，这81个弹性参数并不是完全独立的，其实，完全独立的弹性参数只有21个，即便如此，要求解该方程也是极其复杂的。实际上，在地震勘探中，经常遇到的是一种特殊的各向异性介质，即横向各向同性（TI）介质，这种介质的各向异性主要是由裂隙和薄互层引起的，具有垂直对称轴的TI介质称为VTI介质，具有水平对称轴的TI介质称为HTI介质，HTI介质可以看成是VTI介质的对称轴旋转90o得到的。对于VTI介质，在x-z平面内，（1）式可简化为：

其中，ρ是密度，C11，C13，C33，C44，C66是弹性系数。

上面的方程组可以分解为两部分，一部分为方程（3），它为VTI介质中的SH波

的方程，另一部分为方程（2）和方程（4）组成的方程组，它们为VTI介质中的

准P波和准SV波的波动方程，这两种波是耦合在一起的，即在VTI介质中，P波传播可以引起SV波，SV波传播也可以引起P波。

（一）VTI介质中SH波波场模拟

1．VTI介质中SH波方程差分格式

根据1知，VTI介质中SH波的方程为：

定义C=( C66/C44)1/2，V=( C44/ρ)1/2，则上述方程化为：

式中，V为VTI介质中的SH波沿垂直于层面方向的速度，C是该介质的各向异性参数，表示VTI介质中地震波沿水平方向和沿垂直方向的速度比。对于方程（6），时间导数采用二阶中心差分，空间导数采用2N阶差分精度的交错网格差分，可得到如下的差分格式[10-17]：

式中，[u2(i,j,k)]、[u2(i,j,k)]同上。

2．VTI介质中SH波方程的PML吸收边界条件

通过引入中间变量u2x，u2z，A2x和A2z，可以把式（6）表示为下面的等价方

程组：

很容易证明（8）与（6）是等价的。

由（8）式来构建新的方程组，即引入函数d1(x)和d2(z)，把方程组（8）变为：（9）的解与（8）的解相比是逐渐衰减的，函数d1(x)和d2(z)分别起到了x和z

方向的衰减系数的作用。这样可以利用（9）式采用下面的方式进行数值模拟：在模拟区域的四周引入完全匹配层，在模拟区域内，令d1(x)和d2(z)为0（此时式（9）等价于式（6）），波不受到衰减，在模拟区域外的完全匹配层内，d1(x)和

d2(z)不为0，传播到完全匹配层内的波受到衰减，当波传播到完全匹配层的边界时，能量已经很弱，边界反射几乎为0，这样就可以得到正确的模拟结果。（9）

式即为VTI介质中SH波方程的PML吸收边界条件。

实际计算的时候，衰减系数d1(x) （d2(z)类似）采用下面的形式：

其中，w为PML层的厚度，x为PML层内网格节点到PML内侧边界的距离，从上式中可看出，在PML内侧边界x为0处, 衰减系数d1(x)=0，在PML外侧边界x为w处，衰减系数d1(x)达到最大，为dmax，dmax可根据波传播一个网格的

距离其振幅衰减为原来的十分之一来计算。

3．波场的数值模拟

首先建立宽度和深度均为2000m的VTI介质模型，模型中，密度ρ=2.44g/cm3，弹性系数C44=7.0×109 N/m2，C66=21.0×109 N/m2，震源位于点（1000m，1000m）处，采用垂直于x-z平面的方向的集中力源作为震源[9]，选用30Hz的Ricker子波作为震源子波。模拟波场时，x方向和z方向的步长均取为8m，t的

步长取为1ms，空间导数采用10阶差分精度的差分格式。波场的模拟结果见图1，其中(a)为t 等于250ms时SH波波场的等值线，(b)、(c)和(d)为t =530ms时