非常规油气藏体积压裂数值模拟新进展

唐子春,王 朝,张子珂,胡 蝶,雷征东

(1.中国石油勘探开发研究院，北京 100083；2.中国石油大学(北京))

非常规油气藏体积压裂数值模拟新进展

唐子春1,2,王 朝1,2,张子珂1,2,胡 蝶1,2,雷征东1

(1.中国石油勘探开发研究院，北京 100083；2.中国石油大学(北京))

页岩气藏的高效开发与体积压裂的合理设计密切相关，在调研大量相关文献的基础上，综述了近年来国外非常规油气藏水平井体积压裂数值模拟技术。从天然裂缝与人工缝的相互作用、应力影以及支撑剂模拟等方面叙述了压裂力学模型的发展历程，介绍了非常规压裂模型的原理，为该技术的深入研究提供参考;阐述了连续介质模型和离散裂缝模型在非常规油气藏体积压裂数值模拟中的应用，分析了其原理、数模思路以及适用性;总结了国外研究的主要方向，提出了针对长庆油田致密储层体积压裂数值模拟的发展方向。

长庆油田；非常规油气藏；数值模拟；体积压裂；离散裂缝模型

北美页岩气的勘探与开发至今已有百年的历史，而真正进入商业化生产仅有20多年。水平井开发、体积压裂等技术的突破促进了美国的页岩气革命，使非常规油气在世界能源供给中占有越来越大的比例。从非常规油气藏获得工业性的经济产量有赖于水平井近井地带的改造体积，而在体积压裂过程中产生的水力裂缝会与天然裂缝相互交错，在水平井近井地带形成复杂的裂缝网络。因此，如何合理、准确地描述复杂缝网系统成为非常规油气藏中数值模拟的关键。近年来，国外非常规油气藏体积压裂的数值模拟技术有了飞跃式的发展，而目前国内还没有对此最新研究成果进行系统论述。本文针对该类问题，大量调研了近年来的相关研究，总结了双重介质模型、离散裂缝模型的数值模拟方法，提出了相应的数模思路以及对未来技术的展望。

1 复杂缝网模型

不同于PKN模型、拟三维模型等常规压裂模型，复杂缝网模型需要考虑很多因素，例如天然裂缝和水力裂缝的相互作用、应力影的影响等。Olson等人提出了一个缝网模型来预测水力裂缝的扩展，该模型虽然考虑了人工缝与天然缝的相互作用，但不能模拟压裂液流动和支撑剂输送；Rogers等人[1]提出了一个离散裂缝模型，简化了天然裂缝与人工缝的关系，但没有考虑压裂缝之间的作用；Nagel等人[2,3]提出了一个拟三维缝网模拟器，可模拟天然裂缝与人工缝之间的相互作用，但该模型计算量大，不具有工程适用性；而最为经典的是Xu等人[4,5]提出的线网模型和Weng等人[6]提出的非常规压裂模型。线网模型(Wire-Mesh)将复杂裂缝网络简化为一个椭球体，该椭球体由两组正交、等距分布的垂直裂缝面组成。线网模型应用了岩石力学理论，计算了缝网随着压裂施工的扩展，考虑了天然裂缝与水力裂缝的力学作用，模拟了支撑剂在缝网中的输送。非常规压裂模型(UFM)是Kresse.O[7]于2011年提出的新模型，模型主要考虑以下几个方面：

(1)基本的控制方程。①流动方程：一维泊肃叶圆管层流，假定流体为幂律流体；②质量守恒；③考虑液体滤失效应：限制条件是射孔点流速之和等于总排量。

(2)高度生长模型。主要基于拟三维压裂模型计算，考虑流体压力、地应力、断裂韧性、弹性模量等因素，在裂缝尖端的上部、下部分别计算断裂韧性，判断裂缝高度延伸的可能方位。

(3)裂缝相互作用模型[8]。该方法是最早实现水力裂缝与天然裂缝相互作用的模型之一，模型本身基于线弹性断裂力学，没有考虑液体影响。水力裂缝默认垂直于最小水平主应力方向起裂，单一裂缝时，非常规压裂模型退化为拟三维压裂模型。Warpinski[9]提出，当水力裂缝遇到天然裂缝时，会出现穿过和滑移两种情况，而Blaton[10]的观点是会出现张开、穿过、捕捉三种情况，但二者本质相同，分别由图1b和c所示。

图1 水力裂缝遇到天然裂缝的可能情况

(4)应力影模型。模型使用边界元(BEM)方法，该方法考虑了三维应力修正项，模拟结果与CSIRO十分接近[11]。式(1)和式(2)分别代表第i个裂缝单元在受到应力影作用后的法向与切向应力值。

(1)

(2)

式中n，s分别代表张开、剪切应力；Dn，Ds代表张开、剪切位移；N代表裂缝体内总单元数；Cij是二维平面应变影响因子，该因子是用于作用于计算来自相邻裂缝单元的额外影响应力值，具体表达式参见Crouch[12](1976)；Aij是三维修正系数，修正因子Cij。当距离增加时，三维条件下有限的裂缝高度会削弱裂缝单元间的影响，需要引入该系数进行修正。具体表达式参见Olson[13](2004)。

非常规压裂模型在考虑应力干扰时有两大假定：①裂缝宽度等于法向位移不连续度；②裂缝表面剪切应力等于0。基于这两个假设才能计算切向位移不连续度及法向的应力干扰值。应力影模型对裂缝传递过程的影响体现在两个方面：一是影响井场所有裂缝单元的初始地应力大小与方向，改变了裂缝的压力与宽度；二是影响了裂缝尖端的应力场，可能导致裂缝沿其他方向扩展。多条裂缝的相互干扰模拟结果：当裂缝间距低于10 m时，水力裂缝之间表示出强烈的排斥现象，而间距40 m以上排斥现象变得比较微弱。

(5)支撑剂传输模型。为提高计算效率，水平方向(沿缝长)设定为一维模型，垂向(沿缝高)设定为三层模型，从底部到上部分别为砂堤、混砂液、无砂液，同时考虑支撑剂的沉降与侵蚀作用。Cohen等人[14]对非常规压裂模型的各个参量进行了1800多次的拟合试验，结果表明，支撑剂尺寸为40～70目时，可流动面积比30～50目大了3倍以上，分布也更为均匀。模型计算过程考虑了各个模型的相互作用，对每个时步，该模型基本计算流程如图2所示。

图2 UFM模型计算流程

非常规压裂模型最大的优势在于，除了考虑常规压裂软件的基本功能外，还考虑了多缝间应力干扰，包括单井多段间应力干扰与井间应力干扰，在现有的压裂模型中有独特的优势。该模型还可以模拟与天然裂缝的相互作用，从而进行多井压裂等复杂的压裂设计，为压裂领域较为前沿的压裂设计模型。但是，模型本身的天然裂缝与水力裂缝的扩展准则是基于线弹性断裂力学理论，没有考虑流体的影响，可能会造成较大的偏差；地层天然裂缝识别困难，给建立天然裂缝网络带来一定的不确定性。目前，模型还处于进一步完善之中。

2 体积压裂数值模拟

2.1 连续介质模型

连续介质模型包括双孔模型、双孔双渗模型以及多重孔隙介质模型。传统的双孔模型假设基质到裂缝的流动为拟稳定流，而在非常规油气藏中，流动并不是瞬间发生，需要多级基质网格来模拟致密储层中的微孔隙、裂缝和有机质区域，获取过渡流特征。Pruess[15]提出的MINC模型解决了这个问题，该模型将基质单元剖分成一系列次级单元，这些次级单元形成独立的一维流动系统，并与对应的裂缝单元相连，后续研究又将模型中次级单元与裂缝表面的距离改进为对数变化。Zhang[16]使用ECLIPSE中的多重孔隙结构来描述体积压裂区域基质和裂缝的关系，提出了一个综合模拟方法，研究油藏和水力压裂参数对页岩气藏生产动态的影响，模拟结果表明，该方法很好地获取了基质到裂缝之间的过渡流动特征。

对于水力裂缝的模拟，Kalantari Dahaghi[17]研究认为，对数间距局部网格加密对水平井体积压裂区域有较好的模拟效果。在此基础上，Cipolla和CMG公司的Rubin[18-19]基于双重介质模型提出了一种叫做“双渗模型-对数间距变化-局部网格加密(DK-LS-LGR)”的方法，结合微地震监测结果，对致密储层中的复杂缝网进行整体建模。该方法在体积改造区域内部采用对数间距变化、局部网格加密的双渗模型，在外部采用普通的双渗网格。复杂缝网的建模首先通过微地震监测数据筛选体积改造区域的裂缝网络，之后，通过对数间距局部网格加密模拟水平井的体积改造区域。结果表明，该方法历史拟合效果好，可以较好反映页岩的生产特征。Du等人[20-21]提出了一个综合的工作流对Barnett页岩气藏进行模拟。在地质模型的建立过程中，天然裂缝建模基于井筒成像、测井和微地震监测的解释结果，水力裂缝建模依赖于压裂施工参数、支撑剂的分布以及微地震监测的结果，采用Oda方法将油藏地质模型和离散裂缝模型粗化为双孔模型。现场应用结果表明，这个双孔模型工作流可以较好地模拟非常规油气藏的生产，使人们对致密储层的生产有了更深的认识，同时，也为非常规油气藏一体化模拟工作流的快速发展奠定了基础。

2.2 离散裂缝模型

离散裂缝模型(DFM)是指用一系列指定网格代表裂缝，使裂缝离散分布在模型中，利用网格之间的传导率控制基质与裂缝之间的流动。根据求解方法的不同，可将该方法分为非结构化离散裂缝模型(UDFM)和嵌入式离散列分模型(EDFM)[22-23]。

2.2.1 非结构化离散裂缝模型

常规的离散裂缝模型有赖于非结构化网格以及局部网格加密实现复杂缝网的精细描述，Karimi-Fard[24]、Fu[25]、Marcondes[26]等人应用有限元方法建立了离散裂缝模型，采用非结构化网格对裂缝进行了显性剖分。

近年来，通过将地震、测井、微地震监测以及岩心分析结合到复杂缝网模型，进行不同方式的网格化，实现了离散裂缝模型到油藏数值模型的耦合。Cipolla等人[27-29]提出了地震-数模一体化工作流,该工作流包括以下步骤：①通过地震、测井解释以及岩心分析手段，提供必要的岩石力学性质、地应力分布、断层以及天然裂缝的分布，建立离散裂缝模型和地质力学模型，并运用微地震监测进行约束；②基于已建立的地质模型，使用非常规压裂模型(UFM)输入压裂施工参数和岩石力学性质，通过计算耦合人工裂缝和天然裂缝，生成复杂缝网系统；③通过自动非结构网格化算法进行局部网格加密，生成油藏数值模型，由下一代数值模拟器INTERSECT进行计算。图3为页岩气藏水平井衰竭式开发10年以后压力分布模拟结果，从图中可以看出，该工作流可以精细地描述水平井体积压裂区域的开发特征。

图3 体积压裂水平井衰竭开采10年后的压力分布

在此基础上，Cipolla等人[30]提出了一个完井设计工具，将其耦合到了工作流中。该工具利用应力场、岩石力学特性和成像技术，对页岩孔隙度、黏土含量、有机碳含量(TOC)等参数进行筛选，确定合理的压裂段和射孔位置，进行严格高效的水平井分段压裂完井设计，判别标准如表1所示。Cohen等人[31]也将工作流做了一些改进，包括压裂段之间的应力干扰、压裂液温度和流变性的耦合等因素。目前，该流程已在斯伦贝谢公司的Mangrove软件中实现，工作流如图4所示。

地震-数模一体化工作流实现了压裂模型与油藏数值模拟器之间的耦合，但是没有有效的工具来模拟体积压裂对非常规油气藏单井产能的影响。Cohen等人[32]耦合了一个半解析的非常规产能模型(UPM)，进行多参数敏感性分析，研究压裂施工参数设计和油藏性质与产能之间的关系。在此基础上，Cohen等人[33]耦合了一个温度模型，研究压裂液体系和油藏温度对产能的影响。

图4 地震-数模一体化工作流

油藏性质完井性质总有机碳>=3%低应力值天然气含量>2.83m3/t电阻率>15Ω·m干酪根含量高黏土含量<40%泥岩孔隙度>4%杨氏模量>13.8GPa泊松比<0.2气测渗透率>10-4μm2中子测井孔隙度<35%密度孔隙度>8%

2.2.2 嵌入式离散裂缝模型

双重介质模型建立的缝网形态相对理想，但不能正确描述改造区域的复杂情况；而常规离散裂缝模型依赖于局部网格加密或非结构化网格，可以满足高精度的要求，但网格数目增加，造成计算量呈数量级增长。为了解决这些问题，Li和Lee(2008)[34]首次提出了嵌入式离散裂缝模型(EDFM)，该模型使用结构化网格描述基质系统，通过计算裂缝与基质网格的交会部分，显性地引入裂缝系统控制体作为单独的一部分，将裂缝平面嵌入到基质网格中并进行边界元离散化。

Moinfar[35-36]系统地发展了嵌入式离散裂缝模型，增加了倾斜裂缝的特性，并将其应用于全隐式组分模拟器中。由于两个系统的物质平衡方程之间没有直接关联，需要通过非相邻连接(NNC)实现相互之间的流体流动，即模型中需要计算的四种传导率——裂缝系统中裂缝控制体之间的传导率、基质网格与裂缝控制体之间的传导率、相交裂缝之间的传导率以及裂缝和井筒之间的传导率。通过定义多种类型的传导率，嵌入式离散裂缝模型可以动态地表征裂缝特征，在每个时间步末端更新裂缝的开度和渗透率，重新计算传导率并应用到下一个时间步的流动方程中。Shakiba等人[37-38]利用微地震监测对嵌入式离散裂缝模型进行约束，同时应用组分模拟器UTCOMP更加精准而高效地对体积压裂区域的缝网系统进行了数值模拟。

嵌入式离散裂缝模型继承了双重介质模型的优势，不需要裂缝控制体与基质网格相邻，也不要求基质和裂缝系统具有相同的网格尺寸。同时，该模型离散化表征了裂缝系统，不需要在模型中使用局部网格加密或非结构网格化。上述优势使该模型可以在常规数值模拟器中应用，大大提高了计算效率和精确度。

3 发展趋势

国外页岩气藏开发具有独特的商业化模式，需要在衰竭式开采条件下将采收率最大化，其研究内容主要包括：①体积压裂工艺参数对页岩气产能的影响，例如排量、支撑剂类型、压裂液体系和温度等；②完井工艺对页岩气藏生产动态的影响，例如分段多簇压裂中的压裂段数、段间距以及每段中的射孔簇数等因素，相应的数值模拟技术的突破与创新也着眼于这些方面。Mangrove的一体化工作流主要为单井模拟，而近年来对页岩储层中井间干扰现象的关注度越来越高；同时，随着拉链式压裂和井工厂作业的推广，基于多井模型的非常规油藏数值模拟技术开始在国外发展。

不同于非常规页岩储层，长庆油田主要是致密砂岩储层，在压裂施工和射孔完井工艺上不需要进行有机碳含量、黏土含量等因素的筛选，即"甜点区"的优选；同时，长庆油田主要采用注水开发致密油，在开发方式上与国外存在根本性的差异。因此，体积压裂数值模拟技术的发展，需要考虑超低渗、致密储层中的应力敏感、启动压力梯度等特殊现象，并将研究重点放在注采井网的优化上来。

4 结论

(1)复杂缝网模型逐渐向复杂化、系统化的方向发展。非常规压裂模型可以较好地模拟压裂过程中缝网的形成及缝间应力干扰等问题，能够并入油藏数值模型进行一体化分析，对现场问题提供了很好的解决方案。过于简化的天然裂缝与水力裂缝的扩展准则，使裂缝扩展与实际有一定偏差，且对天然裂缝的描述仍然停留在统计意义上，很难准确描述地层裂缝系统，模拟的裂缝形态有很大不确定性。在以后的研究中，仍然需要大量的现场数据对非常规压裂模型进行完善与修改，使之更加精确与高效。

(2)对于体积压裂水平井的数值模拟，首先需要地质建模，通过力学模型与微地震手段的结合对复杂缝网进行描述；对于流动模型的选择，如果选用双重介质模型，则需要一种合理的粗化技术；如果选择离散裂缝模型，则需要一种高效的非结构化网格技术以及局部网格加密方法对裂缝网络准确刻画。

(3)双重介质模型假设的缝网形态较为均匀，不能准确描述实际情况；常规离散裂缝模型依赖于非结构网格化技术，具有较高的精确度，但是大量的网格数目造成了巨大的计算量；嵌入式离散裂缝模型综合了两种模型的优势，较好地解决了一些复杂的问题。近年来发展较快的一体化工作流，基于连续介质模型和离散裂缝模型，综合应用地质研究、微地震监测和压裂参数等手段对体积压裂水平井中的复杂缝网进行了表征，表现出良好的效果。

[1] STEPHEN R, DAVIDE E, RORY D, et al. Understanding hydraulic fracture geometry and interactions in the horn river basin through DFN and numerical modeling[C].SPE 137488, 2010.

[2] NEAL B N, IVAN G, MARISELA S,et al.Simulating hydraulic fracturing in real fractured rocks-overcoming the limits of pseudo 3D models[C].SPE 140480, 2011.

[3] NEAL B N, IVAN G, MARISELA S, et al. Stress shadowing and microseismic events: A Numerical Evaluation[C].SPE 147363, 2011.

[4]WENYUE X,MARC J, UTPAL G,et al.Wiremesh:A novel shale fracturing simulator[C].SPE 132218, 2010.

[5] WENYUE X, JAMES L, MIKE D. Quick estimate of initial production from stimulated reservoirs with complex hydraulic fracture network[C].SPE 146753, 2011.

[6]XIAOWEI W, OLGA K, CHARLES E C, et al. Modeling of hydraulic fracture network propagation in a naturally fractured formation[C].SPE 140253, 2011.

[7] KRESSE O, COHEN C, WENG X, et al. Numerical modeling of hydraulic fracturing in naturally fractured formations[J].Rock Mechanics,2011,15(5)516-535.

[8] GU H, WENG X, LUND J B, et al. Hydraulic fracture crossing natural fracture at non-orthogonal angles, a criterion, its validation and applications[C].SPE 139984, 2011.

[9] WARPINSKI N R, TEUFEL LW. Influence of geologic discontinuities on hydraulic fracture propagation[C].SPE 13224,1987.

[10] BLANTON T L. An experimental study of interaction between hydraulically induced and pre-existing fractures[C].SPE 10847, 1982.

[11]KRESSE O, WENG X, WU R, et al. Numerical modeling of hydraulic fractures interaction in complex naturally fractured formations[J]. Rock Mechanics and Rock Engieering,2013,46(3):555-568.

[12] CROUCH S L. Solution of plane elasticity problems by the displacement discontinuity method[J]. International Journal for Numerical Methods in Engineering, 1976, 10(2):301-343.

[13] OISON J E. Predicting fracture swarms——the influence of subcritical crack growth and the crack-tip process zone on joint spacing in rock[J]. Geological Society London Special Publications, 2004, 231(1):73-88.

[14] COHEN C E, ABAD C, WENG X, et al. Analysis on the impact of fracturing treatment design and re-servoir properties on production from shale gas reservoirs[C].IPTC 16400, 2013.

[15] KARSTEN P. A practical method for modeling fluid and heat flow in fractured porous media[C].SPE 10509, 1985.

[16] XU Z, CHANGAN D, FRANZ D, et al. Sensitivity studies of horizontal wells with hydraulic fractures in shale gas reservoirs[C]. IPTC 13338, 2009.

[17]AMIRMASOUD K D. Numerical simulation and modeling of enhanced gas recovery and CO2sequestration in shale gas reservoirs: A Feasibility Study[C]. SPE 139701, 2010.

[18] CRAIG L C, ELYEZER L, JAMES C E, et al. Re-servoir modeling in shale-gas reservoirs[C].SPE 125530, 2009.

[19] BARRY R. Accurate simulation of non-darcy flow in stimulated fractured shale reservoirs[C].SPE 132093, 2010.

[20] CHANGAN D, XU Z, BRAD M, et al. A workflow for integrated Barnett shale gas reservoir modeling and simulation[C].SPE 122934, 2009.

[21] CHANGAN M D, XU Z, LANG Z, et al. Modeling hydraulic fracturing induced fracture networks in shale gas reservoirs as a dual porosity system[C].SPE 132180, 2010.

[22] HUI M H， MALLISON B T. System and method for predicting fluid flow characteristics within fractured subsurface reservoirs:US[P],2010.

[23] HUI M H， MALLISON B T. System and method of grid generation for discrete fracture modeling: US[P],2011.

[24] MOHAMMAD K F, ABBAS F. Numerical simulation of water injection in fractured media using the discrete-fracture model and the Galerkin method[C].SPE 83633, 2003.

[25] FU Y, YANG Y K, DEO M. Three-dimensional, three-phase discrete-fracture reservoir simulator based on control volume finite element formulation[C].SPE 93292, 2005.

[26] MARCONDES F. An element-based finite-volume method approach for naturally fractured compositional reservoir simulation.13th Brazilian Thermal Sciences Meeting, Uberlandia, Brazil, December.

[27] CRAIG L C, XIAOWEI W, MARK G M, et al. Integrating microseismic mapping and complex fracture modeling to characterize hydraulic fracture complexity[C]. SPE 140185, 2011.

[28] CRAIG L C, TONY F, MICHAEL J W, et al. Seismic-to-simulation for unconventional reservoir development[C].SPE 146876, 2011.

[29] MAHSANAM M, CRAIG L C. A workflow for modeling and simulation of hydraulic fractures in unconventional gas reservoirs[C].SPE 153022, 2012.

[30] CRAIG L C, XIAOWEI W, HITOSHI O, et al. New algorithms and integrated workflow for tight gas and shale completions[C].SPE 146872, 2011.

[31]HASHEMI R, CAOTHIEN S, ROVIG J. Field measurement of viscosified fluids[C]. SPE 17462, 1988.

[32] COHEN C E, WENYUE X, XIAOWEI W, et al. Production forecast after hydraulic fracturing in naturally fractured reservoirs: coupling a complex fracturing simulator and a semi-analytical production model[C].SPE 152541, 2012.

[33] COHEN C E, ABAD C, WENG X, et al. Analysis on the impact of fracturing treatment design and re-servoir properties on production from shale gas reservoirs[C].IPTC 16400, 2013.

[34] LIYONG L, SEONG H L. Efficient field-scale si-mulation for black oil in a naturally fractured reservoir via discrete fracture networks and homogenized media[C].SPE 103901, 2006.

[35] ALI M, ABDOLJALIL V, KAMY S, et al. Deve-lopment of a novel and computationally-efficient discrete-fracture model to study IOR processes in na-turally fractured reservoirs[C].SPE 154246, 2012.

[36] ALI M, KAMY S, RUSSELL T J, et al. Coupled geomechanics and flow simulation for an embedded discrete fracture model[C].SPE 163666, 2013.

[37] JOSESA C F, MAHMOOD S, ALI M, et al. Implementation of a preprocessor for embedded discrete fracture modeling in an IMPEC compositional reservoir simulator[C].SPE 173289, 2015.

[38] MAHMOOD S, KAMY S. Using embedded discrete fracture model and microseismic monitoring data to characterize the complex hydraulic fracture networks[C]. SPE 175142, 2015.

编辑：王金旗

1673-8217(2017)03-0108-06

2016-10-22

唐子春，1991年出生，2014年毕业于中国石油大学(北京)石油工程专业，现为中国石油大学(北京)石油与天然气工程专业在读硕士研究生，主要从事超低渗透油藏工程与数值模拟研究。

国家科技重大专项“低渗、特低渗油藏水驱扩大波及体积方法与关键技术”(2017ZX05013-002)。

TE357.7

A