基于岩石物理模拟与声学实验识别孔隙—裂隙充填型水合物

2020-12-31景鹏飞胡高伟卜庆涛陈杰万义钊毛佩筱

海洋地质与第四纪地质 2020年6期

景鹏飞，胡高伟，卜庆涛，陈杰，万义钊，毛佩筱

1. 自然资源部天然气水合物重点实验室，中国地质调查局青岛海洋地质研究所，青岛 266071

2. 中国地质科学院，北京 100037

3. 海洋国家实验室海洋矿产资源评价与探测技术功能实验室，青岛 266071

天然气水合物主要赋存于大陆边缘和大陆永久冻土区[1-2]，以层状、块状、分散状、胶结状、瘤状等形式赋存于地层中[3]，可归结为两大类，即孔隙充填型和裂隙充填型。孔隙充填型水合物替代了孔隙中的流体，赋存于沉积物孔隙中[4]；裂隙充填型水合物则以纯水合物形式与沉积物互相接触，一般存在于沉积物裂隙、地层裂隙、断层或泥火山等特殊地质体中。因此，由于赋存状态和岩石物理性质不同，针对孔隙充填型、裂隙充填型水合物往往采用不同的方法进行资源评估。然而，在野外勘探过程中，如何通过简便的方法识别水合物类型，尚无较为成熟的识别方法。

已有钻探结果为认识孔隙充填型和裂隙充填型水合物打下了良好基础[5-6]。中国、美国、印度、韩国等钻探表明，中国南海神狐海域天然气水合物以孔隙充填为主[7]，而美国墨西哥湾[8]、韩国Ulleung盆地[9]和印度K-G 盆地[10-11]的细粒沉积物中均发现了裂隙充填型水合物。2006 年，印度国家天然气水合物计划（NGHP01）在K-G 盆地通过测井和钻探发现大量裂隙充填型水合物富集，测井和红外成像显示水合物多呈结核状、脉状、瘤状，充填于高角度（60°～120°）裂隙以及近水平层理中，其中充填水合物的裂隙层位厚度最大可达135 m，长度达数米以上[1,5,11]。孔隙充填型水合物和裂隙充填型水合物最主要的区别在于非微观尺度研究中，孔隙充填型水合物通常被看作是各向同性介质，而裂隙充填型水合物则是各向异性介质[12]。中国地质调查局广州海洋地质调查局2013 年在中国南海北部海域开展的第二次水合物勘探GMGS2，在13 个站位布设了21 个钻井，并取得了保压岩心[13]，其中5 个站位发现了水合物，这一航次中同时发现了孔隙充填型和裂隙充填型水合物。

由于孔隙充填型和裂隙充填型水合物的形成机理和地球物理响应都存在差异，因此，在无钻探井位的区域，准确识别两类水合物在沉积层中的分布形态，对于水合物资源评价、钻采安全以及环境等方面都至关重要[14]。目前为止，关于孔隙充填型和裂隙充填型水合物识别的文献相对较少。戴丹青等[15]利用中观衰减理论，用斑块饱和模型和双孔隙模型针对水合物在地层中的赋存状态，研究了地震频带内BSR 界面反射系数的变化特征，并将研究结果应用到存在孔隙充填型水合物的南海神狐海域和存在裂隙充填型水合物的布莱克海台，模拟结果与实际情况吻合较好。王吉亮等[3]为识别印度东海岸NGHP01-10D 井的水合物赋存状态，利用等效介质理论（EMT）和层状介质横向各向异性模型分别计算了水合物的饱和度，发现层状介质横向各向异性模型计算的结果与实际情况吻合较好，且存在的裂隙多为高角度裂隙。Lee 等[16]利用各向同性模型和横向各向同性理论计算了韩国东海Ulleung 盆地UBGH2-3B 井中水合物层的电阻率、纵波速度和水合物饱和度，研究结果表明井中局部范围水合物的赋存类型为裂隙充填型，气烟囱上半部分（5～14 mbsf）地层中主要赋存裂隙角度约7°的裂隙充填型水合物，下半部分（45～103 mbsf）中赋存高角度裂隙充填型水合物。Liu 等[17]提出利用纵波速度和电阻率交汇图的斜率来识别两类水合物。Sriram 等[18]对含孔隙充填型/裂隙充填型水合物地层的BSR 进行了AVA（振幅随入射角的变化）模拟，通过模拟将两者的AVA 类别进行归类，从而达到水合物在储层中赋存形态的识别。Liu 等[19]通过岩石物理模拟分别计算了孔隙充填型和裂隙充填型水合物的密度和速度，提出用密度速度属性识别孔隙充填型及裂隙充填型水合物的方法。

迄今为止，对水合物赋存形态识别的研究多以模拟手段为主来计算两种类型水合物储层的电阻率、纵波速度、地震衰减/反射系数及AVA/AVO 属性。截至目前，未发现将模拟结果与实验实测数据进行对比分析的案例。本文针对孔隙充填型水合物选取了 EMT（A、B 模式）、BGTL 和 STPE 模型，针对裂隙充填型水合物（本文仅以裂隙空间发育的纯水合物为研究对象）（图1a）选取了横向各向同性理论（TIT），模拟了两种类型水合物的密度和声波速度与水合物体积分数间的变化特征，并通过声学实验获取了孔隙充填型和裂隙充填型水合物的纵波速度特征。

1 天然气水合物岩石物理模型选取与计算

1.1 孔隙充填型水合物岩石物理模型选取与计算

Liu 等[18]、胡高伟[21]、王吉亮等[3]等关于孔隙充填型水合物岩石物理模拟的研究结果表明，等效介质理论（EMT）、BGTL 理论和简化三相方程（STPE）的计算结果与实验测试波速的吻合程度相比其他岩石物理模型更具优势。因此，本文针对孔隙充填型水合物选取了EMT-A/B、BGTL 和STPE 岩石物理模型，计算了含孔隙充填型水合物的南海沉积物（孔隙度为0.4）的密度和纵波速度。表1 为模拟计算所用的南海沉积物的弹性参数及含量表，模型计算公式见附录。

表1 南海泥质黏土的体积模量（K）、剪切模量（G）、密度（ρ）和含量（C）[22]Table 1 The bulk modulus （K）, shear modulus （G）, density （ρ） and content （C） of the silty clay in the South China Sea [22]

1.2 孔隙充填型水合物模型计算结果

由于EMT-A 和BGTL 的密度计算公式相同，因此，两者模拟的密度曲线重合（图2），而EMT-B 模式中，水合物形成后充当骨架的一部分，同时使得沉积物孔隙度减小，该模式计算密度时对沉积物孔隙度进行了修正，但最终计算结果与STPE 相同。各模式计算的孔隙充填型水合物介质的密度尽管有差异，但最大差值仅为1.25%。

图2 EMT、BGTL、STPE 模型计算的含孔隙充填型水合物的密度Fig.2 The density of pore-filling hydrate calculated by EMT、BGTL and STPE models

3 类岩石物理模型计算的纵波速度均随水合物体积分数的增加而增大（图3），但不同模型的波速计算方式不同，波速初始起跳值也不同。当沉积物孔隙中未生成水合物（体积分数为0）时，EMT-A/B模式计算的沉积物纵波速度大于BGTL 和STPE模型计算的波速，3 类模型计算的纵波速度大小分两种情况：①当水合物体积分数小于17%时，波速大小为：EMT-B＞EMT-A＞STPE＞BGTL；②水合物体积分数大于17%时，EMT-B＞STPE＞EMT-A＞BGTL。EMT-B 模式计算的波速值最大，BGTL 的计算值最小，图3 中横坐标体积分数40%对应水合物饱和度为100%。

图3 EMT、BGTL、STPE 模型计算的含孔隙充填型水合物的纵波速度Fig.3 The P-wave velocity of pore-filling hydrate calculated by EMT、BGTL and STPE models

1.3 裂隙充填型水合物岩石物理模型与计算

裂隙充填型水合物模型用横向各向同性理论（TIT）进行计算。模型由两个端元组成（图4），端元Ⅰ为裂隙，100%由纯水合物充填，端元Ⅱ为饱和水的南海沉积物（泥质黏土，孔隙度为0.4），其弹性参数及含量见表1。端元Ⅱ的体积模量和剪切模量基于 Hill 平均方程（B-5）和 STPE 模型计算。

由横向各向同性理论计算的含裂隙充填型水合物的南海沉积介质的密度和纵波速度结果（图5）发现，随着水合物体积分数增大，沉积介质的总密度减小，且当模型中端元Ⅰ占比为100%时，裂隙中为纯水合物，此时介质密度最小，水合物的密度为0.93 g/cm3；纵波速度曲线最小值代表南海沉积物的模型理论速度（1.84 km/s），最大值代表纯水合物的速度（3.08 km/s），纵波速度随着水合物体积分数增大而增大。

图4 裂隙充填型水合物端元模型[23]Ⅰ为100%纯水合物充填的裂隙，Ⅱ为饱和水泥质黏土。Fig.4 Fracture-filling hydrate end member model[23]Ⅰ: Fracture filled by pure hydrate; Ⅱ: Silty clay with saturated water.

图5 横向各向同性理论计算的含裂隙充填型水合物的密度和纵波速度Fig.5 The density and P-wave velocity of fracture-filling hydrate calculated by layered media model

2 孔隙充填型、裂隙充填型水合物声学实验模拟

2.1 孔隙充填型水合物声学实验模拟

孔隙充填型水合物声学特性模拟实验装置如图6。实验之前需要先对弯曲元超声测速探头的固有声波走时和时域反射技术（TDR）测定沉积物中含水量的误差进行标定。速度测试时需要减去超声探头的固有声波走时，标定结果为10-3-0.01246和y=9.291x×10-4-9.812×10-3，弯曲元换能器中纵波的固有传播时间为5.56 μs，横波固有传播时间为 10.56 μs；TDR 实测值与理论值间的误差范围为0.17%～6.75%，由于TDR 测量的结果为一个场平均值，因此，TDR 系统测试精度满足实验要求。

实验模拟过程为：将饱和度水南海沉积物装入高压反应釜，并加入SDS（十二烷基硫酸钠）溶液（加快水合物的生成），通入纯甲烷气体（纯度99.99%）使釜内压力达到目标压力（约 7 MPa），并放置 24 h使甲烷气溶解到SDS 溶液中，随后开启制冷系统使反应釜内部温度到约1 ℃。随着水合物生成，由TDR 系统测试水合物的饱和度，弯曲元声波系统测试波速。

实验中，孔隙充填型水合物的密度由沉积物中矿物及生成的水合物的密度及其含量计算得来。纵波速度为实验测试结果（图7）。由于南海沉积物中短期很难生成高饱和度水合物，因此，反应釜中生成水合物的最大体积分数为9.6%（对应饱和度为24%）。基于实验波速数据发现，纵波速度随水合物体积分数呈正比趋势，当水合物体积分数为9.6%时，沉积介质纵波速度达到最大值2.06 km/s。

图6 天然气水合物声学模拟实验装置Fig.6 Acoustic experimental apparatus for natural gas hydrate

图7 南海沉积物中孔隙充填型水合物的纵波速度和密度Fig.7 P-wave velocity and density of pore-filling hydrate in deposits of the South China Sea

2.2 裂隙充填型水合物声学实验模拟

先利用实验装置（图8）进行裂隙充填型水合物实验合成，其过程为：利用雾化喷头将纯水雾化后喷洒到液氮中，合成所需量的冰粉，将冰粉与沉积物按照比例装入多孔塑胶管中，反应釜密封后通入甲烷气使釜内压力达到水合物合成所需目标压力，再将反应釜置于循环水域控温装置中进行降温（约1 ℃），直至水合物合成完毕。最后开釜取出样品，连同塑胶管一同置于岩心夹持器（图9）上，利用超声探头发射/接收声波，进行波速测试。

图8 裂隙充填型水合物合成实验装置Fig.8 Experimental apparatus for fracture-filling hydrate synthesis

不加入水合物催化溶液时，南海沉积物中短期很难生成孔隙充填型水合物，且长时间内合成的水合物饱和度也很低[21]。该次裂隙充填型水合物合成实验中未向沉积物中加入SDS 溶液，一个实验周期（约7 天）内沉积物孔隙中几乎不会有水合物生成，实验介质可视为纯裂隙充填型水合物。

含裂隙充填型水合物（裂隙为水平裂隙）的南海沉积介质的密度与纵波速度结果中（图10），密度计算公式为其中和分别为裂隙和沉积物的体积分数，为纯水合物的密度，计算取值为0.93 g/cm3，为南海泥质黏土沉积物的骨架密度（2.66 g/cm3），为沉积物孔隙度，为水的密度（1 g/cm3）；纵波速度为实验超声测试结果。当水合物体积分数为0 时，整个沉积介质的密度即为沉积物的骨架密度（2.66 g/cm3），随着水合物体积分数增大，整个沉积介质（水合物+沉积物）的密度呈减小趋势，而当水合物体积分数为100%时，介质的密度即为纯水合物的密度值（0.93 g/cm3）。水合物体积分数为0 和100%分别对应100%南海沉积物和100%纯水合物的波速，当水合物体积分数增大时，纵波速度也随之增大。

图9 岩心夹持测试波速装置Fig.9 Rock holding device for testing wave velocity

图10 含裂隙充填型水合物的南海沉积介质的密度与纵波速度Fig.10 Density and P-wave velocity of sediment medium of South China Sea with fracture-filling hydrate

3 孔隙充填型、裂隙充填型水合物识别方法

两种类型水合物的岩石物理模拟和实验纵波速度结果中（图11），波速均随水合物体积分数增大而增大。实验中生成的孔隙充填型水合物的最大体积分数为9.6%（对应饱和度为24%）。孔隙充填型水合物的纵波速度在体积分数小于4%（对应饱和度为10%）时低于岩石物理模拟的波速，推断造成该现象的主要原因是当南海沉积物孔隙空间生成的水合物饱和度较低时，水合物呈悬浮状分散于沉积物孔隙内，导致声波信号发生衰减，且这种悬浮状态对沉积物骨架支撑几乎不起作用，故当水合物体积分数（饱和度）较低时，实验模拟的波速低于岩石物理模拟值。而当体积分数大于4%时，实验波速趋向于BGTL 模型计算的理论波速。横向各向同性理论（TIT）计算的裂隙充填型（水平裂隙）水合物理论波速在水合物体积分数较大时（＞40%）与实验测试的纵波速度吻合较好，但当水合物体积分数较小时，两者间的最大差为0.27 km/s，整体吻合较好。

通过将密度与速度参数结合，计算了两种类型水合物的阻抗（ρVp）（图 12）和属性（图 13）。由于裂隙填充型水合物的实验阻抗斜率非严格负值，因此，基于阻抗无法进行两种类型水合物的判别。岩石物理模拟和实验测试的属性中，孔隙充填型水合物和裂隙充填型水合物曲线的斜率符号相反：孔隙充填型水合物的属性曲线斜率为正，而裂隙充填型水合物为负。孔隙充填型和裂隙充填型水合物的这一特征为两种类型水合物的判别提供了有效手段。

4 识别方法野外验证

中国第二次天然气水合物钻探航次（GMGS2）在中国南海珠江口盆地东部第16 井中同时发现了孔隙充填和裂隙充填型水合物，利用该井的密度和测井资料对属性识别水合物类型的方法进行验证。该16 井中孔隙充填和裂隙充填型水合物平均饱和度分别为42%和50%，其中孔隙充填和裂隙充填型水合物分别位于井深（泥线以下：BML）189～198 和 10～20 m[18,24]。

该井上部 10～20 m 验证结果（图 14）表明该范围水合物类型以裂隙充填型水合物为主，夹杂部分孔隙充填型水合物。井深10～11.6、14～18.5 m深度含水合物层与背景层的属性差值小于零，指示该范围内存在裂隙充填型水合物，12～13.6、18.5～20 m 深度属性差值大于零，代表赋存的水合物类型为孔隙充填型。井底部189～198 m 赋存的水合物类型以孔隙充填型为主（图15），在194～196 m 深度夹杂少量裂隙充填型水合物。从验证结果来看，基本与野外观测结果一致。

图11 岩石物理模拟和实验测试的纵波速度Fig.11 The P- wave velocity of petrophysical simulation and experimental test

图12 岩石物理模拟和实验测试的孔隙充填型和裂隙充填型水合物的阻抗Fig.12 The impedance of pore-filling and fracture-filling hydrates obtained by petrophysical simulation and experimental test

5 结论与展望

图13 岩石物理模拟和实验测试的孔隙充填型和裂隙充填型水合物的属性Fig.13 Theproperties of pore-filling and fracture-filling hydrates obtained by petrophysical simulation and experimental test

图14 GMGS2-16 井上部含裂隙充填型水合物储层属性差值结果Fig.14 Result of attribute difference of fracture-filling hydrate reservoir in the upper part range of MGMS2-16 Well

密度和纵波速度间的差异提供了识别孔隙充填型水合物和裂隙充填型水合物在沉积物中赋存形态的基础。通过结合岩石物理模型和实验模拟，获得了南海沉积物中孔隙充填型和裂隙充填型水合物的密度和纵波速度参数，测试了两种类型水合物的阻抗和属性。研究结果表明基于属性可以有效识别孔隙型充填和裂隙充填型水合物：对于含孔隙充填型水合物的南海沉积物，岩石物理模型的属性与实验测试的属性曲线都是正斜率；对于含裂隙充填型水合物的南海沉积物，岩石物理模型和实验测试的属性都为负斜率。对于含孔隙充填型水合物的南海沉积介质，当水合物体积分数大于4%时，岩石物理模型计算的属性与实验数据吻合较好，尤其是STPE 模型与实验结果吻合最佳，因此，STPE 模型在计算南海沉积物中水合物的饱和度时相比其他模型更具优势。对于裂隙充填型水合物，当水合物体积分数较低（＜40%）时，横向各向同性理论计算的属性和实验测试结果存在明显差异。因此，当水合物体积分数低于40%时，计算水合物的属性过程中需要对模型进行适当修正。

图15 GMGS2-16 井下部含孔隙充填型水合物储层属性差值结果Fig.15 Result of attribute difference of pore-filling hydrate reservoir in the lower part of GMGS2-16 Well

附录（岩石物理模型）

A.1 等效介质理论（EMT）

Helgerud 等[25-26]对水合物微观分布形态进行了研究，认为水合物与沉积物颗粒间的共生关系主要是充填于沉积物孔隙中，或胶结了沉积物骨架，与骨架呈支撑形态（图A-1）。对于模式A，水合物的生成替代了沉积物孔隙中的流体并占据了流体空间，因而该模式中，水合物也被认为是孔隙流体的一部分；模式B 中，水合物被当作沉积物骨架的一部分，孔隙空间中只有流体，该模式对沉积物骨架的密度、体积模量和剪切模量产生了影响[21]；模式C 中，孔隙度的降低等同于模式B，另一方面，岩石骨架的体积模量和剪切模量发生了改变。南海沉积物中孔隙充填型水合物模拟选用A（水合物充当流体部分）和B（水合物充当骨架部分）模式进行计算。

图A-1 水合物在沉积物中的 3 种赋存形态Fig.A-1 Three types of hydrate occurrence in sediment

（1）模式 A

沉积介质的纵波速度及密度计算公式如下：

式中，Ksat和Gsat分别是等效介质的体积模量和剪切模量分别为岩石固相和流体相的体积密度。φc为临界孔隙度，一般取 0.36～0.40[27]。当沉积物中充填体积模量为Kf的流体时，可以根据Gassmann 方程，通过如下公式得到沉积物的体积模量Ksat和剪切模量Gsat：

式中，Kma为岩石固相的体积模量，Kdry和Gdry分别为干岩石的体积模量和剪切模量，Kf为流体的体积模量。模式A 中，孔隙中有水合物生成，因此，Kf计算公式为：

式中，Sh、Kh分别代表水合物占孔隙的体积分数和水合物的体积模量，Kw为水的体积模量。

式（A-2）中，Kdry和Gdry的计算公式为：

其中，

式（A-8）中，P为有效压力；Kma、Gma分别是岩石骨架的体积模量、剪切模量，由公式（A-13）计算；ν 是岩石骨架的泊松比，且为临界孔隙度时单位体积内颗粒平均接触的数目，一般取8～9.5[28-29]。

（2）模式 B

模式B 中水合物被认为是岩石骨架的一部分，产生了两个效应：一个是使孔隙度减小，另一个是改变了骨架的体积模量和剪切模量。因此，在模式A 的基础上，需对沉积物孔隙度进行修正，即φr=φ(1-Sh)。同时，应将水合物作为矿物组分代入公式（A-9）来计算岩石的Kma和Gma。此外，沉积物孔隙中只有水，孔隙流体密度和体积模量等直接用水的替代。

A.2 BGTL 理论

BGTL 理论（ Biot-Gassmann Theory modified by Lee 的简称）建立在经典的BGT 理论（Biot-Gassmann Theory）上，在预测速度时不仅考虑了分压的影响，而且还考虑了岩石的孔隙度、固结度等因素的影响，其公式为：

由公式（B-2）和（B-3）得出沉积介质的剪切模量为：

A.3 简化三相方程（STPE）

利用Lee[32]的简化三相方程（STPE）可以计算各向同性的气体水合物填充沉积物孔隙空间的速度。Lee 推导出低频下气体水合物填充沉积物孔隙空间的体积模量和剪切模量，用于测井和地震数据:

式中，α为固结参数，γ是与剪切模量相关参数，可由公式γ=(1+2α)/(1+α)得到，Kma，Kw和Kh是骨架、水和水合物的体积模量，μma是骨架的剪切模量，φ为孔隙度。参数ε是水合物形成使沉积物骨架发生硬化的降低量，Lee 推荐使用ε=0.12 为建模数值。气体水合物填充沉积物孔隙空间的纵横波速度可由下式得：