CN116401897B - 一种基于近似解和能量方程的半解析裂缝扩展模拟方法 - Google Patents
一种基于近似解和能量方程的半解析裂缝扩展模拟方法 Download PDFInfo
- Publication number
- CN116401897B CN116401897B CN202310668991.XA CN202310668991A CN116401897B CN 116401897 B CN116401897 B CN 116401897B CN 202310668991 A CN202310668991 A CN 202310668991A CN 116401897 B CN116401897 B CN 116401897B
- Authority
- CN
- China
- Prior art keywords
- crack
- fracture
- fracturing
- equation
- time
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000004088 simulation Methods 0.000 title claims abstract description 64
- 238000000034 method Methods 0.000 title claims abstract description 58
- 239000000243 solution Substances 0.000 claims abstract description 62
- 239000011435 rock Substances 0.000 claims abstract description 39
- 238000004364 calculation method Methods 0.000 claims abstract description 26
- 238000006073 displacement reaction Methods 0.000 claims abstract description 24
- 230000000704 physical effect Effects 0.000 claims abstract description 18
- 238000002347 injection Methods 0.000 claims abstract description 17
- 239000007924 injection Substances 0.000 claims abstract description 17
- 230000011218 segmentation Effects 0.000 claims abstract description 7
- 239000012530 fluid Substances 0.000 claims description 49
- 230000014509 gene expression Effects 0.000 claims description 39
- 230000003993 interaction Effects 0.000 claims description 15
- 239000007788 liquid Substances 0.000 claims description 13
- 238000005461 lubrication Methods 0.000 claims description 10
- 238000011065 in-situ storage Methods 0.000 claims description 9
- 230000035699 permeability Effects 0.000 claims description 9
- 238000005086 pumping Methods 0.000 claims description 8
- 239000007787 solid Substances 0.000 claims description 8
- 239000000463 material Substances 0.000 claims description 7
- 230000010354 integration Effects 0.000 claims description 5
- 238000005553 drilling Methods 0.000 claims description 4
- 230000004907 flux Effects 0.000 claims description 3
- 230000005484 gravity Effects 0.000 claims description 3
- 239000004576 sand Substances 0.000 claims description 3
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims description 3
- 238000005315 distribution function Methods 0.000 claims 1
- 230000036962 time dependent Effects 0.000 claims 1
- 238000011161 development Methods 0.000 abstract description 8
- 230000018109 developmental process Effects 0.000 abstract description 8
- 238000012795 verification Methods 0.000 description 15
- 230000006870 function Effects 0.000 description 13
- 238000012544 monitoring process Methods 0.000 description 13
- 238000010168 coupling process Methods 0.000 description 11
- 230000008878 coupling Effects 0.000 description 10
- 238000005859 coupling reaction Methods 0.000 description 10
- 230000000694 effects Effects 0.000 description 9
- 230000008569 process Effects 0.000 description 6
- 238000010586 diagram Methods 0.000 description 5
- 238000009826 distribution Methods 0.000 description 5
- 239000000700 radioactive tracer Substances 0.000 description 5
- 238000006243 chemical reaction Methods 0.000 description 4
- 238000004134 energy conservation Methods 0.000 description 4
- 239000012088 reference solution Substances 0.000 description 4
- 230000008901 benefit Effects 0.000 description 3
- 230000015572 biosynthetic process Effects 0.000 description 3
- 238000012937 correction Methods 0.000 description 3
- 238000010790 dilution Methods 0.000 description 3
- 239000012895 dilution Substances 0.000 description 3
- 238000011439 discrete element method Methods 0.000 description 3
- 238000004519 manufacturing process Methods 0.000 description 3
- 239000011159 matrix material Substances 0.000 description 3
- 230000008859 change Effects 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 238000013461 design Methods 0.000 description 2
- 230000005489 elastic deformation Effects 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000007480 spreading Effects 0.000 description 2
- 238000003892 spreading Methods 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000010200 validation analysis Methods 0.000 description 2
- 206010010356 Congenital anomaly Diseases 0.000 description 1
- 101100023016 Methanothermobacter marburgensis (strain ATCC BAA-927 / DSM 2133 / JCM 14651 / NBRC 100331 / OCM 82 / Marburg) mat gene Proteins 0.000 description 1
- 241000009298 Trigla lyra Species 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 230000002596 correlated effect Effects 0.000 description 1
- 238000005336 cracking Methods 0.000 description 1
- 230000006837 decompression Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 229920001971 elastomer Polymers 0.000 description 1
- 239000000806 elastomer Substances 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 101150095438 metK gene Proteins 0.000 description 1
- 238000012821 model calculation Methods 0.000 description 1
- 239000011148 porous material Substances 0.000 description 1
- 238000005381 potential energy Methods 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 238000010008 shearing Methods 0.000 description 1
- 238000003860 storage Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
- Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping
Abstract
本发明公开了一种基于近似解和能量方程的半解析裂缝扩展模拟方法，涉及油气藏开发技术领域；该方法包括以下步骤：S1:获得目标储层的储层物性、岩石力学参数；获得单井的分段分簇参数和每一压裂段的注排量参数；S2:建立多裂缝扩展基础模型；S3:根据储层物性、岩石力学参数、单井的压裂段分段分簇参数和注排量参数求解多裂缝扩展基础模型，得到压裂段中裂缝参数；S4:将裂缝参数带入能量平衡方程和流量平衡方程，进行收敛性判断；若满足收敛条件，得到解值，计算下一步；若不满足收敛条件，则循环步骤S3~S4，直至收敛，获得压裂段的解值；S5:重复步骤S3~S4，得到单井所有压裂段的模拟结果。该方法的精度高、计算效率高。
Description
技术领域
本发明涉及油气藏开发技术领域，具体涉及一种基于近似解和能量方程的半解析裂缝扩展模拟方法。
背景技术
页岩储层孔隙度、渗透率极低，给页岩气的经济高效开发带来了极大的困难和挑战，长水平井段钻井和多段大排量水力压裂施工时页岩气开发的关键和核心技术，能最大程度地增加压裂裂缝地改造体积和表面积，最终达到提高产量和采收率的目的。页岩储层脆性大、天然裂缝和水平层理发育，压裂过程中容易发生剪切滑移和张性破坏，压裂裂缝不再是单一对称的两翼缝，可能形成复杂的网状裂缝，给页岩水力压裂设计、裂缝监测及解释、压后产能预测等带来诸多不便。因此，压裂裂缝的展布和裂缝扩展对页岩压裂设计施工、裂缝监测、产能评价至关重要。
目前常见的裂缝扩展数值模拟方法都在不同程度上存在模拟精度不高或者运算效率较低的缺点，在模拟大尺度多裂缝竞争扩展问题上具有较明显的局限性。如现有的DEM模型，其可以用于非连续面发育的岩体的裂缝展布和扩展，对于在模拟非连续面发育岩体的变形与破坏时具有先天优势，但其仅能间接表征材料的宏观力学参数，存在额外计算误差，因此精度不高，应用范围有局限；BEM模型虽然降低了问题维数，大大降低了计算规模和求解时间，但其难以考虑基质孔隙流体渗流与岩石骨架变形的相互耦合，具有明显的局限性。
发明内容
本发明的目的在于提供一种将无量纲近似解与能量守恒相结合的流固耦合新的基于近似解和能量方程的半解析裂缝扩展模拟方法，该方法在保证计算精度的同时还能提高运算速度。
本发明通过下述技术方案实现：本发明提供一种基于近似解和能量方程的半解析裂缝扩展模拟方法，包括以下步骤：
S1:根据目标储层的测井数据、录井数据获得所述目标储层的储层物性、岩石力学参数；根据单井的钻井资料获得所述单井的分段分簇参数和每一压裂段的注排量参数；
S2:建立多裂缝扩展基础模型；
S3: 根据所述储层物性、岩石力学参数、所述单井的压裂段分段分簇参数和注排量参数求解所述多裂缝扩展基础模型，得到所述压裂段中的包含时空信息的近似解：裂缝开度，裂缝长度和裂缝内流体压力；
S4:将初始值带入能量平衡方程和流量平衡方程求解得到下一迭代包含时空信息的数值解：裂缝开度，裂缝长度和裂缝内流体压力。然后根据迭代前后的流量差进行收敛性判断；若迭代前后的流量绝对差小于设定的阈值，则满足收敛条件并结束迭代。否则循环步骤S3~S4，直至收敛；获得所述压裂段的数值解作为下一时间步的初始值，进入下一时间步的计算。
S5: 重复步骤S3~S4，直至达到所需的总泵送时间t(k)=T；获得所述压裂段的解值；
S6:重复所述步骤S2~S5，得到所述单井的所有压裂段的解值，即为单井全井段模拟结果。
进一步地，在本申请地一些实施例中，所述储层物性包括目标储层的孔隙度、渗透率、流体和固体压缩系数；所述岩石力学参数包括岩石的弹性系数、泊松比和断裂韧性；所述分段分簇参数包括所述单井的段间距、簇间距、射孔长度、孔眼直径、孔密度、单段孔数；所述注排量参数包括每个所述压裂段的泵注总液量、排量、动态粘度和时间；所述泵注总液量包括净液量和混砂液量。
进一步地，在本申请地一些实施例中，所述多裂缝扩展基础模型包括联立的质量平衡方程、局部法向牵引力T和宽度w之间的弹性关系表达式、径向通量q i (r,t)与流体压力梯度的表达式、裂缝扩展准则、体积平衡方程、雷诺润滑方程、各裂缝入口压力平衡和射孔簇的流量总和表达式、各裂缝内的流量总和的表达式。
进一步地，在本申请地一些实施例中，其中，所述质量平衡方程为：
其中：q i 是缝内流速；C L 为滤失系数；k为岩石渗透性；c r 为储层压缩系数；φ为岩石
孔隙度；p 0 为储层压力；为与时间相关的函数，表示裂缝尖端到达坐标为r的点的时
间；
所述局部法向牵引力T和宽度w之间的弹性关系表达式为：
其中F为非局部积分算子；为裂缝在时间点t坐标r处的宽度；为
裂缝在时间点t的半径；E’为弹性模量；
径向通量与流体压力梯度的表达式为：
其中μ为动态粘度；
裂缝扩展准则的表达式为：
其中K I 表示I型裂缝的应力强度因子，K IC 表示I型裂缝的断裂韧性；
体积平衡方程为：
其中R w 为井眼半径；
雷诺润滑方程为：
；
各裂缝入口压力平衡和射孔簇的流量总和表达式为：
各裂缝内的流量总和的表达式为：
。
进一步地，在本申请地一些实施例中，在裂缝扩展准则的表达式中，对于径向裂缝，K I 表示为：
。
进一步地，在本申请地一些实施例中，所述体积平衡方程在缝口处的边界条件、裂
纹尖端的边界条件和裂纹尖端的初始条件下对进行积分得到；
其中在缝口处的边界条件为：
其中R w 为井眼半径；
裂纹尖端的边界条件表达式为：
；
裂纹尖端的初始条件为：
。
进一步地，在本申请地一些实施例中，所述解值包括：裂缝无量纲半径R i (t)/Z、缝内流量q i (R w ,t)、缝口宽度w i (R w ,t)和裂缝面积A(t)。
进一步地，在本申请地一些实施例中，所根据所述储层物性、岩石力学参数对所述
多裂缝扩展基础模型进行赋值中赋值的参数包括：E、v、K IC 、μ、Q、Z、σ min 、R w 、h i,j ，其中，E为杨
氏模量；为泊松比；K IC 为I 型裂缝的断裂韧性；μ为动态粘度；Z为压裂段总长度；Q为单个
分段总注入流量；σ min 为最小原位围压应力；R w 为井眼半径；h i,j 为裂缝i 和j 之间的距离。
进一步地，在本申请地一些实施例中，S3中，所述求解基于第一假设，所述第一假设包括：
（1）裂缝扩展遵循线性弹性断裂力学（LEFM），该力学假设材料在任何地方都遵循线性松弛应力-应变关系；当裂纹尖端的应力强度因子达到材料断裂韧性时，将发生裂纹扩展；
（2）裂缝内为不可压缩牛顿流体，并用润滑方程进行描述；
（3）岩石是不透水的，液体滤失忽略；
（4）所有裂缝均径向平行扩展；
（5）在弹性和流体流动方程中均忽略重力；
（6）流体前沿与裂缝前沿重合；
（7）远场原位应力是均匀且恒定的。
进一步地，在本申请地一些实施例中，所述求解包括以下步骤：
S31：分别设置力学参数以及计算的
初始时间、最终时间和时间步长的值；
S32：设定初始条件为，获得每条裂缝的缝内
压力，裂缝半径，宽度；
其中，；
；
相互作用力为：
；
其中，和通过下式求解：
；
S33：使用牛顿法对方程组进行数值求解。
本发明与现有技术相比，具有以下优点及有益效果：
本发明提供的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，该方法将能量平衡方程和修正逼近函数结合，在考虑流固作用之外，开创性得将能量守恒引入了多裂缝系统加以考虑，从而建立了兼顾应力耦合和能量耦合的多簇裂缝半解析模拟方法。该模拟方法的计算效率相对边界元模型高出103-5倍，且效率优势会随着模拟对象的尺度增加而快速上升。
附图说明
图1为本申请提供的裂缝的几何模型示意图；
图2为本申请中单裂缝验证中近似解和参考解的裂缝半径、井筒宽度(ρ=0)和效率随时间的变化曲线；
图3为本申请中多裂缝验证中段长120m的模拟算例裂缝间距设置示意图；
图4A~图4D为本申请中多裂缝验证中h1=h2=h3=h4=30m时C5Frac与ILSA II的对比曲线；
图5为本申请中多裂缝验证中段长20m的模拟算例裂缝间距设置示意图；
图6A~图6D为本申请中多裂缝验证中均匀间距下模拟输出参数对比图；
图7为本申请中不同模拟时间下裂缝扩展形态对比图；
图8为本申请中利用樊页井进行现场验证时樊页井的全井段各簇裂缝长度和开度；
图9为本申请中利用樊页井进行现场验证时微地震监测缝长折算原理示意图；
图10为本申请中利用樊页井进行现场验证时樊页井的模拟缝长与微地震折算缝长对比结果图；
图11为本申请中利用樊页井进行现场验证时樊页井的压裂第12段模拟缝长与微地震折算缝长对比图；
图12为本申请中利用樊页井进行现场验证时樊页井的模拟结果与示踪剂监测结果对比图；
图13为本申请中利用樊页井进行现场验证时樊页井的井身结构图。
具体实施方式
下面将结合实施例对本申请的技术方案进行清楚、完整地描述，显然，所描述的实施例是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本申请保护的范围。
本技术领域技术人员可以理解，除非另外定义，这里使用的所有术语（包括技术术语和科学术语），具有与本申请所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语，应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样被特定定义，否则不会用理想化或过于正式的含义来解释。
本技术领域技术人员可以理解，除非特意声明，这里使用的单数形式 “一”、“一个”、“所述”和“该”也可包括复数形式。应该进一步理解的是，本申请的说明书中使用的措辞“包括”是指存在所述特征、整数、步骤、操作、元件和/或组件，但是并不排除存在或添加一个或多个其他特征、整数、步骤、操作、元件、组件和/或它们的组。这里使用的措辞“和/或”包括一个或更多个相关联的列出项的全部或任一单元和全部组合。
很多科学问题具有高度的复杂性，无法得到解析解，而数值方法通常会消耗大量的计算资源。比如有限元和离散元方法，保证精度的同时不能兼顾计算效率，难以满足实际生产中对大规模计算的需求。边界元降低了问题的维数，从而降低了计算规模和求解时间。然而采用边界元方法的高精度模型，即使依托高性能工作站，模拟多簇裂缝扩展，仍需要几天至几个星期的运算时间。目前国内外针对复杂多裂缝扩展主要采用的是边界元BEM、离散元DEM、有限元FEM方法、扩展有限元XFEM。相较于其它三种方法，边界元降低了问题的维数，从而降低了计算规模和求解时间。然而即使依托高性能工作站，使用边界元方法在井周尺度模拟多簇裂缝扩展仍需要以天为单位的运算时间以满足精度要求；且边界元公式通常会产生完全填充的矩阵，这意味着存储需求和计算时间会随着模拟对象大小的平方而增长。随着裂缝边界在三维空间的扩展，运算效率会急剧下降，这种特性显然限制了边界元模型在实际生产中的应用。
为了解决这一技术问题，发明人在本申请中提供了一种基于近似解和能量方程的半解析裂缝扩展模拟方法，该方法将能量平衡方程和修正逼近函数结合，在考虑流固作用之外，开创性得将能量守恒引入了多裂缝系统加以考虑，从而建立了兼顾应力耦合和能量耦合的多簇裂缝半解析模拟方法。
其中，在本申请中的逼近函数是由空间函数和时间函数组成的表达式，可以通过正交多项式求解控制方程组得到，该函数不需要迭代运算，因此具有极高的计算效率。此外，发明人在逼近函数的基础上还提出了复合粘度概念和能量等效方法，建立了多裂缝修正逼近函数，综合考虑了应力阴影、地层滤失和非稳定流量对裂缝扩展的影响，改善了逼近函数只能适用于入射流量稳定的单裂缝中使用的应用限制，使其得到模拟方法在保持了计算效率的同时，扩大了应用范围。
该一种基于近似解和能量方程的半解析裂缝扩展模拟方法，具体包括以下步骤：
S1:根据目标储层的测井数据、录井数据获得所述目标储层的储层物性、岩石力学参数；根据单井的钻井资料获得所述单井的分段分簇参数和每一压裂段的注排量参数；
S2:建立多裂缝扩展基础模型；
S3: 根据所述储层物性、岩石力学参数、所述单井的压裂段分段分簇参数和注排量参数求解所述多裂缝扩展基础模型，得到所述压裂段中的包含时空信息的近似解：裂缝开度，裂缝长度和裂缝内流体压力；
S4:将初始值带入能量平衡方程和流量平衡方程求解得到下一迭代包含时空信息的数值解：裂缝开度，裂缝长度和裂缝内流体压力。然后根据迭代前后的流量差进行收敛性判断；若迭代前后的流量绝对差小于设定的阈值，则满足收敛条件并结束迭代。否则循环步骤S3~S4，直至收敛；获得所述压裂段的解值作为下一时间步的初始值，进入下一时间步的计算。
S5: 重复步骤S3~S4，直至达到所需的总泵送时间t(k)=T；获得所述压裂段的解值；
S6:重复所述步骤S2~S5，得到所述单井的所有压裂段的解值，即为单井全井段模拟结果。
其中，根据所述储层物性、岩石力学参数、所述单井的压裂段分段分簇参数和注排量参数求解所述多裂缝扩展基础模型，得到所述压裂段中的裂缝参数，包括：根据所述储层物性、岩石力学参数、所述单井的压裂段分段分簇参数和注排量参数进行赋值，其具体包括：根据所述储层物性、岩石力学参数对所述多裂缝扩展基础模型进行赋值，并根据流量平衡方程设置初始条件；按照所述单井的压裂段分段分簇参数和注排量参数，按照时间步长为t，步长Δt求解所述多裂缝扩展基础模型，得到所述压裂段中的裂缝参数。
其中，所述储层物性包括目标储层的孔隙度、渗透率、流体和固体压缩系数；所述岩石力学参数包括岩石的弹性系数、泊松比和断裂韧性；
所述分段分簇参数包括所述单井的段间距、簇间距、射孔长度、孔眼直径、孔密度、单段孔数；所述注排量参数包括每个所述压裂段的泵注总液量、排量、动态粘度和时间；所述泵注总液量包括净液量和混砂液量。
在单个压裂段中通常存在N条裂缝，其具体模型示意如图1所示。在压裂段中，压裂段总长度为Z，各裂缝的距离依次为h1，h2……hn-1。因此有：
（1）
在向井筒泵注压裂液后，N条裂缝在流体驱动和裂缝之间弹性应力的作用下发生竞争扩展。整个模型要求解的变量包括：
（1）裂缝在时间点t坐标r处的宽度；
（2）裂缝在时间点t的半径；
（3）裂缝在时间点t的缝口流量；
（4）缝内净压力；
（5）来自其他裂缝的弹性相互作用力。
在S3中，所述求解基于第一假设，所述第一假设包括：
（1）裂缝扩展遵循线性弹性断裂力学（LEFM），该力学假设材料在任何地方都遵循线性松弛应力-应变关系；当裂纹尖端的应力强度因子达到材料断裂韧性时，将发生裂纹扩展；
（2）裂缝内为不可压缩牛顿流体，并用润滑方程进行描述；
（3）岩石是不透水的，液体滤失忽略；
（4）所有裂缝均径向平行扩展；
（5）在弹性和流体流动方程中均忽略重力；
（6）流体前沿与裂缝前沿重合；
（7）远场原位应力是均匀且恒定的。
需要说明的是，为了求解上述变量，本申请建立的模型需要将裂缝的弹性变形、流体质量平衡、缝内流动和扩展准则相结合，并添加各裂缝的相互作用应力、井筒内的压力和流体连续性条件联立求解，即上述变量的求解需要联立裂缝的弹性变形表达式、流量质量平衡方程、雷诺润滑方程。首先，基于不可压缩流体和不可渗透岩石的假设，模型的质量平衡方程为：
（2）
其中，q i 是缝内流速；CL为滤失系数；k为岩石渗透性；c r 为储层压缩系数；φ为岩石
孔隙度； p 0 为储层压力；为与时间相关的函数，表示裂缝尖端到达坐标为r的点的
时间。
弹性体被认为是通过作用在每个裂缝表面上的牵引力Ti而发生变形的。在圆形裂缝的假设下，局部法向牵引力T和宽度w之间的弹性关系为：
（3）
其中，F为非局部积分算子，Ti为局部法向牵引力。
此外，由于径向通量与流体压力梯度成比例，因此，可以得到：
（4）
其中，μ为动态粘度。
根据第一假设，裂缝总是在极限平衡状态下扩展，因此裂缝扩展准则可表达为：
K I = K IC （5）
其中K I 表示I 型裂缝的应力强度因子， K IC 表示I 型裂缝的断裂韧性。对于径向裂缝，K I 可以表示为：
（6）
此外，每条裂缝会在各自的缝内流量下发生扩展。因此，从质量平衡角度出发，缝口处的边界条件由下式给出：
（7）
其中R w 为井眼半径。
裂纹尖端的边界条件为：
（8）
初始条件为：
（9）
在初始和边界条件下，可以得到体积平衡方程，即为所述流量质量平衡方程：
（10）
所述雷诺润滑方程表示为：
（11）
由于N 条裂缝中的流量由同一井筒中流进，因此各裂缝具有相同的入口压力：
（12）
同时，由于各条裂缝的流量均来自同一井筒，因此各裂缝内的流量总和等于总泵注排量：
（13）
因此，所述多裂缝扩展基础模型为包括联立的质量平衡方程、局部法向牵引力T和
宽度w之间的弹性关系表达式、径向通量与流体压力梯度的表达式、裂缝扩展准
则、断裂韧性表达式、体积平衡方程、雷诺润滑方程、各裂缝入口压力平衡和射孔簇的流量
总和表达式、各裂缝内的流量总和的表达式，即上述式（2）、式（3）、式（4）、式（5）、式（6）、式
（10）、式（11）、式（12）、式（13）。
其中，所述裂缝参数为裂缝在时间点t的缝口流量Q i (t)，裂缝在时间点t坐标r处的宽度w i (r,t)、裂缝在时间点t的半径R i (t)、裂内净压力p f (r,t)、来自其他裂缝的弹性相互作用力σ i (r,t)。
在一些实施例中，所述解值包括：裂缝无量纲半径R i (t)/Z、缝内流量q i (R w ,t)、缝口宽度w i (R w ,t)和裂缝面积A(t)。
在一些实施例中，所根据所述储层物性、岩石力学参数对所述多裂缝扩展基础模型进行赋值中赋值的参数包括：E、v、K IC 、μ、Q、Z、σ min 、R w 、h i,j ，其中，E为弹性稀释；v为泊松比；K IC 为I 型裂缝的断裂韧性；μ为动态粘度；Z为压裂段总长度；Q为单个分段总注入流量；σ min 为地层最小地应力；R w 为井眼半径；h i,j 为裂缝i 和j 之间的距离。
其中能量平衡方程将裂缝发育的物理过程以能量表征，从而在能量守恒的约束下实现多裂缝系统能量耦合，使用分部积分解决了流体压力在裂缝尖端的奇异问题，能量计算过程采用了不需要离散网格的显示积分，因此计算精度较高且不会因为模拟对象尺度变化而降低。
水力压裂的本质其实能量分化的过程，然而逼近函数像现有多数模拟方法一样，只考虑了流体与固体之间的应力耦合，并没有考虑能量作用。对此，发明人将裂缝发育所需能量被分成两大部分，一部分是通过流体对固体做功消耗掉了，这部分主要由岩石因流体填充的裂缝体积增大而发生的弹性势能变化，岩石破裂消耗能量和为抵消最小地应力及压力阴影做功所消耗能量组成。另一部分是流体耗散能量，由流体在裂缝中流动消耗的能量组成。为了研究地层流体渗流和限流对多裂缝扩展的影响，申请人后续将流体滤失进入地层消耗的能量和流体经过射孔时损失的能量公式化并计入了流体耗散能量，在能量守恒的约束下，与输入总能量建立了能量平衡方程。由于裂缝扩展过程中流体压力、裂缝长度和宽度等复杂时空变量无法用简单表达式准确描述，因此为了解决这一技术问题，发明人提出了地层滤失-压力阴影耦合裂缝扩展逼近函数近似替代裂缝扩展数值解，适用于任意岩石韧度和流体粘度组合。
其方程变量的近似解具体获得方式如下：
（1）压力分布的近似解
通过将流体流动控制方程进行有限差分离散，可消除在每个时间步骤中计算分布的需要，极大降低计算强度。流体压力的形式表示为：
（14）
其中，无量纲压力分布为：
（15）
（2）断裂韧性的近似解
裂缝扩展条件由K I =K IC 给出，其中K I 由牵引力计算：
（16）
方程（15）的前两项对K I 的贡献为零，因此仍然有保持空间均匀的压力，必
须使用这个压力以确保始终满足K I =K IC 。通过替换， 被求解为：
（17）
（3）相互作用应力的近似解
通常，需要根据每个裂缝内部的压力分布计算相互作用力。然而，这种方法不能满足快速近似计算的要求。因此，在该模型中，使用均匀加压裂纹解获取相互作用应力的近似值，则相邻裂缝j在裂缝i上产生的应力的法向分量为：
（18）
其中：
（19）
其中是簇间距h j,i （裂缝i和j之间的距离）与裂缝半径Rj的比值，是径向位
置r与裂缝半径R i 比值，。P j,i 是均匀的净压力。因此，获取裂缝间相互作用力
近似值的关键是选择合适的净压力，以便该的能够最大程度上接近非均匀净压力下裂缝间
的相互作用力。本申请中使用的方法是在每个时间步长为各条裂缝选择某一均匀的净压力
值，在该净压力值下，可以得到非均匀内部净压力下相同的裂缝体积。因此，对于裂缝j，由
均匀内部净压力产生的裂缝体积为：
（20）
相互作用力模型通过对所有相邻裂缝的相互作用应力求和来实现。因此，施加在水力裂缝上的相互作用力近似为：
（21）
（4）全体积守恒
通过对受初始和裂缝尖端边界条件影响的局部体积平衡进行积分，得到整体体积平衡方程：
（22）
为快速计算流体滤失，没有明确地使用演化半径，而是对时间进行近似。但这并不
意味着裂缝被限制遵循幂律增长，相反它们的半径是耦合解的结果。通过对与接触时间相
关的积分的历史依赖性的近似，计算被替换为：
（23）
（5）整体能量平衡
由于井筒压力变化，很难准确估计井筒压力。因此，通过相反的近似来计算井筒的压力，以满足总能量平衡，将与韧性相关的能量添加到能量平衡中，这使得新模型对所有的体系都有效，包括所谓的“韧性主导”体系。更新后的功率平衡为：
（24）
上式左侧是第i个裂缝的能量输入速率，而上式右侧的变量分别为：
Dc是与岩石压裂有关的耗散率：
（25）
Df是与粘性流体流动相关的耗散率：
（26）
DL是与滤失相关的流体损失率：
（27）
是通过变形岩石应变能来增加·应变能的一部分，是可恢复的弹性能：
（28）
是由每个裂缝旁边缝的压缩应力所做的功：
（29）
是原地应力对裂缝的作用。为考虑流体损失，原位应力的功被修正为：
（30）
其中S i 用于量化非均匀的原位应力。
是通过缝口的功率损失，可通过经典压降方程计算：
（31）
通过解压力、宽度和半径对流量的隐性依赖表达式代替未知的和
，得到：
（32）
（33）
（34）
（35）/>
（36）
（37）
基于上述近似解和联立的方程组，该模型的求解步骤如下：
S31：分别设置力学参数以及计算的初
始时间、最终时间和时间步长的值；
S32：设定初始条件为，获得每条裂缝的缝内
压力，裂缝半径，宽度；其中，
；（38）
（39）/>
（40）
相互作用力为：
（41）
其中，和通过下式求解：
（42）；
S33：使用牛顿法对方程组进行数值求解。
然后基于上述计算值，对应力-应变耦合的局部裂缝宽度、净压力和半径进行评估。然后将 和 代入能量平衡方程，即可得到：
（43）
然后，检查初始条件的和之间的相对差异。如果该值低于给定公差
TOL，即：
其中，公差TOL又被作为判断是否收敛的阈值，若迭代前后的流量绝对差小于这个阈值，则判断迭代收敛，否则循环步骤S32~S33，直至收敛，获得所述压裂段的数值解。阈值的设定可根据实际现场作业精度要求变化，阈值越小迭代次数越多，算力消耗越大。 阈值越大，收敛越快，计算耗时越少。在通常 情况下，阈值我们一般设定为1%。
S34将S33步骤中最后一次迭代的流量值作为下一时间步的初始值带入S32得到裂
缝的初次迭代近似解，循环步骤S32~S33，直至达到所需的总泵送时间t(k)=T；得到单个分
段压裂液完全注入后的裂缝长度R(t)，宽度和裂缝内流体压力。
S35:重复所述步骤S32~S35，得到所述单井的所有压裂段的解，即为单井全井段模拟结果。
下面将结合实施例对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
为了验证本申请提供的复杂裂缝扩展模型，发明人分别验证了单裂缝和多裂缝示例。
其中，用作验证的单裂缝的参考解以Dontsov开发的便士形水力裂缝数值近似解计算得出。利用本申请提供的方法获得近似解的给定参数为E’=9.5GPa，KIC=107 Pa·m1/2，μ=0.001 Pa·s，Qo=0.01 m3/s的时间变化，滤失参数的选择为CL={10-5 m/s1/2, 10-4 m/s1/2,10-3 m/s1/2}，其比较结果如图2所示。
从图2可以知道，模型在裂缝半径、井筒宽度和效率方面的计算结果与参考解的匹配度均在1%以内，唯一的例外是CL值较大时的波动表现，但这里的波动表现是在参考解中，而不是在近似解中。因此，结果表明该方法能够准确地考虑滤失的影响。
用作验证的多裂缝示例的参考值利用隐式水平集解析模型（ILSA II）获得。
该模型是一种完全耦合的三维水力裂缝模拟器，弹性方程采用三维位移不连续法求解，流体流动采用有限体积法求解。如图3所示，建立的验证模型段长120m，包含5个均匀间隔的裂缝，簇间距为30m。裂缝的分布情况，裂缝面1至5的Z坐标分别为z1=0，z2=30，z3=60，z4=90和z5=120。下文中，裂缝1和裂缝5称为外部裂缝，裂缝2和裂缝4称为内部裂缝，裂缝3称为中间裂缝。
验证模型涉及的参数设定如表1所示：
表1
| 参数 | 取值 | 单位 |
| 弹性稀释 | 95 | GPa |
| 泊松比 | 0.2 | / |
| 断裂韧性 | 0 | MPa·m1/2 |
| 动态粘度 | 1 | Pa·s |
| 泵注排量 | 0.1 | m3/s |
| 段长 | 120 | m |
| 地层最小地应力 | 70 | MPa |
模拟输出结果包括：裂缝无量纲半径R i (t)/Z、缝内流量q i (R w ,t)、缝口宽度w i (R w ,t)和裂缝面积A(t)，其中裂缝面积定义如下：
其模拟输出结果如图4A~图4D所示，其中图4A为裂缝半径的模拟输出结果，图4B为裂缝流速的模拟输出结果；图4C为井筒宽度的模拟输出结果；图4D为裂缝面积的模拟输出结果，其中图4A~图4C中的外侧对应图3中的第1簇裂缝和第5簇裂缝；内测对应的是图3中的第2簇裂缝和第4簇裂缝；中间对应图3中的第3簇裂缝。
从模拟结果来看，建立的模型得到的裂缝半径、面积和缝口流量与隐式水平集解析模型析模型（ILSA II）非常接近。同时，还观察到明显的非常接近。同时，还观察到明显的应力阴影现象，其中外部裂缝优先应力阴影现象，其中外部裂缝优先扩展，而内部裂缝的内部裂缝的扩展受到抑制。外部裂缝的流量急剧增加至0.013m3/s，在模拟15s后整个泵注排量几乎全部分配给外部裂缝，同时流向其他裂缝的排量降低至接近0，这是因为内部裂缝必须在外部裂缝建立出的诱导应力场中相互竞争扩展，而内部裂缝的任何额外扩展都会增强诱导应力场。随着模拟时间的增加，外部裂缝的扩展变得更加明显，而内部裂缝的尺寸几乎不再发生变化。为了进一步验证模型的准确度，裂缝间距被缩小到5m以检验模型在不同段长下与基准模型ILSA II 的区别。如图5所示，建立的验证模型段长20m，包含5个均匀间隔射孔簇，簇间距为5m。
验证模型涉及的参数设定如表2所示：
表2
| 参数 | 取值 | 单位 |
| 弹性稀释 | 95 | GPa |
| 泊松比 | 0.2 | / |
| 断裂韧性 | 0 | MPa·m1/2 |
| 动态粘度 | 1 | Pa·s |
| 泵注排量 | 0.1 | m3/s |
| 段长 | 20 | m |
| 簇间距 | 5 | m |
模拟输出结果包括：裂缝无量纲半径R i (t)/Z、缝内流量q i (R w ,t)、缝口宽度w i (R w ,t)和裂缝面积A(t)，其中裂缝面积定义如下：
其模拟输出结果如图6A~图6D所示，其中图6A为裂缝半径的模拟输出结果，图6B为裂缝流速的模拟输出结果；图6C为井筒宽度的模拟输出结果；图6D为裂缝面积的模拟输出结果；其中图6A~图6C中的外侧对应图3中的第1簇裂缝和第5簇裂缝；内测对应的是图3中的第2簇裂缝和第4簇裂缝；中间对应图3中的第3簇裂缝。
从模拟结果来看，建立的模型得到的裂缝半径、面积和缝口流量与隐式水平集解析模型非常接近。同时，还观察到明显的应力阴影现象，其中外部裂缝优先扩展，而内部裂缝的扩展受到抑制。
由于总裂缝面积与单井产能高度相关，因此，发明人在已获取参数的基础上，进一步将总裂缝面积作为评价水力压裂效果的重要指标，模拟在不同时间下的裂缝扩展形态，其结果如图7所示。从图7可以知道：在模拟的初期，由于所有裂缝的尺寸都较小，因此它们之间的相互作用力不显著，所有的裂缝都可能以较为接近的扩展速率增大裂缝面积，并且该过程几乎与时间呈线性关系。上述模拟结果可说明该模型的模拟准确性。
为了进一步验证本申请所提供的方法的计算效率，本申请还将该方法（C5Frac）与ILSA II的计算时间进行对比，其计算效率的验证结果如表3所示。
表3
| 模型 | C5Frac | ILSA II |
| 计算时间 | 255秒 | 220612秒 |
| 模拟时间/步数 | 203秒128步 | 203秒128步 |
| 处理器/存储器 | Intel-i74770k4.00GHz32 GB RAM | Intel-XEONE56492.53GHz96 GB RAM |
从表中可以看出，当设置相同的模拟时间时，C5Frac的计算效率要远高于ILSA II的计算效率，大大减少了模型计算所需的时间。
此外，发明人还以胜利油田东营凹陷樊页井为对象进行该井全井段开展压裂模拟，并利用现场实际微地震和示踪剂监测结果进行验证，以验证该基于近似解和能量方程的半解析裂缝扩展模拟方法在实际油藏中的有效性和准确性。该井的井身结构图如图13所示。
根据樊页井的测井、录井数据等资料获得的实际的储层物性、岩石力学参数如表4所示：
表4
| 解释序号 | 井段m | 厚度m | 孔隙度% | 渗透率10-3μ㎡ | 解释序号 | 井段m | 厚度m | 孔隙度% | 渗透率10-3μ㎡ |
| 1 | 3217.5-3225.0 | 7.5 | 4.9 | 0.2 | 41 | 3779.0-3794.0 | 15 | 6.1 | 0.5 |
| 2 | 3246.0-3251.0 | 5 | 5.5 | 0.5 | 42 | 3797.0-3865.5 | 68.5 | 5.1 | 0.3 |
| 3 | 3256.0-3263.0 | 7 | 2.6 | 0 | 43 | 3865.5-4001.0 | 135.5 | 5.9 | 0.4 |
| 4 | 3272.8-3274.3 | 1.5 | 2.8 | 0 | 44 | 4001.0-4007.5 | 6.5 | 5.7 | 0.4 |
| 5 | 3275.8-3283.0 | 7.2 | 2.7 | 0 | 45 | 4007.5-4021.0 | 13.5 | 7.1 | 1.1 |
| 6 | 3296.8-3300.0 | 3.2 | 3.3 | 0 | 46 | 4021.0-4048.0 | 27 | 4.9 | 0.3 |
| 7 | 3311.0-3313.4 | 2.4 | 2.6 | 0 | 47 | 4058.0-4068.0 | 10 | 5.9 | 0.5 |
| 8 | 3335.0-3340.0 | 5 | 4.2 | 0.1 | 48 | 4068.0-4084.5 | 16.5 | 6.4 | 0.6 |
| 9 | 3353.7-3356.6 | 2.9 | 2.2 | 0 | 49 | 4084.5-4135.0 | 50.5 | 6.4 | 0.6 |
| 10 | 3357.8-3360.3 | 2.5 | 2.7 | 0 | 50 | 4201.0-4280.0 | 79 | 5.4 | 0.3 |
| 11 | 3361.8-3365.0 | 3.2 | 2.7 | 0 | 51 | 4289.0-4308.0 | 19 | 5.7 | 0.4 |
| 12 | 3368.6-3370.4 | 1.8 | 2.3 | 0 | 52 | 4308.0-4312.0 | 4 | 3.2 | 0 |
| 13 | 3381.8-3389.0 | 7.2 | 2.4 | 0 | 53 | 4312.0-4319.0 | 7 | 6.5 | 0.6 |
| 14 | 3396.5-3397.5 | 1 | 2.4 | 0 | 54 | 4319.0-4331.0 | 12 | 5.6 | 0.3 |
| 15 | 3422.0-3430.0 | 8 | 6.5 | 0.6 | 55 | 4331.0-4372.5 | 41.5 | 5.8 | 0.4 |
| 16 | 3431.5-3432.8 | 1.3 | 7.1 | 1.7 | 56 | 4372.5-4379.0 | 6.5 | 3.3 | 0 |
| 17 | 3435.0-3442.0 | 7 | 7.1 | 1 | 57 | 4379.0-4402.0 | 23 | 6.4 | 0.5 |
| 18 | 3443.8-3445.0 | 1.2 | 3.1 | 0 | 58 | 4402.0-4446.0 | 44 | 5.3 | 0.3 |
| 19 | 3448.0-3451.5 | 3.5 | 6.5 | 0.6 | 59 | 4463.0-4498.0 | 35 | 4.5 | 0.1 |
| 20 | 3454.0-3457.5 | 3.5 | 6.2 | 0.5 | 60 | 4498.0-4512.0 | 14 | 5.6 | 0.4 |
| 21 | 3459.0-3464.0 | 5 | 6.5 | 0.6 | 61 | 4512.0-4523.0 | 11 | 6.4 | 0.5 |
| 22 | 3466.3-3471.0 | 4.7 | 6.3 | 0.5 | 62 | 4523.0-4527.0 | 4 | 5.4 | 0.3 |
| 23 | 3474.0-3486.0 | 12 | 5.2 | 0.3 | 63 | 4527.0-4538.0 | 11 | 6.2 | 0.5 |
| 24 | 3487.6-3490.0 | 2.4 | 9.8 | 3.5 | 64 | 4538.0-4580.0 | 42 | 5.8 | 0.4 |
| 25 | 3490.0-3494.0 | 4 | 5.9 | 0.4 | 65 | 4599.0-4611.0 | 12 | 5.3 | 0.3 |
| 26 | 3495.7-3500.2 | 4.5 | 5.3 | 0.3 | 66 | 4611.0-4646.0 | 35 | 4.3 | 0.2 |
| 27 | 3505.0-3510.3 | 5.3 | 5.8 | 0.4 | 67 | 4646.0-4714.0 | 68 | 5.4 | 0.3 |
| 28 | 3516.0-3519.0 | 3 | 5.7 | 0.3 | 68 | 4714.0-4717.0 | 3 | 2.3 | 0 |
| 29 | 3522.3-3524.0 | 1.7 | 6.3 | 0.6 | 69 | 4717.0-4730.0 | 13 | 4.5 | 0.2 |
| 30 | 3527.0-3533.8 | 6.8 | 4.4 | 0.1 | 70 | 4737.0-4764.0 | 27 | 4.8 | 0.2 |
| 31 | 3539.0-3541.0 | 2 | 9.1 | 2.6 | 71 | 4772.0-4869.0 | 97 | 4.8 | 0.2 |
| 32 | 3541.0-3545.0 | 4 | 7.5 | 1.3 | 72 | 4882.5-4928.0 | 45.5 | 4.8 | 0.2 |
| 33 | 3554.5-3556.3 | 1.8 | 7.4 | 1.1 | 73 | 4928.0-4931.0 | 3 | 3.8 | 0.1 |
| 34 | 3562.5-3563.6 | 1.1 | 1.4 | 0 | 74 | 4931.0-5000.0 | 69 | 5 | 0.2 |
| 35 | 3567.0-3571.0 | 4 | 5.3 | 0.2 | 75 | 5000.0-5012.0 | 12 | 4.8 | 0.2 |
| 36 | 3576.0-3583.5 | 7.5 | 6.6 | 0.6 | 76 | 5012.0-5142.0 | 130 | 3.7 | 0.1 |
| 37 | 3587.0-3602.0 | 15 | 5.8 | 0.4 | 77 | 5142.0-5179.0 | 37 | 6.6 | 0.7 |
| 38 | 3606.0-3632.8 | 26.8 | 6 | 0.5 | 78 | 5179.0-5189.0 | 10 | 4.5 | 0.1 |
| 39 | 3641.0-3669.0 | 28 | 6.1 | 0.5 | 79 | 5189.0-5199.0 | 10 | 7 | 0.8 |
| 40 | 3673.5-3772.0 | 98.5 | 6.5 | 0.6 | 80 | 5199.0-5342.0 | 143 | 6 | 0.4 |
并按照樊页井真实的分段分簇方案，对段/簇间距、射孔参数进行设置，其分段分簇参数如表5所示。
表5
/>
/>
/>
按照樊页井真实的泵注程序，依次对31个压裂段的排量、粘度进行设置，该井的压裂施工料单如表6所示。
表6
该井的各段的泵注参数统计如表7所示。
表7
然后按照不同簇的空间三维位置，绘制樊页井各段各簇的裂缝展布，如图8所示。分别提取每段每簇裂缝形态，查看模拟结果，发现较明显的应力阴影效应；前半部压裂裂缝改造体积大，裂缝开度大；后半部效果较差，裂缝开度较小；如图8所示。
发明人还利用樊页井微地震监测数据对模拟结果进行验证。其具体方法为：首先在保持裂缝面积不变的前提下，利用微地震裂缝长度和宽度获取折算缝长，再与模拟结果进行对比。这是因为由于微地震监测数据中，每个压裂段包含裂缝长度和宽度两个参数，而模拟结果中由于裂缝为圆形，只有裂缝半径参数。微地震监测缝长折算原理如图9所示。其微地震监测缝长与模拟结果对比如图10所示。从图10可以看出，樊页井模拟缝长与微地震折算缝长的整体趋势基本一致，并且部分段的缝长高度吻合，以压裂第12段为例，该段折算缝长350m，模拟缝长337m，模拟误差仅为3.71%，参阅图11。上述结果说明樊页井的模拟结果准确真实。
发明人还将樊页井各压裂段的缝长模拟结果与示踪剂监测结果进行对比，其对比图如图12所示。从图12发现发现模拟结果与监测数据较好得保持了一致，即裂缝发育呈现“差-好-中-好“趋势，将各压裂段的缝长模拟结果与示踪剂监测结果得到的各段产液贡献率进行对比，发现部分段的模拟结果与示踪剂结果的匹配度较好，且部分段的匹配度还要高于微地震监测结果。
最后应说明的是：以上各实施例仅用以说明本申请的技术方案，而非对其限制；尽管参照前述各实施例对本申请进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本申请各实施例技术方案的范围。
Claims (10)
1.一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，
包括以下步骤：
S1:根据目标储层的测井数据、录井数据获得所述目标储层的储层物性、岩石力学参数；根据单井的钻井资料获得所述单井的分段分簇参数和每一压裂段的注排量参数；
S2:建立多裂缝扩展基础模型；
S3:根据所述储层物性、岩石力学参数、所述单井的压裂段分段分簇参数和注排量参数求解所述多裂缝扩展基础模型，得到所述压裂段中的包含时空信息的近似解：裂缝开度，裂缝长度和裂缝内流体压力；
S4:将初始值带入能量平衡方程和流量平衡方程求解得到下一迭代包含时空信息的数值解：裂缝开度，裂缝长度和裂缝内流体压力；然后根据迭代前后的流量差进行收敛性判断；若迭代前后的流量绝对差小于设定的阈值，则满足收敛条件并结束迭代；否则循环步骤S3～S4，直至收敛；获得所述压裂段的数值解作为下一时间步的初始值，进入下一时间步的计算；
S5:重复步骤S3～S4，直至达到所需的总泵送时间t(k)＝T；获得所述压裂段的解值；
S6:重复所述步骤S2～S5，得到所述单井的所有压裂段的解值，即为单井全井段模拟结果；
其中，t(k)＝T中t为时间，k为迭代次数。
2.根据权利要求1所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，其中，所述储层物性包括目标储层的孔隙度、渗透率、流体和固体压缩系数；所述岩石力学参数包括岩石的弹性系数、泊松比和断裂韧性；所述分段分簇参数包括所述单井的段间距、簇间距、射孔长度、孔眼直径、孔密度、单段孔数；所述注排量参数包括每个所述压裂段的泵注总液量、排量、动态粘度和时间；所述泵注总液量包括净液量和混砂液量。
3.根据权利要求1所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，所述多裂缝扩展基础模型包括联立的质量平衡方程、局部法向牵引力T和宽度w之间的弹性关系表达式、径向通量qi(r,t)与流体压力梯度的表达式、裂缝扩展准则、断裂韧性表达式、体积平衡方程、雷诺润滑方程、各裂缝入口压力平衡和射孔簇的流量总和表达式；
其中，t为时间，r为裂缝上一点的坐标。
4.根据权利要求3所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于：
其中，所述质量平衡方程为：
其中：qi是缝内流速；CL为滤失系数；k为岩石渗透性；cr为储层压缩系数；φ为岩石孔隙度；po为储层压力；ti(r)为与时间相关的函数，表示裂缝尖端到达坐标为r的点的时间；
所述局部法向牵引力T和宽度w之间的弹性关系表达式为：
其中F为非局部积分算子；wi(r,t)为裂缝在时间点t坐标r处的宽度；Ri(t)为裂缝在时间点t的半径；E’为弹性模量；
径向通量qi(r,t)与流体压力梯度的表达式为：
其中μ为动态粘度；
裂缝扩展准则的表达式为：
KI＝KIC
其中KI表示I型裂缝的应力强度因子，KIC表示I型裂缝的断裂韧性；
体积平衡方程为：
其中Rw为井眼半径；
雷诺润滑方程为：
各裂缝入口压力平衡和射孔簇的流量总和表达式为：
pf(1)(Rw,t)＝pf(2)(Rw,t)＝…＝pf(N)(Rw,t)
各裂缝内的流量总和的表达式为：
其中，Qi(t)为裂缝i在时间点t的缝口流量；Q0为单个分段总泵注排量；pf(r,t)为裂缝在时间t、坐标距离井筒长度为r处的缝内流体压力；σo为地应力；Δpc为地应力与储层压力的差值；pf(1)(Rw，t)、pf(2)(Rw，t)……pf(N)(Rw，t),分别表示的各裂缝(裂缝1、2……N)入口在时间t处的入口压力；pf(i)(r，t)为裂缝i在时间为t、坐标r处的流体压力。
5.根据权利要求4所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，在裂缝扩展准则的表达式中，对于径向裂缝，KI表示为：
其中，Ti(ρi，t)中Ti为局部法向牵引力，ρi＝r/Ri(t)。
6.根据权利要求4所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，所述体积平衡方程在缝口处的边界条件、裂纹尖端的边界条件和裂纹尖端的初始条件下对进行积分得到；
其中在缝口处的边界条件为：
其中Rw为井眼半径；
裂纹尖端的边界条件表达式为：
裂纹尖端的初始条件为：
其中，hk为压裂段中各裂缝的距离。
7.根据权利要求1所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，所述解值包括：裂缝无量纲半径Ri(t)/Z、缝内流量qi(Rw,t)、缝口宽度wi(Rw,t)和裂缝面积A(t)。
8.根据权利要求7所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，所根据所述储层物性、岩石力学参数对所述多裂缝扩展基础模型进行赋值中赋值的参数包括：E、v、KIC、μ、Q、Z、σmin、Rw、hi,j，其中，E为杨氏模量；v为泊松比；KIC为I型裂缝的断裂韧性；μ为动态粘度；Z为压裂段总长度；Q为单个分段总注入流量；σmin为最小原位围压应力；Rw为井眼半径；hi,j为裂缝i和j之间的距离。
9.根据权利要求1所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，S3中，所述求解基于第一假设，所述第一假设包括：
(1)裂缝扩展遵循线性弹性断裂力学(LEFM)，该力学假设材料在任何地方都遵循线性松弛应力-应变关系；当裂纹尖端的应力强度因子达到材料断裂韧性时，将发生裂纹扩展；
(2)裂缝内为不可压缩牛顿流体，并用润滑方程进行描述；
(3)岩石是不透水的，液体滤失忽略；
(4)所有裂缝均径向平行扩展；
(5)在弹性和流体流动方程中均忽略重力；
(6)流体前沿与裂缝前沿重合；
(7)远场原位应力是均匀且恒定的。
10.根据权利要求8所述的一种基于近似解和能量方程的半解析裂缝扩展模拟方法，其特征在于，所述求解包括以下步骤：
S31：分别设置力学参数{E,ν,KIC,μ,Q,Z,σmin,Rw,hi,j}以及计算的初始时间、最终时间和时间步长{t0,tf,Δt}的值；
S32：设定初始条件为获得每条裂缝(i＝1,…,N)的缝内压力裂缝半径/>宽度/>
其中，
相互作用力为：
其中，和/>通过下式求解：
S33：使用牛顿法对方程组进行数值求解；
其中，E’为弹性模量，Πi(ρi,t)为无量纲压力分布函数。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310668991.XA CN116401897B (zh) | 2023-06-07 | 2023-06-07 | 一种基于近似解和能量方程的半解析裂缝扩展模拟方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310668991.XA CN116401897B (zh) | 2023-06-07 | 2023-06-07 | 一种基于近似解和能量方程的半解析裂缝扩展模拟方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116401897A CN116401897A (zh) | 2023-07-07 |
| CN116401897B true CN116401897B (zh) | 2023-09-01 |
Family
ID=87016464
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310668991.XA Active CN116401897B (zh) | 2023-06-07 | 2023-06-07 | 一种基于近似解和能量方程的半解析裂缝扩展模拟方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116401897B (zh) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116738622B (zh) * | 2023-08-09 | 2023-10-20 | 成都理工大学 | 一种能够降低能耗的射孔器设计方法 |
| CN117114208A (zh) * | 2023-10-23 | 2023-11-24 | 成都理工大学 | 一种全耦合的井工厂动态调整整体压裂优化方法 |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107545113A (zh) * | 2017-09-08 | 2018-01-05 | 西南石油大学 | 非常规油气藏水力压裂复杂缝网形成过程模拟方法 |
| CN108302833A (zh) * | 2017-08-31 | 2018-07-20 | 环和地能创新科技有限公司 | 封闭式深层地热能采集系统和方法 |
| CN111734380A (zh) * | 2020-07-30 | 2020-10-02 | 西南石油大学 | 一种水平井多级压裂裂缝扩展形态的快速预测方法 |
| CN112347411A (zh) * | 2020-10-30 | 2021-02-09 | 西南石油大学 | 一种考虑闭合应力的裂缝岩样酸蚀缝宽数值计算方法 |
| CN113887045A (zh) * | 2021-09-30 | 2022-01-04 | 西南石油大学 | 一种暂堵裂缝动态压力与扩展轨迹的预测方法 |
| CN114970281A (zh) * | 2022-06-17 | 2022-08-30 | 西南石油大学 | 一种考虑热效应的高温深层页岩地层水力裂缝扩展计算方法 |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US9394778B2 (en) * | 2014-04-24 | 2016-07-19 | Edward O. Anders | Apparatus, systems, and methods for fracturing a geological formation |
| CN112036096B (zh) * | 2020-09-07 | 2021-03-26 | 西南石油大学 | 裂缝性地层井壁强化效果评价的流固耦合数值模拟方法 |
-
2023
- 2023-06-07 CN CN202310668991.XA patent/CN116401897B/zh active Active
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108302833A (zh) * | 2017-08-31 | 2018-07-20 | 环和地能创新科技有限公司 | 封闭式深层地热能采集系统和方法 |
| CN107545113A (zh) * | 2017-09-08 | 2018-01-05 | 西南石油大学 | 非常规油气藏水力压裂复杂缝网形成过程模拟方法 |
| CN111734380A (zh) * | 2020-07-30 | 2020-10-02 | 西南石油大学 | 一种水平井多级压裂裂缝扩展形态的快速预测方法 |
| CN112347411A (zh) * | 2020-10-30 | 2021-02-09 | 西南石油大学 | 一种考虑闭合应力的裂缝岩样酸蚀缝宽数值计算方法 |
| CN113887045A (zh) * | 2021-09-30 | 2022-01-04 | 西南石油大学 | 一种暂堵裂缝动态压力与扩展轨迹的预测方法 |
| CN114970281A (zh) * | 2022-06-17 | 2022-08-30 | 西南石油大学 | 一种考虑热效应的高温深层页岩地层水力裂缝扩展计算方法 |
Non-Patent Citations (1)
| Title |
|---|
| 酸岩反应热对裂缝温度场和酸液有效作用距离的影响;王辽;沈羞月;王荣;;重庆科技学院学报(自然科学版)(第06期);219-224 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116401897A (zh) | 2023-07-07 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN116401897B (zh) | 一种基于近似解和能量方程的半解析裂缝扩展模拟方法 | |
| US11542801B2 (en) | Optimized design method for temporary blocking agent to promote uniform expansion of fractures produced by fracturing in horizontal wells | |
| Lian et al. | A study on casing deformation failure during multi-stage hydraulic fracturing for the stimulated reservoir volume of horizontal shale wells | |
| Li et al. | Numerical simulation of limited-entry multi-cluster fracturing in horizontal well | |
| AU2014404414B2 (en) | Optimizing multistage hydraulic fracturing design based on three-dimensional (3D) continuum damage mechanics | |
| CN113076676B (zh) | 非常规油气藏水平井压裂缝网扩展与生产动态耦合方法 | |
| CN111734380B (zh) | 一种水平井多级压裂裂缝扩展形态的快速预测方法 | |
| CN112949129B (zh) | 一种深层页岩水平井压裂多簇裂缝异步起裂延伸计算方法 | |
| CN113850029B (zh) | 一种页岩气水平井密切割压裂射孔参数优化设计方法 | |
| CN106897544B (zh) | 对来自压裂的非常规地层的油气生产进行建模的方法 | |
| CN110469303B (zh) | 一种基于四类改造体积的体积压裂参数优化设计方法 | |
| CN107130959B (zh) | 一种煤层气产量预测方法 | |
| CN113389534A (zh) | 一种水平井密切割压裂裂缝扩展预测与设计参数优化方法 | |
| US10385659B2 (en) | Evaluation of production performance from a hydraulically fractured well | |
| CN115618759A (zh) | 页岩气地层压裂施工参数优化方法 | |
| CN114239308B (zh) | 一种多尺度高密度压裂参数优化方法 | |
| Huang et al. | A new multi-cluster fracturing simulation model coupled with perforation erosion: based on the continuous–discontinuous method | |
| Cohen et al. | Production Forecast after Hydraulic Fracturing in Naturally Fractured Reservoirs: Coupling a Complex Fracturing Simulator and a Semi-Analytical Production Model | |
| CN114925626A (zh) | 考虑孔眼动态磨蚀的多裂缝扩展预测方法及设备 | |
| CN114943138A (zh) | 一种基于柔性线源解析法的复杂结构井产能预测方法 | |
| CN114169261A (zh) | 一种基于g函数曲线分析的多裂缝参数反演方法和装置 | |
| Michael | Hydraulic fractures from non-uniform perforation cluster spacing: Laboratory study on gelatin and data-driven modeling | |
| Gao et al. | Optimal design of nonuniform multiple transverse fractures for horizontal wells in tight gas reservoirs | |
| Zeng et al. | Production calculation of multi-cluster fractured horizontal well accounting for stress shadow effect | |
| RU2745142C1 (ru) | Способ и система моделирования трещин гидроразрыва пласта бесконечно-конечной проводимости и поперечно-продольного расположения относительно горизонтального ствола скважины |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |