第6章__水力压裂力学

1、 郭建春 副教授、博士 Tel13980785559E-mail: 2003.8 成都 压裂原理、设计与实施：压裂原理、设计与实施： 第第6章：水力压裂力学章：水力压裂力学 第第7章：压裂液化学和支撑剂章：压裂液化学和支撑剂 第第8章：压裂材料的性能章：压裂材料的性能压裂设计（第压裂设计（第10章）章）压裂施工（第压裂施工（第11章）章）碳酸盐岩酸化设计（第碳酸盐岩酸化设计（第17章）章）碳酸盐岩酸化技术（补充）碳酸盐岩酸化技术（补充） 6.1 引言 6.2 早期水力压裂模拟 6.3 三维和拟三维模型 6.4 滤失 6.5 支撑剂铺置 6.

2、6 热传递模型 6.7 缝端效应 6.8 裂缝弯曲以及其它近井筒效应 6.11 泵注程序设计 6.10 多层压裂 6.9 酸压裂 6.12 压裂历史拟合水力压裂力学是对压裂工艺和压裂机理的简单描述。 6.1 引引 言言水力压裂力学流体力学固体力学断裂力学热力学描述单相、两相或三相流体在裂缝中的流动描述由于流体压力变化引起的岩石变形或张开描述与水力裂缝端部附近发生的破坏和裂开描述压裂流体与地层之间的热交换所有的响应是耦合的，相互影响开发和利用水力压裂施工的重要原因 进行经济优化 （确定多大施工规模得到最高回报率） 施工评估 模拟特定的泵注程序得到相应 的裂缝几何形状和支撑剂铺置 泵注程序优化一

3、裂缝延伸二维模型l卡特模型 Carter， 1957年lGDK 模型 Christianovich、Geertsma、Deklerk lPKN 模型 Perkins和Kern 提出, Norgren完善裂缝延伸三维模型l拟三维模型 三维扩展，一维流体流动l真三维模型 三维扩展，至少二维流体流动模型分类： 6.2.1 基本的压裂模拟基本的压裂模拟22)(1)1 (8)(RrEpRrw 椭圆裂缝的体积为：netpERV311632）（ 半径为R的静态扁平裂缝的宽度： 半径为R的裂缝扩展的压力：REpFnet212Sneddon 和 Elliot（1946）（6.1）（6.2）（6.3）对于缝高hf

4、不变和无限大（即平面应变）裂缝其最大宽度为： Ehpwfnet212裂缝的形状为椭圆，平均缝宽 ww4。定义平面应变模量E更为方便：21EE（6.5）（6.4）Perkin 和 Kern （1961） 径向裂缝扩展的压力： 5122233132VEpFnet 泵注排量qi保持不变，裂缝中的流体摩擦阻力不计，没有滤失时： 5122233321323116tqEERqiFi整理得到R：51222211289tqtEqRii（6.6）（6.7）（6.8） 6.2.2 水力压裂二维模拟水力压裂二维模拟o PKN模型 假设每一垂向截面独立作用，即假设截面的压力是由 高度控制的而非由缝长控制的。在缝长远大

5、于缝高的条件下成立没有考虑断裂力学和缝端的影响，而主要考虑了缝内流体的流动以及相应的压力梯度的影响图6.1 PKN裂缝o KGD模型 假设每一水平截面独立作用，即假设裂缝面任一点处裂缝宽度沿垂向变化远比水平方向的变化慢。在缝高远大于缝长或者储积层边界产生完全滑移的条件下成立缝端区域起着很重要的作用，而缝内压力可以估算图6.2 KGD 裂缝6.2.2.1 垂向裂缝的Perkins 和 Ken模型的推导 流动的基本方程：364whqdxdpf将缝宽方程并用注入速度的一半代替q，并假设流速沿缝不变得到：代入上式dxhEqdppfinetnet4334414316LhEqpfinet沿裂缝半长L对上式

6、积分，并利用边界条件pnet=0 得到：（6.9）（6.4）（6.10）（6.11）4/13)(ExLqxwi 实际的裂缝宽度：重要发现： 垂向平面应变特性的假设 断裂韧性可忽略（裂缝延伸所需的能量远比流体沿 缝长方向流动所需的能量最小） 缝中流体滤失和存储或者体积变化可以忽略的假设 固定缝高的假设 没有直接给出作为解的一部分（6.12）6.2.2.2 模型中考虑流体滤失 裂缝任一点处的滤失速度：expttCuLLCL滤失系数 t当前时间texp该点滤失速率 uL持续的时间 质量平衡方程：fLiqqqqL 整个裂缝的滤失速度 qf 缝内流体存储体积流速Cater（1957）（6.13）（6.1

7、4）假设裂缝在空间和时间上都保持恒定，上式变为： tAffLiftAwdAuq02即：tAwdAAtuqffftli02利用拉普拉斯变换得到： 12422SerfceCwqAsLif 压裂设计是通过由Carter的方法得到与时间有关的缝长与由Kern模型确定的缝宽之间反复迭代，直到得到相容解（6.15）（6.16）（6.17）Nordgren（1972） 连续性方程（即质量守恒）：0tAqxqLq 流体通过某一横截面的体积流速A 裂缝的横截面积（对于PKN模型为whf/4）qL单位长度上滤失体积流速 LfLuhq2其中：uL由方程（6.13）得到，横截面面积不是裂缝面的面积Af（6.18）压力

8、用缝宽表示代替，方程（6.18）写为： twxttCxwhELfexp2428128以无量纲形式对该方程数值求解，得到与时间有关的缝宽和缝长。方程解中的无量纲时间定义如下：tqhECtffLD3223564（6.19）（6.20）6.2.2.3 Khristianovich-geerssma-de Klerk 模型的导出 Khristianovich 和 Zheltov（1995）导出了缝高远大于缝长，即离开井筒任意距离时缝宽与垂向位置无关这种水力裂缝延伸的解。 通过假设缝内流速恒定；除缝端没有流体穿透（即没有压力）外，缝中的压力大部分处的压力以定压近似。可用解析法解该问题。 流体滞后的概念一

9、直是缝端力学的中的重要组成部分，已经在现场得到证明（Warpinski，1985）。如果缝端无流体穿透区很小（约为总缝长的百分之几），他们发现裂缝主体中沿整个缝的压力几乎等于井中的压力，只是在靠近缝端剧减。 Geertsma 和 de Klerk（1969）对于缝端区域很小这个问题给出了解。 对于矩形横截面，流动的基本方程为：312whqdxdpf 可以积分形式写为：Lfinetwdxhqp036（6.21）（6.22） 应用Barenblatt缝端条件，意味着应力集中系数为零。 裂缝宽度方程：netwLpEw4 0102LnetLxdxxp图6.3 Barenblatt 的缝端状况（6.23

10、）（6.24） 通过解方程（6.22）至方程（6.24）三个方程，得到Perkins和Kern（1961）给出的表达形式。4132,6421ELhqpfiwnet（6.25）41284fiwhELqw 井壁裂缝宽度：（6.26） 在没有滤失的情况下，解得缝长和缝宽: 32613338. 0thqEtLfi（6.27）31613348. 1thEqwfiw（6.28） 假设流体滤失对裂缝形态或压力分布没有影响，将模型推广到包括流体滤失的情况下：一个两翼KGD裂缝的体积为：wffLwhV2运用体积平衡和与Carter相似的解法，得到： 12642SerfcehCwqLsfLwi其中：wLwtCS8

11、（6.29）（6.30） 为了包括瞬时滤失Sp的影响，应该以ww+(8/)Sp代替ww 。6.2.2.4 PKN 和 KGD 模型的假设平面裂缝（裂缝沿最小主应力垂直方向扩展）流动沿缝长一维流动流体为牛顿流体滤失特性由滤失理论（6.13）得到的简单表达式所控制地层岩石为连续、均匀、各向同性的线弹性体裂缝被认为缝高不变，完全在某一给定的地层中扩展PKN模型假设缝长远大于缝高，忽略了有关断裂力学的影响KGD模型假设缝高远大于缝长，包括了缝端动态过程控制裂缝延伸的假设6.2.2.45 PKN 与GDK 模型比较 1. 模型差异 裂缝几何形态: PKN:垂直剖面为椭园形状;GDK:垂直剖面为矩形 PK

12、N：水平剖面为(2n+2)次抛物线形;GDK:水平剖面为椭园形 PKN：裂缝长而窄;GDK:裂缝短而宽 应变: PKN:平面应变发生于垂直剖面，层间无滑动;GDK:平面应变主要发生于水平剖面，层间有滑动. PKN：裂缝张开在垂直剖面求解;GDK:裂缝张开在水平剖面求解. 压力行为: 井底压力随时间增加而升高，随缝长增加 井底压力随施工时间逐渐降低，随缝长增加而递减 2.选择的一般原则 根据井底压力的变化趋势。 压裂层的埋藏深度：浅层选用GDK，深层选用PKN。 压裂层与上下岩层的地应力差：上下岩层的地应力大于压裂层且差值大于13.8MPa, PKN和GDK都适用。上下岩层的地应力小于压裂层且差

13、值小于13.8MPa,薄层选PKN，块状厚层或射孔段长的井选GDK。 施工参数：PKN用于低粘高滤失压裂液小排量注液施工规模大；GDK用于高粘低滤失压裂液大排量注液施工规模小。3.实例比较 参数地 层 参 数工 程 参 数油层中部深度 2590m渗透率 0.015*10-3um2有效厚度 46m油层压力 29.3Mpa闭合压力 40Mpa油层温度 107 C供油半径 805m裂缝高度 161m压裂液滤失系数 3.0*10-4 m/min1/2总用液量 1336m3前置液百分数 36.3%总加砂量 284m3泵注排量 4.0m3/min平均砂比 33%计 算 结 果PKNGDK动态缝半长，m支撑

14、缝半长,m动态缝宽,mm支撑缝宽,mm井眼附近支撑剂浓度,kg/m2稳定状态下压后日产量,104m36334075.843.304.305.323932449.965.547.324.22 计算结果支撑半长项 目模型229457686PKN56818003785用液量,m3GDK66220804353支撑剂量,m3施工费用,104元导流能力,um2.cm无因次导流能力压后累计产量 （3 年） ,104m3净现值,104元PKNGDKPKNGDKPKNGDKPKNGDKPKNGDKPKNGDK5371507525752541233456645735106230175225295615267434

15、849639443217740835047532691123920410195516563前面简单模型的局限性： l 需要给定缝高或假设产生的是径向缝原因： l 不能断定裂缝是否被限制在某一特定的地层中 l 由井筒（压力最高处）至缝端的过程中缝高是变化的解决办法： l 利用平面三维3D和拟三维（P3D）模型来弥补包括缝高增长的三种主要水力压裂模型o 普通三维模型 l 没有对裂缝方位作假设l 计算量大，需要专人对结果作解释l 模型适合于研究水力裂缝起裂的细节以及 近 井筒的复杂情况，而非裂缝整个延伸过程 在此不对该模型作进一步的讨论。o 平面三维模型 l 假设裂缝是平面的，并且其方向与最小主应力

16、方向垂直，没有考虑由于偏离平面引起的复杂状况l 这种模型的模拟软件也需大量的计算，一般不用于常规压力设计l 模型用于研究裂缝的主体在裂缝起裂地层以外或者压裂液垂向流动比水平流动更强烈的情况这种模型在6.3.1节介绍o 拟三维模型 l 主要类型有块体和单元体两种 l 块体（椭圆）模型中，假设垂向剖面由中心相连的两个半椭圆组成，每一时间步长计算出水平裂缝和井筒中裂缝缝端的垂向延伸，假设的裂缝形态也要拟合到这些位置；采用固有的假设条件，分析得到：流体沿射孔到椭圆边缘的流线流动，而且流线有专门的形状。l 单元体模型将裂缝视为一系列相连的单元对待，不需要对裂缝形态进行假设，但一般假设为平面应变，流体垂向

17、流动计算与裂缝几何形状之间没有做完全耦合。 这种模型在 6.3.2 和 6.3.3 节介绍 6.3.1 平面三维模型平面三维模型定义：缝内流体的二维流动与岩石三维弹性响应耦 合的模型。 任意水力压裂模型求解的复杂性在于不同过程裂缝的几何形状和流体流动的密切耦合。在求解过程中应考虑的问题：已知形态和压力的裂缝的宽度剖面裂缝形态已知形态和宽度（已知几何形状）的裂缝内的流体流动Hirth 和Lother（1968）以及Bui（1977）裂缝中压力和缝宽的关系式：ydxdyxyxpyyxxfyxw,（6.31）式中：应力 f 弹性影响函数，一般情况下只有对于均质线弹性材料，才可以导出该方程的可用的形式

18、（见旁注6E）。在实际应用中，一般假设岩石为各向同性。o 压力和缝宽的关系简单的裂缝形态和缝内压力分布如定压下的椭圆形裂缝，破裂准则：方程（6.3） 复杂的裂缝形态和压力分布，破裂准则用缝端附近的缝宽和临界应力集中系数或断裂韧性KIC表示： xEKxwIC24（6.32）式中：x距离缝端的距离 裂缝的形态随时间不断演化，假设该过程用线弹性断裂力学描述。o 质量守恒方程（描述流体流动）：02Lyxuwtywuxwu（6.32） 可以写为矢量形式：02Luwtuw（6.33） 上式中的前两项与质量流量的矢量的变化有关， 后两项分别表示由宽度增加和滤失引起的流体存储。 方程（6.34）的左边为动量改

19、变速率；右边分别为压力、粘滞力和重力，它可解释为小的流体单元在力的作用下而加速。该方程可以扩展并根据压裂地层的不同形状而简化（见旁注6F）。o 动量守恒方程：gpDtuD（6.34）式中：剪应力 g 重力加速度应力与流速之间的本构方程：对于x分量 方程（6.34）变为：xzxyxxxxgzyxxpDtDuzuxzxxzzuyzyyz方程（6.34可写为）：xxxgzuxpDtDu22（6.35）（6.36）（6.37）假设流动为稳定流动，得到：xxgzuxp22（6.38）6F 水力压裂中的动量守恒水力压裂中的动量守恒方程（6.34）实矢量方程，其分量形式可以写为：iziyixiiigzyxx

20、pdtdu（6F.1）上式的左边为物质导数，它可与偏导数建立关系：zuyuxutdtdzyx（6F.2）因此方程（6F.1）可扩展为：iziyixiiiziyixigzyxxpzuuyuuxuutup（6F.3） 对于非渗透介质，z方向分量可以忽略不计；并假设为稳态流，（6F.3）简化为：iziyixiiiyixgzyxxpyuuxuup（6F.4）式中：i = 1或2 对于渗透介质，也可采用方程（6F.4），在这种情况下，滤失作为插入项包括在质量平衡中，但假设不影响与压力、应力和流体流速有关的方程。牛顿流体对于牛顿流体，仅含有粘度参数，应力分量为：uxuxxx32uyuyyy32uzuzzz

21、32xuyuyxyxxyyuzuzyzyyzzuxuxzxzzx（6F.5） 对于不可压缩流体，正应力方程中的第三项为零 在平行板间一维流动，如不考虑滤失，两个速度分量全部为零质量守恒意味着第三个方程不能随位置而改变。所以全部的正应力方程为零。上面的方程可以简化为剪切流动方程 原因应力分量：zuxzxzuyzy（6F.6）将方程（6F.6）代入（6F.4）得到：ixxyxxzuxpyuuxuup22gzuxpyuuxuupiyyyyx22（6F.7）如为沿缝长一维流动，方程（6F.7）可简化为：xpzux122（6F.8）如不存在滑移（裂缝面处流体流速为零），方程（6F.8）的解为：22221

22、wzxpux（6F.9）对上式积分，得整个通道上的平均流速：xpwux122（6F.10） 单位高度的流速可以通过平均流速乘以缝宽w得到。 在二维流动的情况下，如果惯性可以忽略，方程（6F.7）左边为零。在这种情况下可按方程（6F.10）写出y方向的方程只是多了重力项。6G 非牛顿流体的动量平衡和本构方程非牛顿流体的动量平衡和本构方程牛顿流体的应力和速度之间的关系，用张量表示：ifif（6G.1）式中：形变张量的速率，其分量为：iiiiifxuxu（6G.2） 对于非牛顿流体写出与（6G.1）相似的方程：ijafi（6G.3）式中：a的函数 a可能仅是的函数，而且两者之间存在某种函数关系： 2

23、laa（6G.4）式中：l2 二级张量变量：ijijijl2（6G.5）如：幂律流体的函数a为：122nalK（6G.6） 宾汉塑性的函数a为：2200la（6G.7） 稠度指数K取决于流形，并与流体基本性质综合稠度系数K有关。对于平板流动，关系为：nnnKK312（6G.8） 对于管流关系式为：nnnKK413（6G.9） 平板间一维流动的幂律流体，平均流速：nnnnxwnnxpKu1112211（6G.10） 对当n=1时，上式转化牛顿流体方程，其中K被粘度代替，表6G.1对不同流形下的表达式作了总结。6G 非牛顿流体的动量平衡和本构方程 方程（6.31）到（6.38）总结了牛顿流体的平面

24、模型，非牛顿流体也可得到相似的结果（见旁注6G）。这些方程一般不适合分析解，需要数值模拟这些方程也很难得到高效稳定的数值解 原因 解的不同部分相互耦合得很紧密 （如流体流动与固体变形） 缝宽与压力之间存在非线性关系以及移动边界问题的复杂性o 平面模型的数值模拟 由Clifton 和AbouSayed首次进行的数值模拟方法 将由射孔孔眼起裂的一条裂缝分成数个相同的单元（一般为16方格）然后开始对方程求解。随着边界延伸，单元要变形以符合新的几何形状。这种解法的一种困难就是单元可能变得高宽比大、角度变小，如图6.5所示，这种数值方法一般用来解对于这种几何形状有问题的那些方程。图6.5 平面3D裂缝

25、分成矩形单元 Barree给出的另外一种数值模拟方法 这种方法通过将分层油藏划分为具有相同大小的矩形单元网格的方法避免了网格扭曲的问题，这些单元网格的边界在裂缝可能产生的区域内。 在这种情况下，网格不移动，而是超过破裂准则时，裂缝破裂缝端前面的单元被打开，允许流体流动而成为裂缝的一部分，如图6.6所示。图6.6 表示随裂缝延伸单元打开的固定网格解这种方法的局限性：随着数值模拟的进行，单元数在逐渐增加，所以最初的单元数少，导致误差大。在模拟前要估计裂缝的大致规模以确保采用合理的单元数。专门的计算过程的假设条件： 用一个简化的方法代表模量差异。 用抗张强度准则代替裂缝扩展而不是断裂力学影响。 靠近

26、缝端裂缝诱导的应力随与距缝端距离的平方根而变化，因此，破裂准则取决于网格划分的精度。 该准则用所有边界单元中心的应力与材料抗张强度相比，如果超过了该强度，那么假定单元张开。 6.3.2 以单元为基础的拟三维模型以单元为基础的拟三维模型 在以单元划分的模型中，缝长被分为数个离散的单元。这与平面模型相似，不过是沿一个方向离散而不是两个方向。假设： 流体流动是沿缝长的方向 固体力学简化为任意截面内的平面应变 这个假设对于缝高得到控制的裂缝是合理的，与缝高相比这种裂缝是长裂缝。 这个假设使得可将固体力学和断裂力学解与流体流动分开： 平面应变暗指每一横截面独立作用与其他横截面无关。 一维流动的假设暗指横

27、截面内的压力始终为：gyppcp（6.39）式中：pcp沿射孔中心水平线上的压力 y 到射孔中心的垂直距离 上面的方程仅对裂缝延伸相当缓慢，由垂向流动引起的压力梯度可以忽略不计时才是有效的。这种裂缝垂向缝端基本处于静态的假设称为平衡高度假设。6.3.2.1 固体力学解 Simonson等（1978）导出了对称的具有三个地层的解 Fung等（1978）导出了更具有普遍性非对称多重底层的解。 在高度平衡的假设条件下，固体力学解简化为确定裂缝横截面形状与静压力或pcp的函数关系。 根据Fung等的方法 顶部和底部缝端的应力集中系数KIu和KIl可以分别以射孔中心的压力pcp和各层中的应力i表示ifi

28、fiffniiiffcpfncpfIlhhhhhhhhhhgphK2cos22 4121111ifififfniiiffcpfncpfIuhhhhhhhhhhgphK2cos22 4321111（6.40）（6.41）式中： pf 流体密度 hcp 射孔中心高度 hi 由底部缝端至第i 层底层顶部的高度，如图 6.7所示图6.7 控制缝高增长问题中的变量定义 通过迭代求解该非线性方程，如果pcp取某一值时的解（两个垂向缝端的位置以及压力）已知，假设缝高增量，则两个垂向的缝高就可以计算出来而且满足方程（6.40）和（6.41），裂缝端达到这些位置所需的pcp可计算出来。最后与该解相关的缝宽剖面也

29、可得到。 fiffiififiiniiifncpcphhhyhyhyhhhhhyyyhEyhyyhpEyw2cos 2cosh 4411111（6.42）式中： y 距底部缝端的高度现如考虑如图6.8所示的一个对称的三个地层的情况忽略重力分量的影响，问题就简化为对称的情况，裂缝向两个遮挡层的延伸是相同的。此时方程（6.40）可简化，得到：fpaypaycpfIlIuhhphKK22（6.43）式中： 产层与隔层的应力差 hpay和i产层的厚度和应力图6.8 简单的三层地层中缝高增长问题图6.9为由（6.43）计算得到的缝高和净压力的关系曲线 尽管（6.43）是一个特例，但它表明了两个有实用意义

30、的结论。 在某一临界压力条件下，裂缝会突破边界进入隔层fIccritnethKp2,2（6.44） 净压力不会达到产层与隔层的净压力差值，因为那种情况会产生无限缝高图6.9 上、下遮挡层对称情况的缝高与净压力关系曲线（hs=进入遮挡层）6.3.2.2 流体力学解流体力学解平面3D与P3D模型间的主要区别：流体流动计算；大多数P3D模型中的流体流动模型与Nordgren（1972）的相同（即一维形式用于描述平面三维模型中）。一般的P3D模型不能代表几方面的特征： 垂向缝宽变化对流体速度的影响 局部失水，是用整个缝高同时失水估算的 端部脱砂（TSOs）后的流体滤失，流体流过支撑剂充填层时，忽略了流

31、体的滤失 由于对流或重力引起的支撑剂沉降采用平均流速和宽度（宽度以横截面积除以高度代替）质量守恒方程简化为：iLLhutAxAu2（6.45）式中： u 平均横截面 ul 每层的滤失速度 hl 每层的滤失高度动量守恒方程简化为：zxpxz（6.46）对于具有性质n和K的幂律流体：nxxzzuK（6.47） 结合无滑移边界条件，解方程（6.46）得到通道中的平均流速：nnnxwnnkxpxpu11221sgn（6.48）式中： sgn 量值标志对于牛顿流体，n=1、=K；上式变为：xpwux122（6.49）求沿缝高方向的所有横截面的总流速 将方程（6.45）中换成平均流速 将方程（6.48）从

32、横截面的缝底到缝顶积分 dyyuywqx（6.50）平均流速：nnfhAKxpxpAqu/112sgn（6.51）式中： 通道函数： dyywhnnAhnnhfnnff2121142（6.52）6.3.2.2.1 层流和紊流层流和紊流流体在平板间以低速无滤失流动，除入口较小区域外，任意流体元与流道的壁之间保持固定的距离，即层流流动，雷诺数NRe 2100即为紊流6.3.2.2.2 压裂液的流变性压裂液的流变性压裂液常为幂律流体；幂律模型中的有效参数K和n是剪切速率在以一定范围内通过室内实验导出的；实际中要考虑到剪切速率只会在有限范围内。压裂液的性质随时间和温度而改变。高温下流体粘度降低；交联剂

33、可能引起流体粘度下降前初始粘度增加；模拟软件中通过K和n在不同温度下随时间变化表而考虑进了温度和时间的影响。6.3.2.3 模型的数值解模型的数值解 描述缝高增长的力学（压力缝宽缝高关系），质量守恒和动量守恒（速度压力相互关系）的三个基本解是相互耦合的，需要同时解。解这些耦合方程的方法：网格点随流体移动的显式有限差分法网格点随流体移动的隐式有限差分法 在开始做裂缝演化模拟前，要先按照在第6.3.2节中介绍的“固体力学解”，计算得到压力缝高缝宽的关系（右平衡高度解）表。 显式有限差分方法显式有限差分方法 任意时刻的裂缝中流体被分为n个单元；每个单元的横截面积为Aj；其在xi和xi+1处的两个边界

34、面分别以ui和ui+1。如图6.10所示（网格以数字编号i=1代表端部，如果有必要它可以用于井筒中的新的单元）。A4A3A2A1x5x4x3x2x1u5u4u3u2u1图6.10 裂缝分成位置和速度由网格点定义的单元质量守恒方程：xAuhtAuiiLL2（6.53）将导数以中心有限差分近似代替得到：xtuuAVxhAAiitiLtitti1（6.54）式中：VL 时间步长为内，整个单元流体的滤失体积 速度为网格点处的计算得到的，并假设每个单元的面积是恒定的，这样横截面面积可以由速度和上一时间步长的面积更新。这样压力梯度：2111iiiiixxppxp（6.55） 由方程（6.51）得到新的速度

35、后，网格位置的更新方程： tttutxttxiii（6.56） 这就是拉格朗日动坐标方法，它的局限性在于： l 网格点坐标随流体移动，滤失引起每个单元收缩甚至随着深入裂缝中而消失。 l 井筒必须有新的单元不断加入，这就使得控制任意时刻运用多少单元或单元的大小很困难。另外的方法： 引入平面三维模型中讨论的固定网格；它的优点在于模拟初始阶段不需要精度很高时，所需要的单元数相对较小，而随着模拟的进行单元数要增加。 引入移动网格，网格点以一定的速度移动，如裂缝始终被分成一定数量大小的单元（即采用拉伸坐标，见旁注6H）。6H.1 坐标系拉伸坐标系拉伸如果：（6H.3） x在0L(t)之间变化时，X一直在

36、01之间。网格划分简化了，但差分方程却更复杂了。导数：6H 坐标系拉伸和稳定分析坐标系拉伸和稳定分析 tLxX XLx1（6H.2）XdtdLLXtt（6H.1）方程（6.53）变为：LqXqLXAdtdLLXtA1（6H.4） 假设压力梯度为：（6H.5） 对于PKN模型，缝高是固定的，Cp=hf,，其中定义为：（6H.6） 将方程(6H.5)代入(6.53)，并运用链规则：（6H.7）6H.1 稳定分析稳定分析xACxppwp22211111xAAnxAAxAAAnDtAnnnn 为了进行误差分析，所有变量必须为绝对值，D定义为： 导数扩展为中心差分近似，Ai中的项变为：nnfnpKhCD

37、1112 方程(6H.7)中的最高项为：221111xAAxAAAnDnnnnnAxAxnD131122（6H.10）（6H.9）（6H.8） 为将A用A（1+）代替， A（1+）可近似为（对于小量）：nAxAxnDnn131213112（6H.11） 如取一时间步长，那么增加为：nAxAxnDtEnn131213112（6H.12） 为了减小误差，它必须小于A, 仅当下式成立：vpCACxt22（6H.13）式中：Cv粘度滤失控制系数213nxpvnCv（6H.14） 显式有限差分的一个主要局限：在计算中所用的时间步长不能超过某一临界值，以确保稳定性。隐式有限差分方法没有时间步长的限制显式差

38、分和隐式差分的区别： 显式差分法仅在上一时间步的基础上，求当前步长的解 隐式差分法只用当前值。要得到当前时间步长的所有变量值，就需要建立一组方程并解方程。 隐式差分法 对于线性问题，可用高斯消去法得到解；对于一维流动问题，从隐式有限差分法得到三角方程组。但对非线性问题，这些方法可能比较复杂。 积分或分析元法：既没有时间步长的限制，又能避免形成一组方程。商业时间共享方法PKN模型当在缝端x =的基本方法：Ehpwfnet2364whqdxdpf将方程（6.58）中的p代入方程（6.57）得：（6.57）（6.58）3642whqdxdwhEff（6.59）用 uhfw/4代替q上式变为：uwEh

39、dxdwf232（6.60）或：udxEhdwwf322（6.61）在距离x上积分得： 3/1396 xxxfudxEhxwxxw（6.62）积分项假设为常数，它进一步简化为： 3/1396xuEhxwxxwf（6.63） 对于高度不固定，非牛顿流体，用幂律流变参数写出与此相似的方程： xuEKhnAAnnnfntiti3212式中： dywywnhnnnhff2142（6.64）（6.65）每一时间步长求解的方法：1、估算缝端速度2、对于从缝端至井筒的每一单元 在单元内边确定速度，这样时间步长中的流体滤失和体积变化达到质量平衡。 根据方程（6.64）估算流体流速确定单元内部的横截面积。 根据

40、单元外边界速度和估算的单元内部的速度计算平均流速（第一次迭代中，假设内、外流体流速相同）3、将流入裂缝的实际流速遇上一步迭代计算的井筒流体流速作比较。4、用NewtonRaphson方法对估算的缝端速度进行修正直到达到体积平衡，该过程要24步6.3.2.4 非平衡高度解非平衡高度解 如果裂缝延伸至高渗层或隔层应力不足以遮挡裂缝，裂缝垂向延伸很快的情况下，由于流体垂向流动形成的压力梯度会变得较大，平衡高度假设不再成立。 对于非平衡高度增长，压力梯度要根据缝高增加速度估算。 根据缝端分析解（Lenoach 1994），得到净压力：nnnnnfnetnnnnhEKnnEp2112112sin1cos

41、 2222（6.67）式中：utip 端部速度 = 2 /（2+n） 常压下的裂缝，应力集中系数与净压力的关系：2fIcnethKp（6.68） 将上两式合并， Lenoach 法可求解由于端部速度非零时的表观断裂韧性，这种影响可以纳入实际的岩石断裂韧性中，将两者的和代替实际岩石断裂韧性用于方程（6.40）和（6.41）来确定缝高增长。 在一时间步长中，由一对垂向缝端位置移到另一相应位置的方法：估算单元顶和底部的缝端速度用估算的速度计算该时间步长结束时的新的缝端位置由方程（6.40）和（6.41）计算应力集中系数确定应力集中系数（即计算值减去岩石断裂韧性）用方程（6.67）和（6.68）计算产

42、生该附加应力集中系数所需的速度将该速度与估算的速度对比，继续迭代直到得到正确的速度6.3.2.5 横向耦合横向耦合 在固体力学解中的一个假设条件是：每一横截面独立作用，它隐含在任意一点的压力和宽度是对应的假设条件中。 实际上，任一点的压力不仅与局部缝宽有关，而且它与缝宽在整个裂缝的分布有关。 除非半缝长小于缝高，否则这种横向的耦合一般不是很重要；如果忽略横向耦合，裂缝几何形状不会有很大区别，然而估算的压力偏低。 在泵注过程中横向耦合的影响就是增加井筒与近井筒压力而减小缝端压力。图6.11表示的就是KGD、PKN和横向耦合的PKN模型模拟的缝高控制的压裂施工中的压力变化曲线。图6.11 有、无水

43、平方向耦合情况下的压力由图看出：横向耦合模型预测的压力总比其他两个模型预测的压力高横向耦合模型预测的宽度比其他两个模型小横向耦合模型的最低压力点（KGD和PKN模型此时压力相等）对应于正方形，即缝长等与缝高的一半；此时由横向耦合模型计算得到的压力比其他两个模型计算得到的压力高约40%。 6.3.3 块体拟三维模型块体拟三维模型块体模型的本质： 它是非常简单的模型 该模型能否成功的运用取决于分析问题中系数的恰当选择描述水力压裂过程的方程： 质量守恒方程 方程（6.33） 张开裂缝的分布情况以及净压力分布 xdxwxxIxpLLnet,（6.69） 动量守恒方程 pwqqmnm1241（6.69）

44、式中： 流道系数 m 幂律流体系数， 用于紊流 n 既考虑牛顿流体，又涵盖紊流假设：裂缝形态自相似 由两个水平延伸相同而垂相延伸不同的半圆裂缝组成。利用：空间平均法将上述方程化为：与时间有关的普通微分方程块体方程用于KGD模型质量平衡方程：dtLwdLqwqL3（6.72）式中：LEpwnet1 2225Lwwqnmmm而：1425Errm（6.71）（6.70）（6.73） 这些方程非常简单，但系数值不一定确定，而且还不一定是 常数。模型精度的高低主要取决于确定系数时的工作细致程度。 净压力曲线净压力曲线 是判断裂缝延伸模式的依据：是判断裂缝延伸模式的依据： 线段线段1压力的斜率在压力的斜率

45、在0.1250.2之间。说之间。说明裂缝在预期的缝高和综合滤失系数下向地层深明裂缝在预期的缝高和综合滤失系数下向地层深处延伸。处延伸。 线段线段2压力的斜率不变。可能预示着注入压力的斜率不变。可能预示着注入量等于滤失量，也可能是裂缝在长度上已停止延量等于滤失量，也可能是裂缝在长度上已停止延伸，或裂缝高度即将失控，或缝内即将出现堵塞。伸，或裂缝高度即将失控，或缝内即将出现堵塞。 压裂模型的识别压裂模型的识别 线段线段3压力的斜率为压力的斜率为1。反映了缝内。反映了缝内发生堵塞，裂缝在长度上已停止延伸，注入发生堵塞，裂缝在长度上已停止延伸，注入的液体只能增加裂缝的宽度。的液体只能增加裂缝的宽度。

46、线段线段4压力的斜率为负值。说明裂缝压力的斜率为负值。说明裂缝在高度上已失去控制，延伸到非压裂目的层在高度上已失去控制，延伸到非压裂目的层段，或又压开了新的裂缝，或裂缝在延伸过段，或又压开了新的裂缝，或裂缝在延伸过程中遇到了规模较大的天然微裂隙体系。程中遇到了规模较大的天然微裂隙体系。 压裂施工设计中的一个关键问题是准确了解流体经裂缝滤入地层的速度，没有这方面的知识而要得到一定裂缝几何形状的施工设计是不可能的。小型压裂施工用于估算滤失系数。 Carter、Williams 和 Settari 对滤失的划分提供了三种不同的标准。将滤失划分为三个不同的滤失过程： 油藏流体的移位和压缩 滤液或压裂液

47、的地层侵入 外部滤饼的形成将滤失划分为三个时间段： 初始阶段 未形成滤饼 第二阶段 外部滤饼逐渐形成 最后阶段 动态滤失滤饼与初滤失 滤饼: 压裂液的造壁性和降滤剂共同作用的结果 初滤失: 形成滤饼前的滤失PwPvPcPs 一般滤失模型： 油藏流体的位移和压缩 滤液和压裂液的地层侵入 低渗透岩石外部或高渗岩石内部滤饼的形成 滤饼区的流动滤饼控制过程 侵入区的流动压裂液粘度控制过程 地层流体的压缩地层流体粘度及压缩性控制过程室内定压条件下的滤失实验的假设条件： 沉积的滤饼量与单位面积穿过的流体体积VL成正比。 滤饼渗透率K饼与厚度无关。 流经滤饼遵从达西定律。 通过滤饼的压力降P饼保持不变 (6

48、.74) LfiltLVpKddV饼饼Lu 其中fil为滤液的粘度，假设滤饼厚度L饼与滤液的流体体积成正比，即VL=L饼。对VL积分可以看出： (6.75)由此导出通过井壁滤饼的流体滤失系数为： (6.76) filwPKC2饼饼tCw/uLCarter 认为滤失的流体体积可按下式确定： (6.77) 式中，SP为没有形成滤饼情况下的滤失体积。更恰当的解释结果为： (6.78) 式中，tsp为初损时间。PWLStCV 2pspWLSttCV 2 油藏中的初始带被称为滤失带或侵入带，它假设： 通过该带的压力降 P是常数。 滤液完全取代了地层的流动相，造成了似柱塞式位移以及100%滤液饱和。 流体

49、和岩石是不可压缩的。表达式为： (6.79) 式中，Kfil与滤液有关的渗透率 Lv 侵入带长度。LpKdtdVfilfilLLu 对上式积分并假设： (6.80) 其中 为孔隙度，从而得到： (6.81) 此时粘度控制滤失系数Cv为： (6.82) 任意时刻的滤失体积为： (6.83)dtVLLLu1/tCLufilfilpKC2tCVL2 未受滤液侵入的油藏不含有压裂液，要将油藏流体由裂缝面驱替需要一定的压力，假设： 滤液/油藏界面与远场油藏间存在固定压差Pc 。 可压缩流体具有不变的总体压缩系数Ct。 侵入流体的前缘推进较慢。 无限油藏。 侵入流体前缘可视为无限孔隙介质面，可由分析解（C

50、ollins，1961）得到： (6.84)其中压缩性控制的滤失系数Cc为： (6.85)式中，Kr油藏岩石的渗透率； 油藏流体粘度。CrtrcPCKCtCLur 实际上，三个过程是同时发生的，方程（6.74）（6.79）和（6.84）中的滤失速度必须相同，压力降之和也要与油藏压力和破裂压力之间的总压力差相等： (6.86) (6.87)式中，Ct总滤失系数； 缝中压力与远场油藏压力pr间的压力差。crtrfilfilspfiltpCKtpKttpKtC22饼饼总饼ppppr总p 如果忽略初损体积和初损时间，这些方程综合（Williams等，1979）得到总滤失系数： (6.88)其中系数Cw

51、、Cc和Cv是由总压力差计算得到的。 方程（6.88）仅当滤饼渗透率与压力无关时才成立，如果滤饼压缩性强，其渗透率近似于1/p成正比，Nolte（1988）研究表明流体滤失受滤饼或油藏限制。2222242wvcvWwvwvcwcvtCCCCCCCCCCCC 在这种情况下，流体滤失速率是方程（6.75）的最小值，其中压力降为总压力降，或： (6.89)其中综合滤失系数为： (6.90)其中系数Cc和Cv也是用总压力差计算得到的。tCcvLu2242cvvvccvCCCCCC 如果将侵入流体的压缩性以及油藏流体移动边界等一些更详细的特性考虑在内的话，在滤失模型中可以加入很多复杂因素。当在水力压裂施

52、工中的其他已知参数精度给定的情况下，没有必要再将上述这些影响考虑在内。本节讨论将压力差变化和流体性质改变这两种影响包括在内的模型修正。 Settari（1985）发现侵入带的流体滤失可以用与Cv相当的量 代替Cv来描述： (6.91)式中，Cv用当前所有滤液的平均粘度和平均相 对渗透率计算； VL在这之前滤失到地层的体积。vCLvvVtCC22 Settari还发现以相当值替换当前流体滤饼形成系数，即可解释压裂施工中不同流体滤失特性的区别。在这种情况下，该当量值为 (6.92) 此时先前滤失的流体体积VL也会影响到滤饼的形成。LwwVtCC22 影响滤饼形成的关键流体成分是胶体和/或流体滤失添

53、加剂的浓度。Settari模型的扩展可以在方程（6.74）基础上导出，在该模型 中假设滤饼厚度不再与通过的流体体积成正比而是与沉积胶体体积成正比，这样滤饼不是与时间积分的速度成正比，而是与整个时间段中胶体浓度和流体速度的乘积的积分成正比。这时方程（6.92）由下式代替： (6.93)式中 Cgel 流体中胶体质量浓度； Mgel胶体相对密度。 上式中隐含着CgelCW2为常数这个假设条件，该假设条件与Mayerhofer等（1991）的室内研究结果是一致的。水的粘度随温度的升高而降低，其对滤失系数的影响应包括在压裂模型中。gelWgelWMtCCC22 施工中压力演化 对压力敏感的滤失 孔隙弹

54、性和背应力 水力压裂的目的是在地层中铺置支撑剂从而形成一条油气生产通道，支撑剂在水力裂缝中的出现对流体的性质产生三个主要影响： 支撑剂对流体流变性的影响 对流或重力流 支撑剂输送 一般来讲携砂液的粘度slurry比仅有压裂液时的粘度slurry要高，对于牛顿流体两者之间的关系已通过实验建立，。Nolte（1988）发现对于幂律流体这种相互关系可以用牛顿流体的关系： (6.94)式中， fv 支撑剂体积百分数； fvM 流动砂浆最大百分比。nvMvbaseslurryrff5 . 2/11 流体间密度的差异可能是由于密度较大的流体在密度较小的流体下流动，或密度较小的流体在密度较大的流体上流动造成

55、的。在水力压裂中，当高密度砂浆泵注在前面泵入阶段的流体中或前置液下流动时，以及如Clifton和Wang（1988）考虑的流体流动的其他二维问题就是这种情况。牛顿流体的流动方程可写为： (6.95) (6.96)以及（对于不可压缩流体情况下的方程（6.32）： (6.97)将方程(6.95)、(6.96)代入方程(6.97)得： (6.98)xp12u2x0u2yuxuLyxtygtypxp3L232232121u2121gtp12u 2y 方程（6.98）中右端最后一项是对流项，正如其他两项由于缝宽变化或滤失造成的存储或沉降项一样。Baree和Conway（1994），Vnwin和Hammo

56、nd（1995）以及Smith和Klein（1995）发现如果泵入过量前置液，泵注结束后流体流动可能会引起对流直至前置液滤失掉。方程（6.98）同时表现出了对流对缝宽相当敏感，如果缝宽较大，对流可能更关键。值得庆幸的是低杨氏模量的地层一般是高渗的，这种情况下端部脱砂设计和停泵后滤失快有效地避免了对流。l水力压裂的核心是形成满足导流能力要求的填砂裂缝，支撑剂在裂缝的沉降影响到填砂裂缝几何尺寸和裂缝导流能力。l内容： 支撑剂的沉降特性 沉降布砂设计 悬浮布砂设计l自由沉降： 指单个颗粒在无限流体介质空间内的沉降。l干扰沉降：颗粒群在有限流体空间内的沉降。一、 支撑剂的沉降特性（一）、 单颗粒自由沉

57、降速度固体颗粒在流体中沉降时的受力）固体颗粒的重力）流体对固体颗粒的浮力）固体颗粒在流体中运动时所产生的阻力重力：浮力：阻力：gdFPPg36gdFfPb362228 21PPfDPfDdUdCAUCF213)(4fDfPpPCgdU F=Fg-Fb当FFd时，可得： CD与雷诺数有关，雷诺数与UP有关fPpUdNRe18)(2gdUfPPP 当NRe2时，斯托克斯区域，Re24NCD43.029.044.171.0)(34.20fPfPPdU 当2NRe500时，过渡区，6 .0Re5 .18NCDfPfPPdgU)(74.1 当NRe500时，牛顿区域，44.0DCl支撑剂浓度的影响l壁面

58、的影响l颗粒形状的影响l支撑剂浓度对沉降速度的影响（存在粒间干扰，UHUP） 当NRe2时， UH/UP=5.5 当 2Nre 500时， UH/UP=3.5 当 NRe 500时， UH/UP=2)1(10/ .82.12PHUU)1 ( 9 . 5)/ln( .PHUU)1 (88. 611/ .PHUUl壁面的影响 UH=FWUP FW壁面拉拽系数 当 NRe1时，5430644. 0131. 0147. 06526. 01wdwdwdwdfw当 NRe100时，23)2(1WdPFW当 1NreKIc时，裂缝向前延伸。其中： (6.105)式中，KIc临界应力集中系数的材料临界值； F

59、比裂缝表面能。该理论的优点： 能将复杂的破裂过程用一个参数代表； 假设KIc是与杨氏模量或强度相似的材料常数212FIcEK Irwin（1957）发现二维张性裂缝的缝端附近区域应力分布存在奇异性： （6K.1） （6K.2） （6K.3）式中，在裂缝中心线测得的夹角； r至缝端的距离； KI应力集中系数。) 2/3sin() 2/sin() 2/cos(2) 2/3sin() 2/sin(1) 2/cos(2) 2/3sin() 2/sin(1) 2/cos(2rKrKrKxyyx 通过对由方程（6K.2）和（6K.3）得到的应力场与从-L到L二维裂缝缝端应力场得比较，Rice发现KI可由下

60、式计算： （6K.4）式中，p（x）缝中压力分布。 对于缝中压力为常压的情况，方程（6K.4）简化为方程（6.104）。 对于径向裂缝，相应方程为： （6K.5） 对于常压情况，上式简化为方程（6.105）。dxxLxLxpLKLLI)(1drrRrrpRKRI022)(2 假设KIc为材料常数，不考虑尺寸效应，对方程（6K.4）积分，将断裂韧性加入二维或平面三维模型中得到应力集中系数KI。对于计算得到的KI值： KI大于KIc，裂缝可以继续向前延伸； KI小于或等于KIc，裂缝会维持在其当前位置； KI为负值，裂缝会后退直到KIc0.在水力压裂中，实测的净压力偏高的原因： 缝高延伸控制比预计