地球物理学报  2019, Vol. 62 Issue (11): 4353-4366   PDF    
任意空间取向TI介质qP和qSV波解耦的波动方程
张庆朝1,2, 朱国维1,2, 何登科1,2, 秦良3, 呼邦兵1,2     
1. 中国矿业大学(北京) 煤炭资源与安全开采国家重点实验室, 北京 100083;
2. 中国矿业大学(北京) 地球科学与测绘工程学院, 北京 100083;
3. 中交(天津)生态环保设计研究院有限公司, 天津 300461
摘要:由于构造运动等作用,TI介质对称轴往往沿空间任意方向分布,具有任意空间取向对称轴的TI(ATI)介质更符合实际地质情况.VTI介质与ATI介质的相速度在形式上具有一致性,VTI介质中地震波的相角对应ATI介质对称轴与地震波传播方向的夹角.本文基于Tsvankin的VTI介质精确相速度公式,利用TI介质对称轴和地震波传播方向上单位向量的数量积和向量积来计算ATI介质的精确相速度.根据弱各向异性假设,导出qP波和qSV波的近似相速度,分析了近似公式的误差,讨论总结了ATI介质qP波和qSV波的相速度特征.本文中的单位向量采用观测坐标系表示,通过相角关系,可以较为方便地由ATI介质近似相速度导出频散关系,然后借助傅里叶逆变换推导出时间-波数域qP波和qSV波解耦的波动方程.数值算例表明本文的波动方程是qP波和qSV波解耦的,波场计算结果稳定,未出现明显的数值频散,验证了本文方法的有效性.
关键词: TI介质      qP-qSV波波动方程      相速度      伪谱法     
A decoupled qP-wave and qSV-wave equation in TI media with arbitrary spatial oriented symmetry axis
ZHANG QingChao1,2, ZHU GuoWei1,2, HE DengKe1,2, QIN Liang3, HU BangBing1,2     
1. State Key Laboratory of Coal Resources and Safe Mining, China University of Mining and Technology, Beijing 100083, China;
2. College of Geoscience and Surveying Engineering, China University of Mining and Technology, Beijing 100083, China;
3. CCCC(Tianjin) Eco-Environmental Protection Design & Research Institute Co., Ltd., Tianjin 300461, China
Abstract: Due to the tectonic deformation and other geological factors,the symmetry axis direction of transversely isotropic media is arbitrary. Formally,the phase angle of plane waves travelling in VTI (transversely isotropic media with a vertical symmetry axis) media is consistent with the angle between the symmetry axis and the direction of plane waves travelling in ATI media. With the accurate phase velocity derived by Tsvankin in VTI media,we calculate the accurate phase velocity of quasi-P and quasi-SV waves for ATI media (transversely isotropic media with arbitrary spatial orientation symmetry axis) in terms of the dot product and cross product of the unit vectors of the symmetry axis and the direction of wave propagation. According to the weak anisotropy assumption and Taylor's formula,we derive the approximate formulas of phase velocity of qP- and qSV-waves. The errors of the approximate formulas are analyzed. The phase velocity characteristics of qP-wave and qSV-wave in ATI media are discussed. The unit vectors are represented by the observational coordinate system. Based on the approximate phase velocity formulas of ATI media,the decoupled qP and qSV-wave equations in the time-wavenumber domain are derived with the phase angle relations by using inverse Fourier transform. Numerical examples demonstrate that the qP and qSV wavefields are decoupled and stable without obvious numerical dispersion,showing the validity of the method proposed in this study.
Keywords: TI media    qP-qSV wave equations    Phase velocity    Pseudo-spectral method    
0 引言

波动方程同时包括了波场的运动学特征和动力学特征,准确、全面地描述了各类地震波型在地层中的传播情况.数值求解波动方程解的过程就是数值模拟地震波传播的过程.对于TI介质,如果仍然采用基于各向同性介质假设的正演和偏移成像算法,则会导致成像精度下降以及偏移归位不准确等问题(Zhan et al., 2012).TI介质有5个独立的弹性系数,弹性波动方程复杂,纵横波耦合,波场分离困难,采用近似公式来简化TI介质波动方程是应用较多的一种方法.Alkhalifah通过将对称轴方向的横波速度设置为零,导出了VTI介质声学波动方程.声学近似方程可以较准确地描述纵波运动学特征,避免了复杂的横波波动方程,也减少了弹性波动方程数值模拟所需的计算资源.许多其他的学者(Klie et al., 2001; Zhang et al., 2005; Zhou et al., 2006aFletcher et al., 2008; Fowler et al., 2010)相继对VTI介质声学近似进行了研究.Duveneck等(2008)基于胡克定律和运动学方程,推导了VTI介质一阶速度-应力声波波动方程.Zhang等(2003, 2005)把Alkhalifah(2000)的波动方程扩展到TTI介质,实现了TTI介质声波方程有限差分数值模拟和逆时偏移成像.Zhou等(2006b)推导了二维TTI介质声波方程.Fletcher等(2009)利用旋转波数,推导了三维TTI介质声波方程.

声学近似波动方程主要存在两个问题,一是横波噪声干扰,二是方程解的稳定性条件.首先,把对称轴方向横波速度设置为零,不意味所有方向的横波速度都为零,数值模拟过程中会有剩余伪横波干扰(Grechka et al., 2004),随着迭代次数的增加,还有可能造成严重的数值频散问题;另外在TI介质对称轴倾角变化剧烈的区域,波动方程的数值解也是不稳定的.Duveneck等(2008)通过在震源附近设置椭圆各向异性层(ε=δ)来消除伪横波噪声的影响,但在模型参数变化较大的区域仍会有横波干扰.类似地,Yoon(2010)在倾角快速变化的区域设置椭圆各向异性区来减少方程解的不稳定性.Fowler(2010)等推导了对称轴方向S波速度不为零的VTI介质二阶耦合波动方程.Fletcher(2009)等通过把对称轴方向的横波设置为有限的速度来解决方程解的不稳定性以及数值频散问题.张岩等(2013)通过滤波的方法来消除次生的菱形横波干扰.

推导TI介质qP-qSV波解耦的纯qP波方程是解决声学近似带来的方程不稳定性及qSV波噪声干扰的有效方法.Du等(2005, 2007)基于弱各向异性近似和平方根近似,推导了TTI介质时间-波数域纯qP波解耦方程.Liu等(2009)推导了VTI介质时间-空间域纯qP波动方程.Pestana等(Pestana and Stoffa, 2010; Pestana et al., 2011)基于快速展开法(REM,Rapid Ex-pansion Method),实现了VTI介质时间-波数域纯qP波逆时偏移成像.Chu等(2011)利用泰勒展开公式,推导了TTI介质解耦纯qP波波动方程.Zhan等(2012)利用基于REM的有限差分与伪谱法混合的方法实现了TTI介质纯qP波逆时偏移.Pacheco等(2013)基于Tsvankin(1996)的精确频散关系,利用一阶泰勒和帕德近似导出了纯qP波动方程.黄翼坚等(2011)用待定系数法得到一个近似相速度公式,并在此基础上导出了一个近似的VTI介质纯qP声波方程.黄建平等(2016)基于伪谱法实现了TTI介质一阶qP波方程正演模拟.杨鹏等(2017)利用近似展开的方法,将偏微分算子分解,推导出各向异性介质中二阶纯qP波方程.黄金强和李振春(2017)借助Low-rank分解, 实现了一种间接的纯qP波波场外推方案.郭成锋等(2017)等利用伪谱法和有限差分法联合的方法实现了纯qP波波场延拓,提高了计算效率.

在实际地层中,由于构造运动等作用,TI介质的对称轴往往沿任意空间方向分布,在观测坐标系中表现为既有倾角变化也有方位变化,具有任意空间取向对称轴的TI介质(Transversely isotropic media with arbitrary spatial orientation symmetry axis,ATI)模型更符合实际情况.ATI介质缺少了VTI介质和HTI介质对于对称轴必须垂直或水平的强约束条件,也缺少了TTI介质对称轴局限在铅直入射面内的限制,含义广泛、地质意义明确(郝重涛和姚陈,2007).进行ATI介质波动方程正演模拟计算,有助于更全面地了解TI介质中地震波的传播规律,并且也与偏移成像、反演等问题直接相关.本文从Tsvankin的VTI介质qP波精确相速度出发,利用传播方向与对称轴方向上的两个单位向量来计算ATI介质qP和qSV波相速度.根据弱各向异性近似,推导了三维ATI介质纯qP和qSV波近似相速度及频散关系.基于频散关系,通过傅里叶逆变换构建了时间-波数域qP-qSV波解耦的波动方程.最后对三维均匀ATI介质进行了qP和qSV波数值模拟,说明本文的波动方程中qP与qSV波是解耦的,消除了声学近似方程中的横波干扰.

1 主要原理

各向异性介质中波的相速度、群速度以及偏振方向准确地描述和刻画了其在传播中的特征(张岩和吴国忱,2013).本文从相速度开始,导出ATI介质频散关系,在此基础上建立时域-波数域波动方程.

1.1 ATI介质相速度

将平面简谐波的位移方程代入TI介质波动方程, 可以得到著名的Kelvin-Christoffel方程,求解Kelvin-Christoffel方程的系数行列式便可求解TI介质的相速度和频散关系.用刚度矩阵弹性系数表示的VTI介质相速度为

(1)

根号前面的加号和±号分别对应qP和qSV波相速度.θ是传播方向与TI介质对称轴的夹角,表示平面波的相角.式(1)中的弹性系数无法直接对应物理量,用两个弹性系数VpzVsz和三个各向异性参数εδγ代替5个弹性系数(Thomsen, 1986)来表示TI介质中波的传播速度.其中

(2)

VpzVsz分别是沿对称轴方向(垂直于各向同性面)的纵、横波相速度;εδ表征qP和qSV波的各向异性强弱,取值越大代表各向异性越强;γ是横波的各向异性参数,也是横波分裂强度的度量(王璐琛,2016).Tsvankin(1996)利用式(2),把式(1)改写为

(3)

式中f≡1-Vsz2/Vpz2.为了获得声学近似波动方程,Alkhalifah(2000)令对称轴方向的横波速度为0,Vsz=0,即f=1,则声波方程的相速度为

(4)

把式(4)中的根号移到等号一边取平方运算,便可得到Alkhalifah(2000)推导的VTI介质声学近似频散关系.对方程两边取平方运算会给式(4)中的相速度引入新解,这是Alkhalifah等推导的声波方程除了qP波解以外还有其他干扰解的原因(黄翼坚等,2011).

ATI介质是由VTI介质绕坐标轴旋转得来,引入了倾角和方位角,但并未引入新的物理量.通常讨论的各向异性介质的刚度矩阵是本构坐标系下的刚度矩阵(吴国忱,2009),ATI介质的本构坐标系和观测坐标系不一致,一般通过Bond变换(Winterstein, 1990)将本构坐标系下的刚度矩阵变换为观测坐标系下的刚度矩阵.基于ATI介质刚度矩阵,结合弹性动力学的本构方程、牛顿运动方程和柯西方程,可以得到ATI介质的弹性波波动方程和Christoffel方程(姚振岸等, 2017).吴国忱等(2009)推导的三维TTI介质qP波、qSV波相速度的解析表达式为

(5)

其中

θ是平面波的极化角,即传播方向与观测坐标系z轴的夹角,φ是平面波传播方位角,θ0是TTI介质对称轴与观测坐标系z轴的夹角,x轴是水平测线方向,z轴是垂直方向,本文采用右手直角坐标系,具体几何关系见图 1图 2.类似地,姚振岸等(2017)推导出了ATI介质qP波和qSV波相速度,表达式和式(5)一致,但其中的EF分别换成下列形式

图 1 TTI介质 Fig. 1 Tilted Transversely Isotropic (TTI) medium
图 2 平面波入射角θ和方位角φ Fig. 2 Incident angle θ and azimuth angle φ

增加了φ0项,φ0是ATI介质对称轴在观测坐标系中的方位角,其他参数意义同式(5).

比较式(1)和(5),可以发现VTI介质、TTI介质以及ATI介质的相速度表达式具有相似的特征.ATI介质中体波速度随传播方向变化的速度图案相对TI介质对称轴确定,速度特征只依赖于传播矢量与对称轴的夹角(郝重涛和姚陈,2007).式(5)相对于式(1)显得复杂,这是因为在观测坐标系下,ATI介质相对VTI介质引入了对称轴的倾角和方位角.结合θφθ0φ0的几何意义,式(5)中的F对应式(1)中的sin2θ,而E对应式(1)中cosθ.为了不混淆式(5)和(1)中的参数,这里用θ来改写公式(5):

其中,可以证明θ就是平面波传播方向与ATI介质对称轴的夹角.(sinθcosφ,sinθsinφ,cosθ)是ATI介质中平面波传播方向的单位向量,(-sinθ0cosφ0,-sinθ0sinφ0,cosθ0)是ATI介质对称轴方向的单位向量,通过求取两个单位向量的数量积和矢量积,便可得到它们的余弦和正弦值,即式(5)中的EF.本文中ATI介质对称轴方向的单位向量取(sinθ0cosφ0,sinθ0sinφ0,cosθ0),观测坐标系为右手坐标系,θ0为对称轴与垂直方向夹角,φ0为对称轴方位角.

类似于将公式(1)改写成公式(3),我们将式(5)改写成

(6)

式(6)就是ATI介质中的qP波和qSV波精确相速度表达式,与郝重涛和姚陈(2007)给出的ATI介质qP波和qSV波相速度解析式一致.公式(6)和(3)形式一样,本质含义相同,都表示TI介质相速度随着传播方向与对称轴的夹角而变化.Tsvankin(1996)根据Thomsen(1986)弱各向异性假设,使用泰勒公式展开式(3)中的根号项,舍弃εδ的二次及多次项,得到VTI介质qP波近似相速度,类似地,可以得到ATI介质中qP波和qSV波近似相速度:

(7)

(8)

其中

式(7)和(8)分别代表ATI介质中qP波和qSV波近似相速度.正如郝重涛和姚陈(2007)的研究结果所示,式(6)、(7)、(8)表明ATI介质中的体波速度特征只依赖于传播矢量与对称轴的夹角θθ 由对称轴方向(sinθ0cosφ0,sinθ0sinφ0cosθ0)和地震波传播方向(sinθcosφ,sinθsinφ,cosθ)决定.

1.2 qP和qSV波动方程

在三维观测坐标系内传播的平面波有下列相角关系:

其中kxkykz表示波数,ω是角频率.

对公式(3)进行移项取平方,然后将上面的相角关系代入公式,可得到VTI介质的qP-qSV波频散关系,Fletcher等(2008)利用旋转波数导出了TTI介质的频散关系

(9)

其中kz为观测坐标系下的波数,kxkykz投影到本构坐标系下的波数,有下列关系

Vsz=0,式(9)变为TTI介质的Alkhalifah(2000)声学近似qP波频散关系

(10)

类似于公式(9)的推导过程,也可直接通过式(4)导出式(10).

对公式(10)进行傅里叶逆变换,得到时间-空间域的时间4阶偏微分方程,方程的解较为复杂.通过引入不同的辅助函数(Fletcher et al., 2008),可将式(10)转换为等价的二阶耦合方程.Fletcher等(2008)引入辅助函数:

式(10)两边乘以p(ω, kx, ky, kz),化简后有

(11)

(12)

对式(11)、(12)进行傅里叶逆变换有

(13)

(14)

p为压力波场,q为辅助波场.式中的偏微分算子为

(15)

(16)

(17)

不同的辅助函数产生不同的耦合方程,Fowler等(2010)基于Alkhalifah(2000)的声学近似频散关系,从二阶耦合方程组特征多项式出发,利用待定系数法,构建了不同形式的耦合方程,其中包含了式(13)、(14).选择其中的时间空间二阶且不含混合空间偏导数的方程可大大减小编程复杂度以及数值模拟所需的计算机资源.基于声学近似频散关系导出的波动方程符合qP波运动学特征,但同时也含有伪横波干扰(Grechka et al., 2004).在求解精确频散关系式(9)的过程中,对式(3)取平方以便去掉其所含的根式,该平方运算引入了多余的解.泰勒公式展开根式的过程未引入新的解,本文从ATI介质近似相速度公式(7)、(8)出发,推导了时间-波数域解耦的qP波和qSV波波动方程.

为了计算的方便,引入单位向量,其中 s 表示TI对称轴方向的单位向量,r 表示地震波传播方向上的单位向量,i 表示坐标轴x方向的单位向量,j 表示坐标轴y方向的单位向量,k 表示坐标轴z方向的单位向量,这里的讨论都是基于观测坐标系.根据图 2,地震波传播方向和对称轴方向的单位向量可以写成下列形式

(18)

地震波传播方向 r 与TI对称轴方向 s 的夹角θ可写成

(19)

将式(19)代入式(7)、(8)中,再利用相角关系,得到qP和qSV波频散关系式(20)、(21).频散方程两边乘以时间-波数域波场U(t, kx, ky, kz),同时在频域应用傅里叶逆变换,得到ATI介质时间-波数域qP波和qSV波解耦的波动方程(22)和(23),可以直接利用伪谱法实现数值模拟计算.在波数域对式(22)、(23)进行傅里叶逆变换得到时间空间域波动方程,然后可以利用有限差分实现数值模拟.由于式(22)、(23)含有四次波数项,导致空间时间域的波动方程含有四阶空间偏导数,增加了编程复杂度.式(20—23)的系数见附录A.

(20)

(21)

(22)

(23)

进行波动方程数值模拟计算时,必须保证方程的数值解收敛不发散.附录B给出了qP波和qSV波波动方程的稳定性条件.

2 数值算例 2.1 相速度误差分析

对于椭圆介质(ε=δ),精确相速度表达式(6)和近似相速度表达式(7)、(8)完全相同,此时qP波的波前变为椭圆,而qSV波速变为常值,即qSV波前退化为圆.对于非椭圆各向异性介质(εδ),近似式(7)、(8)和精确式(6)之间存在误差.在给定TI介质参数的条件下,通过比较式(7)、(8)与式(6)的偏差来分析相速度近似公式的误差.表 1是均匀ATI介质的一组Thomsen各向异性参数.

表 1 Thomsen参数 Table 1 Thomsen parameters

图 3是qP波精确相速度公式(6)和近似公式(7)计算的相速度对比图.图 4是qSV波精确相速度公式(6)和近似公式(8)计算的相速度对比图.图 3图 4选取的相速度传播方向位于x-z垂直平面内,即地震波传播方位角φ为0°或180°,其中图 3a4a对应VTI介质,3b、4b对应TTI介质,3c、4c对应ATI介质.

图 3 TI介质qP波相速度 (a) VTI介质,θ0=0°,φ0=0°;(b) TTI介质,θ0=45°,φ0=0°;(c) ATI介质,θ0=45°,φ0=90°. Fig. 3 qP-wave velocities for TI media (a) VTI media, θ0=0°, φ0=0°; (b) TTI media, θ0=45°, φ0=0°; (c) ATI media, θ0=45°, φ0=90°.
图 4 TI介质qSV波相速度 (a) VTI介质,θ0=0°,φ0=0°;(b) TTI介质,θ0=45°,φ0=0°;(c) ATI介质,θ0=45°,φ0=90°. Fig. 4 qSV-wave velocities for TI media (a) VTI media, θ0=0°, φ0=0°; (b) TTI media, θ0=45°, φ0=0°; (c) ATI media, θ0=45°, φ0=90°.

从图中可以看出,qP波和qSV波的相速度变化曲线呈现周期性规律,明显地,公式(6)、(7)、(8)是夹角θ的周期函数.在各向异性参数确定的条件下,相速度仅依赖于地震波传播方向和TI介质对称轴夹角θ.观测坐标系中,VTI介质和TTI介质对称轴同处于方位角为0的x-z垂直平面内,与该平面内的地震波传播方位角相同,对比图 3图 4中的(a)和(b)子图,可以看出TTI介质和VTI介质的相速度曲线类似,只是相位差了θ0.ATI介质对称轴没有位于x-z垂直平面内,所以此平面内的地震波任意传播方向都不会平行于ATI介质对称轴,即θ≠0°或180°,qP波近似公式(7)中的2εsin4θ项无法像在VTI介质和TTI介质中一样取到零值,故图 3c中ATI介质的qP波最小相速度要大于图 3(a)(b)的最小相速度值.qSV波近似公式(8)中没有2εsin4θ项,所以VTI、TTI以及ATI中qSV波相速度的最小值相同.这里只是讨论了εδ大于零时的情况,εδ的相对大小以及正负取值都会导致不同的相速度变化特征.

图 5图 6是极坐标系下TI介质qP波和qSV波的相速度,图中标注的角度是地震波传播方向与观测坐标系z轴的夹角,即传播方向极化角θ.TTI介质和VTI介质的相速度形态类似,只是绕原点逆时针旋转了θ0角度,而ATI介质相速度形态特征相对于VTI和TTI介质发生了改变.从图 35中可以看到,qP波相速度近似公式(7)的值小于精确相速度公式(6),而图 46中的qSV波相速度近似公式(8)的值大于精确相速度式(6),这是因为泰勒展开(6)式中的根号时,只保留了各向异性参数εδ的一次项,而舍去的高次项为正值,当然这是在εδ全大于0的条件下.qP波在各向同性层方向即与对称轴垂直的方向上相速度最大,沿对称轴方向的相速度最小;ATI介质对称轴不在x-z平面内,地震波传播方向与对称轴夹角最小时即为最小相速度方向.对称轴方向和垂直对称轴方向(各向同性层方向)上,即θ=0°(180°)与θ=90°时,式(6)和(8)相同,qSV相速度相等且为最小值.式(8)近似解在夹角为45°时取最小值,而因为包含了(8)式舍弃的泰勒高阶项,(6)式在45°夹角附近取最小值.这里关于qSV波相速度特征的分析,是基于ε大于δ的条件下进行的.qP波相速度特征与εδ正负都有关系,本文不详细讨论;qSV波的相速度特征相对简单,主要与εδ相对大小有关.

图 5 TI介质极坐标系下qP波相速度 (a) VTI介质,θ0=0°,φ0=0°;(b) TTI介质,θ0=45°,φ0=0°;(c) ATI介质,θ0=45°,φ0=90°. Fig. 5 qP-wave velocities in polar coordinates for TI media (a) VTI media, θ0=0°, φ0=0°; (b) TTI media, θ0=45°, φ0=0°; (c) ATI media, θ0=45°, φ0=90°.
图 6 TI介质极坐标系下qSV波相速度 (a) VTI介质,θ0=0°,φ0=0°;(b) TTI介质,θ0=45°,φ0=0°;(c) ATI介质,θ0=45°,φ0=90°. Fig. 6 qSV-wave velocities in polar coordinates for TI media (a) VTI media, θ0=0°, φ0=0°; (b) TTI media, θ0=45°, φ0=0°; (c) ATI media, θ0=45°, φ0=90°.

图 3图 5显示TI介质qP波相速度近似公式(7)较好地吻合了精确公式(6),而qSV波相速度近似公式(8)和精确公式(6)的偏差大于qP波偏差.近似式(7)和(8)的误差主要来自于式(6)中根号项泰勒展开舍弃的εδ高次项,εδ取值不同,误差曲线的形态也会不同;在各向异性参数确定的情况下,舍弃的这些高次项的值随传播方向的变化而变化.图 7是VTI介质和ATI介质中qP波和qSV波近似相速度公式的误差曲线(图 34中近似曲线和精确曲线的偏差),各向异性参数见表 1,ATI介质对称轴的倾角为45°,方位角为90°.从图 7可以看到,相速度近似公式和精确公式之间的偏差随着TI介质中地震波传播方向和介质对称轴的夹角而变化,当夹角为0°或者90°时,误差减小为0.qP波和qSV波相速度误差曲线呈镜像关系,因为二者舍弃的泰勒高阶项符号相反.本文中,ATI介质近似相速度误差曲线形态不同于VTI介质的主要原因是传播方向和对称轴夹角的取值范围少了0~45°(不包括45°角).图 7a的VTI介质qP波近似相速度均方根误差为8.9383 m·s-1图 7b的VTI介质qSV波近似相速度均方根误差为18.2389 m·s-1图 7c的ATI介质qP波近似相速度均方根误差为8.3738 m·s-1图 7d的ATI介质qSV波近似相速度均方根误差为16.4126 m·s-1.近似相速度的均方根误差较小,近似公式(7)和(8)能较为准确地描述qP和qSV波特征.

图 7 近似相速度公式误差 (a) VTI介质qP波;(b) VTI介质qSV波;(c) ATI介质qP波;(d) ATI介质qSV波. Fig. 7 Phase velocity errors of approximate formula (a) qP-wave in VTI; (b) qSV-wave in VTI; (c) qP-wave in ATI; (d) qSV-wave in ATI.
2.2 波场模拟

设计一个均匀介质模型,模型物理参数见表 1,介质对称轴倾角和方位角分别为45°和90°.模型网格数200×200×200,网格尺寸5 m.震源函数采用雷克子波,主频为50 Hz,震源网格坐标(x, y, z)=(100, 100, 100),数值计算的时间步长为0.2 ms,采样时长130 ms.分别进行声学近似方程、qP以及qSV波动方程数值模拟计算.

图 8是130 ms时刻式(13)、(14)声学近似方程的波场快照,波场内部存在明显的qSV干扰波.图 9是解耦的纯qP波场快照切片,没有剩余qSV干扰波存在.图 9b9d的切片波前形态完全一样.本文中对称轴空间取向特殊,对称轴相对于x-z垂直面和x-y水平面的几何关系相同,所以图 9bx-z垂直切片和9d的x-y水平切片波前形态一样,一般情况下二者的形态特征是不同的.图 9cy-z垂直切片的波前特征和图 5b的TTI介质qP波相速度类似,这是因为对称轴位于y-z垂直面内,此时y-z垂直面就是TTI介质.

图 8 ATI介质声学近似方程波场快照切片 (a)正交切片;(b) x-z切片;(c) y-z切片;(d) x-y切片. Fig. 8 Snapshot slices of acoustic wave field by approximate equation for ATI medium (a) Orthogonal slice; (b) x-z slice; (c) y-z slice; (d) x-y slice.
图 9 ATI介质qP波波场快照切片 (a)三维立体正交切片;(b)垂直方向的x-z切片;(c)垂直方向的y-z切片;(d)水平方向的x-y切片. Fig. 9 Snapshot slices of qP-wave field for ATI (a) Orthogonal slice; (b) x-z slice; (c) y-z slice; (d) x-y slice.

图 10是130 ms时的qSV波波场快照切片.本文中,对称轴与x-z垂直面和x-y水平面的几何关系相同,所以图 10b10d的波前特征类似,呈镜像关系,如果从10d的背面来看,则两者的波前特征一致.综合图 8910,本文推导的波动方程,qP波和qSV波解耦,波场快照清晰,没有多余干扰波组.

图 10 ATI介质qSV波波场快照切片 (a)正交切片;(b) x-z切片;(c) y-z切片;(d) x-y切片. Fig. 10 Snapshot slices of qSV-wave field for ATI (a) Orthogonal slice; (b) x-z slice; (c) y-z slice; (d) x-y slice.
3 结论与认识

ATI介质由VTI介质绕坐标轴旋转得到,并未引入新的物理量,但由于包含了对称轴倾角和方位角,导致地震波速度及波动方程比VTI介质时更复杂.VTI介质与ATI介质的相速度在形式上具有一致性,VTI介质中地震波的相角对应ATI介质对称轴与地震波传播方向的夹角.相速度曲线表明,在各向异性参数确定的条件下,ATI介质中地震波的相速度特征与波的传播方向和对称轴之间的夹角有关,近似相速度误差大小也与此夹角相关,误差的主要来源是舍弃的εδ高次泰勒展开项.本文基于ATI介质近似相速度公式,推导了三维解耦的qP波和qSV波频散关系及波动方程,数值模拟计算结果表明,qP波和qSV波波场解耦分离,没有额外的干扰波组;计算结果稳定,未出现明显的数值频散问题,证明本文方法能够较准确地描述ATI介质波场特征,同时消除了声学近似方程的横波干扰.

附录A qP和qSV波动方程(22)、(23)系数

(A1)

(A2)

附录B 波动方程稳定性条件

将平面波解代入式(22)的有限差分递推公式中,有

(B1)

利用倍角公式对式(B1)进行化简

(B2)

波动方程的稳定性与各向异性参数εδ以及介质对称轴的空间取向都有关,式(B2)过于复杂,这里给出VTI介质的稳定性条件.式(B2)等号两边开平方,并且取θ0=0°(VTI),得到

(B3)

上式与对称轴空间取向无关,令α=min[abs(ε), abs(δ)],式(B3)进一步简化为

根据采样定理,信号带宽应该小于奈奎斯特频率(采样频率的1/2),即有

其中Δd为网格尺寸(空间步长).整理得

(B4)

类似qP波推导过程,可以推导整理出qSV波动方程稳定性条件为

(B5)

References
Alkhalifah T. 2000. An acoustic wave equation for anisotropic media. Geophysics, 65(4): 1239-1250. DOI:10.1190/1.1444815
Chu C, Macy B K, Anno P D. 2011. Approximation of pure acoustic seismic wave propagation in TTI media. Geophysics, 76(5): WB97-WB107. DOI:10.1190/geo2011-0092.1
Du X, Bancroft J C, Lines L R. 2005. Reverse-time migration for tilted TI media.//75th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1930-1933.
Du X, Bancroft J C, Lines L R. 2007. Anisotropic reverse-time migration for tilted TI media. Geophysical Prospecting, 55(6): 853-869. DOI:10.1111/j.1365-2478.2007.00652.x
Duveneck E, Milcik P, Bakker P M, et al. 2008. Acoustic VTI wave equations and their application for anisotropic reverse-time migration.//78th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2186-2190, doi: 10.1190/1.3059320.
Fletcher R P, Du X, Fowler P J. 2008. A new pseudo-acoustic wave equation for TI media.//78th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts.
Fletcher R P, Du X, Fowler P J. 2009. Reverse time migration in tilted transversely isotropic (TTI) media. Geophysics, 74(6): WCA179-WCA187. DOI:10.1190/1.3269902
Fowler P J, Du X, Fletcher R P. 2010. Coupled equations for reverse time migration in transversely isotropic media. Geophysics, 75(1): S11-S22. DOI:10.1190/1.3294572
Grechka V, Zhang L B, Rector J W. 2004. Shear waves in acoustic anisotropic media. Geophysics, 69(2): 576-582. DOI:10.1190/1.1707077
Guo C F, Du Q Z, Zhang M Q, et al. 2017. Numerical simulation and reverse time migration using an improved pure P-wave equation in tilted transversely isotropic media. Chinese Journal of Geophysics (in Chinese), 60(1): 258-270. DOI:10.6038/cjg20170121
Hao C T, Yao C. 2007. Analysis of body_wave velocity characteristic for TI media with arbitrary spatial orientation. Chinese Journal of Geophysics (in Chinese), 50(2): 546-555.
Huang J P, Huang J Q, Li Z C, et al. 2016. Staggered grid pseudo-spectral numerical simulation of first-order quasi P-wave equation for TTI media. Oil Geophysical Prospecting (in Chinese), 2016: 487-496. DOI:10.13810/j.cnki.issn.1000-7210.2016.03.010
Huang J Q, Li Z C. 2017. Modeling and reverse time migration of pure quasi-P-waves in complex TI media with a low-rank decomposition. Chinese Journal of Geophysics (in Chinese), 60(2): 704-721. DOI:10.6038/cjg20170223
Huang Y J, Zhu G M, Liu C Y. 2011. An approximate acoustic wave equation for VTI media. Chinese Journal of Geophysics (in Chinese), 54(8): 2117-2123. DOI:10.3969/j.issn.0001-5733.2011.08.019
Klíe H, Toro W, 2001. A new acoustic wave equation for modeling in anisotropic media: 71st Annual International Meeting, SEG, Expanded Abstracts, 1171-1174.
Liu F Q, Morton S A, Jiang S S, et al. 2009. Decoupled wave equations for P and SV waves in an acoustic VTI media.//79th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2844-2848.
Pacheco D F B, Pestana R, Vivas F. 2013. New pseudo-acoustic wave equations for modeling and reverse time migration in TTI media.//75th EAGE Conference & Exhibition.
Pestana R C, Stoffa P L. 2010. Time evolution of the wave equation using rapid expansion method. Geophysics, 75(4): T121-T131. DOI:10.1190/1.3449091
Pestana R C, Stoffa P L, Ursin B. 2011. Separate P- and SV-wave equations for VTI media.//81st Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 163-167.
Tsvankin I. 1996. P-wave signatures and notation for tranversely isotropic media:An overview. Geophysics, 61(2): 467-483. DOI:10.1190/1.1443974
Thomsen L. 1986. Weak elastic anisotropy. Geophysics, 51(10): 1954-1966. DOI:10.1190/1.1442051
Winterstein D F. 1990. Velocity anisotropy terminology for geophysicists. Geophysics, 55(8): 1070-1088. DOI:10.1190/1.1442919
Wu G C, Liang K, Qi Y P. 2009. Phase velocity and group velocity in 3D TTI media. Progress in Geophysics (in Chinese), 24(6): 2097-2105. DOI:10.3969/j.issn.1004-2903.2009.06.023
Wang L C, Chang X, Wang Y B. 2016. Forward modeling of pseudo P waves in TTI medium using staggered grid. Chinese Journal of Geophysics, 59(3): 1046-1058. DOI:10.6038/cjg20160325
Yang P, Li Z C, Gu B L. 2017. Pure quasi-P wave forward modeling method in TI media and its application to RTM. Chinese Journal of Geophysics (in Chinese), 60(11): 4447-4467. DOI:10.6038/cjg20171130
Yao Z A, Sun C Y, Deng X F, et al. 2017. Velocity characteristic analysis of elastic wave for TI media with arbitrary spatial orientation. Oil Geophysical Prospecting (in Chinese), 52(4): 715-733. DOI:10.13810/j.cnki.issn.1000-7210.2017.04.010
Yoon K, Suh S, Ji J, et al. 2010. Stability and speedup issues in TTI RTM implementation.//80th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 3221-3225, doi: 10.1190/1.3513516.
Zhan G, Pestana R C, Stoffa P L. 2012. Decoupled equations for reverse time migration in tilted transversely isotropic media. Geophysics, 77(2): T37-T45. DOI:10.1190/geo2011-0175.1
Zhang L B, Rector III J W, Hoversten G M. 2003. An acoustic wave equation for modeling in tilted TI media.//73rd Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 153-156.
Zhang L B, Rector J W, Hoversten G M. 2005. Finite-difference modeling of wave propagation in acoustic tilted TI media. Geophysical Prospecting, 53(6): 843-852. DOI:10.1111/j.1365-2478.2005.00504.x
Zhang Y, Wu G C. 2013. Methods of removing pseudo SV-wave artifacts in TTI media qP-wave reverse-time migration. Chinese Journal of Geophysics, 56(6): 2065-2076. DOI:10.6038/cjg20130627
Zhou H, Zhang G, Bloor R. 2006a. An anisotropic acoustic waveequation for VTI media.//66th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, H033.
Zhou H B, Bloor R, Zhang G Q. 2006b. An anisotropic acousticwave equation for modeling and migration in 2D TTI media.//76th Ann. Internat Mtg., Soc. Expi. Geophys.. ExpandedAbstracts, 194-198.
郭成锋, 杜启振, 张明强, 等. 2017. 改进的TTI介质纯P波方程正演模拟与逆时偏移. 地球物理学报, 60(1): 258-270. DOI:10.6038/cjg20170121
郝重涛, 姚陈. 2007. 任意空间取向TI介质中体波速度特征. 地球物理学报, 50(2): 546-555. DOI:10.3321/j.issn:0001-5733.2007.02.028
黄建平, 黄金强, 李振春, 等. 2016. TTI介质一阶qP波方程交错网格伪谱法正演模拟. 石油地球物理勘探, 51(3): 487-496. DOI:10.13810/j.cnki.issn.1000-7210.2016.03.010
黄金强, 李振春. 2017. 基于Low-rank分解的复杂TI介质纯qP波正演模拟与逆时偏移. 地球物理学报, 60(2): 704-721. DOI:10.6038/cjg20170223
黄翼坚, 朱光明, 刘池洋. 2011. 一个近似的VTI介质声波方程. 地球物理学报, 54(8): 2117-2123. DOI:10.3969/j.issn.0001-5733.2011.08.019
王璐琛, 常旭, 王一博. 2016. TTI介质的交错网格伪P波正演方法. 地球物理学报, 59(3): 1046-1058. DOI:10.6038/cjg20160325
吴国忱, 梁锴, 戚艳平. 2009. 三维TTI介质相速度和群速度. 地球物理学进展, 24(6): 2097-2105. DOI:10.3969/j.issn.1004-2903.2009.06.023
杨鹏, 李振春, 谷丙洛. 2017. 一种TI介质纯qP波正演方法及其在逆时偏移中的应用. 地球物理学报, 60(11): 4447-4467. DOI:10.6038/cjg20171130
姚振岸, 孙成禹, 邓小凡, 等. 2017. 任意空间取向TI介质弹性波速度特征分析. 石油地球物理勘探, 52(4): 715-733. DOI:10.13810/j.cnki.issn.1000-7210.2017.04.010
张岩, 吴国忱. 2013. TTI介质qP波逆时偏移中伪横波噪声压制方法. 地球物理学报, 56(6): 2065-2076. DOI:10.6038/cjg20130627