青海共和盆地干热岩地应力测量及其储层压裂改造意义分析

In-situ stress measurements in hot dry rock, Qinghai Gonghe Basin and simulation analysis of reservoir fracture modification

摘要

近年来,我国在干热岩资源勘查和钻井技术方面取得一定的进展,但是对于高温干热岩压裂技术的基础研究还较薄弱。对于增强型地热系统,一般采用水力压裂的方式向干热岩储层注入高压流体使储层中的天然裂隙扩展延伸,从而达到增大储层渗透性和热交换面积的目的。因此,可靠的地应力数据对指导干热岩储层改造具有重要意义。本文以我国青海共和干热岩试采工程为例,首先,结合滞弹性应变恢复法(ASR)、直径变形分析法(DCDA)、岩心饼化法和成像测井法获得场区详细的地应力数据;其次,讨论了地应力状态对干热岩储层开发的影响;最后,为了评价现今应力状态下储层压裂改造的效果,利用离散裂隙网格(DFN)方法建立了水力压裂区的三维裂隙地质模型,通过数值模拟计算每个裂隙的滑动趋势(T_s)和膨胀趋势(T_d),并分析了在不同的注水压力下裂隙的活动性。研究表明:(1)共和干热岩场区花岗岩储层范围(3 500~4 000 m)内主要发育逆断层型地应力状态,总体以水平挤压应力为主导;(2)在3 500~4 000 m深度水平最大主应力的平均方向为39.35°±14.23°,整体上呈NE向,这与青藏高原东北缘向NE方向推挤运动有关;(3)适合共和干热岩场区的注水改造压力为46~55 MPa,且裂隙剪切活化后多为张开裂隙(高膨胀趋势),有利于增加热交换的面积,提高干热岩的开采效率;(4)地应力实测结果和微地震监测结果显示,天然裂隙在水力压裂下除了水平方向上的扩展外还有垂直方向上的扩展,最终形成水平—垂直裂隙网络。本文研究成果也可为我国未来的干热岩地应力测量及其在开发评价中的应用提供参考。

关键词

青海共和 / 干热岩 / 地应力测量 / 水力压裂 / 储层评价

Key words

Qinghai Gonghe / hot dry rock / in-situ stress measurement / hydraulic fracturing / reservoir evaluation

0 引言

地热能作为一种极具竞争力的清洁可再生能源,对实现碳达峰、碳中和目标有着重要的现实意义^[1⇓-3]。我国已探明的地热能储量丰富,尤其是可开采的、利用价值高的浅层地热能在全国范围内普遍分布^[1,4-5]。国际上,自20世纪70年代以来,以美国、法国为代表的国家相继开展干热岩开发试验研究^[6]。近年来,我国在干热岩资源勘查和钻井技术方面取得一定的进展,但是对于高温干热岩压裂技术的基础研究还较薄弱^[7-8]。尤其是开发过程中的水力压裂、裂隙改造和注水诱发地震等问题还需要更精细化的研究^[9-10]。其中详细的地应力数据对干热岩开发过程中的水力刺激方案设计、储层改造响应解释和诱发微地震预防有着重要意义^[11]。

增强型地热系统(Enhanced Geothermal System,简称EGS)是干热岩开发的主要方式,即采用水力压裂的方式向干热岩储层注入高压流体使储层中的天然裂隙扩展延伸,从而达到增大储层渗透性和热交换面积的目的^[6,12⇓-14]。不同的地应力类型对裂隙改造的影响不同,通常正断层和走滑型应力易于裂隙激活,而水平应力差较大的逆断层型应力易于激活水平裂隙。地应力量值的增大会使破裂压力增大,影响储层改造效率^[15]。地应力的差异也会影响压裂效果,主应力差越小,裂隙延伸和转向所需的净压力越小,裂隙网络的复杂度越高,进而可以制造足够大的换热面积^[16-17]。此外水力压裂中,人工裂隙面一般与最小主应力方向垂直相交,并沿最大主应力方向扩展延伸,当主应力方向与储层中自然裂隙方向夹角位于优势区间内时,有利于部署生产井和制造规模更大的裂隙网络^[18]。此外,地应力数据还可以辅助评价储层改造效果。在水力压裂过程中,高压流体的注入影响了天然裂隙表面的有效应力,使其最终失稳发生剪切破坏^[19]。因此,在干热岩的开发过程中,详细的地应力数据是不可或缺的。

本文以青海共和干热岩试采工程为例,利用滞弹性应变恢复法(ASR)、直径变形分析法(DCDA)、岩心饼化法和成像测井法获得干热岩场区详细的地应力数据,其次根据地应力数据和压裂参数模拟储层水力压裂改造,评价了该项目的储层改造效果。

1 研究区地质与干热岩开发概况

1.1 研究区地质背景

共和盆地位于青藏高原东北缘,地处中央造山带北西向延伸的构造部位,是东昆仑断裂、阿尔金断裂和南祁连断裂共同左旋走滑断裂作用下控制的新生代陆内盆地^[20]。盆地西侧以瓦洪山走滑断裂与柴达木盆地相邻;盆地东侧以多禾茂断裂为界与巴彦克拉盆地相邻;盆地南侧以阿尼玛卿缝合带为界与松潘—甘孜褶皱带相邻,盆地北侧以青海湖南山断裂为界与祁连造山带相邻。共和盆地整体呈菱形展布,黄河从盆地中间北东向穿越,形成深切峡谷,并将盆地分为东、西两部分。其盆地东侧包括贵南—沙沟—贵德地区,盆地西侧主要包括茶卡—恰卜恰广阔地区(图1)^[21-22]。

共和盆地的基底形成于三叠纪,定型于白垩纪中晚期,总体上经历了盆地形成、湖泊盆地发展、盆地扩展、盆地固定和间隔抬升的复杂过程^[23]。盆地周围发育有东昆仑和西秦岭造山带的火成岩序列和变质沉积岩,以中—晚三叠世变质沉积岩、印支期花岗岩和花岗闪长岩为主,在盆地西南部发育NW向—近NS向构造岩浆岩隆起带,具有走滑性质。共和盆地干热岩赋存岩体形成于中—晚三叠世,中生代经历抬升剥蚀岩体出露至地表^[23]。

1.2 青海共和干热岩工程概况

共和盆地有华力西和印支两期的岩浆活动,其中印支期占主要地位^[24]。火山活动或现代幔源岩浆侵入作用为侵入岩的形成提供了条件。盆地中部有深大隐伏断裂带穿切而过,这种隐伏的深大断裂有利于形成高温对流型地热田^[25]。恰卜恰地热田的成因推测可能是深部熔融体的热量通过隐伏断层传输到较浅的花岗岩体中,上部较厚的第四纪沉积地层阻止了热量逃逸,从而形成埋藏较浅的高品位干热岩资源^[26]。

图2a展示了EGS的工作原理,即通过注入井向干热岩体中注入高压流体使裂隙失稳扩展,并利用注入井与生产井之间的联通性裂隙,将完成热交换的液体从生产井中抽出,达到发电利用的目的。共和EGS工程场地位于共和盆地东北部恰卜恰地区,紧邻海南州共和县城。目前,恰卜恰地热田已经圈定的干热岩面积达246.9 km²,资源总量达13.66×10¹⁸ J,折合标准煤4.66亿t^[27]。场区深度0~<1 360 m为新近纪河湖相沉积盖层,1 360 m深度之下为中三叠世花岗岩体(图2b)。场地蕴含丰富的地热资源,花岗岩地层中的地温梯度接近正常地层的两倍。

2019—2021年,中国地质调查局启动了我国首个大型干热岩试验性开发工程,在共和干热岩场地共钻探了3口4 000 m深井。其中,GH-01井为注入井,GH-02和GH-03井为生产井。本文拟采用水力改造区域深度(3 500~3 900 m)内GH-01井岩心。3口井的3 500 m以上采用套管固井,3 500~3 900 m井身范围为花岗岩裸岩段,也是设计的水力改造区域,注入的液体将会在这个范围内进行热交换从而被生产井回收利用(图2b,c)。

2 干热岩地应力测量方法

主要的地应力测量方法有基于钻孔的水压致裂法和成像测井法,基于岩心的滞弹性应变恢复法(ASR)、直径变形分析法(DCDA)和岩心饼化法,以及基于地质的震源机制解。虽然地应力的测量方法有多种,但每一种都有其局限性^[28]。传统的水压致裂法在非高温深孔测量的最大深度已达到5 300 m,但由于干热岩储层特殊的高温环境难以成功应用。成像测井法是通过识别孔壁的破坏特征,进而约束地应力,只要成像测井探头能够获取孔壁的图像就能够成功应用。基于岩心的地应力测量方法都可以应用于高温干热岩钻孔中^[11]。岩心饼化法适用于高应力环境,共和干热岩储层具有发生岩心饼化的应力条件,因此该方法可以作为一种半定量的地应力评价方法。DCDA法和ASR法都是通过岩心表面的应变来测量地应力,虽然高温能够影响岩心变形,但这种变形是均匀的,对测量带来的误差可以忽略不计。其中,ASR法普遍应用于深部地热及干热岩地应力测量,测量效果理想。

2.1 ASR方法

天然岩石可以被认为是一种黏弹性介质。因此,当岩心从地下被取出时,由于与非均匀地壳应力场分离,一部分弹性应变立即释放,而另一部分应变随着时间的推移逐步释放(滞弹性应变恢复)(图3)。研究表明,滞弹性应变可以对原位应力进行约束^[29]。考虑到3个主应力的方向与3个主滞弹性应变的方向一致,因此可以通过仪器测量岩心至少6个独立方向上的滞弹性应变恢复量来计算主应变及其后的主应力方向。从0到t的时间内,岩心样品在任意方位角的滞弹性应变恢复量(ε_a(t))与下式有关^[29]:

(1)ε_a(t)= \frac{1}{3}[(3l²-1)σ_x+(3m²-1)σ_y+
(3n²-1)σ_z+6lmτ_xy+6mnτ_yz+6nlτ_zx]Jas(t)+
(σ_m-p₀)Jav(t)+α_TΔT(t)

式中:σ_x、σ_y、σ_z、τ_xy、τ_yz和τ_zx是原位应力张量的不同分量;p₀是孔隙压力;σ_m是平均法向应力;σ_T是线性热膨胀系数;ΔT(t)是0到t内温度的变化;Jav(t)和Jas(t)分别是体积和剪切变形模式下滞弹性应变恢复柔度。3个主应力的大小可以利用下式计算^[30]:

(2)σ_i= \frac{e_i\left(t\right)}{jas\left(t\right)}+ \frac{\left[e_m\right(t)-{\alpha }_T\Delta T(t\left)\right]}{\mathrm{j}\mathrm{a}\mathrm{v}\left(t\right)}+p₀

式中:e_m(t)是平均法向应变;e_i(t) (i=1,2,3)是主要的滞弹性应变偏差。当去除地应力时,岩石基质滞弹性应变ε_a(t)主要取决于原位应力张量分量、热膨胀系数和试验过程中温度的变化、孔隙压力和体积与剪切变形模式下滞弹性应变恢复柔度(公式(1))。因此,如果在试验期间环境温度保持恒定,则可以利用公式(2)计算三维主应力的大小。

如上所述,ASR测量结果的准确性受残余应变、岩心采取时间、岩心各向异性、试验环境温度和岩心脱水作用等因素的影响^[30⇓-32]。及时采集充足的应变数据是获得可靠地应力数据的关键,因此普遍在钻井现场就近建立实验室。温度和岩心失水是导致非应力性应变的主要干扰因素,试验过程中须将岩心密封包裹并置于恒温水浴箱内。ASR方法要求岩心样品完整、均质且各向同性。以往测试经验表明,如果没有肉眼可见的微裂隙,花岗岩岩心视为各向同性和均质是合理的^[32]。图4是岩心表面应变片粘贴方式,这样能够测量岩心各个方向上由应力释放引起的滞弹性应变。

2.2 DCDA方法

直径变形分析法(DCDA)的基本原理是,当岩心从高地应力钻孔中取出后,由于应力的释放,会产生瞬时的弹性变形(图5^[33])。理想状态下,取芯钻头所取的岩心应为等值线圆柱形,地应力S_H和S_h不同,导致岩心截面变为椭圆形(图5)^[34]。根据线弹性理论,假设岩石均匀,各向同性且变形较小,则S_H和S_h方向上的拉伸应变ε_max和ε_min分别由下式给出^[33]:

(3)ε_max= \frac{1}{E}[S_H-υ(S_h+S_z)]

(4)ε_min= \frac{1}{E}[S_h-υ(S_H+S_z)]

式中:S_z是平行于钻孔轴线的地应力,对于垂直钻孔,S_z等于S_V;E和υ分别是岩石的杨氏模量和泊松比。此外,ε_max和ε_min还可以通过最大芯径d_max和最小芯径d_min来计算:

(5)ε_max= \frac{d_{\mathrm{m}\mathrm{a}\mathrm{x}}-d_0}{d_0}, ε _min= \frac{d_{\mathrm{m}\mathrm{i}\mathrm{n}}-d_0}{d_0}

式中,d₀为应力释放前岩心直径,根据ε_max和ε_min以及标志线对应圆周角θ处的直径d_θ,计算得到该处的应变ε_θ:

(6)ε_θ= \frac{d_{\theta }-d_0}{d_0}=ε_maxcos²(θ-α)+ε_minsin²(θ-α)

式中,α为最小芯径处的θ,结合公式(3)和(4)得到d_θ的表达式:

(7)d_θ= \frac{d_{\mathrm{m}\mathrm{a}\mathrm{x}}+d_{\mathrm{m}\mathrm{i}\mathrm{n}}}{2}+ \frac{d_{\mathrm{m}\mathrm{a}\mathrm{x}}-d_{\mathrm{m}\mathrm{i}\mathrm{n}}}{2}cos²(θ-α)

由于最小芯径d_min与应力释放前芯径d₀相差很小,用d_min代替d₀,结合公式(3)、(4)和(5),得到表达式:

(8)S_H-S_h= \frac{E}{1+\upsilon } \frac{d_{\mathrm{m}\mathrm{a}\mathrm{x}}-d_{\mathrm{m}\mathrm{i}\mathrm{n}}}{d_0}≈ \frac{E}{1+\upsilon } \frac{d_{\mathrm{m}\mathrm{a}\mathrm{x}}-d_{\mathrm{m}\mathrm{i}\mathrm{n}}}{d_{\mathrm{m}\mathrm{i}\mathrm{n}}}

因此,差应力(S_H-S_h)可以由测定的d_max和d_min确定。

2.3 岩心饼化方法

在高地应力区域,可以用岩心饼化法来估算原地应力状态。岩心饼化的形貌特征一定程度上反映了地应力状态,岩饼外沿轮廓的“马鞍状”凹形轴线可以指示最大水平主应力的方向^[35]。通过构建岩心饼化发生的临界地应力状态与岩饼厚径比的经验和理论关系式,可以估算最大主应力的大小^[36-37]。

对于垂直钻孔,无论岩石类型如何,最大主应力与岩心饼化直径的关系式如下^[36]:

(9)S_H=A+ \frac{B}{exp\left(\frac{t}{C}\right)\mathrm{ }}

式中,A、B和C是曲线拟合参数,对于花岗岩,A=63.70,B=75.11,C=3.20。

3 地应力测量结果

3.1 ASR应变曲线

图6显示了ASR试验样品的原始滞弹性应变曲线和测试环境温度。共和钻井取芯的时间从夏季一直持续到秋季,为了使温度控制器与现场环境温度兼容,将测试系统的温度设置为20~30 ℃。由于工作区海拔较高,现场实验室的昼夜温差较大,有1件样品的最大温度波动约为3 ℃,其余3件样品在近似恒温状态下完成了数据采集。温度波动同步引起了热应变的变化,显示出一些应变“波峰”和“波谷”,但没有改变滞弹性应变恢复的趋势^[38]。因此,其对主应力计算的影响可以忽略不计。滞弹性应变在整个实验过程中以相似的趋势向各个方向延伸。总体上,试验曲线符合滞弹性应变恢复的变化特征。

图7是根据原始滞弹性应变数据计算的主应变和平均应变曲线。在本次研究中,实验的样品为花岗岩,可以认为是均质和各向同性的。因此,主应变和主原位应力方向可以认为是一致的。根据测量的滞弹性应变值,使用最小二乘法将滞弹性应变张量计算为与时间(t)的函数。在图6中,展示了ε₁、ε₂、ε₃和ε_m的计算结果。在测量的时间里,这些花岗岩样品在各个方向上恢复了20~350微应变范围内滞弹性应变值(图6)。与采集系统的测量精度相比,这些应变测量值足够大,可以满足三维主应力分析。

3.2 地应力大小

如图6和图7所示,温度波动并没有影响整体膨胀趋势和应变特征,因此热应变对原位应力计算的影响可以忽略不计^[32]。使用Jas(t)/Jav(t)=2.0来计算原位主应力的大小。Lin等^[30]分析了Jas(t)/Jav(t)=1.0~3.0对主应力的影响,指出σ₁和σ₃的不确定性约为5%,而σ₂的不确定性可以忽略不计。利用4个深度的岩心进行ASR测试,得到如下结果:S_H=104.6~131.1 MPa,S_h=96.6~107.5 MPa,垂直应力(S_V)根据密度测井数据计算,该深度范围内S_V=88.0~99.9 MPa(图8)。利用DCDA方法测量了3 884 m处岩心的直径变形数据,计算得到水平应力差(S_H-S_h)为25.14 MPa,与ASR测量的结果一致。最后,根据共和干热岩场区的岩石力学资料,共和花岗岩的抗拉强度为14.2~21.6 MPa,应力比因子约为6.5,计算得S_H=92.3~140.4 MPa(图8)。

最终的结果如图8所示,在3 500~3 900 m深度范围内,3个主应力的关系为S_H>S_h>S_V,为逆断层应力状态,这说明构造板块运动对区域应力场的影响大于重力作用。

3.3 地应力方向

在地壳各向异性应力、温度差异和井孔泥浆的作用下,钻孔周围的应力集中会导致井壁上发生压缩破坏和拉伸破坏。井壁上的应力分布是不均匀的,当井壁上的挤压应力大于井壁岩石的单轴抗压强度时,会在井壁形成一对相隔180°的条形破碎带,称为钻孔崩落(BBO)(图9b);同样地,当井壁上的拉伸应力大于井壁岩石的抗拉强度时,会在井壁上形成一对相隔180°的线形拉伸裂纹,称为诱发张裂隙(DIF)(图9a)。由于井壁的应力关系,钻孔崩落会沿着水平最小主应力的方向发育,诱发张裂隙会沿着水平最大主应力方向发育^[39-40]。因此,可以利用井壁上这种特殊的现象,结合成像测井资料来确定地应力的方向。

统计了共和干热岩场地钻井3 500 m至井底的成像测井数据,发现钻孔崩落和诱发张裂隙明显发育。最终水平最大主应力的方向如图9c所示,在3 500~4 000 m深度水平最大主应力的平均方向为39.35°±14.23°。这与其他学者根据震源机制解得到的共和区域构造应力场水平最大主应力方向一致,整体上呈NE向^[41-42],这与青藏高原东北缘向NE方向推挤运动有关^[43-44]。

3.4 地应力状态对干热岩储层开发的意义

干热岩储层改造成功的关键是通过水力压裂制造注入井和生产井之间的优势性裂隙网络通道,提高储层岩体的渗透率,进而增强热交换和流体回收的效率^[12]。在不同的地应力状态下,会引起不同方向的水力裂隙活动。通常在正断层或走滑型应力状态下易形成垂直的水力裂隙网络,而在逆断层型应力状态下易形成水平的裂隙网络^[45]。地应力的量值大小也会影响压裂的效率,地应力量值越大,破裂压力越大,压裂的效率越低^[15]。水平主应力的差异控制着水力裂隙的形态,水平主应力的差值越小,水力裂隙延伸和转向所需的压力越小,形成的裂隙网络复杂度越高,进而可以制造足够大的换热面积^[16-17]。此外,岩石基质中形成的水力裂隙较为简单,一般与最小主应力方向垂直相交,并沿着最大主应力方向扩展延伸,但会在与天然裂隙的交界处形成复杂裂隙网络^[18,24]。

由前文可知,干热岩储层的水力压裂区(3 500~3 900 m)为走滑—逆冲应力类型的转换区(图8),S_h与最小主应力S_V的量值十分接近,这可能导致二者对水力裂隙扩展的控制能力相当。因此,本文推断,储层改造时,水力裂隙除了会在水平方向上发展,同时也可能在垂直方向上扩展,裂隙延伸的理论方向应与S_H方位保持一致。天然裂隙的主要产状为NW15°,水平最大主应力的方向为NE39.35°±14.23°。因此,该区域水力裂隙的特征可能是天然裂隙向NW方向扩展,水力压裂产生的裂隙向NE方向扩展,然后,可能在与天然裂隙的交界处形成复杂的裂隙网络。

4 干热岩储层压裂效应评价

4.1 水力压裂效应评价方法

天然裂隙或断层的活动受地应力和孔隙流体压力共同控制,孔隙流体压力的增加会降低裂隙表面的有效应力,从而增加天然裂隙和断层表面重新滑动的可能性。但孔隙流体压力的变化不会影响差应力,也不会影响主应力的方向。因此,可以利用有效正应力来理解裂隙的激活与扩展^[46⇓-48]。有效正应力被定义为主应力减去孔隙流体压力,当流体进入天然裂隙中,孔隙流体压力增加使有效正应力降低,天然裂隙易失稳发生剪切破坏。裂隙表面的正应力(σ_n)和剪应力(τ)可以由以下公式计算:

(10) σ_n=σ₁l²+σ₂m²+σ₃n²

(11) τ=(σ₁-σ₂)²l²m²+(σ₂-σ₃)²m²n²+
(σ₃-σ₁)²l²n²

式中,l、m和n是裂隙面法向方向余弦。

利用滑动趋势(T_s)和膨胀趋势(T_d)可以评价断层或裂隙重新激活的可能性^[49-50]。滑动趋势是作用在断层或裂隙上的剪应力与正应力之比。对于不考虑内聚力的断层或裂隙,如果滑动趋势超过摩擦系数的值,就会发生滑移。在分析孔隙流体压力对断层稳定性的影响时,研究人员通常使用滑动趋势(T_s)来量化给定方向的断层或裂隙重新活动的可能性^[44]。一般情况下,断层或裂隙的临界摩擦系数(滑动趋势)为0.6,即滑动趋势超过0.6时断层或裂隙易发生滑动^[19,51-52]。膨胀趋势是指断层或裂隙在三维应力条件下裂隙面膨胀开启的趋势或可能性,膨胀趋势的范围为0~1,数值越高说明扩展的可能性越高。裂隙的滑动和膨胀趋势可以由以下公式计算:

(12) T_s=τ/σ_n

(13) T_d=(σ₁-σ_n)/(σ₁-σ₃)

在本文中,根据ASR法和成像测井法已经获得了共和干热岩场区详细的地应力实测数据。在这些实测资料的基础上,用滑动趋势(T_s)和膨胀趋势(T_d)来评价压裂段(3 500~3 900 m)天然裂隙在水力压裂下的改造效果。

4.2 储层裂隙建模和压裂模拟

共和干热岩场地的3口钻井均有成像测井资料,能够观察到井壁的裂隙产状,包括方向和倾角。但是这种方法获得的裂隙是局部的,只能反映该钻孔附近的特征,为了能将裂隙推广到三维中,利用离散裂隙网格(DFN)方法。即根据已获得的裂隙信息为基础,添加约束条件,通过随机模拟的方式生成相同概率裂隙。根据该方法建立了共和干热岩压裂区(3 500~3 900 m)的裂隙地质模型。由于本文主要关注在压裂作用下,注入井和生产井之间裂隙改造情况,因此裂隙网格模型的水平设置为以注入井GH-01井为中心,范围100 m的正方形,高度设置为3 500~3 900 m共400 m的范围,最终形成长、宽各200 m,高400 m的长方体裂隙网格。

对于注水压裂的模型,采用Hantush 推荐的方法,即计算孔隙压力随时间和注入位置距离的变化^[53]。假设注入的流体满足扩散方程的解,由点源或条段均匀注入各向异性的渗透介质,由于岩体内部存在各个方向的裂隙,因此注入液在压力的作用下可能渗透到相邻的岩石中,引起不同区域孔隙压力的变化。然后,将这些变化应用于计算裂隙面上的正应力和剪应力大小,进而通过滑动趋势和膨胀趋势分析其活动性。

干热岩储层压裂段(3 500~3 900 m)的花岗岩渗透率非常低。然而,根据成像测井发现原生裂隙发育,这些裂隙往往成为流体扩散的主要通道,因此实际岩体中渗透率往往偏大^[19]。结合现场花岗岩渗透率实测数据,最终将渗透率设置为100 mD。除了岩石基质渗透率外,井下压力同样会对模拟的最终结果产生影响。共和干热岩现场压裂资料显示,实际水力压裂过程中,地面注入泵的压力范围为30~66 MPa。这里,为了评价不同压力对裂隙压裂的影响,模拟过程中将注入压力的范围设置为30~60 MPa。

4.3 储层裂隙活动性结果与评价

由上所述,通过ASR和成像测井方法分别获得了共和干热岩场地的地应力大小和方向,根据DFN方法建立了储层压裂区的裂隙三维地质模型。在3 500~3 900 m处,令σ₁=S_H,方向为NE39°;σ₂=S_h,方向与σ₁垂直;σ₃=S_V,方向是垂直的。基质的渗透率为100 mD,加压泵注水压力分别为30、46、55和60 MPa,最终得到不同压力下裂隙表面的滑动趋势(T_s)和膨胀趋势(T_d)。图10显示了4种不同注水压力下,裂隙表面的滑动趋势变化。随着注水压力的增大,可以看到T_s也在增大,且浅层的裂隙较深层裂隙相比更容易被激活。

Byerlee通过总结不同类型岩石的剪切试验数据,得出岩石的静摩擦系数范围为0.6~1.0,即当滑动趋势T_s大于0.6时可以认为裂隙处于临界失稳状态,即将或已经发生失稳滑动^[54]。

为了更直观地认识到裂隙表面状态,本文利用统计学方法,将模拟结果进行累加分配处理,即将小于某点的数据进行累加后再计算占总体的比例(图11)。如图11a所示,当注水压力为30 MPa时,所有的裂隙T_s都在0.3以下,这说明该压力下天然裂隙不可能滑动,是达不到改造要求的;随着泵注压力继续增大到46 MPa,裂隙面最大T_s达到了0.6,这表明此时压裂区有些天然裂隙已处于临界应力状态;继续增压到55 MPa,约有90%裂隙的T_s超过了0.6,这说明绝大部分裂隙处于剪切破坏应力状态,裂隙可能已经以剪切破坏的方式发生滑动;继续增压到60 MPa,99%的裂隙的T_s超过了1.0,这说明压裂区裂隙绝大部分已经发生剪切破坏。在这样的压力下虽然改造的效果十分理想,但易于诱发有感地震,影响工程安全进行;另一方面,由于储层压裂段裂隙内孔隙压力可能已经超过最小主应力,裂隙将产生Ⅰ型拉张破坏,裂隙张开度过大,进而可能导致水流速度偏大,无法充分吸收岩体热量,甚至有产生热突破的风险^[55-56]。因此,不建议地面泵注压力超过60 MPa。结合模拟结果,建议压裂区优势注水压力范围为46~55 MPa。需要指出的是,本文暂不考虑泵注过程中管柱摩擦阻力的影响,实际泵注压力上限将有所增加。

裂隙的膨胀趋势T_d并不随着注入压力的变化而改变,而只是和应力状态有关,这是因为虽然注入压力增加改变了压裂区的有效应力,但是差应力并没有发生变化。如图11b所示,压裂区裂隙的膨胀趋势处于0.4和1.0之间,这是压裂区处于逆断层应力状态且水平最小主应力(S_h)与垂直应力(S_V)相差较小导致的,这说明裂隙被激活后是处于张开状态的,这有利于构建高渗透性的裂隙联通通道。

图12为注水压裂改造时的诱发微地震云图,反映了从期次1到期次3注水期间,储层微地震事件空间和大小分布情况。由图12a可以看到,3 800 m深度之下微地震事件较少,这与模拟的结果一致(图10c),压裂区浅部裂隙较深部裂隙相比更容易被激活。逆断层型应力一般易于形成水平方向上的水力裂隙,走滑型应力一般易于形成垂直方向上的水力裂隙^[45]。青海共和干热岩场地压裂区水平最小主应力(S_h)与垂直应力(S_V)接近,为走滑—逆冲型应力转换区,应力对水力裂隙在水平和垂直方向上的发展控制效果是接近的,因此天然裂隙会在水力压裂下沿着自身的产状延伸扩展最终形成水平—垂直水力裂隙网络。随着压裂期次的增加,微地震事件由深部向浅部发展,并有部分事件超出了压裂区的上限范围,在3 500 m深度之上出现,这表明在压裂时裂隙除了表现水平方向上的延伸外,还有垂直方向上的扩展,这与应力测量的特征一致(图12a)。由图12b可以看到,微地震事件优势方向为NW,与天然裂隙的走向一致,这说明水力裂隙优先沿着自身产状方向发展,而并未在最大主应力方向,即NE向扩展和延伸,这与已有理论认识不一致。这是因为在实际压裂过程中,期次1注水排量为3.04 m³/min且最大注水压力为48 MPa,期次2、3注水排量为2.4 m³/min且最大注水压力为46 MPa,都只是处于临界应力状态(图10b)。并且在这种小排量的压裂模式下,流体不断注入优先充填天然裂隙,在裂隙尖端形成应力集中,最终导致天然裂缝失稳发生剪切失稳,而很难形成新的水力裂隙,因此水力裂隙的主要扩展方向与天然裂隙一致。未来,应加强高泵注流量下的实验或模拟工作,以分析其对裂隙扩展方向的影响。

5 结论与认识

为了了解青海共和干热岩储层的原位地应力大小和方向,本文在青海共和干热岩场地采用岩心ASR法和成像测井方法,获得了场区详细的地应力实测资料。通过这些资料,评价和对比分析了干热岩水力压裂改造的效果。主要认识如下。

(1)在3 500~4 000 m深度开展了ASR、DCDA和岩心饼化地应力测试,结果表明,青海共和干热岩场区花岗岩储层深度范围内,发育逆断层型应力状态。总体上,以水平挤压应力为主导。

(2)成像测井结果显示,在3 500~4 000 m深度水平最大主应力的平均方向为39.35°±14.23°。结果与共和区域构造应力场水平最大主应力方向一致,整体上呈NE向,这与青藏高原东北缘向NE方向推挤运动有关。

(3)利用滑动趋势和膨胀趋势,通过模拟的方式评价了青海共和干热岩场区压裂改造效果。结果显示,注水压力达到46 MPa时压裂区处于临界应力状态;达到55 MPa时约90%的裂隙处于剪切破坏应力状态;达到60 MPa时99%的裂隙的T_s已经超过了1.0,这说明压裂区裂隙绝大部分已经发生剪切破坏。超过60 MPa压力下不仅易诱发有感地震,还有产生热突破的风险。因此,建议压裂区的地表注水压力范围为46~55 MPa。

(4)地应力实测结果和微地震监测结果显示,天然裂隙在水力压裂下除了水平方向上的扩展外,还有垂直方向上的扩展,最终形成水平—垂直裂隙网络。此外,由于共和干热岩压裂区的小排量压裂工艺,水力裂缝优先沿自身产状方向发展,而并未在最大主应力方向,即NE向扩展和延伸。未来,应加强高泵注流量下的实验或模拟工作,以分析其对裂隙扩展方向的影响。

感谢中国地质科学院水文地质环境地质研究所王贵玲研究员邀请参与本专辑的撰稿,感谢审稿专家对本文的建设性建议,感谢同济大学杜泊潼博士在裂隙建模方面提供的帮助。

参考文献

[1] 王贵玲, 刘彦广, 朱喜, 中国地热资源现状及发展趋势[J]. 地学前缘, 2020, 27(1): 1-9.

[2] 王贵玲, 蔺文静, 刘峰, 地热系统深部热能聚敛理论及勘查实践[J]. 地质学报, 2023, 97(3): 639-660.

[3] LIN W J, WANG G L, GAN H N, et al. Heat source model for Enhanced Geothermal Systems (EGS) under different geological conditions in China[J]. Gondwana Research, 2023, 122: 243-259.

[4] 许天福, 张延军, 曾昭发, 增强型地热系统(干热岩)开发技术进展[J]. 科技导报, 2012, 30(32): 42-45.

[5] WANG G L, GAN H N, LIN W J, et al. Hydrothermal systems characterized by crustal thermally-dominated structures of southeastern China[J]. Acta Geologica Sinica (English Edition), 2023, 97(4): 1003-1013.

[6] 蔺文静, 王贵玲, 邵景力, 我国干热岩资源分布及勘探: 进展与启示[J]. 地质学报, 2021, 95(5): 1366-1381.

[7] 付亚荣, 李明磊, 王树义, 干热岩勘探开发现状及前景[J]. 石油钻采工艺, 2018, 40(4): 526-540.

[8] 陈作, 许国庆, 蒋漫旗. 国内外干热岩压裂技术现状及发展建议[J]. 石油钻探技术, 2019, 47(6): 1-8.

[9] KIM K H, REE J H, KIM Y, et al. Assessing whether the 2017 M_w5.4 Pohang earthquake in South Korea was an induced event[J]. Science, 2018, 360(6392): 1007-1009.

[10] SEITHEL R, GAUCHER E, MUELLER B, et al. Probability of fault reactivation in the Bavarian Molasse Basin[J]. Geothermics. 2019, 82: 81-90.

[11] 许家鼎, 张重远, 缑艳红, 干热岩地应力测量评价方法与前沿挑战[J]. 地球学报, 2023, 44(1): 200-210.

[12] BARIA R, BAUMGÄRTNER J, RUMMEL F, et al. HDR/HWR reservoirs: concepts, understanding andcreation[J]. Geothermics, 1999, 28(4/5): 533-552.

[13] BROWN D W, DUCHANE D V, HEIKEN G, et al. Mining the earth’s heat: hot dry rock geothermal energy[M]. Berlin, Heidelberg: Springer, 2012.

[14] OLASOLO P, JUÁREZ M C, MORALES M P, et al. Enhanced geothermal systems (EGS): a review[J]. Renewable and Sustainable Energy Reviews, 2016, 56: 133-144.

[15] 赵金洲, 任岚, 胡永全. 页岩储层压裂缝成网延伸的受控因素分析[J]. 西南石油大学学报(自然科学版), 2013, 35(1): 1-9.

[16] MAO R B, FENG Z J, LIU Z H, et al. Laboratory hydraulic fracturing test on large-scale pre-cracked granitespecimens[J]. Journal of Natural Gas Science and Engineering, 2017, 44: 278-286.

[17] 沈骋, 赵金洲, 谢军, 海相页岩缝网可压性靶窗空间分布预测: 以川南长宁区块为例[J]. 地质力学学报, 2020, 26(6): 881-891.

[18] JUNG R. EGS: goodbye or back to the future 95[M]//BUNGER A, JEFFRE R. Effective and sustainable hydraulic fracturing. Brisbane, Australia: InTech, 2013: 458-477.

[19] MORRIS A P, FERRILL D A, WALTER G R, et al. Lessons learned from the Youngstown, Ohio induced earthquake sequence from January 2011 to January 2012[J]. Journal of Rock Mechanics and Geotechnical Engineering, 2017, 9(5): 783-796.

[20] 宋博文, 徐亚东, 梁银平, 中国西部新生代沉积盆地演化[J]. 地球科学: 中国地质大学学报, 2014, 39(8): 1035-1051.

[21] 施炜, 马寅生, 吴满路, 青藏高原东北缘共和盆地第四纪磁性地层学研究[J]. 地质力学学报, 2006, 12(3): 317-323.

[22] YUAN D Y, GE W P, CHEN Z W, et al. The growth of northeastern Tibet and its relevance to large-scale continental geodynamics: a review of recent studies[J]. Tectonics, 2013, 32(5): 1358-1370.

[23] 唐显春, 王贵玲, 马岩, 青海共和盆地地热资源热源机制与聚热模式[J]. 地质学报, 2020, 94(7): 2052-2065.

[24] 雷治红. 青海共和盆地干热岩储层特征及压裂试验模型研究[D]. 长春: 吉林大学, 2020.

[25] 蔺文静, 甘浩男, 赵振, 青海共和盆地岩石圈热—流变结构及地热意义[J]. 地球学报, 2023, 44(1): 45-56.

[26] GAO J, ZHANG H J, ZHANG S Q, et al. Three-dimensional magnetotelluric imaging of the geothermal system beneath the Gonghe Basin, Northeast Tibetan Plateau[J]. Geothermics, 2018, 76: 15-25.

[27] 张森琦, 严维德, 黎敦朋, 青海省共和县恰卜恰干热岩体地热地质特征[J]. 中国地质, 2018, 45(6): 1087-1102.

[28] HAIMSON B. The effect of lithology, inhomogeneity, topography, and faults, on in situ stress measurements by hydraulic fracturing, and the importance of correct data interpretation and independent evidence in support of results[M]//XIE F R. Rock stress and earthquakes. Boca Raton: CRC Press, 2010: 11-14.

[29] MATSUKI K. Three-dimensional in-situ stress measurement with anelastic strain recovery of a rock core[C]// Society of geomechanics. 7th ISRM Congress. Aachen, Germany, 1991: 405-414.

[30] LIN W R, KWAŚNIEWSKI M, IMAMURA T, et al. Three-dimensional in-situ stress measurement with anelastic strain recovery of a rock core[C]// Society of geomechanics. 7th ISRM Congress. Aachen, Germany, 1991: 405-414.

[31] TEUFEL L. Determination of in situ stress from paretial anelastic strain recovery measurements of oriented cores from deep boreholes[C]// Society of geomechanics. U.S. Symposium on Rock Mechanics. Madison: University of Wisconsin, 1993: 176-188.

[32] ZHANG C Y, LIN W R, HE M C, et al. Determining in-situ stress state by anelastic strain recovery method beneath Xiamen: implications for the coastal region of Southeastern China[J]. Rock Mechanics and Rock Engineering, 2022, 55(9): 5687-5703.

[33] FUNATO A, ITO T. A new method of diametrical core deformation analysis for in situ stress measurements[J]. International Journal of Rock Mechanics and Mining Sciences, 2017, 91: 112-118.

[34] 杨跃辉, 孙东生, 郑秀华, 岩芯直径变形分析法及其在松科2井深部地应力调查中的应用[J]. 中南大学学报(自然科学版), 2019, 50(12): 3106-3113.

[35] 李彦恒, 谭可可, 冯利. 基于岩饼几何形态测量的原地应力测定方法[J]. 岩土力学, 2012, 33(增刊2): 224-228.

[36] LIM S S, MARTIN C D. Core disking and its relationship with stress magnitude for Lac du Bonnet granite[J]. International Journal of Rock Mechanics and Mining Sciences, 2010, 47(2): 254-264.

[37] 张丰收, 李猛利, 张重远, 高地应力下深部岩芯饼化裂缝发展规律及机制研究[J]. 岩石力学与工程学报, 2022, 41(3): 533-542.

[38] SUN D S, SONE H, LIN W R, et al. Stress state measured at -7 km depth in the Tarim Basin, NW China[J]. Scientific Reports, 2017, 7: 4503.

[39] BELL J S, GOUGH D I. Northeast-southwest compressive stress in Alberta evidence from oil wells[J]. Earth and Planetary Science Letters, 1979, 45(2): 475-482.

[40] PLUMB R A, HICKMAN S H. Stress-induced borehole elongation: a comparison between the four-arm dipmeter and the borehole televiewer in the Auburn Geothermal Well[J]. Journal of Geophysical Research: Solid Earth, 1985, 90(B7): 5513-5521.

[41] 都昌庭. 共和地震震源机制解特征[J]. 高原地震, 2001, 13(4): 1-5.

[42] 郝明, 沈正康, 王庆良. 1990年青海共和7.0级地震震后垂直形变研究[J]. 地震学报, 2010, 32(5): 557-569, 633.

[43] ZHANG C Y, LI B, LI H L, et al. Stress estimation in a 3 km-deep geothermal borehole: a snapshot of stress state in southern Cathaysia Block, China[J]. Tectonophysics, 2023, 864: 230031.

[44] 孟文, 郭祥云, 李永华, 青藏高原东北缘构造应力场及动力学特征[J]. 地球物理学报, 2022, 65(9): 3229-3251.

[45] ZIMMERMANN G, MOECK I, BLÖCHER G. Cyclic waterfrac stimulation to develop an Enhanced Geothermal System (EGS): conceptual design and experimental results[J]. Geothermics, 2010, 39(1): 59-69.

[46] LISLE R J, SRIVASTAVA D C. Test of the frictional reactivation theory for faults and validity of fault-slip analysis[J]. Geology, 2004, 32(7): 569-572.

[47] MOECK I, KWIATEK G, ZIMMERMANN G. Slip tendency analysis, fault reactivation potential and induced seismicity in a deep geothermal reservoir[J]. Journal of Structural Geology, 2009, 31(10): 1174-1182.

[48] SCHWAB D R, BIDGOLI T S, TAYLOR M H. Characterizing the potential for injection-induced fault reactivation through subsurface structural mapping and stress field analysis, Wellington field, sumner county, Kansas[J]. Journal of Geophysical Research: Solid Earth, 2017, 122(12): 10132-10154.

[49] MORRIS A, FERRILL D A, BRENT HENDERSON D B. Slip-tendency analysis and fault reactivation[J]. Geology, 1996, 24(3): 275-278.

[50] WORUM G, VAN WEES J D, BADA G, et al. Slip tendency analysis as a tool to constrain fault reactivation: a numerical approach applied to three-dimensional fault models in the Roer Valley Rift system (Southeast Netherlands)[J]. Journal of Geophysical Research: Solid Earth, 2004, 109(B2): 1-16.

[51] MASSIRONI M, BISTACCHI A, MENEGON L. Misoriented faults in exhumed metamorphic complexes: rule or exception?[J]. Earth and Planetary Science Letters, 2011, 307(1/2): 233-239.

[52] FERRILL D A, MORRIS A P, MCGINNIS R N, et al. Mechanical stratigraphy and normal faulting[J]. Journal of Structural Geology, 2017, 94: 275-302.

[53] HANTUSH M S. Drawdown around Wells of variable discharge[J]. Journal of Geophysical Research, 1964, 69(20): 4221-4235.

[54] BYERLEE J. Friction of rocks[J]. Pure and Applied Geophysics, 1978, 116(4): 615-626.

[55] ZHANG C Y, FAN D J, ELSWORTH D, et al. Mechanisms of stress- and fluid-pressure-driven fault reactivation in Gonghe granite: implications for injection-induced earthquakes[J]. International Journal of Rock Mechanics and Mining Sciences, 2024, 174: 105642.

[56] ZHANG C Y, HU Z J, ELSWORTH D, et al. Frictional stability of laumontite under hydrothermal conditions and implications for injection-induced seismicity in the Gonghe geothermal reservoir, Northwest China[J]. Geophysical Research Letters, 2024, 51(10): e2023GL108103.

基金资助

中国地质调查局地质调查项目(DD20190138)

国家自然科学基金项目(42177175)

中央级科研单位基本科研项目(DZLXJK202204)