地球物理学进展  2015, Vol. 30 Issue (2): 797-804   PDF    
引用本文
刘乐, 林伟军, 张海澜, 王秀明. 2015. 横向各向同性地层倾斜井中的声场模式分析[J]. 地球物理学进展, 30(2): 797-804, doi: 10.6038/pg20150243  
LIU Le, LIN Wei-jun, ZHANG Hai-lan, WANG Xiu-ming. 2015. Mode analyses in a deviated borehole penetrating a transversely isotropic formation. Progress in Geophysics, 30(2): 797-804, doi: 10.6038/pg20150243   
横向各向同性地层倾斜井中的声场模式分析
刘乐, 林伟军 , 张海澜, 王秀明    
中国科学院声学研究所声场与声信息国家重点实验室, 北京 100190
收稿日期: 2014-08-08 修回日期: 2014-11-08 .
基金项目：国家自然科学基金重点项目(11134011,41274134)和中国科学院“十二五”信息化专项超级计算重点应用示范项目(XXH12503-02-02-05)联合资助.
作者简介：刘乐,男,1987年生,在读博士生,主要研究方向为井孔声场理论和声波测井.(E-mail:liule10@mails.ucas.ac.cn)
通讯作者：林伟军,男,1971年生,研究员.(E-mail:linwj@mail.ioa.ac.cn)
摘要：针对横向各向同性地层包围的倾斜井模型, 基于有限元方法求解频率-波数域波动方程, 使用卷积完全匹配层进行数值截断.从计算得到的二维谱中可以直观地看到模式波的分布, 继而得到模式波的相速度.与前人的方法相比, 该方法对地层各向异性的强弱程度没有要求, 可以计算强各向异性介质;此外, 计算不依赖于具体的声源形式, 可以在较宽的频率范围内进行, 且准确度高, 误差基本不超过1%.计算了不同倾斜角度、不同偶极源指向性条件下快、慢地层包围的倾斜井中的弯曲模式波的相速度频散曲线.计算结果表明, 当井孔倾斜时, 偶极源激发的弯曲波分裂成快、慢两条弯曲波, 其低频速度分别趋近于地层的SH和qSV速度.接收器接收到的快、慢弯曲波的相速度由倾斜角决定, 与声源指向角无关.本文中的方法具有较好的拓展性, 可以很方便地修改以分析其他模型.
关键词： 倾斜井     弯曲波模式     频散曲线     有限元法    
Mode analyses in a deviated borehole penetrating a transversely isotropic formation
LIU Le, LIN Wei-jun , ZHANG Hai-lan, WANG Xiu-ming    
State Key Laboratory of Acoustics, Institute of Acoustics, Chinese Academy of Sciences, Beijing 100190, China
Abstract: The wave equation was solved in the frequency-wave number domain based on the finite element method to investigate the wave propagation in a deviated borehole penetrating an transversely isotropic formation. A convolutional perfectly matched layer was implemented to simulate the numerical truncation. The mode distribution was shown directly and clearly on the two dimension spectrum, so the phase velocity of the mode could be obtained afterwards. Compared with the perturbation method, the weak elastic anisotropy is not necessary as the derivation asks for nothing about the anisotropy parameters. And unlike FDTD, the computation was conducted in the frequency domain without carrying a specific source form. Results are obtained in a wider frequency range and with higher accuracy. The relative error of phase velocity is mostly less than 1 percent. Using this method, we calculated the phase velocity dispersion curves of the flexural wave in a deviated borehole surrounded by a fast or slow formation with different deviation angles and dipole orientation. The results show that for such a deviated borehole, the flexural wave splits into fast and slow flexural waves, whose phase velocities are respectively equal to that of the SH and qSV wave of the transversely isotropic formation at low frequencies. Further study show that the phase velocities of the fast and slow flexural waves were independent of the source or receiver orientation. The method used in this paper can be extended to other similar models with little modification.
Key words: deviated borehole     flexural mode     dispersion curve     finite element method    
0 引 言

声波测井是石油勘探及开发的一个重要环节.由于沉积作用，实际地层通常呈现为相互叠加的薄层，当声波波长大于薄层厚度时地层可视为横向各向同性(TI)介质.该介质有一个对称轴，垂直于该轴的各个方向上的介质性质(如弹性波速度等)相同，但沿对称轴方向的介质性质则与垂直于该方向有所不同.在声波测井中，TI介质是最普遍的各向异性模型(张碧星和王克协，1998)，常见的沉积地层如页岩地层即呈现不同程度的各向异性，例如，有些页岩地层的竖直向和水平向的速度各向异性可达20%~40%(Thomsen，1986; 唐晓明和郑传汉，2004).近年来，深海油储的勘探成为热点，而这些储层地层通常呈现较强的各向异性特征，此外，离岸开采通常采用大斜度井以降低钻井成本并增加采收率(Chi and Tang， 2006).因此，针对TI介质地层包围的倾斜井孔中的声场的研究工作具有重要的学术和应用价值.

前人的研究工作表明，当井孔与TI介质对称轴不再平行，系统不再具有对称性，无法采用分离变量法进行解析求解(Ellefsen，1990).针对这样的非对称问题，国内外研究者使用不同的方法进行了大量工作.使用广泛的有时域有限差分(FDTD)方法(Leslie and Randall， 1992; Wang et al., 2002; Sinha et al., 2006; 林伟军等，2006; He et al., 2010; 阎守国等，2011)，该方法直接计算单、偶极子声源激发得到的时域接收波形，通过对波形进行处理(Kimball and Marzetta， 1984; Ekstrom，1995; Nolte et al., 1997)得到频散曲线.Ellefsen等(1991a)，Norris和Sinha(1993)，Sinha 等(1994)，Zhang和Wang(1996)，张碧星和王克协(2000)采用不同的微扰法进行模式波的分析，与FDTD相比，微扰法可以得到近似解析的公式，计算速度快，但限于方法本身的特点，仅适用于各向异性程度较弱的情况.此外，Ellefsen等(1991b)将沿井轴传播的波的特征用指数因子表示，对垂直井轴的截面采用有限元方法，计算得到了弯曲模式波的频散曲线.实验方面，林伟军等(2004)进行了非轴对称井孔声场的模拟实验，朱正亚等(1994)、Zhu 等(2007)及王瑞甲等(2013)实验观察了各向异性条件下弯曲波的分裂.

不久前，张海澜等(2011)采用2.5维方法计算了充液管道中的频率-波数域特征，这种方法可以直接进行波场的模式分析，也可以得到对应的时域波形，为解决TI介质中的倾斜井问题提供了新思路.本文基于该方法，模拟得到了TI介质地层倾斜井中的声场分布，计算了不同倾斜角度的模式波的频散曲线.与FDTD方法相比，该方法在频率域进行计算，可以方便而准确地进行模式波的分析，且计算基于有限元方法，在处理复杂结构和边界问题时更准确；与微扰法相比，该方法没有采取各向异性近似，对TI介质的各向异性程度没有要求，此外，计算模型可以很方便地拓展以进行其他分析.

为了模拟井孔外的无限大介质，本文采用非分裂卷积完全匹配层(C-PML)进行数值截断.该完全匹配层由Roden和Gedney(2000)提出，与传统的PML相比，C-PML对于大角度入射及表面波有更好的吸收效果(Drossaert and Giannopoulos， 2007; Komatitsch and Martin， 2007; 张鲁新等，2012)，Li 等(2010)及李义丰等(2010)使用COMSOL的PDE模块在频率域实现了二维弹性介质的C-PML，本文将其拓展到2.5维. 1 模型及方法

本文以TI介质地层包围的倾斜井孔为研究对象.如图 1，计算在直角坐标系进行，为了研究方便，假设井孔是竖直的，地层是倾斜的，且偏离井孔方向 α角度，不失一般性地，假设地层对称轴的投影与直角坐标系x轴重合.声源位于中心原点处，且偶极源的指向角偏离x轴θ角度，接收器位于井轴，接收广义声压(Kurkjian，1985)，即：单极源激发的声场，井轴上的接收器接收流体中的位移位；偶极源激发的声场，接收器接收位移位的一阶导数，即径向位移，且接收器与声源的指向角度相同.

图 1 TI介质地层倾斜井模型Fig. 1 Model of a deviated borehole penetrating a T Iformation

TI介质的弹性系数矩阵由5个常数确定，在直角坐标系中，如果其对称轴与井轴平行，则其弹性系数矩阵为

对于倾斜井，TI介质的对称轴绕y轴旋转α角度，其弹性系数矩阵为C0= MCMT，M为变换矩阵，MT为 M的转置矩阵，M的表达式为(Auld，1973)

在倾斜条件下，与弹性系数矩阵 C 相比，经变换后的 C 0矩阵多出一些非零项，但矩阵仍是对称的.

当井孔及TI介质地层沿对称轴方向保持不变，声波沿井轴的传播可以用因子eikz表示(Ellefsen et al.，1991b)，其中k是z方向的波数.这样，三维问题的求解可以转变为求解二维方程的频率-波数循环，减小了对计算资源的需求，即：

其中p为广义声压，直接计算结果为p(x，y，ω，k).

与前人的FDTD方法或微扰法不同，本文使用一种新的方法进行计算.该方法利用COMSOL Multiphysics的偏微分方程(PDE)模块进行求解.该软件基于有限元方法，以自定义的方式求解偏微分方程.PDE模块提供系数型通式：

其中式(4a)为关于u的二阶方程的系数型通式，式(4b)和(4c)用以实现边界条件.使用该方法求解偏微分方程的主要步骤是将待求问题的方程、边界条件的表达式整理成通式形式并得到对应系数，之后划分有限元网格从而进行求解及相关后处理. 2 计算模型的PDE模块实现 2.1 2.5维弹性介质C-PML

如图 2，从原点向外依次为圆形流体区域，TI介质地层和C-PML区域，实际计算时图中“—”区域仅考虑x方向吸收，“|”区域仅考虑y方向吸收，“+”区域同时考虑x和y方向的吸收.以x方向为例，引入复坐标转换：

即有：

其中sx为复频移扩展坐标变量(Kuzuoglu and Mittra，1996)，表达式为

上式中的αx(x)、σx(x)及κx(x)均为正实数且κx(x)≥1，选取合适的参数即可实现吸收效果.当σx(x)=0且κx(x)=1时，，为非吸收区域，相当于弹性地层.各自表达式为

本文所取参数为

其中x0为进入C-PML的初始位置，dx为C-PML区域厚度，本文取一个纵波波长，vp为弹性介质中纵波速度，f0为计算频率，R0为理论反射系数.
图 2 计算模型示意图Fig. 2 Schematic diagram of the computational model

对于3维弹性介质，其运动方程的频域形式为

式中下标i，j=1，2，3，对应直角坐标系的三个坐标轴方向，τij表示应力，ρs为弹性介质密度.

此外，应力与应变的本构关系为

其中εkl为应变，表达式为

将应变表达式代入应力，继而代入运动方程，可以得到关于位移u的方程，推导结果见附录A.

由于沿井轴的声波传播使用指数因子eikz表示，有=ik，且无需z方向的吸收，即取sz=1.将附录A公式两端同时乘以sxsy并进行整理以符合式(4a)，得到对应系数，见附录B，没有给出的系数取零.

实际计算时弹性常数矩阵的部分元素为零，但是为了保证计算模型可拓展，这里进行了完整的推导. 2.2 流体区域

计算时，流体中的变量为位移位G，满足频域波动方程

其中vf为流体速度.对比通式(4a)，得到对应系数项为

流-固边界条件为：法向位移连续，法向应力连续，切向应力为零，即

按照整理方程表达式的方法将边界条件整理成式(4b)和式(4c)的形式，得到对应系数，这里不再赘述. 3 算法验证及数值算例 3.1 算法验证

TI介质地层采用慢地层Austin Chalk，其参数见表 1.流体中声波速度vf=1500 m/s，密度ρf=1000 kg/m3，井孔半径为R=0.1 m.声源为正负点源组合而成的偶极源，指向角θ=0°.为保证精度，网格尺寸不超过最小波长的，并在流-固边界及声源处进行加密.

表 1 TI介质参数Table 1 Parameters of the TI formations

图 3为计算得到的频率-波数域的二维谱，可以看到清晰的弯曲模式波，图 4为弯曲波相速度与理论解的对比. 横轴为频率，纵轴为弯曲波的相速度，实线为利用数值外推法求得解析解的频散曲线，“o”为从图 3中提取的结果，与理论解一致，误差大部分在1%之内.

图 3 计算得到的二维谱Fig. 3 The two dimensional spectrum

图 4 2.5维方法计算得到的弯曲波相速度与解析解对比Fig. 4 Comparison of flexural phase velocity by 2.5 dimensional method with the analytic solution

除了准确度高之外，各频率点的计算是完全独立的，因而对于较大的计算规模易于并行，且计算不依赖于声源频谱的具体形式，可以在较宽的频率范围内进行，而采用FDTD计算时域波形再加以算法处理得到频散曲线的方法依赖于声源的具体形式，对于声源频谱外的频率难以进行处理(Leslie and Randall， 1992; 阎守国等，2011).此外，仅需要较少的频率点就足以刻画出一条准确的频散曲线，例如，图 4中的频率间隔是250 Hz，实际上取更宽的频率间隔也是可行的，相比之下，时域计算还需要严格地处理时间步长参数而进行迭代求解，这也是频域计算的一个优势. 3.2 数值算例

在验证算法准确性的基础上，计算了不同倾斜角、不同声源指向性条件下的弯曲波相速度频散曲线.图 5a为偶极源指向角θ=0°时激发的慢弯曲波频散曲线随倾斜角的变化，倾斜角为0°和90°时低频速度一致，且等于TI介质的qSV波速度(√ c44/ρs).此外，相速度随倾斜角的变化并不是单调的，而是逐渐增大至45°时达到最大，之后随着倾斜角的增大，相速度逐渐减小.

图 5 慢地层
(a)慢弯曲波相速度；(b)快弯曲波相速度.Fig. 5 Slow formation
(a)Slow flexural phase velocities；(b)Fast flexural phase velocities.

图 5b为偶极源指向角θ=90°时时激发的快弯曲波频散曲线随倾斜角的变化，与慢弯曲波不同，随着倾斜角的增大，相速度逐渐增大，在倾斜角为90°时相速度低频趋近于TI介质的SH波速度(√ c66/ρs).

采用同样的方法，计算了当TI介质地层为快地层(Cotton Valley shale)时的情况，地层参数见表 1，为了更清晰地进行分析，我们仅将倾斜角为45°和90°时得到的快、慢弯曲波及0°时的结果绘在一张图上，见图 6.随着倾斜角的增大，慢弯曲波速度并不是单调变化的，倾斜角为90°时与倾斜角为0°时低频速度都趋于TI介质qSV波速度，高频段稍有不同；而快弯曲波在倾斜角为90°低频速度趋于TI介质SH波速度，这些都是与软地层的规律一致的.图中同样颜色的实线与虚线的对比直观地说明了井孔倾斜带来的弯曲波分裂现象.值得注意的是，硬地层激发的弯曲波在较高的频率处(图中高于5 kHz)速度几乎一致，各倾斜角的结果区分很不明显，这与软地层的情况稍有不同.

图 6 快地层时弯曲波相速度Fig. 6 Fast formation: flexural phase velocities

由于计算在频率域进行，可以很方便地计算弯曲波的低频速度，图 7中实线为TI介质qSV和SH相速度的理论解(White，1983)，“o”为本方法计算的结果，二者符合很好，相对误差大部分小于1%.图中显示，井孔倾斜时，快、慢弯曲波的低频速度近似为TI介质的qSV和SH波速度.

图 7 弯曲波低频速度Fig. 7 Flexural phase velocities at low frequencies

为了更清楚地从二维谱中观察弯曲波的分裂，我们计算了指向角为60°时不同井孔倾斜角度下的接收二维谱.图 8a是倾斜角为0°时硬地层包围井孔中的接收二维谱，图中可见清晰的弯曲波模式，没有发生分裂现象；图 8b是倾斜角为30°时的结果，弯曲波分裂为快、慢弯曲波；图 8c是倾斜角为90°时的结果，弯曲波的分裂更加明显，但是和图 8b类似，仅在低于5 kHz时可以观察到，频率较高时快、慢弯曲波模式几乎重合.作为对比，图 8d显示了倾斜角为90°时软地层包围井孔中的接收二维谱，与硬地层相比，弯曲波的分裂明显，在整个计算频率范围内都可以观察到.这与图 5和图 6中给出的相速度频散曲线的特征一致.

图 8 偶极源指向角60°时的接收二维谱
(a)倾角0°，快地层；(b)倾角30°，快地层；(c)倾角90°，快地层；(d)90°，慢地层.Fig. 8 2 dimensional spectrum with 60°dipole orientation for
(a)0°(b)30°(c)90°(d)90°borehole deviation，(a)，(b) and (c)for fast formation and (d)for slow formation.

图 9a是TI介质为快地层，倾斜角为75°、偶极源指向角θ=60°时的计算结果，可以看到两条清晰的弯曲波模式，其分别对应快、慢弯曲波.图 9b为相速度比较，改变偶极源指向角计算3次，“+”表示θ=0°时激发的慢弯曲波，“.”表示θ=90°时激发的快弯曲波，“o”表示θ=60°时的计算结果，该结果表明当偶极源指向角不沿着x轴或y轴时，接收器同时接收到快、慢弯曲波，考虑到指向角60°并不具有特殊性，可以认为指向角的改变不会影响接收到的弯曲波的速度，这与时域分析的结果(He et al., 2010)是一致的.

图 9 倾斜角75°，偶极源指向角60°时计算结果(a)二维谱(b)相速度Fig. 9(a)2 dimensional spectrum and (b)phase velocity for 75°borehole deviation and 60°dipole orientation
4 结 论

4.1     本文针对横向各向同性地层包围的倾斜井，采用频率-波数域的2.5维方法，实现非分裂卷积完全匹配层吸收边界条件，从频率-波数域的二维谱中提取模式波得到相速度频散曲线.与前人的方法相比，本文采用的方法对介质的各向异性程度没有要求，可以计算各向异性程度较强的TI介质，且由于计算直接在频率域进行，易于并行，对处理求解频散曲线这样的频率域问题具有较高的精度.

4.2     计算的结果表明：1)当井孔倾斜时，偶极源激发的弯曲模式波发生分裂，其中慢弯曲波的低频速度趋近于TI介质的qSV波速度，其随倾斜角的变化并非单调，但倾角为0°和90°时该速度一致且趋近于 √ c44/ρs ；快弯曲波的低频速度趋近于TI介质的SH波速度，当倾斜角从0°逐渐增大到90°，该速度从√ c44/ρs逐渐变为 √ c66/ρs.2)与慢地层不同，当TI介质为快地层时，随着倾斜角的变化，弯曲波分裂在高频段不明显.3)当偶极源的指向角不再沿着偏振方向时，接收器将同时接收到快、慢弯曲波，其相速度与指向角无关.

4.3     本文中的计算模型具有很好的拓展性，可以很方便地修改以求解类似的2.5维模型.

致 谢 感谢张秀梅、何晓副研究员对本工作的帮助.感谢国家自然科学基金对本工作提供资助. 附录A

附录B

参考文献
[1] Auld B A. 1973. Acoustic fields and waves in solids[M]. New York: Wiley.
[2] Chi S H, Tang X M. 2006. Stoneley-wave speed modeling in general anisotropic formations[J]. Geophysics, 71(4): F67-F77.
[3] Drossaert F H, Giannopoulos A. 2007. Complex frequency shifted convolution PML for FDTD modelling of elastic waves[J]. Wave Motion, 44(7-8): 593-604.
[4] Ekstrom M P. 1995. Dispersion estimation from borehole acoustic arrays using a modified matrix pencil algorithm[C].//Proceedings of the 29th Asilomar Conference on Signals, Systems and Computers. Washington, DC, USA: IEEE Computer Society, 449-453.
[5] Ellefsen K J. 1990. Elastic wave propagation along a borehole in an anisotropic medium[D]. Cambridge, Massachusetts: Massachusetts Institute of Technology.
[6] Ellefsen K J, Cheng C H, Toksäz M N. 1991a. Applications of perturbation theory to acoustic logging[J]. J. Geophys. Res., 96(B1): 537-549.
[7] Ellefsen K J, Cheng C H, Toksäz M N. 1991b. Effects of anisotropy upon the normal modes in a borehole[J]. J. Acoust. Soc. Am., 89(6): 2597-2616.
[8] He X, Hu H S, Guan W. 2010. Fast and slow flexural waves in a deviated borehole in homogeneous and layered anisotropic formations[J]. Geophys. J. Int., 181(1): 417-426.
[9] Kimball C V, Marzetta T L. 1984. Semblance processing of borehole acoustic array data[J]. Geophysics, 49(3): 274-281.
[10] Komatitsch D, Martin R. 2007. An unsplit convolutional perfectly matched layer improved at grazing incidence for the seismic wave equation[J]. Geophysics, 72(5): SM155-SM167.
[11] Kurkjian A L. 1985. Numerical computation of individual far-field arrivals excited by an acoustic source in a borehole[J]. Geophysics, 50(5): 852-866.
[12] Kuzuoglu M, Mittra R. 1996. Frequency dependence of the constitutive parameters of causal perfectly matched anisotropic absorbers[J]. IEEE Microwave and Guided Wave Letters, 6(12): 447-449.
[13] Leslie H D, Randall C J. 1992. Multipole sources in boreholes penetrating anisotropic formations: Numerical and experimental results[J]. J. Acoust. Soc. Am., 91(1): 12-27.
[14] Li Y F, Li G F, Wang Y. 2010. Application of convolution perfectly matched layer in finite element method calculation for 2D acoustic wave[J]. Acta Acustica (in Chinese), 35(6): 601-607.
[15] Li Y F, Matar B O. 2010. Convolutional perfectly matched layer for elastic second-order wave equation[J]. J. Acoust. Soc. Am., 127(3): 1318-1327.
[16] Lin W J, Zhang C Y, Zhang H L, et al. 2004. The experimental study of acoustic field in an axial asymmetric borehole[J]. Acta Acustica (in Chinese), 29(4): 319-322.
[17] Lin W J, Wang X M, Zhang H L. 2006. Acoustic wave propagation in a borehole penetrating an inclined layered formation[J]. Chinese J. Geophys.(in Chinese), 49(1): 284-294.
[18] Nolte B, Rao R, Huang X J. 1997. Dispersion analysis of split flexural waves[R]. Annual report of borehole acoustics and logging reservoir delineation consortia, Massachusetts Institute of Technology.
[19] Norris A N, Sinha B K. 1993. Weak elastic anisotropy and the tube wave[J]. Geophysics, 58(8): 1091-1098.
[20] Roden J A, Gedney S D. 2000. Convolution PML (CPML): An efficient FDTD implementation of the CFS-PML for arbitrary media[J]. Microwave and Optical Technology Letters, 27(5): 334-339.
[21] Sinha B K, Norris A, Chang S. 1994. Borehole flexural modes in anisotropic formations[J]. Geophysics, 59(7): 1037-1052.
[22] Sinha B K, Şimşek E, Liu Q H. 2006. Elastic-wave propagation in deviated wells in anisotropic formations[J]. Geophysics, 71(6): D191-D202.
[23] Tang X M, Zheng C H. 2004. Quantitative borehole acoustic methods[M]. Zhao X M Trans. Beijing: Gulf Professional Publishing.
[24] Thomsen L. 1986. Weak elastic anisotropy[J]. Geophysics, 51(10): 1954-1966.
[25] Wang R J, Qiao W X, Jü X D, et al. 2013. Experimental study of the acoustic field in the borehole surrounded by HTI formations excited by dipole sources with different orientations[J]. Chinese J. Geophys. (in Chinese), 56(2): 707-717, doi: 10.6038/cjg20130235.
[26] White J E. 1983. Underground sound: Application of seismic waves[M]. Amsterdam: Elsevier Science Ltd.
[27] Yan S G, Song R L, Lü W G, et al. 2011. Numerical simulation of acoustic fields excited by cross-dipole source in deviated wells in transversely isotropic formation[J]. Chinese J. Geophys. (in Chinese), 54(9): 2412-2418, doi: 10.3969/j.issn.0001-5733.2011.09.025.
[28] Zhang B X, Wang K X. 1996. Theoretical study of perturbation method for acoustic multipole logging in anisotropic formation[J]. J. Acoust. Soc. Am., 99(5): 2674-2685.
[29] Zhang B X, Wang K X. 1998. Status quo and development of the acoustic multipole logging in anisotropic formation[J]. Progress in Geophysics (in Chinese), 13(2): 1-14.
[30] Zhang B X, Wang K X. 2000. Theoretical study of multipole acoustic logging in anisotropic two-phase medium formation[J]. Chinese J. Geophys. (in Chinese), 43(5): 707-718.
[31] Zhang H L, Lin W J, Wang X M. 2011. Numerical studies on the acoustic field generated by a dipole source in noncircular pipes[J]. AIP Conference Proceedings, 1335: 1424-1431.
[32] Zhang L X, Fu L Y, Pei Z L. 2010. Finite difference modeling of Biot's poroelastic equations with unsplit convolutional PML and rotated staggered grid[J]. Chinese J. Geophys. (in Chinese), 53(10): 2470-2483, doi: 10.3969/j.issn.0001-5733.2010.10.021.
[33] Zhu Z Y, Cheng C H, Toksäz M N. 1994. Experimental study of the propagation of flexural waves in an anisotropic borehole model[J]. Chinese J. Geophys. (in Chinese), 37(6): 811-819.
[34] 李义丰, 李国峰, 王云. 2010. 卷积完全匹配层在两维声波有限元计算中的应用[J]. 声学学报, 35(6): 601-607.
[35] 林伟军, 张澄宇, 张海澜,等. 2004. 非轴对称结构井孔声场的实验研究[J]. 声学学报, 29(4): 319-322.
[36] 林伟军, 王秀明, 张海澜. 2006. 倾斜地层中的井孔声场研究[J]. 地球物理学报, 49(1): 284-294.
[37] 唐晓明, 郑传汉. 2004. 定量测井声学[M]. 赵晓敏译. 北京: 石油工业出版社.
[38] 王瑞甲, 乔文孝, 鞠晓东,等. 2013. 不同偏振方向的偶极子声源在HTI地层井孔中激发声场的实验研究[J]. 地球物理学报, 56(2): 707-717, doi: 10.6038/cjg20130235.
[39] 阎守国, 宋若龙, 吕伟国,等. 2011. 横向各向同性地层斜井中正交偶极子激发声场的数值模拟[J]. 地球物理学报, 54(9): 2412-2418, doi: 10.3969/j.issn.0001-5733.2011.09.025.
[40] 张碧星, 王克协. 1998. 各向异性介质地层中多极源声波测井的现状与发展[J]. 地球物理学进展, 13(2): 1-14.
[41] 张碧星, 王克协. 2000. 各向异性双相介质中多极源声波测井理论研究[J]. 地球物理学报, 43(5): 707-718.
[42] 张鲁新, 符力耘, 裴正林. 2010. 不分裂卷积完全匹配层与旋转交错网格有限差分在孔隙弹性介质模拟中的应用[J]. 地球物理学报, 53(10): 2470-2483, doi: 10.3969/j.issn.0001-5733.2010.10.021.
[43] 朱正亚, 郑传汉 Toksöz M N. 1994. 弯曲模式波在各向异性井孔中传播的实验研究[J]. 地球物理学报, 37(6): 811-819.