干热岩开发断层滑动及诱发地震风险研究:以唐山马头营干热岩场地为例

Fault slip tendency and induced-earthquake risks in hot dry rock development: A case study in the Tangshan Matouying HDR site

摘要

诱发地震是干热岩大规模开发面临的主要问题之一,断层滑动和流体超压是较大地震的两大主要诱因,同时每个场地的诱发地震情况受具体的地质条件控制。本文选取马头营干热岩场地为研究区,在地震地质、构造地质、地热地质条件分析基础上,建立水热力多场耦合数值模型,综合考虑地应力场和多条断裂的影响,通过计算断层滑动趋势评估干热岩开发过程中的断层稳定性,并依据注入扰动体积估算长期注入可能诱发的最大震级。结果显示,10年干热岩开发过程中,马头营干热岩场地断层的孔隙压力变化不超过0.8 MPa,不会引起断层上有效应力的明显变化;注采期间断层滑动趋势均远低于滑动预警线,注采过程安全。优选出两种适用于马头营场地的地震矩计算模型,10年开采过程中可能诱发的最大震级为2.0级至2.1级。研究成果不仅能够为马头营场地干热岩安全规模化开发提供支撑,而且有助于深入理解工程活动影响下的地球动力学过程。

关键词

干热岩 / 诱发地震 / 断层滑动趋势 / 最大震级 / 多场耦合

Key words

hot dry rock / induced earthquakes / fault slip trend / maximum magnitude / multi-field coupling

0 引言

干热岩(HDR)资源是地热资源的一种,通常储存于地下3~10 km处的低孔隙度、低渗透率和高温岩体中^[1⇓-3]。干热岩资源储量丰富、分布广泛^[4⇓-6],占可利用地热能总量的90%以上^[7-8]。

诱发地震是目前干热岩大规模开发面临的主要问题之一^[9⇓-11]。例如巴塞尔干热岩场地注水在2006年11月至2007年3月期间发生了一系列M_w<3.4的诱发地震,最终导致该项目于2009年12月取消^[12]。韩国Pohang地热田于2017年11月发生了M_w 5.5地震,是迄今为止监测到的最大诱发地震^[13]。此外,美国的Geysers、澳大利亚Cooper Basin、意大利Larderello-Travale、法国Soultz-sous-Forets、德国Landau等场地均因水力压裂或后期注水开发引起地震^[14-15]。诱发地震引起了当地群众乃至国际媒体的广泛关注,较大的地震还会导致房屋受到损害,制约了干热岩的进一步开发利用。

较大的震级往往是由两个因素造成的^[16]。其一,断层滑动,根据莫耳—库伦滑动准则,孔隙压力增加会导致断层面上的有效正应力降低,断层面摩擦力降低,当达到临界值时断层发生滑动(断层激活);其二,流体超压,即循环注入过程中诱发地震的震级与流体超压占据的体积相关。许多学者针对干热岩开发诱发地震问题进行了研究,Feng等^[17]采用滑动趋势方法研究了唐山马头营干热岩场地注入井附件断层滑动概率。在注入井15 km范围内,大多数断层的断层滑动概率随流体的注入持续增加。Yue等^[16]综合滑动趋势和流体超压,评价了共和盆地干热岩场地循环过程可能诱发的地震风险。王宁等^[18]研究了马头干热岩场地注水期间的地震活动与注入压力的线性关系,并指出当注水量达到11 821 m³,压力达到23.46 MPa时,存在断层滑动风险。然而,多断层影响下长期循环过程中的断层滑动和诱发地震风险需要进一步研究。

本文选取马头营干热岩场地为研究区,在地热地质条件分析基础上,考虑多条断裂影响,采用断层滑动趋势评估了场地干热岩开发过程中断层稳定性,同时依据超压体积评价了可能诱发的最大震级。研究成果不仅有助于深入理解开采条件下的深部热储动态响应规律,同时为干热岩安全规模化开发提供支撑。

1 马头营干热岩场地地质背景

1.1 区域构造与历史地震

马头营干热岩场地由河北省煤田地质局主持建设,地质构造属于燕山褶皱带中的开滦凹陷褶皱体系的一部分,基底为古生界地层。该地区存在多种断裂系统(图1),包括北北东向(如卢东断裂)、北东东向(如唐山断裂带、野鸡坨断裂、昌黎—宁河断裂、仓东断裂和西南庄断裂)、北西向(蓟运河断裂、滦县—乐亭断裂和柏各庄断裂)、北西西向(冷口断裂)和近东西向(汉沽断裂和海河断裂)。北东向断层是该地区最显著的断层,贯穿整个唐山地区^[17,19-20]。

唐山是华北地区地震活动最活跃的地区之一^[21-22]。唐山7.8级地震是过去100年世界范围内最具破坏性的地震之一。虽然过去40年唐山地区的地震活动有所减少,但是仍有M_L 4.0~5.0级地震发生^[23]。

马头营干热岩场地位于唐山老震区东南方向,黄骅坳陷北部,渤海湾盆地北西向构造和北东向构造交会处。马头营凸起北以马北断层为界与乐亭凹陷相连,西以柏各庄断裂为界与南堡凹陷相接,南以红房子断层为界与石臼坨凹陷相邻,多条断裂纵横交会分布^[18,20,24]。场地附近断裂主要包括柏各庄断层、马北断层、红房子断层和高柳断层。

1.2 区域现今构造应力场

地应力大小和方向对干热岩工程建设至关重要,直接影响断层稳定性^{[16,18,24-25]}。上官拴通等^[24]利用非弹性应变恢复(ASR)地应力测试方法,获取了马头营干热岩场地3~4 km深度范围的地应力状态。3 139~3 934 m水平最小主应力介于 59.0~90.7 MPa,水平最大主应力介于103.7~123.6 MPa。三向主应力总体表现为σ_H>σ_v>σ_h。

2 断层滑动趋势

2.1 评价方法

考虑干热岩场地断层分布,建立多场耦合数值模型,计算断层上的剪应力和有效法向应力,进而计算断层滑动趋势。

断层滑动趋势(T_s)是断层表面上的剪应力与法向应力之比^[26],其计算公式为

T_s= \frac{\tau }{{\sigma }_{n,eff}}≥μ_s

式中:τ为断层面上的剪应力;σ_n,eff为断层面有效法向应力(法向应力,σ_n,减去孔隙压力,P_f);μ_s为断层滑动摩擦系数^[27]。如果T_s等于或超过摩擦滑动阻力μ_s,则断层发生滑动。一般认为μ_s=0.6^[28]。

2.2 概念模型及方案

研究区内主要断层以高角度为主,二维模型能够描述注采条件下断层所处的应力状态,同时考虑到三维模型会消耗大量的计算资源,因此本次研究采用二维模型(图2)。沉积岩中的注入会引起更大、更远的地震,而硅质地层中注入的影响范围较小,一般不会超过15 km^[17]。因此,本文以M1井为中心,15 km为范围建立用于数值模拟的二维结构模型,模型包括注入井M1和开采井M2,以及周围主要断层,如柏各庄断层、高柳断层、马北断层和红房子断层。模型垂向范围为4 300~4 500 m。

据2022年现场循环试验,M1注入量为60 m³/h (取实际的均值),井口压力为4~27 MPa。开采井M2井口压力为0.04 MPa。按照上述注入方案,模拟10年,评价断层滑动趋势。

2.3 数学方程

本次研究基于COMSOL软件平台,使用“多孔弹性”多物理场节点,包含“达西定律”和“固体力学”接口,并在“固体力学”接口中添加流体压力梯度作为应力贡献。各场的方程如下。

2.3.1 达西定律方程

根据Darcy定律,穿过多孔表面的净通量为

u=- \frac{k}{\mu }(∇p+ρg∇D)

式中:u为达西速度或流量矢量,m/s;k为多孔介质的渗透率,m²;μ为流体的动态黏度,Pa·s;p为流体的压力,Pa;ρ为流体密度,kg/m³;g为重力加速度的大小,m/s²;D是重力作用方向上的单位向量矩阵。

当考虑岩石和流体压缩性时,连续性方程可表示为

ρS \frac{\partial p}{\partial t}+∇·ρ \left[-\frac{k}{\mu }(\nabla p+\rho g\nabla D)\right]=Q_m

式中S为储水系数,是多孔材料和孔隙中流体的加权压缩系数,单位为Pa^-1。

2.3.2 固体力学方程

流固耦合多为弹性多孔介质中流体流动与变形之间的相互作用,其表达式为

σ=Cε_σ-α_Bp_fE

式中:σ是Cauchy应力张量;ε_σ是应变张量;α_B是Biot系数^[29];p_f是流体孔隙压力,Pa;C为弹性矩阵,必须在“排水”条件下通过测量恒定孔隙压力下的应力变化引起的应变来确定;E 为单位矩阵。

2.3.3 温度场方程

在固相和液相局部热平衡假设条件下,多孔介质的能量守恒方程为

(ρC_p)_eff \frac{\partial T}{\partial t}+ρ_fC_p,fu·∇T+∇·q=Q

q=-k_eff∇T

式中:ρ_f是流体密度;C_p,f是流体比热容;(ρC_p)_eff是恒压有效体积热容,定义为

(ρC_p)_eff=θ_sρ_sC_p,s+ε_pρ_fC_p,f

式中:ε_p为孔隙度;θ_s为固体基质体积分数;ρ_s是固体基质密度;C_p,s是固体基质恒压热容;k_eff是有效热导率(一个标量,如果热导率是各向异性的,则是一个张量);q是传导热量;u是速度场,可以是解析表达式,也可以从流体流动计算出来;Q是热源。

固体—流体系统的有效热导率k_eff与固体的热导率k_s和流体的热导率k_f有关。如果热传导在固体和液体中同时发生,那么有效导热系数就是导热系数k_f和k_s的加权平均值:

k_eff=θ_sk_s+ε_pk_f

温度对应变的影响表示为

ε_th=α(T-T_ref)

式中:α是热膨胀系数,与温度有关;T为当前温度;T_ref是参考温度。

2.4 初始及边界条件及参数

初始孔隙压力根据静水压力计算,为42.6 MPa。边界条件设定为恒定压力边界,与初始压力相同。

初始温度根据M1井测温资料插值获得,4 500 m为157.12 ℃。温度边界条件设定为恒定温度边界,与初始温度相同。

初始应力条件根据前述地应力大小和方向确定,最大水平主应力为110 MPa,方向为281°;最小水平主应力68 MPa。将模型根据水平最大主应力方向顺时针旋转11°,以保证边界应力与水平应力方向一致。

据测井结果,热储层天然平均渗透率为0.5 mD,经过改造,在注采井范围内形成了改造区,平均渗透率增大为3.5 D。热储岩石平均热导率为2.93 W/(m·K),平均比热容为1.53 J/(K·g),高温高压下的岩石弹性模量约为300 GPa。根据测井结果,岩石泊松比为0.26,根据测试和前人研究结果显示,Biot系数为0.23^[30]。

2.5 评价结果

2.5.1 断层孔隙压力变化

随着流体注入,孔隙压力逐渐升高,并且影响范围逐年扩大(图3)。可以看到在3年后,柏各庄断层孔隙压力开始出现升高,最大升高了0.07 MPa;截至10年后,柏各庄断层孔隙压力升高值最大可达0.8 MPa。马北断层同样受到孔隙压力的影响,在约6年后,可以看到马北断层的孔隙压力升高现象。在10年时,马北断层的孔隙压力升高不超过0.4 MPa。总体来看,注采循环不会引起断层上孔隙压力明显升高。

2.5.2 断层滑动趋势

从马头营干热岩场地主要断层的滑动趋势(图4)可以看出,在天然条件下,马北断层最为稳定,滑动趋势普遍低于0.1。其次为红房子断层,滑动趋势低于0.2。柏各庄断层形态较为复杂,不同断层段滑动趋势差异明显。在与马北断层交汇处,柏各庄断层的滑动趋势较大,约为0.4。在其他断层段滑动趋势较低,小于0.3。造成滑动趋势差异的原因是,断层交汇处应力场比较复杂,会出现应力集中的现象,因此,在断层交汇处某些特定走向的断层滑动趋势会升高。

从时间上来看,10年后的断层滑动趋势与初始状态几乎相同。也就是说,注采循环不会引起断层滑动趋势明显增加。其原因是注采循环影响范围较小,断层上的孔隙压力没有出现明显增高,正如前面模拟结果显示。因此,注采不会引起断层的明显激活,注采过程安全。

3 最大震级预测

3.1 预测方法

3.1.1 方法一

McGarr^[31]统计了非常规储层的地震数据,显示最大地震诱发地震矩与岩石刚度和流体注入扰动范围内的岩石体积有关(M₀=GΔV)。扰动体积(ΔV)在岩石中的变化与岩石孔隙压力升高(+ΔP)的特征一致,与流速和渗透率有关。本次研究使用McGarr模型,将受扰岩石体积大小与诱发地震活动的最大震级联系起来。通过特定的阈值压力计算可能的最大地震矩^[31]:

M_max=GΔV

式中:G为刚度模量,MPa。马头营干热岩场地平均杨氏模量E等于24 GPa,泊松比等于0.26。因此,刚度模量G(G=E/2(1+ν))等于9.52 GPa。

3.1.2 方法二

Galis等^[32]提出了一个基于物理学的定量模型,以描述地震物理学。诱发地震的最大震级( M_0^{\mathrm{m}\mathrm{a}\mathrm{x}})和净注入量(ΔV)之间的关系可以表示为

M_0^{\mathrm{m}\mathrm{a}\mathrm{x}}=γΔV^3/2

γ= \frac{0.425\mathrm{ }5}{\sqrt[ ]{\Delta {\tau }_0}} {\left(\frac{K\mu }{T_{\mathrm{R}}}\right)}^{3/2}

式中:Δτ₀是背景压力降,MPa,这个值通常在0.1~10 MPa^[33];K是储层的体积模量,MPa;μ是摩擦系数;T_R是储层厚度,m。体积模量K(K=E/3(1-2μ))等于16.67 GPa,唐山地区的背景应力降Δτ₀约为6 MPa^[34],摩擦系数μ为0.6^[35],储层厚度h等于200 m(4 300~4 500 m)。

上述方法中,计算最小阈值压力增量(ΔP_crit)所包含的岩石体积,这是能够引起地震的最小压力变化。根据王宁等^[18]的研究,当注入压力达到23.46 MPa时才会发生地震,因此可以认为ΔP_crit=23.46 MPa。

最后,用标准公式把地震矩换算成震级^[36]:

M_w= \frac{\mathrm{l}\mathrm{g}\mathrm{ }{M}_0-9.05}{1.5}

3.1.3 方法三

Van der Elst等^[37]根据已有地震数据统计得出了最大震级M_max的另一种评价方法:

M_max= \frac{1}{b}(Σ+lgΔV)

式中:b是斜率;Σ是孕震指数,并且在地热储层的生命周期内恒定不变^[38]。唐山地区,b的值取为0.87^[39]。由于缺乏马头营场地的Σ,使用澳大利亚Cooper盆地地热系统的孕震指数Σ=-1.2^[38]。

利用上述3种方法,评价了马头营干热岩场地注入情况下诱发地震事件可能的最大震级。

3.2 预测结果

3.2.1 扰动岩石体积(ΔV)

注采循环开始后注入井附近孔隙压力逐渐升高(图5)。1年后孔隙压力最大升高约18 MPa,2年后约25 MPa。在6年后趋于稳定,最高升高至35 MPa,之后最大压力基本保持不变。孔隙压力扰动体积也逐年扩大。

3.2.2 最大震级预测

图6显示了3种方法预测的最大震级分布。McGarr模型预测的10年内最大震级由1.0级升高至2.5级。Galis模型预测的震级由-0.25级升高至2.11级。Van der Elst模型预测的最大震级由-0.7级升高至2级。

总体来看,McGarr模型预测的最大震级较大,Galis和Van der Elst模型预测的最大震级略小。McGarr模型是根据全球诱发地震数据统计得出,对于特定的场地可能预测不准。Galis和Van der Elst模型综合考虑了场地地质特征,对场地条件的描述更加准确,因此认为Galis和Van der Elst模型预测的结果更加准确。总之,马头营干热岩场地按照当前的注入方案(50.7 m³/h),连续循环1年,可能引起的最大震级低于0级;连续循环10年,可能引起的最大震级为2级至2.1级。

以上评价结果是针对当前的储层条件和注入方案进行的,可以采取一定措施,如提高储层渗透性能等来降低可能诱发的最大震级;同时在注采循环期间加强井筒压力和微震等监测,及时改变注采方案,可有效避免产生较大地震。

4 结论

本文选取马头营干热岩场地为研究区,采用断层滑动趋势评估了场地干热岩开发过程中断层稳定性,同时依据扰动体积估算了注入可能诱发的最大震级,为干热岩规模化安全开发提供支撑。

预测了注采循环引起的孔隙压力变化。10年内马头营干热岩场地断层的孔隙压力变化不超过0.8 MPa。与天然地应力相比,干热岩的注采过程不会引起断层上有效应力的明显变化。

明确了天然条件和注采循环条件下的断层滑动趋势。天然条件下,柏各庄断层与马北断层交汇处的滑动趋势最大,达到0.4,其余断层滑动趋势均低于0.3。10年的注采循环不会引起断层滑动趋势的明显升高,注采过程安全。

定量预测了可能诱发的最大地震震级。优选出Galis和Van der Elst模型,对长期注入可能诱发的最大震级进行预测,10年的连续循环可能诱发的最大震级为2级至2.1级。

建议在干热岩开发过程中采取提高储层渗透性、优化注入方案等措施,以降低诱发地震的风险。建议加强井筒压力和微震活动的监测,以便在监测到异常信号时及时调整注采方案,确保开发活动的安全性。

需要注意的是,在风险评价过程中,采用的参数均基于目前对场地地质条件的认识。然而,有些参数可能随着注入的进行发生改变,如储层的力学参数和断层滑动摩擦系数等,本次研究并未全面考虑这些参数的变化,因为这往往很复杂,目前的研究不足以支撑这些描述。在未来的工作中,会进一步细化模型、深入研究参数的影响因素,使得预测结果更加准确。

感谢中国地震台网中心国家地震科学数据中心(http://data.earthquake.cn)提供数据支撑。

参考文献

[1] Panel M L. The future of geothermal energy: impact of enhanced geothermal systems (EGS) on the United States in the 21st century[R]. Cambridge, USA: Massachusetts Institute of Technology, 2006.

[2] 许天福, 胡子旭, 李胜涛, 增强型地热系统: 国际研究进展与我国研究现状[J]. 地质学报, 2018, 92(9): 1936-1947.

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

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

[5] 王贵玲, 刘峰, 蔺文静, 我国陆区地壳生热率分布与壳幔热流特征研究[J]. 地球物理学报, 2023, 66(12): 5041-5056.

[6] 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.

[7] LEI Z H, ZHANG Y J, YU Z W, et al. Exploratory research into the enhanced geothermal system power generation project: the Qiabuqia geothermal field, Northwest China[J]. Renewable Energy, 2019, 139: 52-70.

[8] ZHANG C, JIANG G Z, JIA X F, et al. Parametric study of the production performance of an enhanced geothermal system: a case study at the Qiabuqia geothermal area, northeast Tibetan Plateau[J]. Renewable Energy, 2019, 132: 959-978.

[9] 程钰翔. EGS诱发地震特征及风险评价研究[D]. 长春: 吉林大学, 2021.

[10] RATHNAWEERA T D, WU W, JI Y L, et al. Understanding injection-induced seismicity in enhanced geothermal systems: from the coupled thermo-hydro-mechanical-chemical process to anthropogenic earthquake prediction[J]. Earth-Science Reviews, 2020, 205: 103182.

[11] 甘浩男, 王贵玲, 蔺文静, 增强型地热系统环境地质影响现状研究与对策建议[J]. 地质力学学报, 2020, 26(2): 211-220.

[12] KRAFT T, DEICHMANN N. High-precision relocation and focal mechanism of the injection-induced seismicity at the Basel EGS[J]. Geothermics, 2014, 52: 59-73.

[13] GRIGOLI F, CESCA S, RINALDI A P, et al. The November 2017 M_w 5.5 Pohang earthquake: a possible case of induced seismicity in South Korea[J]. Science, 2018, 360(6392): 1003-1006.

[14] BREEDE K, DZEBISASHVILI K, LIU X L, et al. A systematic review of enhanced (or engineered) geothermal systems: past, present and future[J]. Geothermal Energy, 2013, 1(1): 4.

[15] ZANG A, OYE V, JOUSSET P, et al. Analysis of induced seismicity in geothermal reservoirs-an overview[J]. Geothermics, 2014, 52: 6-21.

[16] YUE G F, LI X Y, ZHANG W. Risk assessment of earthquakes induced during hdr development: a case study in the Gonghe Basin, Qinghai Province, China[J]. Geothermics, 2023, 111: 102721.

[17] FENG C J, GAO G L, ZHANG S H, et al. Fault slip potential induced by fluid injection in the Matouying EGS field, Tangshan seismic region, North China[J]. Natural Hazards and Earth System Sciences, 2022, 22(7): 2257-2287.

[18] 王宁, 王亚玲, 张晓刚, 马头营干热岩开采试验场地人工注水诱发地震探讨[J]. 大地测量与地球动力学, 2024: 1-8.

[19] 上官拴通. 马头营区干热岩地热资源赋存分布特征及开发利用前景[J]. 能源与环保, 2017, 39(5): 155-159, 165.

[20] 齐晓飞, 上官拴通, 张国斌, 河北省乐亭县马头营区干热岩资源孔位选址及开发前景分析[J]. 地学前缘, 2020, 27(1): 94-102.

[21] 江娃利. 有关1976年唐山地震发震断层的讨论[J]. 地震地质, 2006, 28(2): 312-318.

[22] 丰成君, 戚帮申, 王晓山, 基于原地应力实测数据探讨华北典型强震区断裂活动危险性及其对雄安新区的影响[J]. 地学前缘, 2019, 26(4): 170-190.

[23] 张素欣, 王晓山, 陈婷, 唐山老震区40年地震时空演化特征分析[J]. 华北地震科学, 2017, 35(1): 32-37.

[24] 上官拴通, 孙东生, 张国斌, 唐山地区3-4 km深部地应力测量及断层稳定性分析[J]. 地质学报, 2021, 95(12): 3915-3925.

[25] 王丹, 文冬光, 杨用彪, 干热岩开发循环试验的研究进展和发展建议[J]. 地质科技通报, 2024: 1-15.

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

[27] BYERLEE J. Friction of rocks[M]. Basel: Birkhäuser Basel, 1978.

[28] ZOBACK M D. Reservoir geomechanics[M]. Cambridge: Cambridge University Press, 2007.

[29] BIOT M A. General theory of three-dimensional consolidation[J]. Journal of Applied Physics, 1941, 12(2): 155-164.

[30] LI Q, AGUILERA R, CINCO-LEY H. A correlation for estimating the Biot coefficient[J]. SPE Drilling and Completion, 2020, 35(2): 151-163.

[31] MCGARR A. Maximum magnitude earthquakes induced by fluid injection[J]. Journal of Geophysical Research: Solid Earth, 2014, 119(2): 1008-1019.

[32] GALIS M, AMPUERO J P, MAI P M, et al. Induced seismicity provides insight into why earthquake ruptures stop[J]. Science Advances, 2017, 3(12): eaap7528.

[33] ABERCROMBIE R. Earthquake source scaling relationships from -1 to 5 ML using seismograms recorded at 2.5 km depth[J]. Journal of Geophysical Research, 1995, 100(24): 15-24,36.

[34] XIAO W W, HUANG Q, WANG R J. Distribution and fractal characteristics of stress drop in different earthquake regions[J]. Acta Seismology Sinica, 1992, 14: 281-288.

[35] TROISE C, DE NATALE G, PINGUE F, et al. Evidence for static stress interaction among earthquakes in the south-central Apennines (Italy)[J]. Geophysical Journal International, 1998, 134(3): 809-817.

[36] HANKS T C, KANAMORI H. A moment magnitude scale[J]. Journal of Geophysical Research: Solid Earth, 1979, 84(B5): 2348-2350.

[37] VAN DER ELST N J, PAGE M T, WEISER D A, et al. Induced earthquake magnitudes are as large as (statistically) expected[J]. Journal of Geophysical Research: Solid Earth, 2016, 121(6): 4575-4590.

[38] SHAPIRO S A, DINSKE C, LANGENBRUCH C, et al. Seismogenic index and magnitude probability of earthquakes induced during reservoir fluid stimulations[J]. The Leading Edge, 2010, 29(3): 304-309.

[39] 杜航, 王俊, 朱音杰, 2020年河北唐山M_S 5.1地震b值时空变化特征[J]. 华北地震科学, 2021, 39(3): 92-98.

基金资助

河北省唐山市马头营凸起区干热岩开发关键技术集成研究与示范项目(20374102D)

中国地质科学院基本科研业务费专项(YK202309)