气藏型储气库圈闭应力场扰动规律及影响因素分析

张 玥，贾善坡,2，温曹轩，张品金，李 健，张鲁钢，施少军

( 1. 东北石油大学 地球科学学院，黑龙江 大庆 163318；2. 东北石油大学 油气藏及地下储库完整性评价重点实验室，黑龙江 大庆 163318；3. 中国石化胜利油田有限公司 勘探开发研究院，山东 东营 257000；4.中国石油吉林油田有限公司 勘探开发研究院，吉林 松原 138000 )

0 引言

地下储气库作为天然气产业链中的关键一环，建设与运行程度越来越高[1-2]。中国天然气地下储气库主要有4种类型:枯竭油气藏型、盐穴型、废弃矿井型及含水层型[3-4]。截止2021年底，中国建成28座地下储气库，其中25座为气藏型储气库[5]。枯竭气藏型储气库具有地质资料齐全、圈闭条件良好、储气空间大等优点，同时具有完整配套的天然气地面系统工程设施可供选择[6]。相比其他3种类型储气库，枯竭气藏型储气库具有更高的可靠性、安全性和经济性[7-8]。地下储气库季节性周期注采引起储层孔隙压力变化，间接引起储层、盖层及其他地层中的应力场交替变化，应力路径的改变可能激活盖层中存在的天然裂缝及断层，影响盖层密封性及力学完整性。作为阻止油气向上逸散的保护层，盖层封隔储集层中的气体，应力场变化影响储气库选址、设计和建设运行。在大流量高速注采时，储气库储集层局部高压改变盖层应力场，盖层发生宏观、微观力学破坏[9]，出现隆起而产生裂缝，导致天然气泄露[10]。有必要研究枯竭气藏型储气库在注采过程中盖层应力场变化，减少储气库泄露危险，确保地下枯竭气藏型储气库安全运行。

人们从不同角度对储气库圈闭进行研究。基于流体渗流弹性介质的Biot本构理论，SEGALL P[11]研究流体注采时区域位移场的变化，表明沉降量随流体抽采率线性增加，但未研究流体注采时圈闭应力场的变化。CHAMANI A等[12]建立三维数值概念模型，研究流体注入储气库时应力场变化，分析流体注入时弹性模量变化对应力场的影响，但未考虑流体采出时应力场的变化。ZHOU X等[13]研究储层与盖层交界面的剪应力对不同盖层渗透率的响应，表明最大剪应力发生在近井区域，当盖层渗透率变大时，最大剪应力转移到远井区，且近井区的剪应力可能在界面发生逆转。LI C等[14]采用有限元数值模拟方法，研究不同温度下流体注入储气库时应力场的变化，表明在注入井附近，随流体压力的增加，平均有效应力不断减小，但未指明流体采出时的应力场变化。FAVERO V等[15]研究CO2注入含水层后盖层的沉降变形，但未考虑圈闭应力场的变化。王迪晋等[16]应用GPS系统研究枯竭型储气库在注气过程中地表盖层的变形响应，表明储气库注采过程引起地表盖层形变，在水平方向上呈明显的“呼吸效应”，垂直方向上地表盖层的运动方向与注采过程呈负相关关系。李宏等[17]建立数值模型研究储气库注气后地表的变形，但忽略圈闭应力场的变化。根据多孔介质弹性理论，王成虎等[18]研究储气库在周期注采下地层应力的变化，表明水平总应力增量随储气库注采过程出现线性正相关关系的变化。张广权等[19]建立三维地应力模型与数值模型，研究储气库注采时圈闭的动态密封性，分析孔隙压力与三轴主应力的分布。

人们对储气库圈闭应力场变化及应力路径因数的研究相对较少，忽略对盖层物理力学参数的敏感性分析。盖层物理力学参数对盖层应力场的变化产生不同程度的影响，盖层的渗透率是决定注入气体能否被密封的重要性质之一[20]。在分析圈闭应力场变化时需要对盖层物理力学参数进行分析。基于渗流—力学耦合理论，笔者建立天然气注采过程的储气库力学模型，探讨储气库在注气和采气阶段圈闭扰动应力场的变化特征及应力路径的发展趋势，分析盖层物理力学参数变化对盖层应力场的影响，为气藏型储气库建库设计及安全运行分析提供依据。

1 渗流—力学耦合理论

1.1 渗流模型

假设储气库中的天然气注采符合达西定律，注采过程视为等温过程，建立地下储气库圈闭渗流模型。基于达西定律的渗流连续方程为

(1)

式中:ρf为流体密度；φp为岩石孔隙度；d为渗流速度；Qm为源汇项；t为时间。

考虑流体及岩石骨架的可压缩性，分析地下储气库的地层受力，则

(2)

式中:S为弹性储存系数；p为孔隙流体压力。

储气库注入与采出的流体质量表现为孔隙产生的体积应变对时间的导数，则质量守恒方程表示为

(3)

式中:εv为体积应变，εv=εx+εy+εz，εx、εy、εz为体积应变分量；α为Biot系数。

1.2 力学模型

在场地尺度下，储气库围岩产生的变形是弹性的、小变形的[16-21]。具有小变形的弹性三维控制方程[22]表示为

(4)

式中:ε和σ分别为应变和应力；u为位移；f为体力；C为四阶弹性系数张量。

假设储气库地下岩土体为各向同性的材料，则式(4)中的本构方程为

(5)

式中:E为弹性模量；ν为泊松比；εxx、εyy、εzz和γxy、γxz、γyz分别为应变在纵向和剪切方向上的分量。

(6)

式中:σij为总应力；δij为Kroneker张量，当i=j，δij=1，i≠j，δij=0。

基于式(5-6)，应力平衡方程的张量形式为

σij-αδijp+Fi=0，

(7)

式中:Fi为不同方向的体载荷。

将式(4-6)代入式(7)，得到有效应力下具有小变形的弹性三维控制方程为

，

(8)

1.3 耦合方程

储气库场地尺度下，渗流—力学耦合方程的控制方程由渗流连续性方程与力学平衡方程组成，结合式(3)与式(8)，分析渗流—力学耦合问题[23-24]。适用于储气库圈闭应力场的渗流—力学耦合方程为

(9)

储气库注采过程的力学场向渗流场的耦合程度取决于岩石骨架和孔隙流体的性质，围岩中流体质量或流体压力的变化对岩石内部的孔隙结构产生体应变，在注采气压力作用下，岩石骨架产生体应变和弹性应力。力学场建立孔隙压力、渗流、体应变与有效应力之间的联系，压力作用产生渗流，引起孔隙压力消散和扩散；渗流场产生体积变化，根据达西定律，孔隙压力使流体流动随时间发生变化，产生随时间变化的弹性应力和体应变，进而耦合回渗流场[25]。

2 储气库力学数值模型

枯竭气藏型储气库有限元计算模型采用二维平面模型，对储层与盖层进行渗流—力学耦合计算。根据地质构造资料，区域构造平缓，属于低幅度背斜构造，地层近水平状，且无断层发育。在不影响整体数值模拟结果的情况下，对实际地质构造及工况进行简化，建立储气库力学数值模型(见图1(a))。该模型长为400 m，高为600 m，由上覆层、盖层、储层、围岩及底板组成。储层长为80 m，高为40 m；盖层长为400 m，高为80 m。上覆层上部有深度为2 000 m的地质层，以等效覆盖层表示。根据储气库“夏注冬采”季节性调峰需求等特点，将夏季与秋季作为注入阶段，周期为180 d；冬季与春季作为采出阶段，周期为180 d。数值模拟方案选择储气库运行阶段，地层压力曲线见图1(b)。

根据储气库地质资料与岩石力学实验资料，盖层弹性模量为5.0 GPa，泊松比为0.25，渗透率为10-6μm2，孔隙度为0.01。储层弹性模量为2.5 GPa，泊松比为0.30，渗透率为890×10-3μm2，孔隙度为0.34。由于该区域勘探程度较低，岩心测试资料较少，上覆层、底板及围岩的地层计算参数采取邻区的平均值，平均弹性模量为5.0 GPa，平均泊松比为0.25，平均渗透率为10-6μm2，平均孔隙度为0.01。

假设储气库注采作业的流体为单相流体，地层中的流体及岩石骨架可压缩，上覆层施加的荷载考虑上覆岩体和流体的质量，根据垂向应力计算理论[26]，施加在上覆岩层的荷载约为54.31 MPa。区域构造活动并不强烈，储气库力学数值模型不考虑构造应力场。模型初始水平应力为各向异性，地应力因数为0.33。储气库上限压力因数为1.2，下限压力因数的选取应考虑储气库的工作气量与经济效益，为了突出注采过程应力场扰动变化及应力路径因数的差异，将下限压力因数选取为0.8。

数值模拟采用超孔压方法，力学边界:上边界施加上覆压力；左、右边界为沿水平方向的指定位移约束；下边界为沿垂直方向的指定位移约束。渗流边界:所有边界为无流动边界。数值模型网格划分见图1(c)，将储层与盖层的网格部分进行加密，网格加密采用映射分布法，单元数为16，单元大小比为1，增长公式采用等差数列，模型求解的自由度数为6.432 3×104。

图1 储气库力学数值模型及网格划分Fig.1 Mechanical numerical model and grid division of gas storage

3 数值模拟分析

3.1 初始应力场平衡

为了准确还原岩石实际受力状态，对地层进行地应力平衡。地应力平衡后的水平总应力σH与垂直总应力σV分布见图2(a-b)，垂向位移场的最大量级为10-15m，最小量级为10-17m(见图2(c))，储层及盖层的水平总应力为-23.00～-22.00 MPa，垂直总应力为-69.00～-66.00 MPa。

图2 储气库地应力平衡后总应力及垂向位移分布Fig.2 Distribution of total stress and vertical displacement after ground stress balance

3.2 扰动应力场特征

在注采气作业后，储气库孔隙压力的变化对盖层、储层和底板的应力场扰动较大。储气库注气和采气的总应力增量见图3-4。由图3-4可知，注气阶段，储气库圈闭应力场的增量为正值，水平总应力增幅范围为0～1.90 MPa，垂直总应力增幅范围为0～3.82 MPa；采气阶段，圈闭应力场的增量为负值，水平总应力增幅范围为-2.90～0 MPa，垂直总应力增幅范围为-6.10～0 MPa。储气库水平总应力增量小于垂直总应力增量，采气阶段的总应力增幅范围大于注气阶段的。

受孔隙压力扰动影响的圈闭应力场不均匀分布，为进一步分析储气库在注采期间总应力场的变化，将总应力场按扰动程度进行分区(见图5)。由图3和图5可知，注气阶段，B区与D区的水平总应力与垂直总应力大于其他区域的，D区最易形成应力集中区；A区作为B区的上覆岩层，孔隙压力增大，有效应力减小，A区形成压实区；C区作为储层外部的围岩区域，应力场扰动程度比储层的小，总应力场增量为正值，说明C区处于压缩状态，属于压实区。由图4-5可知，采气阶段，圈闭处于拉伸状态，发生拉伸变形，D区垂直总应力增量的绝对值最大，大于采出的孔隙流体压力，其他区域的小于采出的孔隙流体压力，越靠近B区与D区，应力场受扰动程度越大；水平总应力增量与垂直总应力增量相比，变化较小，但不能忽视D区发生的应力集中现象。在储气库季节性周期注采作业时，应考虑储层与盖层边界处及底板区域内的应力场变化及扰动程度。

图3 储气库注气阶段总应力增量分布 Fig.3 Total stress increment distribution during gas injection in gas storage

图4 储气库采气阶段总应力增量分布 Fig.4 Total stress increment distribution during gas production in gas storage

图5 储气库注采气阶段圈闭总应力场分区特征Fig.5 Division characteristics of trap total stress field during gas injection and production in gas storage

3.3 储层与盖层应力路径效应

储气库注采过程中，储层与盖层的总应力场受扰动程度较大，应力路径可能发生突变。在储层内，根据(0 m，-500 m)、(400 m，-500 m)两点设置储层中心水平测线，根据(200 m，-400 m)、(200 m，-600 m)两点设置储层中心竖向测线；在盖层区域内，在垂向距离-480、-440和-400 m处分别设置储层与盖层交界面、盖层中部和盖层顶部3条水平测线，分析5条测线位置的储层和盖层应力路径效应(见图6)。

图6 储气库储层和盖层应力测线分布Fig.6 Stress line distribution of reservoir cap in gas storage

在储气库周期性季节注采作业时，假设储层满足无水平应变的单轴应变条件，结合储气库圈闭中岩石实际受力状态，由孔隙压力变化导致总应力变化的计算公式为

(10)

受注采气作业影响后，储气库应力场孔隙压力呈不同程度的增长(见图7)。由图7可知，储层的孔隙压力增量最大，其次是盖层和围岩的。在垂向范围内，储层注气阶段的孔隙压力变化程度稍大；在水平范围内，由储层中心到盖层顶部，孔隙压力增量越来越小，盖层顶部的最大孔隙压力增量是储层中心的0.68倍。根据不同测线的注采总应力增量(见图8-12)，其中βH为水平应力路径因数，βV为垂直应力路径因数，距离注采区越远，总应力场受扰动程度越小；储层中心到盖层中部，注气阶段的总应力场小于采气阶段的，且垂直总应力增量比水平总应力增量大；盖层顶部在注气阶段的总应力增量比采气阶段的大，水平总应力增量比垂直总应力增量大。

由式(10)计算储层区域内的应力路径因数为0.57。由图8-12可知，在水平距离内，应力路径因数易在储层与盖层交界面发生突变，在垂向距离内，应力路径因数易在储层发生突变。距离注采区越远，应力路径因数发生突变的位置从储层与盖层交界面变为储层中心，采气阶段的水平应力路径因数和垂直应力路径因数达到的峰值和谷值变大。垂直应力路径因数整体发展趋势为先下降后上升，水平应力路径因数发展趋势为先上升后下降。采气阶段的垂直应力路径因数大于解析解求解的(0.57)，注气、采气阶段的水平应力路径因数在储层内小于0.57，注气阶段的垂直应力路径因数由储层中心到储层与盖层交界面，逐渐接近解析解求解的值，并在储层边界处交会。应力路径因数整体表现为垂直应力路径因数大于水平应力路径因数，采气阶段的应力路径因数大于注气阶段的。

图7 储气库注气和采气阶段不同范围的孔隙压力增量Fig.7 Pore pressure increment in different range during gas injection and production in gas storage

图8 储气库储层中心水平测线总应力增量及应力路径因数曲线Fig.8 Total stress increment and stress path coefficient curves of horizontal survey line in reservoir center of gas storage

图9 储气库储层中心竖向测线总应力增量及应力路径因数曲线Fig.9 Total stress increment and stress path coefficient curves of vertical survey line in reservoir center of gas storage

图10 储气库储层与盖层交界面测线总应力增量及应力路径因数曲线Fig.10 Total stress increment and stress path coefficient curves of survey line in reservoir and cap interface of gas storage

图12 储气库盖层顶部测线总应力增量及应力路径因数曲线Fig.12 Total stress increment and stress path coefficient curves of survey line in top cap rock of gas storage

3.4 剪切应力场

储气库注气、采气阶段圈闭剪切应力及增量分布见图13-14。由图13-14可知，储层边界更易受剪切应力的影响，储层中心的剪切应力为0 MPa，盖层距离储层较远，剪切应力场受扰动程度较小。注气阶段的剪切应力增量为-1.00～0.60 MPa，采气阶段的剪切应力增量为-1.20～0.80 MPa，当储气库由注气阶段转为采气阶段时，剪切应力场发生扭转，注气阶段剪切应力增量的正值变为负值，负值变为正值。储气库注气、采气阶段不同深度下剪切应力增量曲线见图15。由图15可知，储层与盖层交界面的剪切应力增量更易出现尖点，注气阶段的剪切应力增量最大约为0.80 MPa，采气阶段的剪切应力增量最大约为1.20 MPa，是注气阶段剪切应力增量的1.5倍，距离注采区越远，剪切应力增量越不容易达到峰值。根据储气库注采气阶段剪切应力路径因数分布(见图16)，与水平应力路径因数和垂直应力路径因数相比，剪切应力路径因数变化量极小。注气阶段，距离储层越远，应力路径因数越平滑，起伏越小。与注气阶段相比，采气阶段的剪切应力增量略大，距离储层较远，应力路径因数不再沿储层中心分布，逐渐接近垂直分布。

图13 储气库注气、采气阶段圈闭剪切应力分布Fig.13 Trap shear stress distribution during gas injection and production in gas storage

图14 储气库注气、采气阶段圈闭剪切应力增量分布Fig.14 Incremental distribution of trap shear stress during gas injection and production in gas storage

图15 储气库注气、采气阶段不同深度下剪切应力增量分布Fig.15 Shear stress increment distribution at different depths during gas injection and production in gas storage

图16 储气库注气、采气阶段不同深度下剪切应力路径因数分布Fig.16 Distribution of shear stress path coefficients at different depths during gas injection and production in gas storage

4 影响因素分析

储气库封盖应力场的变化主要取决于两个要素:一是盖层内压力的变化幅值；二是封盖岩石的变形(弯曲梁效应)。若盖层内孔隙压力变化较大，则有效应力变化较大，主要与盖层渗透率有关。盖层可视为球形弯曲梁，注气、采气阶段盖层变形使盖层内部产生附加应力场。对盖层物理力学参数进行敏感性分析，研究天然气注气、采气过程的盖层参数对盖层应力场变化幅度的影响。由于储层与盖层交界面的应力场受扰动程度更大，选择储层与盖层交界面进行应力场变化分析。

4.1 盖层渗透率

储气库注气、采气阶段不同盖层渗透率的应力增量曲线见图17。由图17可知，盖层渗透率K越大，注气阶段的水平总应力、垂直总应力和剪切应力的增量峰值越大；垂直总应力增量变化最大，水平总应力增量次之，剪切应力增量最小；盖层渗透率越小，在储层与盖层交界面处达到峰值的速度越快，储层上方区域的应力增量变化越大。

图17 储气库注气、采气阶段不同盖层渗透率应力增量曲线Fig.17 Stress increment curves of different caprock permeability during gas injection and production in gas storage

4.2 盖层弹性模量

储气库注气、采气阶段不同盖层弹性模量的应力增量曲线见图18。由图18可知，采气阶段，储层与盖层交界面的应力扰动程度大于注气阶段的，盖层弹性模量E越大，储层与盖层交界面的应力增量增长速度越快，且应力增量的谷值与峰值的差距越大。总应力场增量特征总体表现为:水平总应力增量增幅大于垂直总应力增量增幅大于剪切应力增量增幅。储层上方区域的垂直总应力增加发展趋势比水平总应力增加的更平缓。当盖层弹性模量较小时，剪切应力场可能发生应力扭转。

图18 储气库注气、采气阶段不同盖层弹性模量应力增量曲线Fig.18 Stress increment curves of different cap modulus during gas injection and production in gas storage

4.3 盖层泊松比

储气库注气、采气阶段不同盖层泊松比的应力增量曲线见图19。由图19可知，采气阶段，储层与盖层交界面受应力扰动程度大于注气阶段的，随盖层泊松比的增大，水平总应力增量、垂直总应力增量及剪切应力增量达到的峰值更大，谷值更小。注气阶段，盖层泊松比的改变，对盖层区域内的垂直总应力影响较小，储层区域内对峰值有影响，谷值影响不大。随盖层泊松比的增大，注气阶段应力增量增幅由大到小依次为水平总应力增量、垂直总应力增量、剪切应力增量。

图19 储气库注气、采气阶段不同盖层泊松比应力增量曲线Fig.19 Stress increment curves of Poisson's ratio of different cap layers during gas injection and production in gas storage

盖层物理力学参数对储层与盖层交界面的应力场扰动程度由大到小依次为:盖层弹性模量、盖层渗透率、盖层泊松比。采气阶段应力增量增幅比注气阶段的大，盖层物理力学参数的改变对水平总应力增量、垂直总应力增量及剪切应力增量的影响幅度不一致:当盖层弹性模量与泊松比变大时，水平总应力增量增幅大于垂直总应力增量增幅大于剪切应力增量增幅；当盖层渗透率变大时，垂直应力增量增幅大于水平总应力增量增幅大于剪切应力增量。

5 结论

(1)基于岩石渗流—力学耦合理论，建立储气库力学数值模型，研究注气、采气阶段储气库圈闭应力场变化，对盖层物理力学参数进行敏感性分析。

(2)储气库圈闭应力场经过季节性注采作业后不均匀分布，距离储层越近，应力场受扰动程度越大，表现为储层与盖层交界面的应力场受扰动程度最大，其次是盖层、底板的，最小是上覆层的。储层与盖层交界面的应力场受扰动程度由大到小依次为盖层弹性模量、盖层渗透率和盖层泊松比。

(3)采气阶段的应力场受扰动程度大于注气阶段的，垂直总应力增量大于水平总应力增量大于剪切应力增量；采气阶段的应力路径因数大于注气阶段的，垂直应力路径因数大于水平应力路径因数，应力路径因数易在储层内发生突变，在储气库圈闭水平范围内，垂直应力路径因数总体发展趋势为先下降后上升，水平应力路径因数发展趋势为先上升后下降。