谭 琳 刘 芳

(同济大学土木工程防灾国家重点实验室， 上海 200092， 中国)

0 引 言

我国南海北陆坡分布着丰富的天然气水合物资源。水合物是一种新型绿色能源，其产业化开采对我国能源结构优化具有非常重要的意义。然而，南海水合物富集区发育海底滑坡以及与之相关的地质灾害因素，如浅层气、浅水流、气烟囱、岩浆底辟与火山等(何健等， 2018)，水合物本身也是海底滑坡的重要诱因之一，开采方法不当可能会诱发海底滑坡，甚至造成海啸等地质灾害，构成重大安全隐患(Sun et al.， 2020； Tan et al.，2021)。

针对南海水合物富集区，学者已开展大量关于海底滑坡的地震勘测、室内实验及数值模拟研究(孙运宝等，2008； Wang et al.，2018； 霍沿东等，2019； Liu et al.，2020； Nian et al.，2020)，在海底滑坡的几何形态与变形特征分析方面已取得有建设意义的研究成果(He et al.，2014； Li et al.，2014；Yin et al.，2015)，对于海底滑坡的成因机制和灾害评估也有了初步探索(马云等， 2012； 雷亚妮等， 2018； Nian et al.，2019； Guo et al.，2020)。研究表明，地震、活动断层和水合物分解是触发海底滑坡的3大主要因素，其中约10%的海底滑坡由水合物分解诱发(Hance， 2003)。自然条件下，水合物成藏与海底滑坡处于伴生动态发展过程。一方面，滑移过程中，根部沉积物受张拉作用，一般会出现犁式断裂，犁式断裂自浅至深倾角变缓，总体呈上陡下缓的形态(孙运宝等， 2008)； 发生滑移后，在滑坡主体底部形成明显的滑移面，地层中的断裂与滑移面为气体提供了运移通道与储集空间，促进水合物成藏(苏丕波等， 2020)。另一方面，水合物藏在温压环境受扰动后再次分解，对沉积物的胶结作用消失，同时，产生的气体侵入古滑移面，造成孔压积聚并沿横向扩展，使沉积物抗剪强度降低，可能触发新一期滑坡(Dugan， 2012；Sun et al.，2017)。

水合物开采过程中，逃逸的气体可能迁移至古滑移面(图 1)，在其中横向扩展，会加速沉积物的弱化，从而触发大规模滑坡(Tan et al.，2021)。目前对于水合物开采诱发滑坡的机制研究尚处于起步阶段，大多针对均质储层边坡(Moridis et al.，2018；Song et al.，2019)，鲜有考虑实际储层中海底古滑坡的影响。海底古滑坡的存在使开采过程中的气液运移路径非常复杂，基于均质储层边坡的数值模拟忽略了古滑移面对气液运移路径的影响，限制了孔压积聚和材料强度弱化的影响范围，无法全面地反映水合物开采对边坡稳定性的影响机制。

古滑坡的几何形态与空间分布非常复杂(孙运宝等， 2008)、滑坡土体特性具有高度不确定性(Yoneda et al.，2019； Lall et al.， 2022)，增加了地质建模的难度。初步探究水合物开采过程中古滑坡发育边坡的稳定性，宜从简化模型入手。按照古滑坡与水合物藏在空间上的相对位置，水合物藏可分为下卧型和伴生型两种。如图 1所示，当气体运移通道未到达古滑坡体或者气体扩散范围在古滑坡下方地层时，水合物成藏于古滑坡以下，形成下卧型储层(Wu et al.，2018)； 当气体运移通道延伸至古滑坡体或气体扩散至古滑坡体以内时，水合物成藏于古滑坡体内，形成伴生型储层(Li et al.，2013)。

图 1 下卧型与伴生型水合物储层成藏机理和开采逃逸气体被阻滞于古滑移面的示意图Fig. 1 Illustration of the formation mechanism of the under-burden- and associated-type hydrate reservoirs and the escaping gas trapped in the ancient slip surface from hydrate production

不同类型的储层开采，可能触发不同类型的滑坡(Wan et al.，2016)。目前对海底滑坡分类尚未有统一的认识，从不同的角度对海底滑坡分类，术语存在一定差异(朱超祁等， 2015)。对于南海北陆坡的海底滑坡，我国学者从滑坡体的几何形态与变形特征出发，分为滑动型、滑塌型和蠕动型滑坡(He et al.，2014； Chen et al.，2016； Wan et al.，2016)。滑动型滑坡体发生平动，沉积物几乎不发生变形，滑移面较长，几乎与顺坡向平行； 滑塌型滑坡体发生转动，其中沉积物发生明显变形，滑移面为相对较短的曲面； 蠕动型滑坡是指正在发生缓慢顺坡运动的滑坡(Moscardelli et al.， 2008)。

本文针对下卧型与伴生型水合物储层，考虑古滑移面渗透率的空间变异性建立边坡模型，通过数值模拟计算瞬态孔压与土体强度，建立水合物开采二维边坡极限平衡分析法，研究水合物开采过程中边坡稳定性的演变，分析边坡失稳的滑坡类型，初步探究水合物开采过程中深海古滑坡体的再启滑机制。

1 物理模型

神狐海域位于南海珠江口盆地南部，水深400～2000 m，地形平缓，坡度不超过10°，分布着17条峡谷，走向为北北西—南南东向，谷壁较陡，坡度可达15°～22°。周庆杰(2015)利用峡谷区的多波束水深数据结合地震剖面分析识别出145个海底滑坡，这些滑坡分布于峡谷头部和侧壁，发育时期较晚，多为现代滑坡。统计显示，滑坡陡坎的高度8～125 m，坡度3.4°～26.9°，水平滑移距离0.35～3.6 km，面积0.04～15.88 km2，厚度十几米到几十米，埋深几米到几十米(周庆杰， 2015)。勘探表明，水合物藏位于海底以下300 m范围以内(Yang et al.，2017)，古滑坡体内以及下方地层均有水合物成藏，即存在下卧型和伴生型水合物藏(Wu et al.，2018)。

根据神狐海域峡谷区的海床坡度与古滑坡特征，建立下卧型和伴生型水合物藏二维边坡模型，如图 2所示，计算域宽3000 m，其中边坡段宽1800 m，长高比3︰1(坡度约20°)，坡肩和坡趾的水深分别为700 m和1300 m。储层沉积物以泥质粉砂为主，古滑坡体位于30～70 mbsf(meters below seafloor，海底以下…米)，古滑移面位于70～80 mbsf。将下卧型水合物藏设置在古滑移面以下，厚度为40 m(80～120 mbsf)，并设置竖向裂隙区域，连通古滑坡与水合物藏，促使水合物分解产生的气体迁移至古滑移面； 伴生型水合物藏与古滑坡沉积物重叠，即位于30～80 mbsf。

2 数值分析方法

2.1 水合物开采二维边坡动态极限平衡分析模型

在边坡稳定极限平衡分析框架内引入水合物开采过程中的瞬态孔隙压力和土体抗剪强度，建立水合物开采二维边坡动态极限平衡分析模型。边坡稳定极限平衡分析法的基本概念是利用安全系数来量化边坡的稳定性。安全系数定义为滑动力与抗滑力之间的比值，利用条分法求解。水合物开采过程采用热-流-化学耦合数值方法来模拟，从而获取边坡中瞬态孔隙压力和水合物饱和度的空间分布。假设土体黏聚力随水合物饱和度线性减小直至消失，通过莫尔-库仑准则计算边坡中土体抗剪强度的空间分布。

为实现以上方法，借助TOUGH+HYDRATE程序(Moridis， 2014)与SLOP/W软件(GEO-SLOPE InternationalLtd.， 2017)单向耦合来计算边坡稳定安全系数，确定相应的潜在滑移面。首先，采用TOUGH+HYDRATE程序考虑热-流-化学耦合机制模拟水合物开采过程，获取储层边坡中孔压分布场与水合物分解锋面的扩展。然后，根据水合物饱和度确定土体强度参数在储层边坡中的时空分布，基于极限平衡分析框架，利用SLOPE/W软件计算边坡稳定安全系数与潜在滑移面的时空演变。分析方法细节请参考我们的前期研究(谭琳等， 2020)。

2.2 物理模型参数

在热-流-化学耦合数值模拟中，将边坡模型离散为8688个单元(图 2)。设置孔压为静水压，压强梯度为10.5 MPa·km-1，海床表面温度为4 ℃，地温梯度为45 ℃ · km-1。将古滑坡主体作为致密盖层，渗透率设置为10-5mD。考虑古滑移面渗透率的空间变异性，生成渗透率随机场，平均值为38.3 mD，标准差为56.7 mD，顺坡向变程为500 m。表 1和表 2分别列出计算模型与地层的物理力学参数。

图 2 计算边坡模型及不同类型边坡的地层设置Fig. 2 The calculation slope model and the stratum settings of different types of slopes

表 1 计算模型

表 2 地层的物理力学参数Table 2 Physical and mechanical parameters of the strata

2.3 开采场景

开采方法采用单水平井热吞吐法，即通过周期性交替的降压与注热实现水合物开采(Li et al.，2011)。一般地，热吞吐法分为3个阶段： (1)注水阶段，向储层注入高温热水、蒸汽或者其他高温介质； (2)焖井阶段，将井口关闭一段时间，目的在于增大储层升温范围，促使水合物分解； (3)生产阶段，重启开采井，抽出气液混合流体。由于储层沉积物的热导率较低，焖井阶段对升温范围的扩大效果并不明显，数值模拟一般忽略该阶段(Li et al.，2011)，仅考虑注水阶段和生产阶段，本文也采用同样的设置。降压与注热共用同一口井，下卧型和伴生型水合物藏的开采井分别设置在97.5 mbsf和72.5 mbsf的位置。为了保证井周有足够的注水空间，在热吞吐循环之前，首先进行30 d的降压开采，促使井周水合物分解，降压压差为原位压强的50%。然后，开始热吞吐循环开采。注热阶段，以0.10 kg · s-1的速率向储层注入90 ℃热盐水(盐度3.5%)，持续时长2 d； 开采阶段，以0.12 kg · s-1的速率自储层抽出气液流体，持续时长2 d。

图 3 水合物分解锋面的扩展：(a～c)下卧型储层； (d～f)伴生型储层Fig. 3 The propagation of hydrate dissociation front in：(a～c) underburden-type reservoir； (d～f) associated-type reservoira. 435 d; b. 1200 d; c. 3019 d; d. 155 d; e. 1139 d; f. 3439 d

3 结果分析

3.1 水合物分解过程

开采注热过程中，注入的热水在储层中扩散，储层温度升高，促使水合物分解； 降压过程中，井周地层压力降低，促使水合物分解，并将注热与降压阶段产生的气体抽出。

图 3a～图 3c给出了下卧型储层开采过程中水合物分解锋面的扩展。可以看出，下卧型储层中开采井上方和下方均出现分解锋面的扩展，其中上方分解区域较下方的范围大。在注热过程中，热水沿开采井与古滑移面之间的裂隙区向上迁移，裂隙周围水合物分解。在降压过程中，井周地层压力降低，同时将下伏层中温热流体抽至井周，使井周水合物分解。由于上方分解区同时受降压与注热的影响，而下方分解区仅受降压作用的影响，因此上方分解区范围较大。

随着水合物的分解，沉积物黏聚强度降低，直至消失，形成强度弱化区(图 4a～图 4c)。同时，水合物分解后产生的气体发生膨胀，导致超孔压积聚。超孔压随着气液流体在古滑移面中传播，形成高压区(图 5a～图 5c)。高压区沿古滑移面横向扩展，促使逃逸其中的气体形成二次水合物(图 3a～图 3c)。

图 3d～图 3f给出了伴生型储层开采过程中水合物分解锋面的扩展。伴生型储层水合物赋存于古滑坡，古滑坡主体渗透率低，古滑移面渗透率高，所注热水在古滑移面中横向流动，古滑移面中水合物首先分解，锋面逐渐向古滑坡主体扩展，强度弱化区也随之自下而上扩展(图 4d～图 4f)。同时，分解产生的气体发生膨胀，造成高压区的发展。与下卧型储层相比，伴生型储层开采过程中超孔压较高，高压区在古滑移面的横向扩展范围较大。

3.2 海床稳定性演变

初始条件下，由于下卧型储层边坡中古滑移面沉积物的抗剪强度低，潜在滑移面贯穿其中(图 6a)； 伴生型储层边坡古滑移面中发育水合物藏，沉积物剪切强度较高，潜在滑移面位于浅层(图 6b)。

图 4与图 5给出了两种类型边坡开采过程中潜在滑移面的演变(红色网格所示)，水合物分解后，潜在滑移面向开采区移动。由图 4a～图 4c与图 5a～图 5c可以看出，下卧型储层的潜在滑移面贯穿古滑移面，与高压区相交，但未与弱化区相交，说明该情况下储层黏聚力降低对边坡稳定性几乎没有影响，而超孔压积聚对边坡稳定性的影响较大。

图 4 有效黏聚力的演变： (a～c)下卧型储层； (d～f)伴生型储层Fig. 4 The evolution of effective cohesion in：(a～c) underburden-type reservoir； (d～f) associated-type reservoira. 435 d; b. 1200 d; c. 3019 d; d. 155 d; e. 1139 d; f. 3439 d

图 5 超孔压的演变：(a～c)下卧型储层； (d～f)伴生型储层Fig. 5 The evolution of excess pore pressure in：(a～c) underburden-type reservoir； (d～f) associated-type reservoira. 435 d; b. 1200 d; c. 3019 d; d. 155 d; e. 1139 d; f. 3439 d

图 6 初始潜在滑移面：(a)下卧型储层； (b)伴生型储层Fig. 6 The initial potential slip surface of：(a) underburden-type reservoir； (b) associated-type reservoir

图 7 边坡稳定安全系数演变：(a)下卧型储层； (b)伴生型储层Fig. 7 The evolution of the slope stability factor of safety of：(a) underburden-type reservoir； (b) associated-type reservoir

图 8 滑移面有效黏聚力演变：(a)下卧型储层； (b)伴生型储层Fig. 8 The evolution of the effective cohesion on the potential slip surface of：(a) underburden-type reservoir； (b) associated-type reservoir

如图 7a所示，下卧型储层开采过程中边坡稳定安全系数先降低后升高。开采初期，高压区的发展使边坡稳定安全系数逐渐降低。然而，开采中后期，古滑移面中生成二次水合物(图 3a～图 3c)，提高了潜在滑移面上沉积物的黏聚强度(图 8)，使边坡稳定安全系数升高。开采期间，边坡稳定安全系数始终高于1.35，说明边坡稳定性虽有所下降但不会造成古滑坡复活。此外，潜在滑移面始终贯穿古滑移面，与边坡方向平行，表现为滑动型的潜在滑坡模式。

由图 4d～图 4f与图 5d～图 5f可以看出，伴生型储层边坡的潜在滑移面同时穿过高压区与弱化区，边坡稳定性受超孔压积聚与储层黏聚力降低的共同影响。开采过程中，边坡稳定安全系数降至小于1，说明开采可能造成古滑坡复活。潜在滑移面集中在开采区附近，呈现近圆弧形状，表现为滑塌型滑坡模式。

图 9给出了下卧型与伴生型储层边坡潜在滑移面上的超孔压演变。伴生型边坡高压区的超孔压远远高于下卧型边坡，且高压区范围较下卧型边坡大。这是由于伴生型储层开采产生的气体直接在古滑移面中扩展，同时，古滑坡体中水合物的存在进一步降低了渗透性，对气体的阻滞作用更明显，加剧了超孔压的积聚。因此，伴生型储层开采边坡稳定性降低得更加显著。

图 9 滑移面超孔压演变：(a)下卧型储层； (b)伴生型储层Fig. 9 The evolution of the excess pore pressure on the potential slip surface of：(a) underburden-type slope； (b) associated-type slope

4 结 论

针对深海古滑坡发育区下卧型与伴生型水合物藏，考虑古滑移面渗透率的空间变异性，引入瞬态孔压与土体抗剪强度，建立水合物开采二维极限平衡分析模型，针对水合物藏与古滑坡不同的相对位置，初步研究水合物开采过程中古滑坡再启滑机制，分析了不同储层边坡的失稳模式，得到如下结论：

(1)水合物开采过程中，分解区形成强度弱化区，同时产生的气体迁移至古滑移面，造成超孔压积聚，形成沿古滑移面横向扩展的高压区，土体抗剪强度降低，可能导致古滑坡复活。

(2)下卧型储层边坡稳定性主要受超孔压积聚的影响，潜在滑移面始终贯穿古滑移面，边坡稳定安全系数先降低后升高，导致古滑坡复活的可能性较低，潜在的失稳模式为滑动型滑坡。

(3)伴生型储层边坡稳定性同时受黏聚力降低与超孔压积聚的影响，整个潜在滑移面集中在开采区域周围，水合物开采可能会造成古滑坡复活，表现为滑塌型滑坡模式。

此外，水合物开采中古滑坡的再启滑机制受到多重因素的影响，例如古滑坡沉积物的渗透性，该特性具有高度不确定性，这种不确定性对开采古滑坡再启滑机制的影响有待进一步研究。