VTI介质井间地震高斯束正演数值模拟

张 林， 杨飞龙， 史朝阳， 江 桂， 陈继川

(1.陕西省地质调查中心，西安 710068；2. 西安石油大学 地球科学与工程，西安 710065；3. 陕西省地质调查院，西安 710054)

VTI介质井间地震高斯束正演数值模拟

张 林1， 杨飞龙2， 史朝阳1， 江 桂1， 陈继川3

(1.陕西省地质调查中心，西安 710068；2. 西安石油大学 地球科学与工程，西安 710065；3. 陕西省地质调查院，西安 710054)

为了研究垂直横向各向同性(VTI)介质井间地震中,地震波的传播方向及波场特征，采用高斯束方法对复杂构造进行波场正演模拟。分析VTI介质中地震波的群相关系，并提出群相数据库思想来解决地震波在VTI介质中传播的方向和速度问题，改善了波场正演模拟的算法。同时高斯束正演方法不仅能够解决复杂构造的盲区问题，且其动力学特征更能反映各向异性介质对地震波传播的影响。较VTI介质频率-空间域弹性波正演模拟方法，正演结果一致且计算效率高，能有效分辨地震波场特征，为地震数据的处理与解释提供了借鉴。

井间地震； 群相数据库； 垂直横向各向同性介质； 高斯束； 正演模拟

0 引言

井间地震[1]因其高精度和高分辨率以及贴近目标探测的特点和优势，将油气勘探开发带入了精细勘探时代。井间地震勘探是将激发系统与接收系统分别放置于不同井中进行数据采集，可以获得多种类型的地震波场信息。丰富的波场信息为成像提供了大量的信息基础，但同时多种波型也会产生相互干扰，使得地震波场难以分辨。地震波在各向异性介质中传播时，地震波场特征还会受到介质弹性参数的变化而发生变化，对有效波场地识别与分离造成更大困难。地震波场的识别与分离方法有很多，其中井间地震正演数值模拟方法是最有效、最直观的方法之一，不仅可以帮助井间地震任务的设计、制定适合的野外采集方案、指导野外采集工作的进行，还可以识别井间复杂波场特征，指导地震数据处理工作的顺利开展。

井间地震正演数值模拟方法主要包括波动方程类和射线类正演方法。波动方程类正演方法的优点在于能够获得较完整的波场信息，但是计算效率低，尤其是面对大规模地震勘探及三维地质模型。张文波等[2]使用交错网格有限差分方法对井中地震弹性波进行正演模拟，以射线理论为基础的正演方法在计算效率上都优于波动方程类正演方法，但是面对物性变化较大的复杂构造，波场信息不完整，出现盲区。cerveny等[3]提出了使用高斯射线束方法进行地震波场正演数值模拟，高斯束是对波动方程的高频近似，它是将地震波场分解到一定频率范围内的射线束上来实现波场的数值模拟，不仅包含波场传播的运动学特征，还具有地震波传播的动力学变化。并且具有高效、精确等特点，适用于复杂地质构造，对复杂构造的盲区、临界区等具有较好效果。杨飞龙等[4]将高斯束正演方法应用在各向同性介质的井间地震波场正演数值模拟当中。针对射线类正演方法计算效率快、波动方程类正演方法精度高的优点，笔者在研究VTI介质井间地震正演时，结合二者优点采用高斯射线束方法进行正演数值模拟。

1 方法技术

1.1 群相数据库

在地震勘探中，我们所研究的地球介质因地质作用影响出现各向异性，介质的各向异性表现在地震波的传播速度与衰减方向上。各向异性介质中，地震波沿着同一点出发的不同方向上介质物理性质不同，主要体现在群速度和相速度上，由于各向异性的作用使得群速度与相速度分离，导致群角与相角也产生分离。

各向异性介质中，群速度代表了地震波能量的传播速度，体现在波场动力学特征，也与地震波的射线路径有关[5]；相速度表示的是波矢量的传播速度，方向与波前面垂直，它是求解简谐平面波解时假定的速度，通过求解Christoffel方程获取[6]。群角是指介质的对称轴与射线传播方向之间的夹角；而相角指的是介质的对称轴方向和波矢量方向的夹角，群速度的方向就是群角的方向，相速度的方向与相角的方向一致(图1)。从图1中可以看到，波矢量的方向垂直于波前面，相速度的方向即为波前面的传播方向，群速度的方向与波射线的方向一致。波前面在各向异性介质下为一椭球面，当相角不等于群角时，相速度与群速度也不相等。

图1 各向异性介质的群相关系Fig.1 Group-phase relations of anisotropic medium

Daley等[7]根据VTI介质的应力与应变的关系，推导出地震波在横向各向同性介质中传播的相速度，使用Thomsen参数来表示弹性参数并将介质近似为弱各向异性介质，其地震波传播的相速度可表示为：

(1)

式中：θ为相角；ε为P波的各向异性参数，ε越大，介质的各向异性越强；δ为影响VTI介质中垂直对称轴附近纵波的速度的参数；γ表示的是横波的各向异性参数，γ与各向异性的强度成正比，当γ=0时，表示介质中不存在横波的各向异性。

Berryman[8]提出了使用相速度计算群速度的计算公式：

(2)

当入射角为0°或者90°时，群速度与相速度相等。

地震波在各向异性介质中传播时，波矢量方向与地震波传播方向不一致，因此群角和相角也不相等。Thomsen[8]给出了弱各向异性介质中地震波传播的群角和相角之间的关系为：

(3)

将式(3)分别写成P波、SV波及SH波的形式：

(4)

式中：vp、vsv、δ、ε、γ为Thomsen参数；θ是相角；φ是群角。

地震波在各向异性介质中传播时，群角为射线角，当射线遇到地层界面发生反射或者透射时地震波的入射角用相角来表示。式(4)为相角计算群角的计算公式，可是射线追踪过程中我们常遇到已知群角计算相角的问题。Byun等[9]均对该问题进行了研究，提出了采用相角计算群角的近似公式。赵爱华等[10]将弱各向异性介质中相角与群角的关系表示为θ=φ-g(φ)，由群角计算相角需要进行近似取舍，并且运算过程比较复杂，计算效率低。基于此，本次研究中仍以式(4)为基础，首先计算出任意相角对应的群角大小，然后计算出所对应的相速度和群速度，将它们保存在一个数据库中，在此称为群相数据库。如图2所示，为群相数据库中群角、相角、群速度和相速度之间的关系。当进行地震波场正演模拟时，射线以某一群角入射，遇到地层界面发生反射或透射时读取群数据库中与之对应的相角，进行SNELL定理计算。高斯束正演模拟时每个检波器的能量是由多条射线能量加权而来，在VTI介质高斯束正演时,需要进行多次SNELL运算，若提前计算相角对应的群角大小，便可在运算时直接调用数据库中的信息，避免运算过程中产生的误差，同时会提高运算的效率。

图2 群相关系示意图Fig.2 Group-phase relations diagram

图3 VTI介质地震波传播示意图Fig.3 Seismic wave propagation diagram in VTI medium

1.2 反射与透射

如图3所示，地震波在VTI介质中传播时在界面处发生透射、反射，仍然遵循斯奈尔定律[11]。但是对于VTI介质，地震波传播的速度和入射角度都存在群相分离现象，因此较各向同性介质来说斯奈尔定律更加复杂。

在VTI介质中，介质的对称轴方向与垂直方向一致，斯奈尔定律可以表示为：

(5)

式中：P为射线参数；α为入射角；θ为相角；v(θ)为相速度。

在各向异性介质中，可以将反射系数与透射系数分为两部分：①各向同性介质性质的反射、透射系数项;②各向异性介质性质的反射、透射系数项。可以表示为式(6)。

(6)

Ripp、Tipp、Rips和Tips分别为P波的反射系数、透射系数，SV波的反射系数和透射系数。

以PP波(P波入射P波反射)为例，反射系数的各向同性项可以表示为式(7)。

(7)

反射系数的各向异性项可以表示为：

(8)

其中：vp01和vs01为界面以上的P波速度和SV波速度；vp02和vs02为界面以下的P波速度和SV波速度；ρ1和ρ2分别为界面上下两种介质的密度。

1.3 高斯束正演

高斯束正演是将波场分解到具有一定频率范围的射线束上实现地震波场的数值模拟，它是将动力学方程集中在射线附近的高频渐近时间调和解。高斯束正演包括运动学射线追踪、动力学射线追踪和波场叠加三个阶段。

高斯射线束运动学追踪即求解从震源发出的所有射线经过地层反射、透射后最终到达接收井的射线路径、旅行时和中心射线的振幅。运动学射线追踪方法很多，常用的是求解程函方程的射线追踪方法。

(9)

式中：τ为位移；v为速度；x、z为水平和垂直坐标，其中：x、y、z为直角坐标系下空间坐标，τ为旅行时，v为地震波传播速度。已知各反射界面的分布函数，炮点的坐标和出射角，即射线的函数表达式是已知的，当按一定角度步长打出一条条射线后，根据界面函数和射线函数关系可以求出两者交点的坐标，然后判断交点是否在我们所设计模型的有效范围内。如果在就记录此交点坐标，然后利用Snell定律求反射或透射角，进一步求出生成射线方程，再利用上述方法求与其他界面的交点，直至到达接收井上，求出整条射线路径的所有控制点坐标，最后把各控制点坐标输出到存储器中；如果不在就直接进行下一个角度的射线追踪。

经过运动学追踪，就有了中心射线，在此基础上可以进行这高斯射线束的动力学射线追踪，即求得p、q的函数值。函数p(s)、q(s)在高斯射线束中起着非常重要的作用，它们决定了高斯射线束能量的分布状态，也表征沿射线传播方向的高频地震波场动力学特征。

用U(R,t)表示地震波传播到接收点R的波场，根据傅里叶变化得到时间域波场：

(10)

写成离散形式为：

(11)

式中：φ为从震源发出射线的入射角；φ0和φN分别为起始入射角和终止入射角；g为高斯波包；Δφ为入射角的间隔。

根据文献[3]，可知波包g的近似解析表达式为式(12)。

g(R,φ)= (2πfm)1/2|Aφ|exp{-[2πfm(t-θ)/

γ]2+(2πfmG/γ)-2πfmG}·

(12)

式中:f*=fm·(1-4πfmG/γ2)，为高斯波包的主频。

2 数值模拟

图4 VTI介质井间地震高斯束正演地质模型Fig.4 Geological model of cross-well seismic Gaussian beam forward in VTI medium

表1 VTI介质井间地震地质模型参数

表2 井间地震观测系统参数

图5 VTI介质井间地震高斯束正演射线路径Fig.5 The ray path of cross-well seismic Gauss beam forward in VTI medium

如图4所示，为VTI介质井间地震高斯束正演地质模型，模型参数如表1所示。井间地震高斯束正演的观测系统如表2所示。

对该VTI介质的井间地质模型使用上述观测系统进行高斯束正演，得到如图5所示的射线路径。井间地震的波场较地面地震和VSP更加复杂，在VTI介质中，地震波传播的射线路径也不像均匀介质中那么有序，地震波传播时遵循VTI介质下SNELL定律。图5中可以看到井间地震波场正演包含直达P波、直达S波，以及上行反射P波、上行反射S波与下行反射P波、下行反射S波。在实际的井间地震数据处理中，识别上下行反射波场是在波场分离中重要的工作之一。透过正演研究，能够给井间地震数据处理提供一个依据，指导数据处理进一步开展。

图6 VTI介质井间地震高斯束正演射线路径Fig.6 The ray path of cross-well seismic Gauss beam forward in VTI medium

图7 VTI介质井间地震高斯束正演波场记录Fig.7 The wave record of cross-well seismic Gauss beam forward in VTI medium

图6为第二炮井间地震高斯束正演的射线路径，图7为其相应的波场记录。从图6中可以看到，在VTI介质中地震波传播时产生的射线分布不均匀，尤其面对复杂地质构造的时候会产生射线的扭转，在地震记录上也表现为反转的现象。图7的波场记录中直达P波记录中有空缺的地方，也出现一小段斜率变化较大的地方，这是因为在VTI介质中射线传播时遇到界面后反射透视所遵循的VTI介质SNELL定律。从图7中的波场记录上可以看到，VTI介质的井间地震波场更加复杂，在复杂地质构造中地震波场表现的没有那么整齐，出现了一些小的错段。在实际井间地震数据中，由于井间地震数据的复杂性，再加上复杂构造，对于识别分辨波场十分困难，透过VTI介质井间高斯束正演数值模拟，可以有效模拟实际地下构造井间复杂波场特征，帮助我们分辨所需波场，进行下一步地震数据处理。

在实际井间地震数据处理时，波场分离工作是处理的重要任务。通过不同波场类型的正演，可以帮助我们有效地分辨井间地震波场传播的特点，为进一步井间地震数据处理提供依据。图8(d)中的下行反射P波的时距曲线上出现一段斜率不同的下行反射P波波场特征，在接收井深度410 m～510 m范围里。造成此种现象的原因是，地震波在各向异性介质中传播时受到不同方向上速度差异的影响，在实际的井间地震资料上有很多类似这种由于地层的各向异性和复杂构造造成的波场特点。

图8 VTI介质井间地震上下行反射波场特征Fig.8 The up and down reflection wave field characteristics of cross-well seismic Gauss beam forward in VTI medium(a)上行S波的射线路径；(b)对应的正演记录；(c)下行P波的射线路径；(d)对应的正演记录

3 不同井间地震正演方法对比研究

图9 VTI介质地质模型Fig.9 Geological model in VTI medium

为了验证VTI介质井间地震高斯束正演算法的准确性，选择简单的层状地质模型将其与VTI介质频率-空间域弹性波正演模拟方法进行对比。建立如图9所示的VSP勘探地质模型(VSP勘探是一种特殊的井间地震勘探)，模型参数如表3所示。炮点放置于地面300 m处，偏移距为300 m。检波点安置在接收井(直井)上，道间距为5 m，共120道。分别使用高斯束方法和波动方程方法进行正演模拟，得到如图10所示的波场特征。既包含纵波震源产生的波场，也有横波震源产生的波场。对比其纵波震源产生的波场特征，可以看出，VTI介质井间地震高斯束正演的波场形态、旅行时与高阶交错网格弹性波数值模拟结果一致，验证了该方法的准确性。波动方程正演模拟方法不仅具有地震波传播的动力学特征，也能清晰反映地震波在VTI介质中传播的运动学特点。将两种方法对同一地质模型正演效率进行对比，如表4所示，可见在相同电脑配置下，高斯束正演方法的计算速度要比波动方程法正演效率高近乎50倍。在此将两种方法进行对比，更加说明本文所提出的群相数据库思想应用在VTI介质井间地震高斯束正演模拟中，能够改善VTI介质正演模拟算法。

图10 不同正演方法波场对比示意图Fig.10 The contrast diagram of wave field between two different forward methods(a)VTI介质高斯束正演的波场特征；(b)高阶交错网格弹性波数值模拟方法正演的结果

表3 VTI介质井间地震地质模型参数

表4 单炮正演时间统计表

4 结论

笔者研究了VTI介质井间地震高斯束正演方法，提出了使用群相数据库思想来解决群角计算相角的问题，有效地改善了VTI介质的射线正演计算方法。高斯束正演方法不仅能解决复杂地质构造产生的盲区问题，且作为射线类正演方法具有较高的运算效率。将其与井间地震高阶交错网格弹性波数值模拟方法进行对比，运算结果一致，证明笔者所研究的VTI介质井间地震高斯束正演方法的正确性，为在井间地震复杂构造采集设计与波场分离工作中提供了借鉴。

[1] 何惺华.井间地震[M].北京：石油工业出版社，2008.

HE X H. Cross-well seismic[M]. Beijing: Petroleum Industry Press，2008.(In Chinese)

[2] 张文波. 井间地震交错网格高阶差分数值模拟及逆时偏移成像研究[D].西安:长安大学，2005.

ZHANG W B. The studies of staggered grid high-order finite difference numerical simulation and reverse-time migration imaging in cross-well seismic[D].Xi an: Chang’an University，2005.(In Chinese)

[4] 杨飞龙, 孙渊, 李绪宣,等. 基于高斯射线束的斜井VSP正演方法[J]. 地球物理学进展, 2014, 29(6)：2791-2799.

YANG F L, SUN Y, LI X X, et al. Deviated hole VSP forward method based on the Gaussian beam[J]. Progress in Geophysics, 2014,29(6)：2791-2799.(In Chinese)

[5] 杜丽英, 刘国明, 杜丽娟,等. VTI介质中地震波反射波合成记录的方法研究[J]. 地球物理学进展, 2001, 2 (16): 58-64.

DU L Y, LIU G M, DU L J, et al. The method of seismic and reflection wave synthetic seismo gram in VTI media[J]. Progress in Geophysics,2001, 2 (16): 58-64. (In Chinese)

[6] DALEY P F, HRON F. Reflection and transmission coefficients for transversely isotropic media [J]. Bulletin of the Seismological Society of America, 1977, 67 (3): 661-675.

[7] BERRYMAN J G. Long-wave elastic anisotropy in transversely isotropic media [J]. Geophysics, 1979, 44 (5): 896-917.

[8] THOMSEN,L.Weadly elostic anisotropy[J].Geophysics,1986,51(10):1954-1966.

[9] BYUN B,CORRIGAN D,GAISER J.Anisotropic velocity analysis for lithology discrimination [J].Geophysics,1989,54 (12):1564-1574.

[10] 赵爱华, 丁志峰. 一种弱各向异性介质地震波群速度的近似表示新方法[J]. 地球物理学进展, 2005, 20 (4): 916-919.

ZHAO A H, DING Z F. New approximate expressions of seismic group velocities for weakly anisotropic media[J]. Progress in Geophysics,2005, 20 (4): 916-919. (In Chinese)

[11] 李建国, 李彦鹏, 郭晓玲. VTI 介质试射射线追踪 [J]. 石油地球物理勘探, 2010, (4): 491-496.

LI J G, LI Y P, GUO X L.VTI medium test-firing ray tracing[J].OGP,2010(4): 491-496. (In Chinese)

Gaussianbeamforwardnumericalsimulationofcross-wellseismicinVTImedium

ZHANG Lin1, YANG Feilong2, SHI Zhaoyang1, JIANG Gui1, CHEN Jichuan3

(1.Shaanxi Center of Geological Survey, Xi'an 710068， China；2.School of Earth Sciences and Engineering, Xi'an Shiyou University, Xi'an 710065， China；3.Shaanxi Geological Survey Institute, Xi'an 710054， China)

The Gaussian beam method is used in wave field forward modeling in complex structure, in order to study the direction of seismic wave and the wave field characteristics in vertical transversely isotropic(VTI) medium in cross-well seismic. The group-phase relations of seismic wave are analyzed in VTI media, the group-phase database thought is put forward to solve the problems of direction and velocity of seismic wave in VTI media, and the algorithm of wave field forward modeling is improved. At the same time, Gaussian beam forward modeling method is not only solving the problem of blind area of complex structure, and its dynamic characteristics can reflect the anisotropy medium effects on seismic wave propagation. Compared with elastic wave forward modeling method in frequency-spatial domain in VTI media, the forward modeling results are consistent and the computational efficiency is higher. The characteristics of seismic wave field can be effective separation. The method provides accurate basis for seismic data processing and interpretation.

cross-well seismic; group-phase database; vertical transversely isotropic medium; gaussian beam; forward modeling

2017-03-16 改回日期： 2017-06-09

张林(1978-)，男，高级工程师，主要从事综合地球物理勘探，E-mail：136751912@qq.com。

陈继川(1970-)，男，高级工程师，主要从事物探、遥感、测绘等方面研究，E-mail：527379546@qq.com。

1001-1749(2017)06-0791-08

P 631.4

A

10.3969/j.issn.1001-1749.2017.06.12