10 各向异性介质技术

10 各向异性介质技术
利用地震资料研究裂隙裂缝发育的方向和密度意义重大。

对于油气勘探而言，碳酸盐岩是一个有利的高产油气层，世界上约有60％左右的油气来自碳酸盐岩储层，而碳酸盐岩储层与裂隙裂缝的关系极为密切。

对于煤矿开采而言，研究裂隙裂缝的作用更为重要，主要表现在煤层底板突水和瓦斯突出两个方面。

华北大部分矿区的煤系地层基底为奥陶系灰岩，区内张裂性、张剪性断裂及陷落柱非常发育，奥灰水往往借助于小断层或岩溶陷落柱等导水通道突破煤层底板涌入工作面，造成矿井涌水量的增加甚至淹井的煤矿灾害，简称“水害”。

瓦斯突出是指煤矿生产过程中，从煤层、岩层及采空区放出的各种有害气体在工作面上富集并涌出，从而引起瓦斯爆炸的煤矿灾害，简称“火灾”。

无论是“水害”还是“火灾”，其罪魁祸首是岩层中的裂隙裂缝。

由于裂隙裂缝是水及瓦斯富集、存储、运移的场所，因此查明采区内断层、裂隙裂缝的分布有利于预防煤层底板突水和瓦斯突出，直接涉及到煤矿的安全生产。

大量的研究工作和观测数据表明，含裂隙裂缝介质的性质可以用各向异性介质理论进行解释，而传统的地震理论仅研究各向同性介质。

本章讨论各向异性介质中弹性波传播理论的意义也在于此。

10.1 各向异性介质
各向异性现象在地球介质中是普遍存在的。

广义上讲，当介质的特性在同一点处随方向发生变化时，则认为介质是各向异性介质。

10.1.1 地震各向异性
在地震勘探中，各向异性是指在地震波长的尺度下介质弹性特征随方向发生变化。

图10—1给出各向同性介质与各向异性介质的地震波速度变化。

(a)各向同性介质 (b)各向异性介质
图10—1 各向同性介质与各向异性介质的地震波速度变化一般地说，引起地震各向异性的主要因素为：
(1) 结构各向异性(如薄的层状结构，见图10—2)；
(2) 地层中方向应力导致的各向异性；
(3) 岩性各向异性(如颗粒的定向排列)；
(4) 地层中岩石晶体定向排列导致的各向异性；
(5) 岩石定向裂隙裂缝导致的各向异性，见图10—3。

图10—2 薄层状结构引起的各向异性
图10—3 岩石定向裂隙裂缝导致的各向异性
10.1.2 各向异性介质的类型
各向异性介质按其弹性性质变化的程度进行分类。

(1) 极端各向异性介质
如果介质中任一点处沿任意方向的弹性性质都是不同的，则这种介质称为极端各向异性介质，具有21个独立弹性参数。

(2) 正交各向异性介质
如果介质中存在一个平面，在平面对称的方向上弹性性质是相同的，则该平面称为弹性对称面，垂直弹性对称面的方向称为弹性主方向。

如果介质中有三个相互正交的弹性对称面，且它们的弹性主方向上的弹性性质互不相同，则这种介质称为正交各向异性介质，具有9个独立弹性参数。

(3) 横向各向同性介质
如果介质中存在一个弹性对称面，在平面内沿所有方向的弹性性质都是相同的，而垂直平面各点的轴向都是平行的，则称该平面为各向同性面，垂直各向同性面的轴为对称轴。

具有各向同性面的介质称为横向各向同性介质，简称TI(Transverse Isotropy)介质。

当TI介质的对称轴垂直时，称其为VTI(Transverse Isotropy with a Vertical axis of symmetry)介质，即具有垂直对称轴的横向各向同性介质。

它近似地表示水平层
状介质周期性沉积的薄互层各向异性介质，因此VTI介质也称为PTL(Periodic Thin-Layer)各向异性介质。

图10—4为PTL介质示意图。

当TI介质的对称轴水平时，称其为HTI(Transverse Isotropy with a Horizontal axis of symmetry)介质，即具有水平对称轴的横向各向同性介质。

HTI介质近似地表示空间排列垂直裂隙而引起的各向异性，也称为扩容各向异性介质，简记EDA (Extensive Dilatancy Anisotropy)介质。

EDA介质是典型的方位各向异性介质，图10—5为EDA介质示意图。

图10—4 PTL介质示意图图10—5 EDA介质示意图根据各向异性介质的对称特性，各向异性介质可分为8类(Crampin,1981)。

(1) 三斜对称各向异性介质
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
=
66
65
64
63
62
61
56
55
54
53
52
51
46
45
44
43
42
41
36
35
34
33
32
31
26
25
24
23
22
21
16
15
14
13
12
11
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
C (10-1)
(2) 单斜对称各向异性介质
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
=
66
64
55
53
52
51
46
44
35
33
32
31
25
23
22
21
15
13
12
11
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
c
C (10-2)
⎥⎥⎥⎥⎥⎥
⎥⎥⎦
⎤
⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡=665544333231
2322211312110
0000000000000000000c c c c c c c c c c c c C (10-3) (4) 四方对称各向异性介质
⎥⎥⎥⎥⎥⎥
⎥⎥⎦
⎤
⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡--=6616
16
4444331313
1613111216131211
00000000000000000c c c c c c c c c c c c c c c c C (10-4) ⎥⎥⎥⎥⎥⎥
⎥⎥⎦
⎤
⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=664444331313
131112131211
0000000000000000000c c c c c c c c c c c c C (10-5) (5) 三角对称各向异性介质
⎥⎥
⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡-----=)(5.000000000
0000
1211142514
44252525
441414331313
25141311122514131211
c c c c c c c c c c c c c c c c c c c c c c c c c C (10-6)
⎥⎥
⎥⎥
⎥⎥⎥⎥⎦
⎤
⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡---=)(5.00
000
00000
00000
0121114
14
44441414331313
1413111214131211
c c c c c c c c c c c c c c c c c c c C (10-7)
⎥⎥
⎥⎥
⎥⎥⎥⎥⎦
⎤⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡-=)(5.00
000
0000000
000000
0012114444331313
131112131211c c c c c c c c c c c c c C (10-8) (7) 立方对称各向异性介质
⎥⎥⎥⎥⎥⎥
⎥⎥⎦
⎤⎢⎢⎢⎢⎢
⎢⎢⎢⎣⎡=444444111212
121112121211
0000
000000000000000c c c c c c c c c c c c C (10-9) (8) 各向同性介质
⎥
⎥⎥⎥⎥⎥
⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣
⎡=444444111212
121112121211
00000
000000000000000c c c c c c c c c c c c C (10-10) 10.2 各向异性介质中的弹性理论
弹性波在各向异性介质中与在各向同性介质中遵循不同的传播规律，满足不同的波动方程，具有不同的波型、极化、相速度、群速度等波动特征。

10.2.1 各向异性介质中地震波的传播
对于各向异性介质，通常利用广义虎克定律来描述应力与应变之间的关系，即介质的本构方程
∑∑====313
1
)3,2,1,,,(k l kl
ijkl ij l k j i c εσ (10-11)
其中，ij σ为应力张量，kl ε为应变张量，ijkl c 为弹性劲度(elastic stiffness)常数，简称弹性系数，也称刚度张量或刚度矩阵。

因为式(10-11)中包含九个方程(下标ij 所有可能的组合)，每个方程中有九个应变变量，故有81个弹性系数。

当应力满足对称(ij σ=ji σ)时，可以增加条件
jilk ijlk jikl ijkl c c c c ===
这样81个弹性系数减少到36个，广义虎克定律有如下形式
⎥⎥⎥⎥⎥⎥
⎥⎥⎦⎤
⎢⎢⎢⎢⎢⎢⎢⎢⎣
⎡⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎣
⎡xy zx yz zz yy xx xy zx yz zz yy xx c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c c εεεεεεσσσσσσ6665
64
63
62
61565554535251464544434241363534333231262524232221161514131211 (10-12) 其中，ij c (6,5,4,3,2,1,=j i )弹性系数。

由于弹性系数是应变的单值函数，即ji ij c c =。

因此，描述一个复杂的弹性介质需要21个弹性系数。

10.2.2 TI 介质中的波动方程
TI 介质的弹性沿横向是各向同性的，沿纵向则是变化的，相当于薄互层(薄层厚度<<地震波长)或水平裂隙情况。

对于VTI(横向各向同性)介质，yoz xoz xoy ,,面都是对称的，所以有
04645363526251615========c c c c c c c c ，056342414====c c c c
于是弹性系数矩阵为
⎥⎥⎥⎥⎥⎥
⎥⎥⎦
⎤
⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡665544333231232221131211
0000000000000000000c c c c c c c c c c c c (10-13) 其中322331132112,,c c c c c c ===。

因为它在y x ,方向完全相同，又有
⎪⎪⎩⎪⎪⎨
⎧====1211231366
552211c c c c c c c c 66
2c + 所以，对于横向各向同性介质，只有6644331311,,,,c c c c c 这5个独立的弹性系数。

将式(10-13)代入应力运动方程
)3,2,1(3
122=∂∂=-∂∂∑=i t u f x k i i k
ik ρσ (10-14)
则得
)3,2,1(3
1
,,2
22=∂∂=-∂∂∂∑
=i t
u f x x u c q p k i i k p q
ikpq
ρ (10-15)
其中，i f 为外界作用力，321,,u u u 为对应质点z y x ,,三个方向的位移分量。

在二维介质中，设0=i f ，用式(10-13)对式(10-15)简化，即可得二维横向各向同性介质中的弹性波动方程
z x u c c z
u c x u c t u ∂∂∂++∂∂+∂∂=∂∂32441321
24421211212)(ρ (10-16)
2
2
24422266222z u c x u c t u ∂∂+∂∂=∂∂ρ (10-17) z x u c c z
u c x u c t u ∂∂∂++∂∂+∂∂=∂∂1
244132323323244232)(ρ (10-18)
上述方程组可分解为一个独立的波动方程(10-17)及一个波动方程组(10-16)和(10-18)。

方程(10-17)是横向各向同性SH 波波动方程；方程组(10-16)和(10-18)则是准P 波和SV 波波动方程，它们是偶合在一起的，即在横向各向同性介质中，P 波传播可引起SV 波，SV 波传播也可引起P 波，它们联合称为横向各向同性准P —SV 波动方程。

图10—6是点震源S 激发，地震波分别在各向同性和TI 介质无限空间中传播的数值模拟结果。

从图10—6中可以看出，波场有以下特征：
(1) 与均匀各向同性介质中只产生P波不同，在TI介质中，P波和SV波是偶合的，即在无限均匀TI介质中，不仅产生P波，而且产生SV波。

(2) P波的波前面不再是一个圆(二维情况)，而是一个椭圆。

说明在TI介质中的同一点，波沿不同方向的传播速度是不同的。

波前面不总是与波的传播方向(射线)垂直，这种波又称准纵波。

(3) SV波的波前面更加复杂，其偏振方向也不再与传播方向(射线)垂直，但与P波的偏振方向仍然正交，这种波又称准SV波。

图10—6 TI介质的数值模拟结果
(a) 各向同性；(b) 各向异性
10.3 地震资料裂隙裂缝检测技术
目前，国内外油气勘探中利用地震资料研究各向异性、检测裂隙裂缝系统的主要方法有三类。

(1) 多波多分量裂隙裂缝检测技术
理论上，含有裂缝的介质是各向异性介质，地震多波多分量是研究各向异性介质中波传播规律的理想方法。

但是，由于多波多分量地震采集的高成本和转换波处理、解释的复杂化，使得这种技术的推广应用受到极大限制。

(2) S波裂隙裂缝检测技术
S波法裂隙裂缝检测主要有质点振动矢量分析、旋转分析和协方差矩阵分析等。

S波分裂是地下各向异性存在的最好证据，利用S波分裂检测裂隙裂缝已经成为研究裂缝性储层的一种可靠技术手段。

但是，S波地震采集的成本高，而且难以得到高质量的S波数据，这就限制了S波方法的应用。

(3) P波裂隙裂缝检测技术
由于P波地震勘探成本低，从20世纪90年代起，地球物理学家把目光转向P 波勘探，用P波代替S波/转换波检测裂隙裂缝成为一个重要研究课题。

目前，国内外学者在研究利用反射P波进行裂隙裂缝检测方面做了大量工作，其研究重点是利用叠前P波数据的振幅、速度和旅行时差随炮检距或方位角的变化规律来识别裂缝，其研究方法包括裂缝介质中地震P波波场特征的物理模型研究和数学模拟计算。

研究结果表明，反射P波对裂缝性地层所表现出的方位各向异性特征很敏感。

本节的重点讨论第三种方法。

10.3.1 垂直裂缝介质方位各向异性特征的物理模拟
所用的垂直裂缝模型是由一组平行排列的有机玻璃片叠合而成，利用有机玻璃片之间的缝隙来模拟平行排列的裂缝。

有机玻璃片四周用螺杆固定并将其放置于水槽中，有机玻璃片之间的缝隙充满了水，超声波观测在其上方进行。

将偏移距固定，并在裂缝模型上部过中心点使测线方位进行
360旋
180或 转，则能够通过物理模型实验记录直观地观测裂缝介质的P波方位各向异性特征，见图10—7。

图10—7 (a)裂缝物理模型；(b)固定偏移距观测系统示意图
以图10—7所示的观测系统对垂直裂缝模型进行测线方位 360旋转采集数据。

测线方位每旋转 15采集一道数据，共采集25道数据，实验记录见图10—8。

图10—8 垂直裂缝介质方位各向异性特征物理模型实验结果
图10—8中A 同相轴是不同方位测线所接收的裂缝层顶界面的反射，B 同相轴为射线经过裂缝层后，不同方位测线所接收的裂缝层底界面的反射。

A 同相轴的反射时间与振幅基本一致，而B 同相轴的反射振幅、反射时间和传播速度曲线均表明，反射波在通过裂缝体后呈现出方位各向异性特征，具体表现为：
(1) 当测线方位与裂缝走向平行时(夹角为 0)，反射波振幅和速度最大；随着测线方位与裂缝走向之间的夹角增大，反射波振幅和速度逐渐减小，当夹角为 90时达到最小；此后反射波振幅和速度随着夹角增大而逐渐增大，当夹角为 180时又达到最大，变化周期为 180。

(2) 当测线方位与裂缝走向平行时，反射时间最小；随着测线方位与裂缝走向之间的夹角增大，反射时间逐渐增大，当夹角为 90时达到最大；此后反射时间随着夹角增大又逐渐减小，当夹角为 180时达到最小，变化周期仍为 180。

分析这些地震属性的变化规律，反射波振幅和速度可简化表达为
ϕϕ2cos )(B A F += (10-19)
式中，ϕ是炮检方向与裂缝走向的夹角，A 是与偏移距有关的偏置因子，B 是与偏移距和裂缝特征有关的调制因子。

式(10-19)给出了反射P 波振幅、速度等的方位各向异性特征，可用图10—9来近似表示。

图10—9 垂直裂缝介质的P 波方位各向异性示意图
10.3.2 利用P 波方位属性确定岩溶裂隙带的空间分布
假设从某个CMP 道集中，抽取三个方位角ϕ、αϕ+、βϕ+的资料，对应的振幅分别为)(ϕF 、)(αϕ+F 、)(βϕ+F ，注意α和β是已知的。

根据式(10-19)，可以得到
ϕϕ2cos )(B A F += (10-19)
)(2cos )(αϕαϕ++=+B A F (10-20)
)(2cos )(βϕβϕ++=+B A F (10-21)
现在的任务是根据)(ϕF 、)(αϕ+F 、)(βϕ+F 、α和β求取ϕ。

用式(10-19)减去式(10-20)，并利用三角判别式
)2sin()2sin(2cos cos 1
22
121θθθθθθ-+=-
可以得到
ααϕαϕϕsin )2sin(2)()(+=+-B F F (10-22)
类似地用式(10-19)减去式(10-21)，可以得到
ββϕβϕϕsin )2sin(2)()(+=+-B F F (10-23)
解式(10-22)和式(10-23)，得到
πα
αβϕϕββαϕϕβαϕϕαβϕϕϕn F F F F F F F F ±+--+-+--+-=]cos sin ))()((cos sin ))()((sin ))()((sin ))()((arctan[2122 (10-24)
其中， ,2,1,0=n 。

利用式(10-24)，可以通过CMP 道集中三个方位角的P 波资料计算裂缝方向ϕ。

在求取ϕ后，也可以从式(10-19)、式(10-20)和式(10-21)中求得偏置因子A 和调制因子B 。

地震模拟研究结果表明，裂缝密度(每单位长度的裂缝数)和B/A 的值成正比，即裂缝密度越高，B/A 的值越高。

因此，B/A 的计算提供了一种定量测量裂缝密度的方法。

对于叠前资料，可以对每个偏移距使用式(10-24)求取所有偏移距的裂缝方向，然后进行加权平均得到总的裂缝方向。

对于叠加资料，可以利用三个方位测线的叠加数据来计算裂缝方向，值得注意的是叠加数据已经把所有偏移距进行了校正。

综上所述，利用多方位的P 波叠前或叠加资料可以确定裂缝的方位角和密度。

因此，利用P 波方位属性可以确定岩溶裂隙带的空间分布。

利用地震P 波方位属性研究岩溶裂隙带的基本步骤是：
(1) 为了增加用于裂缝属性计算的有效覆盖次数，将4个面元(m 1010⨯)的数据形成超级CMP 道集；
(2) 然后按30°的方位角增量抽取6个方位角道集；
(3) 对6个方位角道集进行速度分析、NMO 校正、叠加和偏移，得到6个方位偏移数据体；
(4) 对6个方位偏移数据体进行波阻抗反演，得到6个方位波阻抗数据体；
(5) 从6个方位偏移数据体和方位波阻抗数据体中提取与岩溶裂隙密度有关的地震属性参数；
(6) 利用从方位偏移数据体、方位相干/方差数据体、方位波阻抗数据体中提取的地震属性(主要包括振幅、频率、分形参数、相干/方差系数、波阻抗、速度)对下组煤层到太原组灰岩、奥陶系灰岩之间的岩溶裂隙带进行解释。

研究不同方位角P波地震振幅、速度、波阻抗随入射角变化关系(AVA、VAV、IPAV)，利用AVA 、VAV、IPAV方法定量计算裂缝属性。

把上述观点进行延伸，研究多个地震属性(Seismic Attributes)随入射角变化的规律，利用地震属性参数随方位角变化的特征提取裂隙属性，从而确定岩溶裂隙带的空间分布，这种方法称为SAVA(地震属性随入射角变化)方法。

10.4 应用实例
利用淮南张集煤矿西三采区三维地震P波资料进行处理和解释，观测系统为束状8线8炮中点激发，8
m10
10⨯，
48⨯道接收，覆盖次数为24次，CMP网格为m
在偏移前通过插值加密为m
5⨯。

该区的三维地震资料采集过程中存在问题，由
m5
于采用的观测系统方位角较窄，导致了各个方位角上的炮检距和覆盖次数分布不十分均匀，因此在抽取方位道集时，要尽量使各个方位上的炮检距和覆盖次数分布均匀，减小其对处理结果的影响。

在地震资料处理过程中，先抽取方位道集，并对其进行相应的处理，得到6个方位偏移数据体；然后利用地震属性技术，针对某一目标层位，从时间属性、振幅属性、频率属性和相位属性中分别提取了波峰相位时间、时域平均能量、峰值频率和平均频率相位四种属性；其次用这四种地震属性分别计算出裂隙的发育方向和密度，得到以下主要成果。

图10—10(a)为利用波峰相位时间属性识别裂隙带，图10—10 (b)为利用时域平均能量属性识别裂隙带，图10—10 (c)为利用峰值频率属性识别裂隙带，图10—10 (d)为利用平均频率相位属性识别裂隙带。

图10—10中，直线方向表示裂隙的方向，直线长度表示裂隙的密度。

综合各种地震属性预测结果得到最终的预测成果图，见图10—11。

将图10—11的预测成果与图10—12的地震解释成果比较后，发现两者之间非常吻合，这就验证了利用地震P波方位属性预测裂隙发育带的正确性。

结果还表明，利用平均频率相位属性识别裂隙发育带的效果最为明显，利用峰值频率的效果次之，而利用波峰相位时间和时域平均能量的效果最弱。

(a) (b)
(c) (d)
图10—10 利用地震P波属性识别裂隙发育带
(a)波峰相位时间；(b) 时域平均能量；(c)峰值频率；(d)平均频率相位
图10—11 最终预测成果图图10—12 地震解释成果图。