CN110348031A - 水平井压裂近井筒裂缝扭曲形态数值模拟方法 - Google Patents

Abstract

公开了一种水平井压裂近井筒裂缝扭曲形态数值模拟方法。该方法可以包括：根据地层参数与施工参数，建立多层非均质地层模型；根据流体渗流方程、岩石变形方程、裂缝扭曲准则与裂缝面流体流动方程，建立流固耦合数值方程；根据流固耦合数值方程，在多层非均质地层模型内进行扩展有限元数值模拟，确定水平井压裂近井筒裂缝每一个时刻的扭曲形态。本发明通过计算裂缝宽度变化，分析施工参数等因素对裂缝扭曲变化的影响，根据裂缝扭曲变化程度和变化时间，提高水平井压裂优化设计的针对性，为现场水平井压裂优化设计提供理论依据。

Description

水平井压裂近井筒裂缝扭曲形态数值模拟方法

技术领域

本发明涉及油气田开发领域，更具体地，涉及一种水平井压裂近井筒裂缝扭曲形态数值模拟方法。

背景技术

目前对于页岩气、致密气等非常规油气藏的开发大多采用水平井分段压裂技术，水平井压裂与直井压裂不同之处在于近井筒区域可形成裂缝扭曲，理论计算和物模实验验证了这种裂缝扭曲现象可能性。Zhang et al.(2005)应用拉格朗日数值计算方法模拟了近井筒附近三维裂缝裂缝转向问题，并与实验结果进行了对比，但他们的模型没有考虑裂缝及孔隙内流体流动与岩石变形间的耦合关系；Lecampion et al.(2013)分析了水平井压裂近井筒附近横向和纵向裂缝间竞争关系，分析结果主要与地层特征、施工参数、射孔参数及应力场有关，该模型的主要局限是无法实时显示裂缝扩展及扭曲形态，不能直观观察到裂缝变化过程；Sherman et al.(2015)应用三维有限元模型模拟了水平井近井筒附近裂缝扩展与转向，结果显示了近井筒附近裂缝复杂程度，记录了与这种复杂程度相一致的压力关系，但模型建立在单层均质储层基础上，且没有考虑施工参数、应力差等因素对裂缝复杂程度影响，没有量化裂缝从横向裂缝扭曲为纵向裂缝时的宽度变化。因此，有必要开发一种水平井压裂近井筒裂缝扭曲形态数值模拟方法。

公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

发明内容

本发明提出了一种水平井压裂近井筒裂缝扭曲形态数值模拟方法，其能够通过计算裂缝宽度变化，分析施工参数等因素对裂缝扭曲变化的影响，根据裂缝扭曲变化程度和变化时间，提高水平井压裂优化设计的针对性，为现场水平井压裂优化设计提供理论依据。

所述方法可以包括：根据地层参数与施工参数，建立多层非均质地层模型；根据流体渗流方程、岩石变形方程、裂缝扭曲准则与裂缝面流体流动方程，建立流固耦合数值方程；根据所述流固耦合数值方程，在所述多层非均质地层模型内进行扩展有限元数值模拟，确定水平井压裂近井筒裂缝每一个时刻的扭曲形态。

优选地，所述地层参数包括：初始地应力场、初始渗流场、初始孔隙度与裂缝面滤失系数、层理数量与层理厚度。

优选地，所述施工参数包括施工流体排量与施工流体粘度。

优选地，所述岩石变形方程为：

其中，σ_ij,e为ij平面的弹性应力，σ_ij，e0为ij平面的初始弹性应力，ε_ij,e为ij平面的弹性应变，Δ_ij表示在ij平面前一个参数的变化量，ε_kk,e为垂直于K平面方向的弹性应变，G为弹性剪切模量，K为弹性体积模量，i、j代表i、j坐标方向，e表示弹性。

优选地，所述流体渗流方程为：

其中，β和M为Biot系数，k为岩石渗透率，γ为孔隙流体比重。

优选地，根据所述裂缝面流体流动方程，获得裂缝顶面法向滤失速率与裂缝底面法向滤失速率，进而建立流固耦合数值方程。

优选地，所述裂缝面流体流动方程为：

其中，v_t为裂缝顶面法向滤失速率，v_b为裂缝底面法向滤失速率，g_f为裂缝距，q为单位缝宽压裂液流体积。

优选地，所述裂缝顶面法向滤失速率为：

v_t＝l_t(p_f-p_t) (9)

其中，v_t为裂缝顶面法向滤失速率，l_t为滤失系数，p_t为裂缝顶面孔隙流体压力，p_f为压裂液压力；

所述裂缝底面法向滤失速率为：

v_b＝l_b(p_f-p_b) (10)

其中，v_b为裂缝底面法向滤失速率，l_b为滤失系数，p_b为裂缝底面孔隙流体压力，p_f为压裂液压力。

优选地，根据所述裂缝扭曲准则为裂缝扩展临界能量释放率准则。

优选地，所述裂缝扩展临界能量释放率准则为：

其中，G_S＝G_s+G_t，G_T＝G_n+G_S,为法向断裂临界应变能释放率；为两切向断裂临界能量释放率，B-K准则认为η为与材料本身特性有关的常数；G^C为复合型裂缝临界断裂能量释放率。

本发明的有益效果在于：分析水平井压裂近井筒附近裂缝空间扭曲面的形成过程，实时显示不同时间段裂缝扩展、扭曲形态，量化裂缝扭曲形成时间，计算裂缝宽度变化，分析施工参数等因素对裂缝扭曲变化的影响，根据裂缝扭曲变化程度和变化时间，提高水平井压裂优化设计的针对性，为现场水平井压裂优化设计提供理论依据。

本发明具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。

附图说明

通过结合附图对本发明示例性实施例进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施例中，相同的参考标号通常代表相同部件。

图1示出了根据本发明的水平井压裂近井筒裂缝扭曲形态数值模拟方法的步骤的流程图。

图2a与图2b分别示出了根据本发明的一个实施例的水平井压裂近井筒裂缝扭曲初始模型中水平井筒与初始纵向裂缝的示意图。

图3示出了根据本发明的一个实施例的裂缝扭曲开始形成的模拟图。

图4示出了根据本发明的一个实施例的裂缝扭曲完全形成的模拟图。

具体实施方式

下面将参照附图更详细地描述本发明。虽然附图中显示了本发明的优选实施例，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施例所限制。相反，提供这些实施例是为了使本发明更加透彻和完整，并且能够将本发明的范围完整地传达给本领域的技术人员。

图1示出了根据本发明的水平井压裂近井筒裂缝扭曲形态数值模拟方法的步骤的流程图。

在该实施例中，根据本发明的水平井压裂近井筒裂缝扭曲模拟方法可以包括：步骤101，根据地层参数与施工参数，建立多层非均质地层模型；步骤102，根据流体渗流方程、岩石变形方程、裂缝扭曲准则与裂缝面流体流动方程，建立流固耦合数值方程；步骤103，根据流固耦合数值方程，在多层非均质地层模型内进行扩展有限元数值模拟，确定水平井压裂近井筒裂缝每一个时刻的扭曲形态。

在一个示例中，地层参数包括：初始地应力场、初始渗流场、初始孔隙度与裂缝面滤失系数、层理数量与层理厚度。

在一个示例中，施工参数包括施工流体排量与施工流体粘度。

在一个示例中，岩石变形方程为：

其中，σ_ij,e为ij平面的弹性应力，σ_ij,e0为ij平面的初始弹性应力，ε_ij,e为ij平面的弹性应变，Δ_ij表示在ij平面前一个参数的变化量，ε_kk,e为垂直于K平面方向的弹性应变，G为弹性剪切模量，K为弹性体积模量，i、j代表i、j坐标方向，e表示弹性。

在一个示例中，流体渗流方程为：

其中，β和M为Biot系数，k为岩石渗透率，γ为孔隙流体比重。

在一个示例中，根据裂缝面流体流动方程，获得裂缝顶面法向滤失速率与裂缝底面法向滤失速率，进而建立流固耦合数值方程。

在一个示例中，裂缝面流体流动方程为：

其中，v_t为裂缝顶面法向滤失速率，v_b为裂缝底面法向滤失速率，g_f为裂缝距，q为单位缝宽压裂液流体积。

在一个示例中，裂缝顶面法向滤失速率为：

v_t＝l_t(p_f-p_t) (9)

其中，v_t为裂缝顶面法向滤失速率，l_t为滤失系数，p_t为裂缝顶面孔隙流体压力，p_f为压裂液压力；

裂缝底面法向滤失速率为：

v_b＝l_b(p_f-p_b) (10)

其中，v_b为裂缝底面法向滤失速率，l_b为滤失系数，p_b为裂缝底面孔隙流体压力，p_f为压裂液压力。

在一个示例中，裂缝扭曲准则为裂缝扩展临界能量释放率准则。

在一个示例中，裂缝扩展临界能量释放率准则为：

其中，G_S＝G_s+G_t，G_T＝G_n+G_S,为法向断裂临界应变能释放率；为两切向断裂临界能量释放率，B-K准则认为η为与材料本身特性有关的常数；G^C为复合型裂缝临界断裂能量释放率。

具体地，根据地层参数与施工参数，建立页岩地层模型，其中，地层参数包括初始地应力场、初始渗流场、初始孔隙度与裂缝面滤失系数、层理数量与层理厚度，其中，施工参数包括施工流体排量与施工流体粘度。构建一个内嵌有水平井筒的孔隙弹性地层几何模型，从水平井筒起裂出纵向裂缝。初始化后，裂缝可在模型域中的任意位置扩展。模型足够大，可最大限度地减小边界效应。模型长250m，宽180m，高60m，储层厚度为25m。

图2a与图2b分别示出了根据本发明的一个实施例的水平井压裂近井筒裂缝扭曲初始模型中水平井筒与初始纵向裂缝的示意图。

模拟地层与线性耦合的孔隙流体扩散/应力元素相结合，并为XFEM模型富集了一组连续介质，裂缝转化将可能在模拟过程中发生。为提高建模精度，对水平井筒附近单元格进行加密。三维几何模型如图2所示。模拟过程中，压裂液为牛顿流体，为粘度/存储占主导的传播介质。假设初始孔隙介质饱和应力和孔隙流体压力。位移与所有边界和对称表面垂直，并在内部角虚节点限制移动。孔隙流体压力在模型边界和孔隙介质内部保持一致。以恒定速率注入压裂液，使压裂液流进入已存在的纵向裂缝中，到达富集元素伪边界节点处，进一步求解压裂液流动方程。

根据流体渗流方程、岩石变形方程、裂缝扭曲准则与裂缝面流体流动方程，建立流固耦合数值方程。

(1)岩石形变

假设孔隙岩石为均质、各向同性孔隙弹性材料，开展准静态压裂，低应力条件下的初始本构方程可表示为：

式中，σ_ij为ij平面的应力，σ_ij,0为ij平面的初始应力，ε_ij为ij平面的应变，ε_kk为垂直于K平面方向应变，p为孔隙压力，p₀为初始孔隙压力，β为Biot系数。对于全饱和多孔介质，Terzaghi有效应力可表示为：

σ_ij,e＝σ_ij+pΔ_ij (2)，

有效应变可表示为：

其中，α为Biot系数，因此，本构方程可简化为公式(4)，即为岩石形变方程。

(2)孔隙流体流动

对于小体积应变，孔隙流体连续性方程可表示为：

其中，v_kk为流体渗流速率，β和M为Biot系数。根据达西定律，流体通过相互连接的孔隙网络时，渗流速率v_k可表示为：

其中，k为岩石渗透率，μ为孔隙流体粘度，γ为孔隙流体比重，p'_k为k平面孔隙压力张量。由此得到孔隙扩散方程即为流体渗流方程为公式(7)。

(3)压裂液流动

根据雷诺润滑理论，裂缝面流体流动方程为裂缝内压裂液连续性方程可表示为(8)，根据裂缝面流体流动方程，获得裂缝顶面法向滤失速率与裂缝底面法向滤失速率，分别为公式(9)与公式(10)，p_t、p_b与p_f之间的压差导致压裂液滤失。

(4)裂缝起裂和扩展

裂缝扭曲准则是裂缝起裂和扩展准则，即作用于裂缝面上的最大拉伸有效应力大于岩石的抗拉强度时，岩石起裂，随之的扩展方向沿着目前的最大主应力方向。应用复合型裂缝扩展B-K准则，即裂缝扩展临界能量释放率准则为公式(11)。

流固耦合数值方程即包括上述公式(4)、公式(7)、公式(8)、公式(9)、公式(10)与公式(11)，根据流固耦合数值方程，在多层非均质地层模型内进行扩展有限元数值模拟，确定水平井压裂近井筒裂缝每一个时刻的扭曲形态。依据线性软化粘性牵引分离定律，将压裂过程模拟为零厚度界面中强度的逐渐丧失，其方向和范围不是预先定义，而是在模拟过程中自动计算。在牵引分离定律中，内聚能量E等于线性软化面积和内聚强度的乘积。裂缝开启前，内聚能量与初始刚度呈线性关系。从裂缝开启到扩展，界面强度从I₀降为0，界面自由张开，裂缝距为g_f0。如果在完全破坏前卸载，随着损坏刚度Mp的下降，界面牵引呈线性下降。由于水力裂缝扭曲过程中同时受到张力和剪切力，应用一种基于能量的混合模式损伤演化定律研究裂缝扩展过程，即B-K裂缝准则，建立水平井压裂裂缝扭曲模型，根据设定的控制方程、初始和边界条件，运行计算有限元数值程序，随着流体的不断注入，逐渐达到裂缝起裂和扩展条件，这时，水平井附近岩石多个不同位置处点破裂，会形成多个分支缝起裂与止裂，它们的转向半径及方向不同，造成主裂缝扭曲，横向裂缝逐渐形成，这一过程会根据时间的不同在软件计算工具中实时呈现。

本发明通过计算裂缝宽度变化，分析施工参数等因素对裂缝扭曲变化的影响，根据裂缝扭曲变化程度和变化时间，提高水平井压裂优化设计的针对性，为现场水平井压裂优化设计提供理论依据。

应用示例

为便于理解本发明实施例的方案及其效果，以下给出一个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。

XX区块水平井压裂近井筒裂缝扭曲现象有待评价，模拟过程中储层及边界孔隙弹性性能参数如表1所示。

表1

初始、边界条件参数如表2所示。

表2

初始条件	储层	边界层
			初始最大水平主应力(MPa)	35	38
初始最小水平主应力(MPa)	25	28
			孔隙流体初始饱和度	1	1
孔隙流体初始压力(MPa)	20	20

图3示出了根据本发明的一个实施例的裂缝扭曲开始形成的模拟图。

图4示出了根据本发明的一个实施例的裂缝扭曲完全形成的模拟图。

利用本发明进行模拟，模拟结果显示，水平井纵向裂缝扩展过程中，随着不同位置处点破裂，会形成多个分支缝起裂与止裂，它们的转向半径及方向不同，造成主裂缝扭曲，横向裂缝逐渐形成，不断形成的不同位置处的破裂点、分支缝造成压力曲线的波动，裂缝扩展压力随排量和粘度的增加而增加；在缝高不同位置处，因扭曲形成了缝宽的交替变化，裂缝扭曲开始时间为60s-90s，如图3所示，水平井近井筒纵向裂缝转变为横向缝时间大概为240s，如图4所示；在裂缝扭曲的形成过程，裂缝宽度波动程度较大，平均0.5mm，增加了施工摩阻；降低粘度和排量利于横向缝的形成，根据扭曲发生时间，可调整段塞浓度和数量。

综上所述，本发明通过计算裂缝宽度变化，分析施工参数等因素对裂缝扭曲变化的影响，根据裂缝扭曲变化程度和变化时间，提高水平井压裂优化设计的针对性，为现场水平井压裂优化设计提供理论依据。

本领域技术人员应理解，上面对本发明的实施例的描述的目的仅为了示例性地说明本发明的实施例的有益效果，并不意在将本发明的实施例限制于所给出的任何示例。

以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。

Claims (10)

1.一种水平井压裂近井筒裂缝扭曲形态数值模拟方法，包括：

根据地层参数与施工参数，建立多层非均质地层模型；

根据流体渗流方程、岩石变形方程、裂缝扭曲准则与裂缝面流体流动方程，建立流固耦合数值方程；

根据所述流固耦合数值方程，在所述多层非均质地层模型内进行扩展有限元数值模拟，确定水平井压裂近井筒裂缝每一个时刻的扭曲形态。

2.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述地层参数包括：初始地应力场、初始渗流场、初始孔隙度与裂缝面滤失系数、层理数量与层理厚度。

3.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述施工参数包括施工流体排量与施工流体粘度。

4.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述岩石变形方程为：

其中，σ_ij,e为ij平面的弹性应力，σ_ij,e0为ij平面的初始弹性应力，ε_ij,e为ij平面的弹性应变，Δ_ij表示在ij平面前一个参数的变化量，ε_kk,e为垂直于K平面方向的弹性应变，G为弹性剪切模量，K为弹性体积模量，i、j代表i、j坐标方向，e表示弹性。

5.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述流体渗流方程为：

其中，β和M为Biot系数，k为岩石渗透率，γ为孔隙流体比重，p'_kk为垂直于K平面的孔隙压力张量，p'为孔隙压力张量，ε'_kk为垂直于K平面方向应变张量。

6.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，根据所述裂缝面流体流动方程，获得裂缝顶面法向滤失速率与裂缝底面法向滤失速率，进而建立流固耦合数值方程。

7.根据权利要求6所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述裂缝面流体流动方程为：

其中，v_t为裂缝顶面法向滤失速率，v_b为裂缝底面法向滤失速率，g_f为裂缝距，q为单位缝宽压裂液流体积。

8.根据权利要求6所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述裂缝顶面法向滤失速率为：

v_t＝l_t(p_f-p_t) (9)

其中，v_t为裂缝顶面法向滤失速率，l_t为滤失系数，p_t为裂缝顶面孔隙流体压力，p_f为压裂液压力；

所述裂缝底面法向滤失速率为：

v_b＝l_b(p_f-p_b) (10)

其中，v_b为裂缝底面法向滤失速率，l_b为滤失系数，p_b为裂缝底面孔隙流体压力，p_f为压裂液压力。

9.根据权利要求1所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述裂缝扭曲准则为裂缝扩展临界能量释放率准则。

10.根据权利要求9所述的水平井压裂近井筒裂缝扭曲形态数值模拟方法，其中，所述裂缝扩展临界能量释放率准则为：

其中，G_S＝G_s+G_t，G_T＝G_n+G_S,为法向断裂临界应变能释放率；为两切向断裂临界能量释放率，B-K准则认为η为与材料本身特性有关的常数；G^C为复合型裂缝临界断裂能量释放率。