﻿ TTI介质qP波伪谱法正演模拟
 石油地球物理勘探  2019, Vol. 54 Issue (2): 302-311  DOI: 10.13810/j.cnki.issn.1000-7210.2019.02.008 0
 文章快速检索 高级检索

### 引用本文

ZHANG Qingchao, ZHU Guowei, ZHOU Junjie, ZHU Congcong, LIU Weigang. qP-wave numerical simulation in TTI media with pseudo-spectral method. Oil Geophysical Prospecting, 2019, 54(2): 302-311. DOI: 10.13810/j.cnki.issn.1000-7210.2019.02.008.

### 文章历史

TTI介质qP波伪谱法正演模拟

qP-wave numerical simulation in TTI media with pseudo-spectral method
ZHANG Qingchao , ZHU Guowei , ZHOU Junjie , ZHU Congcong , LIU Weigang
State Key Laboratory of Coal Resources and Safe Mining, China University of Mining and Technology(Beijing), Beijing 100083, China
Abstract: Based on Tsvankin's exact phase velocity for VTI media, we obtain the exact phase velocity of qP-wave for 3D TTI media with the coordinate transformation. Then following Thomsen anisotropic parameters, we gain the approximate formula of phase velocity for 3D TTI media with Taylor series. The phase velocity analysis shows that the approximate formula could well fit the exact formula. From the approximate formula of phase velocity, we derive the qP-wave dispersion relation and wave equation in the time-wavenumber domain for 3D TTI media. The stability of the wave equation is analyzed. Numerical examples demonstrate that the solution of the wave equation proposed in this paper is stable and free of S-wave artifacts.
Keywords: TTI media    phase velocity    pure qP-wave    finite difference method    pseudo-spectral method
0 引言

1 方法原理 1.1 TTI介质的相速度

 $\begin{array}{l} \frac{{V_{\rm{P}}^2\left( \theta \right)}}{{V_{{\rm{P}}z}^2}} = 1 + \varepsilon {\sin ^2}\theta - \frac{f}{2} + \\ \frac{f}{2}\sqrt {1 + \frac{{4{{\sin }^2}\theta }}{f}\left( {2\delta {{\cos }^2}\theta - \varepsilon \cos 2\theta } \right) + \frac{{4{\varepsilon ^2}{{\sin }^4}\theta }}{{{f^2}}}} \end{array}$ (1)

 $\frac{{V_{\rm{P}}^2\left( \theta \right)}}{{V_{{\rm{P}}z}^2}} = 1 + 2\delta {\sin ^2}\theta {\cos ^2}\theta + 2\varepsilon {\sin ^4}\theta$ (2)

 $\begin{array}{l} \frac{{V_{\rm{P}}^2\left( \gamma \right)}}{{V_{{\rm{P}}z}^2}} = 1 + \varepsilon {\sin ^2}\gamma - \frac{f}{2} + \\ \;\;\;\;\frac{f}{2}\sqrt {1 + \frac{{4{{\sin }^2}\gamma }}{f}\left( {2\delta {{\cos }^2}\gamma - \varepsilon {{\cos }2}\gamma } \right) + \frac{{4{\varepsilon ^2}{{\sin }^4}\theta }}{{{f^2}}}} \end{array}$ (3)

 $\frac{{V_{\rm{P}}^2\left( \gamma \right)}}{{V_{{\rm{P}}z}^2}} = 1 + 2\delta {\sin ^2}\gamma {\cos ^2}\gamma + 2\varepsilon {\sin ^4}\gamma$ (4)

 $\left\{ \begin{array}{l} \mathit{\boldsymbol{s}} = \sin {\theta _0}\mathit{\boldsymbol{i}} + \cos {\theta _0}\mathit{\boldsymbol{k}}\\ \mathit{\boldsymbol{r}} = \sin \theta \cos \varphi \mathit{\boldsymbol{i}} + \sin \theta \sin \varphi \mathit{\boldsymbol{j}} + \cos \theta \mathit{\boldsymbol{k}} \end{array} \right.$ (5)
 图 1 地震波传播方向与坐标轴的关系

sr的夹角为γ，通过sr的矢量积及数量积得到

 $\cos \gamma = \sin \theta \cos \varphi \sin {\theta _0} + \cos \theta \cos {\theta _0}$ (6)
 $\begin{array}{l} {\sin ^2}\gamma = {\sin ^2}\theta {\sin ^2}\varphi + \\ \;\;\;\;\;\;\;\;\;\;\;\;{\left( {\sin \theta \cos \varphi \cos {\theta _0} - \cos \theta \sin {\theta _0}} \right)^2} \end{array}$ (7)

 图 2 由式(3)、式(4)得到的相速度 (a)VTI介质，θ0=0°，φ=0°，ΔVRMS=13.4154m/s; (b)VTI介质，θ0=0°，φ=90°，ΔVRMS=13.4154m/s; (c)TTI介质，θ0=45°，φ=0°，ΔVRMS=13.4605m/s; (d)TTI介质，θ0=45°，φ=90°，ΔVRMS=18.5240m/s
1.2 时间—波数域波动方程的推导

 $\left\{ \begin{array}{l} \sin \theta \cos \varphi = \frac{{{V_{\rm{P}}}\left( {\theta ,\varphi } \right){k_x}}}{\omega }\\ \sin \theta \sin \varphi = \frac{{{V_{\rm{P}}}\left( {\theta ,\varphi } \right){k_y}}}{\omega }\\ \cos \theta = \frac{{{V_{\rm{P}}}\left( {\theta ,\varphi } \right){k_z}}}{\omega } \end{array} \right.$ (8)

 $\begin{array}{l} {\omega ^4} = \left[ {\left( {V_{{\rm{P}}x}^2 + V_{{\rm{S}}z}^2} \right)\left( {\hat k_x^2 + \hat k_y^2} \right) + \left( {V_{{\rm{P}}z}^2 + V_{{\rm{S}}z}^2} \right)\hat k_z^2} \right]{\omega ^2} - \\ \;\;\;\;\;\;\;\;V_{{\rm{P}}x}^2V_{{\rm{S}}z}^2{\left( {\hat k_x^2 + \hat k_z^2} \right)^2} - V_{{\rm{P}}z}^2V_{{\rm{S}}z}^2\hat k_z^4 + \\ \;\;\;\;\;\;\;\;\left[ {V_{{\rm{P}}z}^2\left( {V_{{\rm{P}}n}^2 - V_{{\rm{P}}x}^2} \right) - V_{{\rm{S}}z}^2\left( {V_{{\rm{P}}n}^2 - V_{{\rm{P}}z}^2} \right)} \right] \times \\ \;\;\;\;\;\;\;\;\left( {\hat k_x^2 + \hat k_y^2} \right)\hat k_z^2 \end{array}$ (9)

 $\left\{ \begin{array}{l} {{\hat k}_x} = {k_x}\cos {\theta _0}\cos {\varphi _0} + {k_y}\cos {\theta _0}\sin {\varphi _0} - {k_z}\sin {\theta _0}\\ {{\hat k}_y} = - {k_x}\sin {\varphi _0} + {k_y}\cos {\varphi _0}\\ {{\hat k}_z} = {k_x}\sin {\theta _0}\cos {\varphi _0} + {k_y}\sin {\theta _0}\sin {\varphi _0} + {k_z}\cos {\theta _0}\\ \hat k_x^2 + \hat k_y^2 + \hat k_z^2 = k_x^2 + k_y^2 + k_z^2 \end{array} \right.$ (10)

 $\begin{array}{l} {\omega ^4} = \left[ {V_{{\rm{P}}x}^2\left( {\hat k_x^2 + \hat k_y^2} \right) + V_{{\rm{P}}z}^2\hat k_z^2} \right]{\omega ^2} + \\ \;\;\;\;\;\;\;V_{{\rm{P}}z}^2\left( {V_{{\rm{P}}n}^2 - V_{{\rm{P}}x}^2} \right)\left( {\hat k_x^2 + \hat k_y^2} \right)\hat k_z^2 \end{array}$ (11)

f=1，将式(6)、式(7)代入式(3)也可导出式(11)。对式(11)进行傅里叶逆变换，得到时间—空间域的4阶偏微分方程，方程的解较为复杂。通过引入辅助函数

 $\begin{array}{l} q\left( {{k_x},{k_y},{k_z},\omega } \right) = \frac{{{\omega ^2} + \left( {V_{{\rm{P}}n}^2 - V_{{\rm{P}}x}^2} \right)\left( {\hat k_x^2 + \hat k_y^2} \right)}}{{{\omega ^2}}} \times \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;p\left( {{k_x},{k_y},{k_z},\omega } \right) \end{array}$ (12)

Fletcher等[15]将四阶偏微分方程转换为等价的2阶耦合方程。式(11)两边乘以p(kx, ky, kz, ω)，化简后有

 $\begin{array}{l} {\omega ^2}p\left( {{k_x},{k_y},{k_z},\omega } \right) = V_{{\rm{P}}x}^2\left( {\hat k_x^2 + \hat k_y^2} \right)p\left( {{k_x},{k_y},{k_z},\omega } \right) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;V_{{\rm{P}}z}^2\hat k_z^2q\left( {{k_x},{k_y},{k_z},\omega } \right) \end{array}$ (13)
 $\begin{array}{l} {\omega ^2}q\left( {{k_x},{k_y},{k_z},\omega } \right) = V_{{\rm{P}}n}^2\left( {\hat k_x^2 + \hat k_y^2} \right)p\left( {{k_x},{k_y},{k_z},\omega } \right) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;V_{{\rm{P}}z}^2\hat k_z^2q\left( {{k_x},{k_y},{k_z},\omega } \right) \end{array}$ (14)

 $\frac{{{\partial ^2}p}}{{\partial {t^2}}} = V_{{\rm{P}}x}^2\left( {\frac{{{\partial ^2}p}}{{\partial {{\hat x}^2}}} + \frac{{{\partial ^2}p}}{{\partial {{\hat y}^2}}}} \right) + V_{{\rm{P}}z}^2\frac{{{\partial ^2}q}}{{\partial {{\hat z}^2}}}$ (15)
 $\frac{{{\partial ^2}q}}{{\partial {t^2}}} = V_{{\rm{P}}n}^2\left( {\frac{{{\partial ^2}p}}{{\partial {{\hat x}^2}}} + \frac{{{\partial ^2}p}}{{\partial {{\hat y}^2}}}} \right) + V_{{\rm{P}}z}^2\frac{{{\partial ^2}q}}{{\partial {{\hat z}^2}}}$ (16)

 $\begin{array}{l} \frac{{{\partial ^2}}}{{\partial {{\hat x}^2}}} = {\cos ^2}{\theta _0}{\cos ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {x^2}}} + {\cos ^2}{\theta _0}{\sin ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {y^2}}} + \\ \;\;\;\;\;\;{\sin ^2}{\theta _0}\frac{{{\partial ^2}}}{{\partial {z^2}}} + {\cos ^2}{\theta _0}\sin 2{\varphi _0}\frac{{{\partial ^2}}}{{\partial x\partial y}} - \\ \;\;\;\;\;\;\sin 2{\theta _0}\cos {\varphi _0}\frac{{{\partial ^2}}}{{\partial x\partial z}} - \sin 2{\theta _0}\sin {\varphi _0}\frac{{{\partial ^2}}}{{\partial y\partial z}} \end{array}$ (17)
 $\frac{{{\partial ^2}}}{{\partial {{\hat y}^2}}} = {\sin ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {x^2}}} + {\cos ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {y^2}}} - \sin 2{\varphi _0}\frac{{{\partial ^2}}}{{\partial x\partial y}}$ (18)
 $\begin{array}{l} \frac{{{\partial ^2}}}{{\partial {{\hat z}^2}}} = {\sin ^2}{\theta _0}{\cos ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {x^2}}} + {\sin ^2}{\theta _0}{\sin ^2}{\varphi _0}\frac{{{\partial ^2}}}{{\partial {y^2}}} + \\ \;\;\;\;\;{\cos ^2}{\theta _0}\frac{{{\partial ^2}}}{{\partial {z^2}}} + {\sin ^2}{\theta _0}\sin 2{\varphi _0}\frac{{{\partial ^2}}}{{\partial x\partial y}} + \\ \;\;\;\;\;\sin 2{\theta _0}\cos {\varphi _0}\frac{{{\partial ^2}}}{{\partial x\partial z}} + \sin 2{\theta _0}\sin {\varphi _0}\frac{{{\partial ^2}}}{{\partial y\partial z}} \end{array}$ (19)

 $\begin{array}{l} \frac{{{\partial ^2}{U_P}\left( {{k_x},{k_y},{k_z},t} \right)}}{{\partial {t^2}}} = - V_{{\rm{P}}z}^2\left[ {k_x^2 + k_y^2 + k_z^2 + } \right.\\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\cos }^4}{\theta _0}} \right)\frac{{k_x^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {\delta \sin 4{\theta _0} - 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_x^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}2{\theta _0} - \sigma {{\sin }^2}2{\theta _0} + 3\varepsilon {{\sin }^2}2{\theta _0}} \right)\frac{{k_x^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} - \\ \;\;\;\left( {\delta \sin 4{\theta _0} + 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_z^3{k_x}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\sin }^4}{\theta _0}} \right)\frac{{k_z^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0} + 4\varepsilon {{\cos }^2}{\theta _0}} \right)\frac{{k_x^2k_y^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta \sin 2{\theta _0} - 4\varepsilon \sin 2{\theta _0}} \right)\frac{{{k_x}k_y^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}{\theta _0} + 4\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_y^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left. {2\varepsilon \frac{{k_y^4}}{{k_x^2 + k_y^2 + k_z^2}}} \right]{U_{\rm{P}}}\left( {{k_x},{k_y},{k_z},t} \right) \end{array}$ (20)

 $\frac{{{V_{{\rm{P}}z}}\Delta t}}{{\Delta d}} < \frac{2}{{{\rm{ \mathsf{ π} }}\sqrt {4\alpha + 3} }}$ (21)

2 数值计算 2.1 均匀模型

 图 3 由式(15)、式(16)(上)及式(20)(下)得到的VTI均匀介质模型在139ms的波场快照 (a)三维正交切片；(b)xoz切片；(c)yoz切片；(d)xoy切片 模型网格数为201×201×201，网格间距为dx=dy=dz=5m，各向异性参数：VPz=3000m/s；ε=0.25；δ=0.1，VTI介质对称轴倾斜角θ0=0°。震源函数采用雷克子波，主频为100Hz，震源坐标为(x, y, z)=(500m, 500m, 500m)。数值计算的时间步长为0.1ms，总采样时长为150ms

 图 4 由式(15)、式(16)(上)及式(20)(下)得到的TTI均匀介质模型在139ms的波场快照 (a)三维正交切片；(b)xoz切片；(c)yoz切片；(d)xoy切片模型对称轴倾角θ0= 45°，其他参数同图 3
2.2 Hess VTI模型

 图 5 Hess VTI模型 (a)VPz；(b)ε；(c)δ

 图 6 Hess VTI模型在1.5s的波场快照 (a)式(15)、式(16)；(b)式(20)
2.3 BP 2007 TTI模型

 图 7 BP 2007 TTI模型 (a)VPz；(b)ε；(c)δ；(d)对称轴倾角θ0

 图 8 BP 2007 TTI模型在1.5s的波场快照 (a)式(15)、式(16)；(b)式(20)网格间距为dx=dz=12.5m, 震源坐标为(x, z)=(68750m, 2625m), 其他参数同图 6

3 结束语

(1) 在各向异性参数确定的条件下，任意对称轴取向的TI介质中qP波的传播速度取决于传播方向与TI对称轴的夹角。

(2) 把TTI介质对称轴的方位角取为0°，简化了计算，文中的各向异性模型为2.5维，但不同方向的切片仍然可以表征TTI介质的方位各向异性。

(3) 本文的纯qP波波动方程可以较好地应用于TTI介质正演模拟，得到的波场快照没有伪横波干扰，对复杂TTI模型仍然保持较好的波场稳定性。伪谱法计算量较大，可以和有限差分法联合使用以减少计算时间。

 $\begin{array}{l} {U_\text{P}}\left( {{k_x},{k_y},{k_z},t + \Delta t} \right) - 2{U_\text{P}}\left( {{k_x},{k_y},{k_z},t} \right) + \\ \;\;\;{U_P}\left( {{k_x},{k_y},{k_z},t - \Delta t} \right) = - \Delta {t^2}V_{{\rm{P}}z}^2\left[ {k_x^2 + k_y^2 + k_z^2 + } \right.\\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\cos }^4}{\theta _0}} \right)\frac{{k_x^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {\delta \sin 4{\theta _0} - 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_x^3{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}2{\theta _0} - \delta {{\sin }^2}2{\theta _0} + 3\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_x^2{k_z^2}}}{{k_x^2 + k_y^2 + k_z^2}} - \\ \;\;\;\left( {\delta \sin 4{\theta _0} + 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_z^3{k_x}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\sin }^4}{\theta _0}} \right)\frac{{k_z^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0} + 4\varepsilon {{\cos }^2}{\theta _0}} \right)\frac{{k_x^2k_y^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta \sin 2{\theta _0} - 4\varepsilon \sin 2{\theta _0}} \right)\frac{{{k_x}k_y^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}{\theta _0} + 4\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_y^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left. {2\varepsilon \frac{{k_y^4}}{{k_x^2 + k_y^2 + k_z^2}}} \right]{U_{\rm{P}}}\left( {{k_x},{k_y},{k_z},t} \right) \end{array}$ (A-1)

 $\begin{array}{l} {{\rm{e}}^{ - {\rm{i}}\omega \Delta t}} - 2 + {{\rm{e}}^{{\rm{i}}\omega \Delta t}} = - \Delta {t^2}V_{{\rm{P}}z}^2\left[ {k_x^2 + k_y^2 + k_z^2 + } \right.\\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\cos }^4}{\theta _0}} \right)\frac{{k_x^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {\delta \sin 4{\theta _0} - 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_x^3{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}2{\theta _0} - \delta {{\sin }^2}2{\theta _0} + 3\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_x^2{k_z^2}}}{{k_x^2 + k_y^2 + k_z^2}} - \\ \;\;\;\left( {\delta \sin 4{\theta _0} + 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_z^3{k_x}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\sin }^4}{\theta _0}} \right)\frac{{k_z^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0} + 4\varepsilon {{\cos }^2}{\theta _0}} \right)\frac{{k_x^2k_y^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta \sin 2{\theta _0} - 4\varepsilon \sin 2{\theta _0}} \right)\frac{{{k_x}k_y^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}{\theta _0} + 4\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_y^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left. {2\varepsilon \frac{{k_y^4}}{{k_x^2 + k_y^2 + k_z^2}}} \right] \end{array}$ (A-2)

 $\begin{array}{l} 4{\sin ^2}\frac{{\omega \Delta t}}{2} = \Delta {t^2}V_{{\rm{P}}z}^2\left[ {k_x^2 + k_y^2 + k_z^2 + } \right.\\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\cos }^4}{\theta _0}} \right)\frac{{k_x^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {\delta \sin 4{\theta _0} - 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_x^3{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}2{\theta _0} - \delta {{\sin }^2}2{\theta _0} + 3\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_x^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} - \\ \;\;\;\left( {\delta \sin 4{\theta _0} + 8\varepsilon \sin {\theta _0}{{\cos }^3}{\theta _0}} \right)\frac{{k_z^3{k_x}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0}{{\cos }^2}{\theta _0} + 2\varepsilon {{\sin }^4}{\theta _0}} \right)\frac{{k_z^4}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\sin }^2}{\theta _0} + 4\varepsilon {{\cos }^2}{\theta _0}} \right)\frac{{k_x^2k_y^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta \sin 2{\theta _0} - 4\varepsilon \sin 2{\theta _0}} \right)\frac{{{k_x}k_y^2{k_z}}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left( {2\delta {{\cos }^2}{\theta _0} + 4\varepsilon {{\sin }^2}{\theta _0}} \right)\frac{{k_y^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + \\ \;\;\;\left. {2\varepsilon \frac{{k_y^4}}{{k_x^2 + k_y^2 + k_z^2}}} \right] \end{array}$ (A-3)

 $\begin{array}{l} 2 > \Delta t{V_{{\rm{P}}z}} \times \\ \;\;\;\sqrt {k_x^2 + k_y^2 + k_z^2 + 2\varepsilon \frac{{k_x^4}}{{k_x^2 + k_y^2 + k_z^2}} + 2\delta \frac{{k_x^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + 4\varepsilon \frac{{k_x^2k_y^2}}{{k_x^2 + k_y^2 + k_z^2}} + 2\delta \frac{{k_y^2k_z^2}}{{k_x^2 + k_y^2 + k_z^2}} + 2\varepsilon \frac{{k_y^4}}{{k_x^2 + k_y^2 + k_z^2}}} \end{array}$ (A-4)

α=min[abs(ε), abs(δ)]，则式(A-4)简化为

 $2 > \Delta t{V_{{\rm{P}}z}}\sqrt {k_x^2 + k_y^2 + k_z^2 + 2\alpha \left( {k_x^2 + k_y^2} \right)}$ (A-5)

 ${k_x},{k_y},{k_z} < \frac{1}{2}\frac{{2{\rm{ \mathsf{ π} }}}}{{\Delta d}} = \frac{{\rm{ \mathsf{ π} }}}{{\Delta d}}$

 $\frac{{{V_{{\rm{P}}z}}\Delta t}}{{\Delta d}} < \frac{2}{{{\rm{ \mathsf{ π} }}\sqrt {4\alpha + 3} }}$ (A-6)

 [1] 黄金强, 李振春. 基于Low-rank分解的复杂TI介质纯qP波正演模拟与逆时偏移[J]. 地球物理学报, 2017, 60(2): 704-721. HUANG Jinqiang, LI Zhenchun. Modeling and reverse time migration of pure quasi-P-waves in complex TI media with a low-rank decomposition[J]. Chinese Journal of Geophysics, 2017, 60(2): 704-721. [2] 王娟, 李振春, 孙小东, 等. TTI介质逆时偏移成像[J]. 石油地球物理勘探, 2012, 47(4): 573-577. WANG Juan, LI Zhenchun, SUN Xiaodong, et al. Reverse time migration in tilt transversely isotropic(TTI) media[J]. Oil Geophysical Prospecting, 2012, 47(4): 573-577. [3] Alkhalifah T. Acoustic approximations for processing in transversely isotropic media[J]. Geophysics, 1998, 63(2): 623-631. [4] Alkhalifah T. An acoustic wave equation for aniso-tropic media[J]. Geophysics, 2000, 65(4): 1239-1250. [5] Zhou H, Zhang G, Bloor R.An anisotropic acoustic wave equation for VTI media[C]. Extended Abstracts of 68th EAGE Conference & Exhibition, 2006, H033. [6] Zhou H, Zhang G, Bloor R.An anisotropic acoustic wave equation for modeling and migration in 2D TTI media[C]. SEG Technical Program Expanded Abstracts, 2006, 25: 194-198. [7] Zhang L, Rector J W, Hoversten G M. Finite-differ-ence modelling of wave propagation in acoustic tilted TI media[J]. Geophysical Prospecting, 2005, 53(6): 843-852. DOI:10.1111/gpr.2005.53.issue-6 [8] Zhang H, Zhang Y.Reverse time migration in 3D heterogeneous TTI media[C]. SEG Technical Program Expanded Abstracts, 2008, 27: 2196-2200. https://www.researchgate.net/publication/240736933_Reverse_time_migration_in_3D_heterogeneous_TTI_media [9] Zhang Y, Zhang H.A stable TTI reverse time migration and its implementation[C]. SEG Technical Program Expanded Abstracts, 2009, 28: 2794-2798. https://www.researchgate.net/publication/301434048_A_Stable_TTI_Reverse_Time_Migration?ev=auth_pub [10] Fletcher R P, Du X, Fowler P J.A new pseudo-acoustic wave equation for TI media[C]. SEG Technical Program Expanded Abstracts, 2008, 27: 2082-2086. https://www.researchgate.net/publication/249862440_A_New_Pseudo-acoustic_Wave_Equation_for_VTI_Media [11] 梁锴, 吴国忱, 印兴耀, 等. 三维TTI介质波动方程分解[J]. 石油地球物理勘探, 2009, 44(1): 19-27. LIANG Kai, WU Guochen, YING Xinyao, et al. Wave equation decomposition in 3-D TTI medium[J]. Oil Geophysical Prospecting, 2009, 44(1): 19-27. DOI:10.3321/j.issn:1000-7210.2009.01.005 [12] Yoon K, Suh S, Ji J, et al. Stability and speedup issues in TTI RTM implementation[C]. SEG Technical Program Expanded Abstracts, 2010, 29: 3221-3225. https://www.researchgate.net/publication/269069126_Stability_and_speedup_issues_in_TTI_RTM_implementation?ev=auth_pub [13] 张岩, 吴国忱. TTI介质qP波逆时偏移中伪横波噪声压制方法[J]. 地球物理学报, 2013, 56(6): 2065-2076. ZHANG Yan, WU Guochen. Methods of removing pseudo SV-wave artifacts in TTI media qP-wave reverse-time migration[J]. Chinese Journal of Geophysics, 2013, 56(6): 2065-2076. [14] Fowler P J, Du X, Fletcher R P. Coupled equations for reverse time migration in transversely isotropic media[J]. Geophysics, 2010, 75(1): S11-S22. [15] Fletcher R P, Du X, Fowler P J. Reverse time migration in tilted transversely isotropic (TTI) media[J]. Geophysics, 2009, 74(6): WCA179-WCA187. DOI:10.1190/1.3269902 [16] Du X, Bancroft J C, Lines L R.Reverse-time migration for tilted TI media[C]. SEG Technical Program Expanded Abstracts, 2005, 24: 1930-1934. [17] Du X, Bancroft J C, Lines L R. Anisotropic reverse-time migration for tilted TI media[J]. Geophysical Prospecting, 2007, 55(6): 853-869. DOI:10.1111/gpr.2007.55.issue-6 [18] Liu F Q, Morton S A, Jiang S S, et al. Decoupled wave equations for P and SV waves in an acoustic VTI media[C]. SEG Technical Program Expanded Abstracts, 2009, 28: 2844-2848. [19] Pestana R C, Stoffa P L. Time evolution of the wave equation using rapid expansion method[J]. Geophy-sics, 2010, 75(4): 121-131. [20] Pestana R C, Ursin B, Stoffa P L.Separate P-and SV-wave equations for VTI media[C]. SEG Technical Program Expanded Abstracts, 2011, 30: 163-167. [21] Chu C, Macy B K, Anno P D. Approximation of pure acoustic seismic wave propagation in TTI media[J]. Geophysics, 2011, 76(5): WB97-WB107. DOI:10.1190/geo2011-0092.1 [22] Zhan G, Pestana R C, Stoffa P L. Decoupled equations for reverse time migration in tilted transversely isotropic media[J]. Geophysics, 2012, 77(2): 37-45. DOI:10.1190/geo2011-0175.1 [23] 王伟国, 熊水金, 徐华宁, 等. TTI介质各向异性伪谱法逆时偏移[J]. 石油地球物理勘探, 2012, 47(4): 566-572. WANG Weiguo, XIONG Shuijin, XU Huaning, et al. Reverse-time migration using pseudo-spectral method for tilted TI media[J]. Oil Geophysical Prospecting, 2012, 47(4): 566-572. [24] 黄建平, 黄金强, 李振春, 等. TTI介质一阶qP波方程交错网格伪谱法正演模拟[J]. 石油地球物理勘探, 2016, 51(3): 487-496. HUANG Jianping, HUANG Jinqiang, LI Zhenchun, et al. Staggered grid pseudo-spectral numerical simulation of first-order quasi P-wave equation for TTI media[J]. Oil Geophysical Prospecting, 2016, 51(3): 487-496. [25] 郭成锋, 杜启振, 张明强, 等. 改进的TTI介质纯P波方程正演模拟与逆时偏移[J]. 地球物理学报, 2017, 60(1): 258-270. GUO Chengfeng, DU Qizhen, ZHANG Mingqiang, et al. Numerical simulation and reverse time migration using an improved pure P-wave equation in tilted transversely isotropic media[J]. Chinese Journal of Geophysics, 2017, 60(1): 258-270. [26] Tsvankin I. P-wave signatures and notation for tranversely isotropic media:An overview[J]. Geophysics, 1996, 61(2): 467-483. DOI:10.1190/1.1443974 [27] Thomsen L. Weak elastic anisotropy[J]. Geophysics, 1986, 51(10): 1954-1966. DOI:10.1190/1.1442051 [28] 郝重涛, 姚陈. 任意空间取向TI介质中体波速度特征[J]. 地球物理学报, 2007, 50(2): 546-555. HAO Chongtao, YAO Chen. Analysis of body-wave velocity characteristic for TI medium with arbitrary spatial orientation[J]. Chinese Journal of Geophysics, 2007, 50(2): 546-555. DOI:10.3321/j.issn:0001-5733.2007.02.028 [29] 姚振岸, 孙成禹, 邓小凡, 等. 任意空间取向TI介质弹性波速度特征分析[J]. 石油地球物理勘探, 2017, 52(4): 715-733. YAO Zhen'an, SUN Chengyu, DENG Xiaofan, et al. Velocity characteristic analysis of elastic wave for TI media with arbitrary spatial orientation[J]. Oil Geophysical Prospecting, 2017, 52(4): 715-733. [30] Cerjan C, Kosloff D, Kosloff R, et al. A nonreflecting boundary condition for discrete acoustic and elastic wave equation[J]. Geophysics, 1985, 50(4): 705-708. DOI:10.1190/1.1441945