赵海洋 刘志远 唐旭海 李新勇 耿宇迪 张俊江

1.中国石化西北油田分公司石油工程技术研究院；2.中国石化缝洞型油藏提高采收率重点实验室；3.武汉大学土木建筑工程学院

缝洞型碳酸盐岩油气藏储量丰富，天然裂缝普遍存在于储层中，在压裂增产过程中天然裂缝对水力裂缝延伸方向和扩展形态具有重要影响，能否利用天然裂缝沟通溶洞储集体是储层改造成功与否的关键［1-4］。

许多专家学者针对储层天然裂缝对油气田开发的影响问题进行了理论和实验研究。在油气藏开采方面，樊建明等［5］采用高斯函数对地应力、天然裂缝、人工裂缝各自优势方向进行耦合，将利用自然能量开发井网作为研究对象建立数值模型，研究不同优势方向下油藏的开采方法。姚军等［6］研究了油藏中大裂缝分布、长度、方位对水驱油的影响。苏玉亮等［7］通过对天然裂缝导流能力、长度、密度、方位如何影响油藏开采展开分析，得出天然裂缝的密度、方位对采收率有较大影响。李博［8］利用嵌入式离散裂缝模型对缝洞型油藏进行模拟研究，模拟结果表明天然裂缝数量越多，储层连通性越好，充分连通水力裂缝和天然裂缝对增产具有重要意义。孙若凡［9］在水平井压裂产能模型中考虑了天然裂缝，发现天然裂缝所占比例相对于它的渗透率对产能影响更大。朱圣举等［10］结合现场工程实际，对影响天然裂缝开启压力的因素进行了研究，结果表明对注水压力的控制可以优化开发效果。刘熙远［11］建立了三维地质模型，在储层注水开发过程中考虑天然裂缝，认为沿着裂缝方向注采井连通性更强。吴忠宝等［12］基于等值渗流阻力法原理和现场施工资料，建立微裂缝发育的地质模型和模拟方程，模拟结果和现场结果吻合较好。王宁［13］对缝洞型油藏中天然微裂缝发育程度进行量化分析，发现裂缝发育程度对油藏压后产能有较大影响。在储层改造方面，曾青冬等［14］建立井筒和裂缝流体流动模型，讨论天然裂缝逼近角与水力裂缝扩展的关系。马耕等［15］进行了存在天然裂缝的水力压裂实验，实验结果表明天然裂缝形态对水力裂缝扩展有不同影响。宋晨鹏等［16］建立二维模型对水力裂缝和天然裂缝相交扩展进行模拟，发现小尺寸的天然裂缝对水力裂缝影响不大。林鹤等［17］利用微震检测和施工数据研究发现活动性天然裂缝相对于非活动性天然裂缝对压裂改造效果的影响更明显。赵金洲等［18］研究了水力裂缝沟通天然裂缝后不同因素对裂缝网络的影响，发现天然裂缝发育、地应力差小的环境有利于形成缝网。赵海洋等［19］、荣元帅等［20］对塔河缝洞型碳酸盐岩高效酸压改造进行了大量的现场论证和分析。

上述研究主要集中在天然裂缝连通对储层产能的影响，对于缝洞型碳酸盐岩油藏中天然裂缝如何影响压裂裂缝扩展及沟通溶洞的机理尚不明确，没有提出明确的技术思想。本文基于缝洞型碳酸盐岩储层特点，提出了循缝找洞改造储层思想，同时通过建立天然裂缝和溶洞储集体同时存在的压裂模型，研究单条天然裂缝及天然裂缝网络对压裂裂缝扩展轨迹的影响，对缝洞型油藏中压裂裂缝如何沟通溶洞进行模拟研究，并对研究成果进行推广应用。

1 TOUGH2-AiFrac耦合求解算法

1.1 固相控制方程

考虑一个包含裂缝面Γf的 三维体Ω，拉力F施加在外边界Γh上，边界的法向量为ν，加载在准静态情况下，控制方程的有限元强形式为

式中，σ为应力，Pa；ρ为材料密度，kg/m3；b为三维体的体力，N；t为边界单位法向载荷，Pa。

在水力压裂过程中，在裂缝表面 Γf上，拉力满足在裂缝面两侧连续条件

式中，f为裂缝面局部的单位法向量，无量纲；p为作用在裂缝面表面的水压，Pa。

将总位移u分解为连续部分uc和非连续部分ud，可以得到总位移为

将裂缝面从正面到负面的位移差定义为

根据变形体虚功原理，得到控制方程的有限元弱形式最终表示为

式中，δ为开度，m；ε为应变；u为位移，m。

1.2 TOUGH2流体模拟原理

TOUGH2流体模拟采用基于积分有限差分的空间离散化技术(IFD)和Newton-Raphson方法，对流体各分量的质量平衡方程和各时刻的能量平衡方程进行线性化。在TOUGH2中，每一流体组分的总质量平衡由每一阶段的相关贡献进行累积。对于任意组分得到

式中，MK为 组分K的单位体积质量，kg/m3；QK为组分K的能量通量，kg/(m3· s )；qKΨ为Ψ相中K组分的流量(液体Ψ 用l表示，气体 Ψ用g表示；空气K用a表示，水K用w表示)，kg/(m3· s)。

MK的计算公式为

式中，ϕ为孔隙度；S为饱和度；ρΨ为 Ψ相密度，kg/m3；为在Ψ 相中K组分的质量分数。

在Ψ相中，每个组分的质量流量为达西流量和扩散流量之和

式中，ka为绝对渗透率，m2；krΨ为Ψ相相对渗透率，%；µΨ为 Ψ相黏度，Pa· s ；pΨ为Ψ相压力，Pa；a为重力加速度，m/s2；为Ψ相中K组分的扩散通量，kg/(m3· s )。

用Fick定律计算扩散通量

式中，SΨ为 Ψ相饱和度；DΨ为介质中依赖于温度、气压、介质曲折率和气体饱和度的有效分子扩散系数，在模拟过程中假设温度是恒定的；τ为介质曲折率。

利用上述方程及时间的一阶有限差分离散化，在时间ek+1=ek+∆e，平衡方程重新定义为

式中，R为残差，kg/m3；M为单位体积质量，kg/m3；V为单元体积，m3；A为单元面面积，m2；q为通量，kg/(m2· s)；Q为源项，kg/(m3· s)；k为迭代步数；m为离散单元号。

对于非线性，采用牛顿-拉弗逊迭代法生成线性方程，在每个时间步中迭代，直到满足某个收敛条件。

1.3 AiFrac模拟原理

AiFrac模拟器基于有限无网格法(FEMM)，有限无网格法主要用来模拟非连续模型［21-22］。如图1所示的横截面存在3种不同的单元类型：红色的裂缝单元、黄色的过渡单元和绿色的FEMM单元。裂缝穿过的单元称为裂缝单元，与裂缝单元相邻的单元称为过渡单元，剩下的单元为FEMM单元。存在两种结点：黑色实心的FE结点和黑色空心的PU结点，裂缝单元包含的所有结点为PU结点，FEMM单元包含的所有结点为FE结点，过渡单元同时包含PU结点和FE结点。对于任意点e={x,y,z}，在单元 Ω上，位移场u(e)可以由节点形函数和节点位移插值近似为

式中，u(e)为 FEMM中的全局近似，ωi(e)为非负权重函数且权重函数总和等于1，ui(e)为 结点i的局部近似。与FE结点相关的局部近似设为1，用最小二乘法构造在节点上的局部近似。

对于断裂单元，应正确反映断裂表面的不连续位移场。根据Shepard公式建立分形单元的权函数，对于四面体单元Ω 的节点检索集为 ζΩ={1,2,3,4}，按照可视化准则，假设裂缝面是不可见的，已知单元内计算点e，定义节点对应的一组函数为φ(e)={φ1(e),φ2(e),φ3(e),φ4(e)}， 对于计算点e对应的四面体单元每个节点的有限元形函数 φi(e)可以由体积坐标得到

图1 裂缝单元、过渡单元和FEMM单元的定义Fig.1 Definition of fracture unit, transition unit and FEMN unit

式中，vol(E1E2E3E4)为 四面体单元体积，vol(E(e)EiEjEk)为任意点E(e)和其他3个顶点{Ei,Ej,Ek}构成的四面体单元体积。当计算点e与单元节点的连线与裂缝面相交时，将其φi(e)设为零函数。

定义 ζΩ为FEMM中与单元域Ω 相关的一个节点检索集，定义可见区为计算点对应的所有非零加权函数的集合。对于包含裂缝面的任意单元Ω中的计算点e，φi(e)可以表示为

利用下式得出包含裂缝面与结点相关的有限元法单元权函数为

过渡单元的体积权函数与标准有限元公式相同，传统的有限元单元形函数用标准有限元形函数计算。

1.4 耦合求解器

TOUGH2-AiFrac耦合求解算法，将TOUGH2流体求解器［23-24］与AiFrac固体求解器［24］相结合，耦合过程见图2。流体流动和流体压力计算采用TOUGH2求解器，岩石变形和裂缝扩展计算采用AiFrac固体求解器。

图2 TOUGH2-FEMM仿真原理Fig.2 TOUGH2-FEMM simulation principle

TOUGH2和AiFrac建模是在不同的网格上进行的，但是几何结构和节点编号是相同的。如图3所示，TOUGH2分析每个单元的中心点，AiFrac需要在每个单元的8个节点上进行三维分析。

图3 TOUGH和FEMM计算节点示意图Fig.3 Sketch of TOUGH and FEMM calculate nodes

建模数据将TOUGH2单元的中心点插入到AiFrac单元的顶点，生成TOUGH2和AiFrac网格的过程见图4。首先，用TOUGH2建立了一套用于模拟裂缝和围岩基质中流体流动的六面体单元，然后，在之前六面体单元的基础上，使用AiFrac生成了一组四面体单元，用于模拟岩石变形和压裂，。

图4 TOUGH网格构造出AiFrac网格Fig.4 AiFrac network constructed on the basis of TOUGH network

2 单条天然裂缝对缝洞沟通的影响

研究缝洞型储层中单条天然裂缝对缝洞沟通的影响，建立如图5所示的模型。在35 m×35 m×1.5 m的假三维平板岩体模型中存在长10 m的水力裂缝，x方向为水平最大地应力方向，假设模拟岩体各向同性、线弹性，模型参数设置如下：岩石弹性模量30 GPa，泊松比0.2，孔隙度0.05，岩石渗透系数1×10−18，裂缝孔隙度0.25，裂缝渗透系数 3×10−10。在水力裂缝表面施加20 MPa的均匀压力。在假三维平板岩体模型中，将天然溶洞几何形状简化为圆柱体，圆柱体上下面圆心坐标为(25,17,0.75)，圆柱体厚度为1.5 m，溶洞内充满流体，压力为20 MPa。

图5 单条天然裂缝-溶洞的沟通模拟模型Fig.5 Simulation model of single natural fracture-vug connection

考虑到岩体网格尺寸引起的裂缝扩展敏感性，当岩体网格尺寸越小时，裂缝扩展模拟结果越接近于解析解，但同时会造成计算量大幅度增加。通过多组网格尺寸的模拟对比，选取h=0.5 m作为模拟网格尺寸。岩体四面体网格总数117 600，将圆柱体上下面圆心坐标(25,17,0.75)，半径4 m内框选的所有网格视为溶洞网格。

2.1 无天然裂缝时缝洞沟通模拟计算

建立无天然裂缝的缝洞模型(图6)，基于TOUGH2-AiFrac耦合求解器计算水力裂缝扩展及缝洞沟通情况。模拟设置6个计算步，结果表明无天然裂缝时水力裂缝在压力作用下水平扩展，未与天然溶洞连通。

图6 无天然裂缝的缝洞模型Fig.6 Fracture-vug model without natural fracture

2.2 天然裂缝偏离溶洞时缝洞沟通模拟计算

基于无天然裂缝的缝洞模型，在压裂裂缝扩展路径上增加一条与水平最大地应力方向夹角50°的天然裂缝，但裂缝走向偏离溶洞(图7)，研究天然裂缝对缝洞沟通的影响。假设天然裂缝未与水力裂缝连通时，缝内压力忽略不计。模拟设置3个计算步，结果表明压裂裂缝与天然裂缝连通，但未与溶洞连通(图8)。

图7 天然裂缝偏离溶洞的缝洞模型Fig.7 Fracture-vug model with natural fractures deviating from vugs

图8 天然裂缝偏离溶洞时缝洞沟通结果Fig.8 Fracture-vug connection result of natural fractures deviating from vugs

2.3 天然裂缝指向溶洞时缝洞沟通模拟计算

对比于上述天然裂缝偏离溶洞的情况，在压裂裂缝扩展路径上增加一条与水平最大地应力方向夹角43°的天然裂缝，裂缝走向指向溶洞(图9)。模拟设置3个计算步，结果表明压裂裂缝与天然裂缝连通后，与附近的溶洞储集体连通(图10)。

图9 天然裂缝指向溶洞的缝洞模型Fig.9 Fracture-vug model of natural fractures deviating towards vugs

图10 天然裂缝指向溶洞时缝洞沟通结果Fig.10 Fracture-vug connection result of natural fractures deviating towards vugs

3 天然裂缝网络对缝洞沟通的影响

天然裂缝作为缝洞型油藏中沟通溶洞储集体的重要通道，在储层中广泛存在。研究天然裂缝网络影响水力裂缝扩展及缝洞沟通，对于提高缝洞型油藏油气产量具有实际意义。储层中天然裂缝的特征表述包括长度L、面密度α、走向和倾角，考虑到研究区储层的特性，设置一组长度5±1 m，面密度2%，走向平行于z轴，倾角45±5°的天然裂缝网络。如图11所示在50 m×50 m×1.5 m的假三维平板岩体模型中存在长7 m的水力裂缝，x方向为水平最大地应力方向，在水力裂缝表面施加20 MPa的均匀压力，将天然溶洞几何形状简化为圆柱体，圆柱体圆心坐标(40,13,0.75)，厚度为1.5 m，溶洞内部充满流体，压力20 MPa，模型边界位移设置为0，模型参数取值与第2节相同。模型选取h=0.5 m作为岩体模拟网格尺寸，与第2节模型模拟网格相同。

图11 天然裂缝网络-溶洞的缝洞沟通模型Fig.11 Fracture-vug connection model of natural fracture network-vug

基于无天然裂缝网络时缝洞沟通模型，加入指向溶洞的天然裂缝网络(图12)。假设每条天然裂缝是独立存在的，将多条已经贯通的天然裂缝视为一条裂缝，不考虑天然裂缝中的流体压力。

图12 天然裂缝网络指向溶洞时缝洞沟通结果Fig.12 Fracture-vug connection result of natural fracture network deviating towards vugs

存在天然裂缝网络时，水力裂缝在压力作用下先沿着水平向扩展，之后与天然裂缝沟通，并将流体注入天然裂缝内，天然裂缝在压力作用下扩展，连通更多的天然裂缝，形成相互连通的裂缝网络，最终连通远端的溶洞储集体。

4 缝洞型油藏三维缝洞沟通模拟计算

建立如图13所示的三维缝洞型油藏模型，模拟岩体尺寸50 m×50 m×50 m，位于岩体中心的圆形水力裂缝长5 m，x方向为水平最大地应力方向。设置两组长度5±2 m，走向平行于z轴，倾角分别为45°和135°的天然裂缝网络，裂缝面密度2%。在水力裂缝表面施加20 MPa的均匀压力，天然溶洞圆心坐标分别为(10,30,25)、(38,14,25)，溶洞内部充满流体，压力20 MPa，模型边界位移设置为0，模型参数设置与第2节相同。

图13 三维天然裂缝网络-溶洞的缝洞沟通模型Fig.13 Fracture-vug connection model of 3D natural fracture network-vug

计算结果显示，水力裂缝在压力作用下向四周扩展，连通上、下两条倾斜的天然裂缝，并注入流体，上、下两条天然裂缝在压力作用下扩展，连通更多的天然裂缝，形成相互连通的裂缝导流网络，最终沟通远端的溶洞储集体(图14)。

图14 三维天然裂缝网络-溶洞的缝洞沟通结果Fig.14 Fracture-vug connection result of 3D natural fracture network-vug

上述模拟研究表明，对于不在水平最大地应力方向上的溶洞，人工主裂缝难以对其进行有效沟通。为增加溶洞储集的沟通几率和数量，需充分利用储层中的天然裂缝，在不同方向形成相互连通的裂缝通道，实现对不同方向的多个溶洞进行沟通，即通过控制施工参数，使注入流体遵循原生天然裂缝的展布形态进行流动、连通，利用缝洞型碳酸盐岩储层缝-洞伴生发育的特点，沟通井周不同方向的溶洞储集体，即循缝找洞。

5 现场应用

以循缝找洞思想为指导，制定储层改造技术方案，现场应用85井次，有效率100%，与常规酸压邻井相比，2020年增加原油产量16.12万t，增油效果显著。以W1井为例，该井区域油气富集，但生产无产量，拟通过储层改造实现增产，改造层段为奥陶系6 413.98~6 543.00 m，层位为O2yj+O1−2y，岩性为石灰岩，缝洞储层发育，且分布在井筒附近。施工过程中，通过在压裂初期注入压裂液形成一定长度的人工主裂缝，之后注入酸液充分激活人工主裂缝侧面及未端的天然裂缝，形成相互连通的裂缝通道，以增加沟通溶洞的几率。注入普通胶凝酸后在注入排量不变的条件下，压裂压力呈现较大幅度下降，显示人工裂缝与溶洞储集体存在沟通。

W1井在压裂前基本无产量，压裂改造后产量呈现大幅度上升，平均日产量37.7 t，且持续稳产。

6 结论

基于TOUGH2-AiFrac耦合求解算法建立了缝洞型油藏压裂裂缝扩展模型，研究了天然裂缝对缝洞沟通的影响机制,研究表明：缝洞型碳酸盐岩循缝找洞的储层改造思想是可行的；天然裂缝网络对压裂裂缝扩展轨迹具有重要影响，可充分利用天然裂缝的导向作用沟通溶洞；基于缝洞型碳酸盐岩循缝找洞储层改造思想，可实现不同方位溶洞储集体的沟通；现场应用证实缝洞型碳酸盐岩循缝找洞储层改造可大幅度提升井周油气动用程度，提高采收率。