地球物理学报  2021, Vol. 64 Issue (3): 864-875   PDF    
流体-地质力学耦合建模表征水力压裂诱发地震: 以加拿大Fox Creek地区为例
惠钢1, 陈胜男1, 顾斐2     
1. 加拿大卡尔加里大学化学与石油工程系, 卡尔加里 T2N1N4;
2. 中国石油勘探开发研究院, 北京 100083
摘要:随着水力压裂技术在页岩气开发中的广泛应用,加拿大西部盆地的诱发地震活动显著增加.目前对于诱发地震的综合表征方法还不成熟.本文采用一种综合地质、岩石力学及流体力学的研究方法,对Fox Creek地区2015年2月8日发生的M3.0诱发地震事件进行了综合表征.首先,利用高分辨率三维反射地震资料,采用蚂蚁体追踪技术识别潜在断层.其次,利用测井曲线和压裂施工数据等资料定量求取岩石力学及地应力参数,建立三维地质力学模型,明确水力压裂缝的空间扩展规律.最后,建立流体-地质力学耦合模型,计算水力压裂过程中断层附近的孔隙压力及局部应力变化,利用摩尔-库仑破裂准则判定断层激活的时间与空间位置,揭示本次诱发地震事件的触发机制并提出风险控制对策.结果表明,三条由Precambrian基底向上延伸至Duvernay地层的近垂直断层在水平井压裂过程中被激活.由于水平井的部分压裂缝与断层沟通,注入流体沿断层的高渗透破裂带向下迅速扩散,在基底位置激活断层并诱发M3.0地震事件.其中孔隙压力增加是本例中断层活化的主要因素.现场措施表明,增大压裂水平井与已知断层之间的距离被可以有效地降低地震风险.因此在进行水平井钻井及压裂作业之前,明确地下断层的分布位置至关重要.
关键词: 诱发地震      水力压裂      断层活化      流体-地质力学耦合模型      摩尔-库仑破裂准则      风险控制策略     
Coupled fluid-geomechanics modeling to characterize hydraulic fracturing-induced earthquakes: Case study in Fox Creek, Canada
HUI Gang1, CHEN ShengNan1, GU Fei2     
1. Department of Chemical and Petroleum Engineering, Calgary T2N1N4, Canada;
2. Research Institute of Petroleum Exploration and Development, Beijing 100083, China
Abstract: Along with the wide application of hydraulic fracturing techniques in the development of shale gas reservoirs, a salient increase in induced earthquakes has been reported recently in the Western Canada Sedimentary Basin. However, the comprehensive charaterzation of fracturing-induced seisicmicity has not been well investigated. In this paper, an integrated method combining geology, geomechanics, and hydrodynamics is utilized to investigate the M3.0 induced earthquake that occurred in Fox Creek on February 8, 2015. First, the ant tracking approach is utilized to identify the pre-existing fault through the high-resolution 3D seismic survey. The mechanical parameters and in-situ stress tensors are derived from the well logging and treatment data. Additionally, the 3D geomechanical model is built up and utilized to simulate the propagation of hydraulic fractures. The fluid-geomechanical coupling model is then established to calculate the pore pressure and local stress variation in the vicinity of the fault during and after fracturing operations. Finally, the Mohr-Coulomb failure criterion is employed to determine the spatiotemporal activation of related faults and reveal the triggering mechanisms of induced earthquakes. The results show that three near-vertical faults were activated during hydraulic fracturing, all of which extended upward from the Precambrian basement to the Duvernay Formation. Because of the hydrology communication between inferred faults and the stimulated well, some injected fracturing fluids diffused downward along the high-permeability fault damage zone, activating the fault in the basement and triggering the M3.0 earthquake. The increase in pore pressure is the predominant controlling factor in the process of fault activation. Field measurements show that enlarging the distance between the horizontal wellbore and the known fault had been demonstrated effective in mitigating the potential risks of earthquakes. Therefore, it is essential to determine the distribution of the subsurface faults prior to the drilling and hydraulic fracturing operations of horizontal wells in the unconventional reservoirs.
Keywords: Induced earthquakes    Hydraulic fracturing    Fault activation    Coupled fluid-geomechanics model    Mohr-Coulomb failure criterion    Mitigation strategy    
0 引言

近年来,随着水平井多级压裂技术在页岩气开发中的广泛应用,水力压裂引发的人为地震活动显著增加(张东晓和杨婷云,2015Atkinson et al., 2016Lei et al., 2019Schultz et al., 2017).全球一些较大震级的诱发地震事件,例如在中国四川盆地发生的2017-01-28筠连4.9级地震,2018-12-18兴文5.7级地震,2019-09-07威远5.4级地震和加拿大西部盆地发生的2019-03-04 Red Deer 4.2地震事件,都被证实与附近的水平井压裂作业有时空相关性(Lei et al., 2017, 2019, 2020Schultz and Wang, 2020).

水平井多级压裂被认为是一项有效开发页岩气资源的前沿技术.在水力压裂过程中,压裂液被加压并泵入目标地层,造成岩石的拉伸破坏,产生特定方向的压裂缝(一般平行于最大水平主应力).然后将配伍的支撑剂注入到新生压裂缝中,使其在后续油气开采过程中保持张开状态,从而使页岩气被高效快速采出.一般而言,该压裂过程会产生不连续的、瞬时的微地震事件,而且其震级通常小于0,多数在-1级以下(Chen等,2018).然而,当水力压裂缝在地下扩展并沟通潜在的断层或者大裂缝时,就会诱发较大震级的地震活动.一些中等强度(0.5 < M < 7)的诱发地震事件已被证实与压裂缝沟通潜在断层有关(Bao and Eaton, 2016Eyre et al., 2019Lei等, 2017, 2019, 2020).为减少诱发地震带来的财产及人身安全,各国相继设立了信号灯预警系统,即当水平井压裂过程中产生的诱发地震事件达到警戒震级(红灯)时,所有的压裂作业需立即停止.英国将警戒震级设为0.5, 而加拿大阿尔伯塔省设为4.0(Eaton, 2018).警戒震级的不同显示了各国能源管理部门对于工业活动诱发地震事件的管控程度,并与区域背景下的地震活动频率密切相关.

加拿大西部盆地具有丰富的页岩气资源,水平井多级压裂技术在该盆地取得了巨大成功(Atkinson et al., 2016).近年来,该盆地发生了很多中等强度的地震事件.这些事件都被证实与附近的水力压裂作业密切相关(Bao and Eaton, 2016Eaton et al., 2018Schultz and Wang, 2020).图 1显示了加拿大西部盆地记录的M>2.5地震事件.其中,震级2.5被认为是加拿大西部盆地诱发地震事件的长期监测阈值(Schultz et al., 2017).图 1还显示出诱发地震事件群具有明显的聚集特征,即沿落基山脉变形边缘带(图 1虚线)分布.其中的震源球指示地震事件的震源机制(Wang et al., 2018).图中显示Fox Creek地区是该盆地诱发地震事件的频发区域之一.前人研究表明,该地区诱发地震与地层孔隙压力(Eaton and Schultz, 2018)、纵向接近基底的距离(Pawley et al., 2018)、平面接近生物礁边界(Schultz et al., 2017)及断层应力状态(Zhang et al., 2019)等地质因素密切相关.此外,Fox Creek地区单口水平井压裂液注入总体积达到32944 m3.根据压裂液注入体积与诱发地震震级的经验线(McGarr,2014),该注入体积所对应的地震事件最大震级可以达到M4.0.综上可知,Fox Creek地区水力压裂诱发的地震事件受到地质和工程因素的共同制约(Pawley et al., 2018Schultz et al., 2018).

图 1 加拿大西部盆地震级大于2.5的地震事件分布(截止至2019-04-30) 圆圈代表诱发地震事件.圆球显示工业活动诱发地震的震源机制(Wang等, 2018).圆球上方数字代表事件发生事件.虚线代表落基山脉变形边缘带. Fig. 1 Map of historical seismicity of ML≥2.5 up to 2019-04-30 in Western Canada The cycles denote the induced earthquakes. The magnitude-scaled "beach balls" show the focal mechanisms of induced earthquakes (Wang et al., 2018). Four digits above the beach ball represent the nucleated time for associated earthquakes. The dashed line indicates the deformation edge of the Rocky Mountains.

近年来,国内外学者采用不同的研究方法对诱发地震进行了表征.国外方面,Bao和Eaton(2016)研究了诱发地震活动的触发模式,提出压裂过程中的应力变化能够激活断层.Lele等(2017)建立了压裂诱发地震的概念模型,认为水力压裂缝与附近断层沟通、并引起断层内部压力上升是导致断层活化并诱发地震事件的主要原因.国内方面,张致伟等(2012)研究了四川地震活动与注水井压力与注水量的关系, 认为诱发地震频次、强度与注水量呈正相关.何登发等(2019)利用高精度三维地震资料探索诱发地震的形成机制,认为四川长宁地区的地震活动与现今青藏高原东南缘的大地构造背景密切相关.本文在Lei等(2017)研究的基础上,综合分析地质、岩石力学及流体力学等特征,来表征加拿大西部盆地水力压裂诱发的地震事件.

本文采用一种综合地质、岩石力学及流体力学的研究方法,对加拿大西部盆地Fox Creek地区2015年2月8日发生的M3.0诱发地震事件进行精细表征.首先,建立区域地质、完井、测井、三维地震及压裂施工等资料库.其次,利用该区高分辨率三维地震勘探资料, 采用蚂蚁体追踪技术识别潜在断层.然后,根据测井曲线和施工数据等定量求取岩石力学参数及地应力,建立三维地质力学模型,计算出水力压裂缝的几何形状和缝内流体压力的分布特征.最后,建立流体-地质力学耦合模型,计算水力压裂过程中断层附近的孔隙压力及局部应力变化,利用摩尔-库仑破裂准则判定断层激活的时间与空间位置,从而揭示本次诱发地震事件的触发机制并提出相应的风险控制对策.

1 诱发地震触发机理及断层失稳准则

一般认为,诱发地震的触发机制主要包括两方面(见图 2).一种是注入流体通过高渗透单元(如压裂缝)直接沟通断层,造成断层内流体压力增加;另一种是在流体注入或者开采过程中,弹性应力从目标地层传递到断层,引起断层面应力变化(Ellsworth,2013).上述两种情况均导致断层失稳并发生移动,从而诱发地震事件(Healy et al., 1968Raleigh et al., 1976).

图 2 流体注入或采出诱发地震示意图(Ellsworth, 2013) Fig. 2 Schematic diagram of fluid injection or extraction-induced seismicity(Ellworth, 2013)

具体而言,当断层内部孔隙压力或者弹性应力足够高、超过断层的破坏强度时,断层将发生滑动进而诱发地震活动.为表征断层的失稳特征,国内外研究者们提出了一系列判别标准,其中最常用的是摩尔-库仑破坏准则.该准则是确定断层失稳和预测强震活动后余震分布的有效工具(Catalli et al., 2013).摩尔-库伦破坏准则定义如下:

(1)

其中,CFS为断层面上的库仑破坏应力,MPa;μ为摩擦系数,无因次;τ为剪切应力(滑动方向为正),MPa;σn为正应力(拉伸方向为正),MPa;Pp为孔隙压力,MPa.在三维应力场中,断层面(l, m, n)上的剪切应力τ及正应力σn计可以通过最大,垂向及最小主应力(σ1σ2σ3,MPa)计算,其计算公式为(Lei et al., 2020)

(2)

(3)

为方便运算,本文对上述公式进行了二维近似.在二维摩尔圆中,断层面左右侧的剪切应力τ及有效正应力σn的计算公式为

(4)

(5)

(6)

其中,τlτr分别为断层左、右侧的剪切应力,MPa;σ1σ3分别为断层面承受的最大、最小主应力,MPa;β为断层走向与最大水平主应力的夹角.

在断层滑动过程中,断层周围的岩石会发生弹性变形,导致断层附近局部应力场发生变化,产生库仑破坏应力(CFS)正变化和负变化区域.研究表明,在CFS正变化区域,断层将处于不稳定状态,断层发生失稳、滑移的几率显著增加(Stiros and Kontogianni, 2009King and Devès,2015).其判别方程基于(1)式得出,Δ为各参数变化量.

(7)

摩尔圆可以直观显示上述判别准则.通常根据断层面的剪切应力及有效应力来绘制摩尔圆及库仑破坏线.具体来说,在压裂过程中,压裂缝扩展并沟通地层的潜在断层时,注入流体会增加断层内的孔隙压力,从而降低断层面上的正应力,使断层向摩尔-库仑失稳向移动,直至到达库仑破坏线引起断层滑动(见图 3a).在这种情况下,注入流体通常会以超过1 km·d-1的速率在已沟通的断层内扩散,因此会很快触发地震事件(Tadokoro et al., 2000).

图 3 (a-b)压裂过程中孔隙压力变化和应力变化引起的库仑应力失稳 Fig. 3 (a-b)The Coulomb Failure tiggered by changes of pore pressure and in-situ stress

另一种情况是压裂缝没有直接沟通断层.但压裂液的注入引起地层孔隙压力增大,并导致区域岩石弹性变形,从而改变断层附近的局部应力场.在这种情形下,摩尔圆也将接近库仑破坏线并激活断层(见图 3b).其中,垂直主应力变化是由负载压力变化和Biot系数决定,而水平主应力变化则与岩石的泊松比有关,并且通常小于垂直应力变化(Haddad et al., 2020).

综上所述,我们可以计算出压裂过程中断层附近的孔隙压力及局部应力变化,利用摩尔-库仑破裂准则来判定断层激活的时间及空间位置,从而揭示诱发地震的触发机制.

2 区域地质背景及M3.0地震事件

Fox Creek地区位于加拿大西部盆地中部(见图 1).该地区地层发育较全,自上至下依次为泥盆系、寒武系及前寒武系花岗岩基底,上述地层的总厚度超过2000 m(图 4a).其中Duvernay地层的页岩气资源为该地区主要的油气开发对象.该地层广泛发育富有机质页岩,上覆泥岩为主的Ireton地层,下覆灰岩为主的Swan Hills地层,地层之间均为整合接触.Duvernay地层厚度约为40 m,埋藏中深低于海平面约2520m.地层平均有效孔隙度为6.6%,平均渗透率为392 nD(3.92×10-20m2),平均总有机质(TOC)含量为4.5%(Creaney and Allan, 1990Weir et al., 2019).

图 4 Fox Creek地区发育的地层及岩性统计 Fig. 4 Statistics of stratigraphy and its lithology developed in the Fox Creek region

Fox Creek地区作为诱发地震事件的高发区域,其地震活动最早可追溯到2013年12月(Schultz et al., 2017).此后,随着水平井多级压裂技术的广泛应用,该区诱发地震活动频繁发生.其中有超过10次震级大于3.0的地震事件被证实与水力压裂作业密切相关(Bao and Eaton, 2016Eyre et al., 2019).这些地震事件主要归因于该地区特定的地质因素和工程因素,如Duvernay地层高孔隙压力、与Precambrian基底的距离、与Swan Hills生物礁边缘的距离、被激活断层的初始应力状态等地质因素及压裂液的累计注入体积等工程因素(Eaton et al., 2018Pawley et al., 2018).截止至2018年12月,本地区有606口水平井进行了水力压裂作业.统计表明,平均每口水平井压裂级数为28.7段,水平段长度为2135 m,压裂液注入速率为9.5 m3·min-1,累计注入体积达到32944 m3.如前文所述,该注入体积可诱发震级超过M4.0的地震事件(Mcjarr,2014).

本文以2015年2月8日发生的M3.0地震事件群为例,对诱发地震的触发机制进行综合研究.该地震群发生在Fox Creek地区中部的一口压裂水平井附近.图 5a表明,研究区与本区发育的大型断层都有一定的距离(Pawley et al., 2018).其中与最近断层的距离已超过5 km,表明本次事件与大型断层无直接关系.因此,需要利用高精度的三维地震勘探资料,识别发育规模相对较小的断层,从而揭示构造与诱发地震的关系.

图 5 (a) 研究区位置及大断层分布, 底图为Duvernay顶部构造等值图;(b) 地震事件及水平井空间分布特征;(c) 压裂施工与诱发地震事件时间顺序 Fig. 5 (a) Location of the study area and distribution of large faults.The base map dnotes the elevation of top Duvernay Formation; (b) Spatial distribution of the induced earthquakes and horizontal well; (c) Time sequence of fracturing operation and induced events

2015年2月5日-14日,该水平井从南向北实施了水力压裂作业.根据压裂参数的统计结果,该水平井共压裂级数20段,水平段长度1.5 km.平均施工压力55.8 MPa,平均压裂液注入速率为8.7 m3·min-1,平均每个压裂段注入1904 m3,累计注入量为38078 m3.区域地震台站监测结果显示,压裂过程中共有109起地震事件发生(Bao and Eaton, 2016).图 5b显示了该诱发地震群与水平井压裂段的空间分布特征.可以看出,大部分地震事件发生在水平井的东侧,以及Duvernay地层下方200~1300 m范围内.

现场资料表明,压裂开始40 h后,南部地震群在水平井筒东侧约200 m处被激活.其中包括本次研究的M3.0事件.该事件的震源分析表明(图 5b震源球),被激活的断层呈近N-S走向展布(Schultz et al., 2017).南部地震群的发生持续了40 h,直到2月8日结束(图 5c中第6级压裂段结束).随后,中部地震群被激活,结束后,北部地震群在距井筒东侧300 m处被激活,其震级范围为M1.5~M2.9(图 5b-c).

利用地震震级与地震事件频率关系拟合曲线的斜率值(b值),可识别地震事件群的触发类型.一般而言,对于天然断层激活所诱发的地震事件群,其b值一般小于1(Burridge and Knopoff, 1967).而水力压裂活动所诱发的地震事件群,通常b值要大于1.5(Eaton et al., 2014).本次地震事件群的b值为0.84,为典型的天然断层活化特征.根据水平井压裂缝的位置推断,本次地震群可能是由于压裂缝扩展并沟通了断层的基底部分,从而激活断层并诱发一系列地震事件.为进一步明确其触发机制,需要结合测井、三维地震、岩石力学及应力场、施工数据等资料对本次诱发地震进行综合表征.

3 蚂蚁体追踪技术识别断层

该地区进行过三维地震勘探,因此可用来探寻构造特征与诱发地震事件的关系.在三维反射地震解释中,蚂蚁体跟踪技术通常用于识别小型断层,其工作流程为:(1)井震结合建立关键井合成地震记录,提取工区相关层位;(2)对三维地震数据进行预处理,包括中值滤波、结构平滑、高斯空间滤波和带通滤波处理;(3)从地震属性中提取Chaos、方差和倾角偏差,进行地震反射边缘检测;(4)将滤波方法应用于地震解释后处理,开展蚂蚁体自动跟踪,刻画出断层的空间分布(Pedersen et al., 2002; Hui et al., 2021).

图 6a-6c分别为Duvernay、Cambrian、Precambrian地层的蚂蚁体属性切片及地震事件叠合图.可以看出,大部分地震事件集中在基底,并与基底复杂的断层分布特征大致相同.基于蚂蚁体属性特征,分析断层与地震事件之间的空间关系,共识别出3条近垂直的天然断层(图 6d).三条断层均由Precambrian基底向上延伸至Duvernay地层顶部.其中,断层1的走向和倾角分别为175°和89°,与前人的震源反演结果一致(Schultz et al., 2017).断层1高度约为1300 m,延伸距离为810 m.其发育规模表明,该断层被激活后具有诱发最大震级4.8级地震的潜力(Zoback and Gorelick, 2012).断层2和断层3走向大致相同,为北东-南西向,延伸长度分别约为820 m和590 m,倾角分别为89°和88°.根据前人的研究结果(Yehya et al., 2018Fan et al., 2019),本次研究将三条断层带的厚度设置为45 m,包括5 m的断层核及两侧各20 m的断层破裂带.图 6d还显示出压裂诱发的地震事件在纵向上沿三条断层两侧分布.结合图 5b-5c, 这些地震事件与水平井的某些压裂段具有时间与空间相关性.这表明,压裂过程中水平井压裂缝扩展并沟通了三条基底断层,引起断层失稳并诱发了M3.0地震事件.

图 6 (a-c) 三个地层的蚂蚁体属性切片及地震事件(黑点)叠合图; (d) A-A′地震剖面断层解释结果及断层两侧的地震事件.A-A′剖面位置见图 5b Fig. 6 (a-c) Horizontal cross-section view of ant tracking property and seismic event (black spot) in three formations; (d) Interpretation of faults bounded by related seismic events in the A-A′ seismic profile.A-A′ section position is shown in Fig. 5b
4 三维地质力学模型与压裂缝扩展 4.1 地层异常高压与局部地应力场

研究表明,Fox Creek地区的高地层压力及局部应力场对诱发地震有重要影响(Eaton et al., 2018Pawley et al., 2018).本文综合测井和压裂施工资料,对本区的地层孔隙压力及地应力参数进行了定量求取.

首先,利用压裂过程中的瞬时关井压力来估算孔隙压力.本例水平井20级压裂段的平均瞬时关井压力及测试深度分别为56.8 MPa和3387 m.因此孔隙压力Pp梯度估算为16.77 kPa·m-1,接近前人研究结果的16.8 kPa·m-1(Eaton et al., 2018).相对于其他地区的地层压力,Fox Creek地区的地层异常高压成为诱发地震事件在该区频繁发生的重要地质因素(Pawley et al., 2018).

根据井眼垮塌和钻井诱导缝等资料,结合Shen等(2019)的研究结果,确定本地区自Ireton至基底地层的最大水平主应力SHmax方向均为NE45°.同时,为便于耦合模拟,本文对最大及最小应力轴进行了水平近似处理.垂直主应力Sv梯度则由附近直井的密度测井数据求得,为24.6±0.25 kPa·m-1.其中0.25为分析误差,设定为计算值的1%(Zoback,2007).根据施工过程中的闭合压力,估算本区最小水平主应力Shmin梯度为20.92±0.56 kPa·m-1.最后,利用经验公式SHmax= 3Shmin - 2Pp(Zoback, 2007),估算SHmax梯度为29.16±1.68 kPa·m-1.主应力计算结果与前人研究结果基本一致(Eyre et al., 2019Zhang et al., 2019Shen et al., 2019).

通过比较三应力张量的大小(SHmaxSvShmin),可以确定本区发育的断层类型为走滑断层.根据公式(4)-(6),绘制出三条被激活断层的摩尔圆,以表征其初始应力状态.如图 7所示,如果忽略压裂过程中的局部应力变化,则激活断层1、断层2和断层3分别需要增加3.5、2.0、0.5 MPa的断层内部孔隙压力.

图 7 三条断层的初始应力状态 斜线为库仑破坏线.圆点代表各断层相关地震事件的震源机制.摩尔圆内的等值线表示激活断层所需要增加的孔隙压力(MPa). Fig. 7 Original stress state of three faults The pink line represents the Coulomb failure line.The dots denote focal mechanisms of related seismic events.The contours within the Mole circle show the increased pore pressure (MPa) required for the fault activation.
4.2 岩石力学及物性参数的确定

利用岩石力学参数可以定量表征压裂过程中岩石弹性变形对于断层激活的影响.研究区附近的一口直井测量了从Duvernay地层到Muskeg顶部地层的纵波(VP,m·s-1)及横波(VS,m·s-1)测井数据.从Muskeg底部到Precambrian地层的速度数据则采用前人的研究成果(Ronald et al., 2019).可据此计算各地层的动态泊松比(υ, 无因次)和杨氏模量(E,GPa),再结合三轴应力实验得到的静态参数对其进行校正.此外,根据附近取心井的岩心化验数据,对各地层的孔隙度和渗透率进行了统计(Weir et al., 2019).表 1统计了从Duvernay地层到Precambrian基底各地层的弹性参数及物性参数.其中静态泊松比范围为0.19~0.30, 静态杨氏模量范围为48~70 GPa, 孔隙度和渗透率分别为0.02~0.07和5.2×10-20~6.46×10-18 m2.根据以上结果,建立了整合应力场、岩石力学参数和物性参数的三维地质力学模型.

表 1 各地层弹性参数及物性参数统计 Table 1 Statistics of elastic and physical parameters of each formations
4.3 人工压裂缝扩展模拟

当对水平井实施水力压裂作业时,压裂缝扩展方向通常与水平最大主应力SHmax平行.对于水力压裂缝扩展的表征,通常是利用施工数据和地质参数,进行压裂过程中的净压力(施工压力减去最小主应力及摩擦阻力等)的历史拟合,从而模拟压裂缝的空间扩展.本文采用PKN模型模拟压裂缝的三维实时扩展.在PKN模型中,断裂面处于平面应变状态,断裂界面为椭圆形,裂缝高度为固定值.假定忽略滤失效应,则在某一时刻压裂缝半长和宽度的计算公式为(Yew and Weng, 1997)

(8)

(9)

(10)

其中,L(t)为t时刻压裂缝半长,m;W(t)为t时刻压裂缝宽度,m;P(t)为t时刻压裂缝内初始压力, MPa;G为剪切模量, GPa;qo为注入速率,m3·min-1μ为注入流体黏度,cp;h为压裂缝缝高,m.

模拟结果表明,压裂缝半长为104~166 m,平均为131 m(图 5b).裂缝宽度为0.0135~0.0173 m,平均为0. 0166 m.根据裂缝渗透率与裂缝宽度的经验公式,计算裂缝渗透率约为7.5×10-12 m2,远超过Duvernay地层基质渗透率3.94×10-19m2.上述压裂缝属性特征将用于流体-地质力学耦合模型.

5 耦合流体-地质力学有限元模拟 5.1 流体-地质力学耦合方程

本文采用线性孔隙弹性理论来表征压裂过程中的孔隙压力扩散和应力扰动.基于多孔介质连续性方程及固体弹性变形方程,可推导出流体-地质力学耦合方程(Wang,2000):

(11)

(12)

其中,f(x, t)为单位体积应力,N·m-3q(x, t)为流体注入速率,m3·s-1u为平移向量,m;ν为泊松比,无因次;ε为体积应变,无因次;α为Biot系数,无因次;M为Biot模量,无因次;k为渗透率,μm-2η为动态流体黏度,cp.

上述方程可以利用有限元方法进行求解.本文使用COMSOL Multiphysics软件进行流体-地质力学耦合建模, 定量表征压裂过程中的孔隙压力及应力变化,再利用摩尔-库仑失稳准则来确定断层被激活的位置及时间.

5.2 建立有限元耦合模型及初始化设置

三维模型在xyz方向的尺寸分别设为2.5 km×3.5 km×1.8 km.z底部设为-4.1 km.根据实际地层特征,自Ireton至Precambrian地层厚度依次设置为160、40、100、100、200、1200 m.水平井深度设为-2520 m.压裂缝沿NE45°扩展,其几何参数及属性采用前文计算结果.三条解释断层嵌入到模型中,其断层带厚度设为45 m,包括中间5 m的致密断层核部分及两侧各20 m的断层破裂带(Fan et al., 2019).其中断层核及破裂带渗透率分别设为经验值1×10-20 m2、1×10-14 m2(Yehya et al., 2018).同时,模型中各地层物性及岩石力学参数设定为表 1对应数值.压裂液密度设定为1200 kg·m-3,压缩系数设定为4.6×10-10/Pa,动态黏度设定为0.4 mPa·s.图 8a-8b为耦合模型的三维视图及网格划分.

图 8 (a) 三维耦合模型初始化示意图; (b) 模型网格划分.在压裂缝及断层附近加密网格 Fig. 8 (a) 3D view of poroelastic model initialization; (b) 3D mesh using triangular elements. The mesh surrounding the fractures and faults are refined

在耦合模型中还设置了初始条件和边界条件.假设模型的初始流体系统处于流体静力平衡状态,则可以设置初始孔隙压力p(t=0)=0和初始应力张量σ(t=0)=0.因此, 耦合模拟结果即为孔隙压力和应力张量的相应变化值.将各压裂段的压裂液注入体积转化为各段净注入压力,并作为其耦合模拟过程中的驱动力,其计算公式如(10)所示.压裂缝延伸长度与宽度如公式(8)和(9)所示.净注入压力的作用时间同各压裂段实际施工时间保持一致.前一段压裂结束后,压裂缝末端延伸净压力接近于0(接近于地层孔隙压力),后一段压裂开启,如此循环至最后一级压裂结束.模型边界设定为固定边界,即顶部处于无牵引状态,侧向和底面边界的位移保持为零.为便于耦合模拟,本文对最大及最小应力轴进行了水平近似处理.即垂直上覆压力为变量,最大及最小主应力取决于上覆压力,因此三轴应力简化为单轴应力状态.设定耦合模型受控于线性弹性力学和达西定律,并在压裂缝及断层附近加密网格(图 8b).通过固体力学与多孔介质流体流动的耦合模拟,计算出压裂过程中的应力扰动和孔隙压力变化.

5.3 模拟结果及M3.0地震事件的诱发机理

图 9a-9c显示压裂80h后压裂层位(Duvernay)及震源位置所在基底(Basement)的孔隙压力平面变化情况.可以看出,断层1两侧破裂带区域的孔隙压力变化明显.表明致密断层核对流体流动的阻碍作用不明显.图 9dMW3.0事件震源位置处(图 9b十字)ΔPp及ΔCFS随时间变化的情况.MW3.0地震时间发生时,ΔPp及ΔCFS瞬时值分别达到3.6 MPa及2.2 MPa.基于有效正应力变化值(ΔσnPp=-3.12 MPa)和剪切应力变化值(Δτ=0.14 MPa),图 7所示的断层1将达到失稳状态,可以断定断层1被激活.同样,压裂过程中孔隙压力及局部应力的变化使得断层2和3也被激活.值得注意的是,孔隙压力增加仍然是影响MW3.0诱发地震事件的主要因素(图 9d).

图 9 (a-b) 压裂80 h后压裂层位及基底层位ΔPp分布图; (c) 地震事件叠合于(b)图上; (d) 原边界及大边界模型震源处(图 9b十字)的ΔPp及ΔCFS随时间变化情况.压裂段及地震事件也在图中绘出 Fig. 9 (a-b) ΔPp at the Duvernay Formation and basement, respectively, at t=80 hours after the onset of fracturing operations; (c) Induced events overlapped Fig. 9b; (d) Temporal ΔPp and ΔCFS at the nucleation position of MW3.0 earthquake. Fracturing stages and induced events are also depicted

另外,为检测耦合模型是否受到设定边界的影响.将模型尺寸扩大为5 km×7 km×1.8 km(原模型体积的4倍).模拟结果如图 9d所示.可以看出,新模型震源位置的ΔCFS及ΔPp与原模型差别不大(变化值7%左右).因此在原模型尺度大小及边界条件的设定下,边界大小对于本次模型的结果输出影响不大.

6 讨论

目前,针对Fox Creek地区诱发地震的风险控制对策主要体现在施工控制方面,如优化新钻水平井与已知断层的空间位置、控制压裂液注入速率及累计注入体积等(Schultz et al., 2018).其中,增大水平井筒与已知断层的距离被证实可以有效地降低地震风险.

图 10a所示,根据蚂蚁体反演属性,本例中的断层1向南延伸3 km,呈SE-NW走向.2015年11月,即北部水平井(本例)完井9个月之后,在断层西侧1.5 km左右完钻一口NS向、水平井段长度为3 km的新水平井.自11月11日至29日,对该水平井实施压裂施工作业,共压裂级数31段.压裂液注入速率及累计注入体积分别为9.5 m3·min-1和38303 m3,与北部水平井压裂施工参数8.7 m3·min-1和38078 m3大致相同.

图 10 (a) Duvernay地层顶部的蚂蚁体属性分布.下部圆球显示南部水平井压裂过程中的最大震级M1.25事件的震源机制.黑色虚线为D-D′剖面位置; (b) D-D′剖面的地震事件空间分布特征.圆球大小和颜色分别代表地震事件的震级大小和发生时间 Fig. 10 (a) The horizontal cross-section of ant tracking attributes showed the simulated pore pressure changes at the Basement of Fault 1 at t=80 hours after the onset of HF operations; (b) Vertical F-F′ section view of earthquakes distributions. The balls denoted induced earthquakes, colored by time and scaled by magnitude

地震监测结果表明,压裂过程中在该井以东600~1500 m处诱发了一系列震级较小的地震事件,其中最大为M1.25事件(图 10a).图 10b比较了南、北水平井压裂过程中的两次主要地震事件的空间分布特征.可以看出,南部水平井诱发的事件震级明显变小,而且深度主要集中于压裂地层(Duvernay)上下50 m范围内,这与北部地震事件主要集中在基底的情况形成鲜明的对比.此外,南部最大地震事件M1.25的震源机制表明(图 10a下部震源球),被激活断层的走向为NE 28.2°±8.3°(Zhang et al., 2019).这与北部被激活断层(断层1)的近南北走向显著不同.因此可以推断,由于南部水平井与断层1位置较远,其压裂施工并未激活断层1.相反,在南部水平井和断层1之间或存在一个小规模断层(图 10a断层4).南部地震群的产生可能由于该断层与南部井的压裂缝连接而被激活.此外,模拟结果表明,南部井压裂过程中的CFS变化仅为0.15 MPa,远小于北部井的1.53 MPa.由此可见,优化水平井与已知断层的距离可以有效降低诱发地震活动的风险.该现场措施为优选水平井井位、降低地震风险提供了可靠的现场依据.

7 结论

本文以Fox Creek地区2015年2月8日发生的M3.0地震事件为例,基于蚂蚁体识别的潜在断层及三维地质力学模型,通过流体-地质力学的耦合模拟,明确了水力压裂施工过程中孔隙压力及局部应力的变化特征,揭示出本次地震事件的触发机制,并提出相应的风险控制对策.具体包括:(1)利用蚂蚁体追踪方法识别出三条由Precambrian基底向上延伸至Duvernay顶部的近垂直断层;(2)基于三维地质力学模型,模拟出水平井压裂缝半长约为131 m;(3)部分压裂缝与解释断层相沟通,导致部分注入流体沿断层高渗透破裂带向下扩散,在基底地层激活断层并诱发MW3.0地震事件.其中孔隙压力的增加是断层失稳的主要因素;(4)增大水平井与已知断层之间的距离,被现场证实可以有效地降低地震风险.因此在致密储层进行水平井钻井及压裂施工作业之前,明确地下断层的分布状态至关重要.

致谢  图 1地震活动目录来自加拿大Alberta地震活动记录(www.inducedseismicity.ca/catalogues/).MW3.0诱发地震事件及压裂施工数据来自Bao和Eaton(2016).本次研究使用的三维地震数据由加拿大TGS公司提供.该研究得到Canada First Research Excellence Fund (CFREF)的资助.感谢卡尔加里大学地质系David Eaton教授的帮助.感谢两位审稿人的宝贵意见.
References
Atkinson G M, Eaton D W, Ghofrani H, et al. 2016. Hydraulic fracturing and seismicity in the western Canada sedimentary basin. Seismological Research Letters, 87(3): 631-647. DOI:10.1785/0220150263
Bao X W, Eaton D W. 2016. Fault activation by hydraulic fracturing in western Canada. Science, 354(6318): 1406-1409. DOI:10.1126/science.aag2583
Burridge R, Knopoff L. 1967. Model and theoretical seismicity. Bulletin of the Seismological Society of America, 57(3): 341-371.
Catalli F, Meier M A, Wiemer S. 2013. The role of Coulomb stress changes for injection-induced Seismicity: The Basel enhanced geothermal system. Geophysical Research Letter, 40(1): 72-77. DOI:10.1029/2012GL054147
Chen H C, Meng F B, Niu Y L, et al. 2018. Microseismic monitoring of stimulating shale gas reservoir in SW China: 2. spatial clustering controlled by the preexisting faults and fractures. Journal of Geophysical Research: Solid Earth, 123(B2): 1659-1672. DOI:10.1002/2017JB014491
Creaney S, Allan J. 1990. Hydrocarbon generation and migration in the Western Canada Sedimentary basin. Geological Society, London, Special Publications, 50(1): 189-202. DOI:10.1144/GSL.SP.1990.050.01.9
Eaton D W, Davidsen J, Pedersen P K, et al. 2014. Breakdown of the gutenberg-richter relation for microearthquakes induced by hydraulic fracturing: Influence of stratabound fractures. Geophysical Prospecting, 62(4): 806-818. DOI:10.1111/1365-2478.12128
Eaton D W. 2018. Passive Seismic Monitoring of Induced Seismicity: Fundamental Principles and Application to Energy Technologies. Cambridge: Cambridge University Press.
Eaton D W, Igonin N, Poulin A, et al. 2018. Induced seismicity characterization during hydraulic-fracture monitoring with a shallow-wellbore geophone array and broadband sensors. Seismological Research Letters, 89(5): 1641-1651. DOI:10.1785/0220180055
Eaton D W, Schultz R. 2018. Increased likelihood of induced seismicity in highly overpressured shale formations. Geophysical Journal International, 214(1): 751-757. DOI:10.1093/gji/ggy167
Ellsworth W L. 2013. Injection-induced earthquakes. Science, 341(6142): 1225942. DOI:10.1126/science.1225942
Eyre T S, Eaton D W, Garagash D I, et al. 2019. The role of aseismic slip in hydraulic fracturing-induced seismicity. Science, 5(8): eeav7172.
Fan Z Q, Eichhubl P, Newell P. 2019. Basement fault reactivation by fluid injection into sedimentary reservoirs: Poroelastic effects. Journal of Geophysical Research: Solid Earth, 124(7): 7354-7369. DOI:10.1029/2018JB017062
Haddad M, Eichhubl P. 2020. Poroelastic models for fault reactivation in response to concurrent injection and production in stacked reservoirs. Geomechanics for Energy and the Environment, 24: 100181. DOI:10.1016/j.gete.2020.100181
He D F, Lu R Q, Huang H Y, et al. 2019. Tectonic and geological background of the earthquake hazards in Changning shale gas development zone, Sichuan Basin, SW China. Petroleum Exploration & Development (in Chinese), 46(5): 993-1006.
Healy J H, Rubey W W, Griggs D T, et al. 1968. The Denver earthquakes. Science, 161(3848): 1301-1310. DOI:10.1126/science.161.3848.1301
Hui G, Chen S G, Chen Z X, et al. 2021. An integrated approach to characterize hydraulic fracturing-induced seismicity in shale reservoirs. Journal of Petroleum Science and Engineering, 196: 107624. DOI:10.1016/j.petrol.2020.107624
King G C P, Devès M H. 2015. Fault interaction, earthquake stress changes, and the evolution of seismicity.//Schubert G ed. Treatise on Geophysics. 2nd ed. Oxford: Elsevier, 4: 243-271.
Lei X L, Huang D J, Su J R, et al. 2017. Fault reactivation and earthquakes with magnitudes of up to MW4.7 induced by shale-gas hydraulic fracturing in Sichuan Basin, China. Science Report, 7(1): 7971. DOI:10.1038/s41598-017-08557-y
Lei X L, Wang Z W, Su J R. 2019. The December 2018 ML5.7 and January 2019 ML5.3 earthquakes in south sichuan basin induced by shale gas hydraulic fracturing. Seismological Research Letters, 90(3): 1099-1110. DOI:10.1785/0220190029
Lei X L, Su J R, Wang Z W. 2020. Growing seismicity in the Sichuan Basin and its association with industrial activities. Science China Earth Sciences, 63(11): 1633-1660. DOI:10.1007/s11430-020-9646-x
Lele S P, Tyrrell T, Dasari G R. 2017. Geomechanical analysis of fault reactivation due to hydraulic fracturing. //87th Ann. Internat Mtg., Soc. Expi. Geophys. . Expanded Abstracts.
McGarr A F. 2014. Maximum magnitude earthquakes induced by fluid injection. Journal of Geophysical Research: Solid Earth, 119(2): 1008-1019. DOI:10.1002/2013JB010597
Pawley S, Schultz R, Playter T, et al. 2018. The geological susceptibility of induced earthquakes in the Duvernay play. Geophysical Research Letters, 45(4): 1786-1793. DOI:10.1002/2017GL076100
Pedersen S I, Randen T, Sonneland L, et al. 2002. Automatic 3D fault interpretation by artificial ants.//72nd Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, G037.
Raleigh C B, Healy J H, Bredehoeft J D. 1976. An experiment in earthquake control at Rangely, Colorado. Science, 191(4233): 1230-1237. DOI:10.1126/science.191.4233.1230
Schultz R, Wang R J, Gu Y J, et al. 2017. A seismological overview of the induced earthquakes in the Duvernay play near Fox Creek, Alberta. Journal of Geophysical Research: Solid Earth, 122(1): 492-505. DOI:10.1002/2016JB013570
Schultz R, Atkinson G, Eaton D W, et al. 2018. Hydraulic fracturing volume is associated with induced earthquake productivity in the Duvernay play. Science, 359(6373): 304-308. DOI:10.1126/science.aao0159
Schultz R, Wang R J. 2020. Newly emerging cases of hydraulic fracturing induced seismicity in the Duvernay East Shale Basin. Tectonophysics, 779: 228393. DOI:10.1016/j.tecto.2020.228393
Shen L W, Schmitt D R, Haug K. 2019. Quantitative constraints to the complete state of stress from the combined borehole and focal mechanism inversions: Fox Creek, Alberta. Tectonophysics, 764: 110-123. DOI:10.1016/j.tecto.2019.04.023
Stiros S C, Kontogianni V A. 2009. Coulomb stress changes: From earthquakes to underground excavation failures. International Journal of Rock Mechanics and Mining Science, 46(1): 182-187. DOI:10.1016/j.ijrmms.2008.09.013
Tadokoro K, Ando M, Nishigami K. 2000. Induced earthquakes accompanying the water injection experiment at the Nojima fault zone, Japan: Seismicity and its migration. Journal of Geophysical Research: Solid Earth, 105(B3): 6089-6104. DOI:10.1029/1999JB900416
Wang H F. 2000. Theory of Linear Poroelasticity with Applications to Geomechanics and Hydrogeology. Princeton, NJ: Princeton University Press.
Wang R J, Gu Y J, Schultz R, et al. 2018. Faults and non-double-couple components for induced earthquakes. Geophysical Research Letters, 45(17): 8966-8975. DOI:10.1029/2018GL079027
Weir R, Lawton D, Lines L, et al. 2019. Application of structural interpretation and simultaneous inversion to reservoir characterization of the duvernay formation, Fox Creek, Alberta, Canada. The Leading Edge, 38(2): 151-160. DOI:10.1190/tle38020151.1
Yehya A, Yang Z, Rice J R. 2018. Effect of fault architecture and permeability evolution on response to fluid injection. Journal of Geophysical Research: Solid Earth, 123(11): 9982-9997. DOI:10.1029/2018JB016550
Yew H C, Weng X W. 1997. Mechanics of Hydraulic Fracturing. Oxford: Gulf Professional Publishing.
Zhang D X, Yang T Y. 2015. Environmental impacts of hydraulic fracturing in shale gas development in the United States. Petroleum Exploration & Development (in Chinese), 42(6): 801-807.
Zhang H L, Eaton W D, Rodriguez G, et al. 2019. Source-Mechanism analysis and stress inversion for hydraulic-fracturing-induced event sequences near Fox Creek, Alberta. Bulletin of the Seismological Society of America, 109(2): 636-651. DOI:10.1785/0120180275
Zoback M D. 2007. Reservoir Geomechanics. Cambridge: Cambridge University Press.
Zoback M D, Gorelick S M. 2012. Earthquake triggering and large-scale geologic storage of carbon dioxide. Proceedings of the National Academy of Sciences of the United States of America, 109(26): 10164-10168. DOI:10.1073/pnas.1202473109
何登发, 鲁人齐, 黄涵宇, 等. 2019. 长宁页岩气开发区地震的构造地质背景. 石油勘探与开发, 46(5): 993-1006.
张东晓, 杨婷云. 2015. 美国页岩气水力压裂开发对环境的影响. 石油勘探与开发, 42(6): 801-807.
张致伟, 程万正, 梁明剑, 等. 2012. 四川自贡-隆昌地区注水诱发地震研究. 地球物理学报, 55(5): 1635-1645.