循环应力对储气库岩石微观结构及力学性质影响数值模拟

隋义勇, 林堂茂, 刘 翔, 董长银, 程 威, 张广明, 廖正毅

(1.中国石油大学(华东)石油工程学院,山东青岛 266580; 2.中国石油勘探开发研究院,北京 100083)

中国越来越多的枯竭型油气藏储气库投入建设和运行,而对储气库建设仍处于初级技术阶段[1-4]。针对循环应力对储气库岩石微观结构和力学性质的影响,马新华等[5]通过室内岩心三轴加卸载试验研究岩石变形破坏特征,分析了累积塑性变形的变化规律。郑得文等[6]通过室内试验研究储气库生产过程中的剪切破坏、拉伸破坏及断层剪切滑移失稳风险等问题。孙岩等[7]根据流体物性参数和相应生产动态数据,优化了多周期注采水侵量计算模型。Wang等[8]利用有限元法研究储气库循环注采过程中的疲劳损伤、弹塑性变形及孔隙度变化。Guo等[9]研究循环应力作用下黏土弹性模量变化规律及永久累积变形特征。Wang等[10]研究循环注采过程中产生的微地震与应力历史对岩石内部损伤及破坏的影响。Sun等[11]建立了一种储气库综合地质力学模型,评估储气库长期运行过程中盖层完整性破坏和故障泄漏风险。陈金平等[12]使用有限元法研究储气库地应力场和内压作用下,不同长短轴比椭球形储气库腔体极限压力。Dai等[13]使用离散元法研究花岗岩在加载卸载过程中应变和损伤演变情况。前人研究主要集中在宏观层次,缺少对微观结构变化和机制的研究,岩石微观结构变化是影响岩石宏观力学性质变化的根源[14-16]。笔者采用颗粒离散元数值模拟方法研究循环应力对储气库岩石微观结构及力学性质的影响。

1 真三轴颗粒离散元模型

1.1 建立离散元模型

离散元法将岩体看作是各离散单元的黏结集合体,假设两离散单元之间存在若干弹簧连接,通过两者重叠量以及弹簧刚度确定接触力,离散单元的运动规律符合宏观牛顿定律,离散单元间碰撞符合胡克定律以及弹塑性力学规律[17-21],离散元法将非连续介质与连续介质力学问题有机结合起来。

首先,根据国际岩石力学学会对真三轴岩心尺寸要求,确定数值模型尺寸为2.5 cm×2.5 cm×10 cm,再根据岩心室内试验测得岩心粒度组成及力学性质,数值模型粒度组成既要确保颗粒体系能表征岩心的力学特性,又要确保模型中颗粒数量适中提高计算效率,最终确定模型颗粒半径为0.14～3.62 mm,生成了一个包含27 000多颗粒的纯颗粒模型;然后,考虑到枯竭型油气藏储气库建库之前经过多年的油气开采,储层岩石中存在大量微裂缝,通过PFC3D中缝网生成功能建立随机原生裂缝模型;最后,将纯颗粒模型和裂缝模型叠加,得到离散元真三轴岩心数值模型,如图1所示。

图1 岩心离散元模型建立过程Fig.1 Modeling process of core discrete element model

1.2 模型接触参数及应力加载参数

根据岩体的非连续特性及颗粒间力学接触特性,选用PFC3D中自带的线性平行黏结接触模型模拟岩体颗粒间接触摩擦特性与接触黏结特性。选用PFC3D中自带的光滑节理接触模型模拟岩体原生裂缝内颗粒间的接触特性。根据岩心的室内力学试验测量结果标定模型接触参数,最终确定的模型接触参数如表1所示。室内试验曲线与数值模拟曲线对比如图2所示。由图2可以看出,数值模拟获得的应力-应变曲线与岩心室内试验测得的应力-应变曲线吻合度较好,建立的岩心真三轴离散元模型以及标定的模型参数能够表征储气库储层岩石力学特性,可以用来进行后续的多周期等幅值轴向循环应力加载数值模拟研究。

表1 岩心颗粒接触模型参数

根据储气库储层岩石的原位地应力状态确定数值模型的初始应力状态,已知新疆某储气库储层3个方向的原位地应力梯度,结合储层埋深,确定储气库储层3个方向的原位地应力状态。利用伺服原理控制侧向无摩擦刚性墙体对数值模型施加水平方向的围压,其中水平最小主应力方向围压为30 MPa,水平最大主应力方向围压为48 MPa。数值模型上、下墙体为加载板,通过控制上、下墙体缓慢匀速往复运动对数值模型施加轴向循环加载应力,结合储气库生产上、下限库内压力和储层上覆地层压力,确定数值模型轴向循环加载应力上限为65.23 MPa,轴向循环加载应力下限为49.23 MPa,储气库通常要求安全平稳生产50～60 a,因此决定对数值模型进行等应力幅值轴向循环加载60个周期,研究循环应力加载对储气库岩石微观结构和力学性质的影响。循环应力加载示意图如图3所示。

图2 室内试验曲线与数值模拟曲线对比Fig.2 Comparison of laboratory experimental curves and numerical simulation curves

图3 循环应力加载示意图Fig.3 Cyclic stress loading diagram

2 循环应力对储层岩石内部结构的影响

2.1 循环应力对模型内部黏结状态的影响

不同加载周期模型内生成的无黏结接触分布如图4所示,红色部分显示的是循环应力加载1、30、60个周期时模型内颗粒间黏结接触断裂后生成的无黏结接触分布情况。由图4可以看出,随着循环应力加载周期增加,模型内颗粒间生成的无黏结接触数量逐渐增加,对比原始裂缝分布与无黏结接触分布位置,发现循环应力加载过程中模型内颗粒间生成的无黏结接触主要集中在原始裂缝分布密集的区域。

多周期循环应力加载对无黏结接触数量影响如图5所示。由图5可以看出,循环应力加载过程中模型内颗粒间生成的无黏结接触数量增加速率呈递减趋势,循环应力加载前30个周期模型内颗粒间无黏结接触数量增加速度较快,循环应力加载后30个周期模型内颗粒间无黏结接触数量增加速度减缓。

图4 不同加载周期模型内生成的无黏结接触分布Fig.4 Non-bonded contacts distribution in model with different loading cycles

图5 多周期循环应力加载对无黏结接触数量影响Fig.5 Effect of cyclic stress loading on number of non-bonded contacts

循环应力加载对模型内颗粒间黏结状态的影响与模型内颗粒间接触强度分布和裂缝分布有关。初始时模型内颗粒间黏结接触强度服从Weibull分布,模型内各区域颗粒间黏结接触强度分布不均匀,循环应力加载时应力上限虽未达到模型破坏强度,但已经超过模型内颗粒间弱黏结接触的破坏强度,导致模型内颗粒间黏结接触发生断裂。模型内颗粒在加载外力、不平衡力和颗粒间接触力的共同作用下发生分离运移并产生新接触,力和位移共同作用导致模型内颗粒间强黏结接触逐渐弱化并发生断裂,因此随着循环应力加载周期增加,模型内颗粒间产生的无黏结接触数量逐渐增加。另外,模型内原始裂缝分布密集区域存在软弱界面,在加载外力作用下应力会在界面处集中导致颗粒沿裂缝面发生滑移,因此模型内颗粒间产生的无黏结接触集中在原始裂缝分布密集区域。循环应力加载前30个周期模型内颗粒的位移量较大,颗粒间黏结接触断裂产生的无黏结接触数量较多。循环应力加载后30个周期模型内颗粒逐渐被压密而分布较均匀,模型内部应力分布逐渐趋向均匀,模型内颗粒位移量减小,颗粒间黏结接触断裂产生的无黏结接触数量减少,因此随循环应力加载周期增加模型内颗粒间产生的无黏结接触数量增加速率递减。

2.2 循环应力对模型内部微裂缝发育的影响

图6为不同应力加载周期时模型内生成的微裂缝分布与原始裂缝分布对比,其中红色部分代表循环应力加载后模型内生成的微裂缝分布,绿色部分代表原始裂缝分布。由图6可以看出,随着循环应力加载周期增加,模型内生成的微裂缝数量和微裂缝分布范围逐渐增加,且新生成的微裂缝分布受原始裂缝分布的控制明显,新生成的微裂缝主要分布在原始裂缝分布密集区域,在原始裂缝的边缘和尖端萌生后逐渐延展并相互贯穿连通。

图6 不同加载周期模型内微裂缝分布与原始裂缝分布对比Fig.6 Comparison of micro-cracks distribution with original crack distribution with different loading cycles

循环应力加载对模型内微裂缝发育的影响如图7所示。由图7可以看出,随着循环应力加载周期增加,模型内生成的微裂缝数量增加速率呈递减趋势,循环应力加载前30个周期模型内生成的微裂缝数量增加速度较快,循环应力加载后30个周期模型内生成的微裂缝数量增加速度减缓。

图7 循环应力加载对模型内微裂缝发育的影响Fig.7 Effect of cyclic stress loading on development of micro-cracks in the model

循环应力加载过程中,模型内原始裂缝在张应力的作用下克服模型内颗粒间摩擦力往复开合,导致原始裂缝分布区域及其周边区域的颗粒承受往复的张应力以及剪应力作用,因此原始裂缝分布区域及其周边区域更容易产生微裂缝。同时在加载外力、不平衡力和颗粒间接触力的共同作用下,原始裂缝周围颗粒易沿原始裂缝面产生较大位移量的滑移。这有利于微裂缝的延展和贯穿连通,因此模型内生成的微裂缝主要集中在原始裂缝分布区域及其周边区域。循环应力加载前30个周期,模型内新生成的微裂缝主要受原始裂缝控制,因此生成的微裂缝数量增加速率较快。循环应力加载后30个周期,模型内颗粒体系重新分布导致模型内应力趋向均衡,模型内新生成的微裂缝受原始裂缝控制减弱,模型内生成的微裂缝数量增加速率减缓,但模型内生成的微裂缝分布更加广泛。

2.3 循环应力对模型孔隙度的影响

图8为模型归一化孔隙度φn/φ1(φn和φ1分别为第n和第1循环应力加载周期模型孔隙度，%)随循环应力加载周期增加的变化。由图8可以看出,随着循环应力加载周期增加,轴向加载应力处于循环应力下限时模型归一化孔隙度逐渐降低,且模型归一化孔隙度降低速率呈递减趋势,循环应力加载前30个周期模型归一化孔隙度降低速率较快,循环应力加载后30个周期模型归一化孔隙度降低速率减缓。其中前10个循环应力加载周期内归一化孔隙度降低量约占整个循环应力加载过程中归一化孔隙度总降低量的30%。

图8 循环应力加载对模型孔隙度的影响Fig.8 Effect of cyclic stress loading on model porosity

循环应力加载过程中模型的变形主要有可恢复的弹性变形和不可逆的塑性变形,模型孔隙度降低主要是不可逆塑性变形引起的。随着循环应力加载周期增加,模型内生成的无黏结接触数量和微裂缝数量逐渐增加,模型内产生大量游离颗粒。随着循环应力加载的进行模型内游离颗粒运移至大孔洞和裂缝中,导致模型中大孔洞逐渐被填充,裂缝逐渐被压密甚至闭合,模型在相同应力水平下颗粒分布更加均匀,模型的孔隙度逐渐减小。随着循环应力加载周期增加,模型内生成的无黏结接触数量和微裂缝数量增加速率逐渐减缓,因此模型内不可逆塑性变形量逐渐减小,模型孔隙度降低量也随之减小,模型孔隙度降低速率逐渐减缓。

3 循环应力对储层岩石力学性质的影响

3.1 循环应力对储层岩石应变的影响

图9为模型轴向和侧向应变随着环应力加载周期的变化。由图9可以看出,随着循环应力加载周期增加,在轴向加载应力达到循环应力上限时,模型轴向应变与最小主应力方向侧向应变均逐渐增大,但两者应变增加速率差异明显,最小主应力方向侧向应变增加速率约是轴向应变增加速率的两倍。这表明在多周期等幅值轴向循环应力加载的疲劳损伤过程中,模型侧向膨胀的影响远大于轴向压缩的影响,模型侧向膨胀占主导作用。

图9 循环应力加载对轴向应变和侧向应变影响Fig.9 Effects of cycle stress loading on axial strains and lateral strains

3.2 循环应力对储层岩石弹性模量的影响

图10为模型弹性模量随循环应力加载周期的变化。

图10 循环应力加载对弹性模量影响Fig.10 Effect of cyclic stress loading on elastic modulus

由图10可以看出,随着循环应力加载周期增加,模型弹性模量逐渐减小,模型弹性模量降低速率呈递减趋势,循环应力加载前30个周期模型弹性模量降低速率较快,后30个周期模型弹性模量降低速率减缓,循环应力加载60个周期模型弹性模量较原始弹性模量降低0.33 GPa,降幅为22%。

3.3 循环应力对储层岩石泊松比的影响

图11为模型泊松比随循环应力加载周期的变化。由图11可以看出,随着循环应力加载周期增加,模型泊松比逐渐增大,模型泊松比增加速率呈递减趋势。循环应力加载前30个周期模型泊松比增加速率较快,后30个周期模型泊松比增加速率减缓,循环应力加载60个周期模型泊松比较原始泊松比增加0.157,增幅为78%。

图11 循环应力加载对泊松比影响Fig.11 Effect of cycle stress loading on Poissons ratio

3.4 循环应力对储层岩石损伤量的影响

图12为模型损伤量随循环应力加载周期的变化。由图12可以看出,随着循环应力加载周期增加,模型损伤量逐渐增大,则模型抵抗变形和破坏的能力降低,模型损伤量增加速率呈递减趋势。循环应力加载前30个周期模型损伤量增加速率较快,后30个周期模型损伤量增加速率减缓,循环应力加载60个周期模型损伤量增加22%。

图12 循环应力加载对模型损伤量影响Fig.12 Effect of cycle stress loading on model damage amount

循环应力加载过程中模型内新生成的大量无黏结接触无法抵抗张应力,只能抵抗压应力,同时模型内新生成的微裂缝数量和模型损伤量逐渐增加,模型整体强度降低,这导致模型抵抗变形的能力降低。因此随着循环应力加载周期增加,相同应力水平下模型侧向应变及轴向应变量均逐渐增加,模型弹性模量逐渐降低。根据岩石力学理论可知,岩石抵抗张应力的能力远低于抵抗压应力的能力,而且模型内新生成的无法抵抗张应力的无黏结接触数量逐渐增多。因此随着循环应力加载周期增加,相同应力水平下模型侧向应变增加速率大于轴向应变增加速率,循环应力加载过程中模型侧向变形占主导地位,模型泊松比逐渐增加。在循环应力加载过程中,模型内新生成的无黏结接触数量和微裂缝数量增加速率呈递减趋势,因此模型损伤量增加速率同样呈递减趋势。

4 结 论

(1)随着循环应力加载周期增加,新生成的无黏结接触数量和微裂缝数量均逐渐增加,且新生成的无黏结接触分布和微裂缝分布均受原始裂缝分布控制明显,而孔隙度逐渐降低,且无黏结接触数量、微裂缝数量和孔隙度变化速率均呈递减趋势。

(2)随着循环应力加载周期增加,相同应力水平下侧向应变和轴向应变均增大,但侧向应变占主导地位,弹性模量逐渐减小,而泊松比和损伤量逐渐增加,且弹性模量、泊松比和损伤量变化速率均呈递减趋势。