复杂裂缝中低密度支撑剂铺置数值模拟

任 岚 林 辰 林 然 刘 军 胡哲瑜 王 潘

（西南石油大学油气藏地质及开发工程国家重点实验室， 四川成都 610500）

0 引 言

近年来， 全球非常规油气勘探开发取得了一系列重大突破， 非常规油气产量快速增长， 在全球能源供应中的地位日益凸显［1］。 缝网压裂作为非常规储层开发的一种重要技术手段［2］， 其首要目标是经济有效地建立油气从地层到井底的高导流能力缝网， 而支撑剂在缝网中的有效铺置决定了压裂增产的最终效果。 滑溜水作为非常规储层开采最常用的压裂液， 携砂能力本身较弱， 加上复杂裂缝中分支缝间的阻力效应强， 常规支撑剂迅速沉降， 导致其在复杂裂缝内运移铺置效果并不理想， 从而影响导流能力。 低密度支撑剂沉降速率慢的特点使得在复杂裂缝中运移性更好， 在裂缝闭合前， 支撑剂处于悬浮状态， 从而优化支撑剂剖面， 达到更好的铺置效果。 目前低密度支撑剂在缝网压裂中的应用日益增多， 而低密度支撑剂在复杂裂缝中运移铺置规律研究较少， 数值模拟相对于物理实验方法成本更低， 效率更高。 因此， 开展低密度支撑剂在复杂缝网中的铺置模拟研究非常重要。

目前国内外针对支撑剂铺置规律的数值模拟研究多基于在单缝中的常规支撑剂。 徐暖筑等［3］建立二维模型研究了2 600 kg／m3的支撑剂在单直缝内的铺砂剖面与输送过程； 张涛等［4］建立欧拉两相流模型研究支撑剂在清水压裂中的输送沉降行为， 对支撑剂密度为2 000 kg／m3与2 500 kg／m3运移进行了模拟； Sahai 等［5］研究了普通陶粒在不同的裂缝复杂度、 排量、 浓度、 粒径等因素中的运移情况； Kong 等［6］基于双欧拉模型， 模拟了常规支撑剂在不同类型流体及T 形缝网和交叉缝网裂缝模型中的沉降过程； 刘平礼等［7］运用CFD 数值模拟软件， 研究了单缝中不同影响因素下常规支撑剂（2 800 kg／m3） 铺置运移规律； 梁莹等［8］通过可视化平板物理实验研究了低密度支撑剂在不同施工参数下的砂堤铺置形态； 韩琦［9］通过CFD-DEM方法分析了超低密度支撑剂（1 055～1 750 kg／m3）在裂缝中的输送过程， 但该文只针对单缝模型， 未能研究复杂裂缝中的传输规律。

鉴于此， 基于计算流体力学（CFD）， 建立了欧拉—欧拉液固两相流模型， 运用有限体积法进行求解， 通过物理实验验证该数值模型的准确性， 对比常规支撑剂， 分析了低密度支撑剂在复杂裂缝中沉降运移规律， 以及低密度支撑剂在不同泵注排量、 砂比、 压裂液黏度和裂缝夹角下的铺置效果，从而达到优化施工参数目的， 为现场压裂施工设计提供帮助。

1 低密度支撑剂液固两相流模型

低密度支撑剂在裂缝中的传输过程， 实际为复杂的液固两相流动状态， 而欧拉法（双流体模型）将固相颗粒看作拟流体， 将液体与固相颗粒都放在欧拉坐标下处理的方法， 大大减小了数值计算时间。 该方法是目前解决工程多相流问题常用的方法。

1.1 流体控制方程

连续性方程为：

式中：α——体积分数；

ρ——密度， kg／m3；

∇——哈密顿算子；

t——时间， s；

v——速度， m／s；

下标s、 l 分别代表固相、 液相。

动量守恒方程为：

式中：g——重力加速度， m／s2；

p——分压， Pa；

τ——剪切应力张量， Pa；

β——相间动量交换系数， kg／ （m3·s）。

湍流方程。 湍流影响通过k—ε湍流模型加以描述， 该模型的控制方程包括湍动耗散率ε方程、湍动能k方程， 即：

式中：k——连续相的湍动能， m2／s2；

ε——湍动能k的耗散率， m2／s3；

μt——连续相黏性系数， Pa·s；

σk， σε——湍动能和耗散率对应的普朗特数， 分别取1 和1.3；

Gk，l——湍动能的产生项， kg／ （m·s3）；

Πk，Πε——两 相 间 湍 流 交 换项， kg／ （m·s3）；

C1ε，C2ε——经验常数， 分别取 1.44和1.92。

1.2 相间作用力方程

在液固两相流模型中， 液相会传递给固相多种作用力， 包括曳力、 压力梯度力、 浮力等。 本文主要考虑压裂液对低密度支撑剂作用较大的浮力和曳力以及颗粒相因加速运动对流体相产生的虚拟质量力。

相间动量传递项为

其中：

式中：CL——升力系数， 取0.25；

CVM——虚拟质量力系数， 取0.5；

μl——流体运动黏度， m／s2；

d——颗粒相直径， m；

αs、αl——固相、 液相体积分数，；

ρl——密度， kg／m3；

vl，vs——液相、 固相速度， m／s；

曳力基于固液两相之间的动量交换， 采用Gidaspow 模型， 此模型适用于稠密流化床模拟中，同时支撑剂铺置过程也使用此模型［10-16］。 液固两相间的动量交换系数β为

其中

式中：CD——固相颗粒曳力系数；

Res——两相滑移速度的雷诺数。

2 裂缝模型及模拟方案

在非常规储层中通常使用体积压裂的方式形成复杂缝网， 本裂缝物理模型基于常见的3T 缝网形态， 利用相似准则建立主缝缝长2 m、 缝高0.6 m、缝宽0.01 m 的三维缝网模型［17-20］（图1）。 在几何模型左端边界设3 个等距离进口， 作为压裂液与低密度支撑剂颗粒入口， 右端及分支缝边界设置为排液出口， 模型采用全6 面体网格划分。 设置左进口为速度入口边界， 主缝右端出口及各分支缝出口为压力边界， 固壁边界为光滑静止无滑移边界。 初始化采用进口边界条件， 裂缝模型中压裂液体积分数为100%， 控制方程离散格式采用一阶迎风格式，采用相耦合的SIMPLE 算法进行迭代计算， 收敛标准为各残差小于10-4［21-23］。

通过数值计算模拟不同情况下砂堤形成的过程， 采用控制变量法针对不同的支撑剂类型、 泵注排量、 砂比、 压裂液黏度和二级缝夹角对比分析支撑剂在复杂裂缝中的铺置情况［24-26］。 具体模拟方案见表1。

表1 模拟方案Table 1 Schemes of the numerical simulation

3 低密度支撑剂铺置过程模拟结果

3.1 实验对比验证

通过大型可视化平板单缝模型物理实验进行对比验证， 实验平板与裂缝模型主缝的参数一致。 实验选择支撑剂为30／50 目自悬浮低密度支撑剂， 平均密度为1 300 kg／m3， 采用黏度为2 mPa·s 压裂液， 以10%的砂比， 0.02 m3／min 的排量注入， 根据雷诺数相似准则， 对应现场泵注排量5 m3／min。采用表1 中序号5 的模拟结果对比， 参数一致得到两者平衡时的铺砂剖面如图2 和图3 所示。 实验的砂堤长度为1.8 m， 平衡高度为0.18 m， 平衡时间105 s， 数值模拟得到砂堤长度为2 m， 平衡高度为0.2 m， 平衡时间79 s。 在实际运移过程中， 由于存在颗粒间作用力， 壁面效应等因素， 使得物理模型与数值模型模拟运移时间节点上存在差异， 但两者模拟的剖面和平衡参数基本一致， 表明建立的欧拉—欧拉两相流数值模型能够较为准确地描述低密度支撑剂在裂缝中的铺置情况。

3.2 与常规支撑剂铺置过程对比

对方案表1 中序号2、 9 进行模拟， 得到平衡参数表2， 在相同砂量下， 低密度支撑剂在复杂缝网中铺置更均匀， 能运移到裂缝深处， 传输能力更好。 图4 和图5 为低密度支撑剂（ρ=1 300 kg／m3）与常规支撑剂（ρ=2 600 kg／m3） 在30 s 时铺置情况。 在常规低排量的铺置过程中， 低密度支撑剂相对于常规支撑剂， 其上方形成大面积的悬砂区域，只有部分支撑剂达到饱和状态沉降形成砂丘， 单从常规的平衡高度来评价， 不能说明支撑剂铺置效果的好坏。 因此， 本文采用支撑剂传输前缘浓度梯度角（αc） 与平衡时悬砂浓度梯度高（hc） 来评价支撑剂铺置效果。 支撑剂传输前缘浓度梯度角越大， 支撑剂传输越接近活塞流， 流动能量越大， 压裂液的携砂性能就越好。 为了更好对比， 取用支撑剂前缘刚到缝口时的角度进行评价。 由于低密度支撑剂悬砂部分在裂缝闭合后也具有导流能力， 可以通过平衡悬砂浓度梯度高评价支撑效果。

表2 不同支撑剂在复杂裂缝中平衡参数Table 2 Balance parameters of different types of the proppant in complex fractures

从图4、 图5 中看到， 两者在相同砂量情况下， 常规支撑剂以沉降铺置为主， 分支缝对支撑剂运移阻力大， 其沉降位置主要集中在分支缝前的缝口处以及二级缝缝口处。 低密度支撑剂以悬浮状态铺置为主， 受到的阻力效应小， 在裂缝中沉降较少， 有更多的支撑剂被压裂液携带到裂缝深处中沉降， 在复杂裂缝中铺置更均匀。 同时低密度支撑剂体系前缘浓度梯度角更大， 说明该体系在裂缝末端仍具有较大的流动性， 携砂性能更好。 但由于自身重力低， 在压裂液带动下形成悬浮状态， 沉降速度较慢， 平衡时间较长。 因此， 在复杂缝网中， 采用低密度支撑剂进行铺砂效果更佳， 能有效地支撑深处裂缝。

3.3 泵注排量的影响

选取模拟方案表1 中序号1、 2、 3 进行研究，分析入口速度分别为0.5、 1、 2 m／s 时对低密度支撑剂铺置的影响， 得到平衡时参数见表3。 随着排量的增加， 主缝中hc有所降低， 更多的支撑剂被携带到分支缝中， 特别是在三级缝中填充增加明显。 图6、 图7 为30 s 时不同速度下主缝以及分支缝铺砂情况， 入口速度为0.5 m／s 时， 低密度支撑剂在缝口端堆积较多， 到达分支缝时， 因阻力作用， 悬浮支撑剂迅速沉积， 分支缝中支撑剂多以滑落进入； 入口速度为1 m／s 时， 主缝铺砂较为均匀， 分支缝口的阻力效应明显减弱， 悬浮支撑剂部分沉积； 入口速度为2 m／s 时， 低密度支撑剂大量悬浮在主缝和二级缝中， 三级缝也形成较好的砂堤， 但当入口速度提升后， 湍流效应增强， 支撑剂被卷起， 导致入口处砂堤高度减小， 沉降不稳定。

表3 不同泵注排量下的平衡参数Table 3 Balance parameters at different pumping rates

综上， 提高排量可减缓分支缝的阻力效应， 能使支撑剂更好地铺置分支缝。 但大排量的湍流效应会导致主缝入口处砂堤高度减小， 铺置不均匀。 在现场施工时， 建议初期采用大排量、 高黏压裂液携带低密度支撑剂铺置缝网远端， 后期可以采用大排量、 低黏度尾追中—高密度支撑剂铺置裂缝近端。

3.4 砂比的影响

选取模拟方案表1 中序号2、 4、 5 进行研究，分析5%、 10%、 20%的砂比对低密度支撑剂铺置的影响， 得到平衡参数见表4。 初期提高砂比， 支撑剂在裂缝中hc增加明显， 随着砂比的继续增加，支撑剂扰动效应使得hc增加速率下降，αc逐渐降低， 传输动能减少使得平衡时间逐渐减少。 图8、图9 为20 s 时的复杂裂缝砂堤展布。 对比砂比5%和10%， 主缝中已沉降砂堤厚度与悬浮态支撑剂明显增多， 分支缝中砂堤含量也更多， 表明随着砂比增加， 支撑剂铺置效果更好。 对比砂比10%和20%， 主缝与分支缝中已沉降的砂堤变化不大， 由于缝口支撑剂颗粒间扰动效应增强，αc更小， 支撑剂推进动能降低， 集中铺置裂缝前段， 而低密度支撑剂更易卷起的特性， 导致缝口处悬浮态支撑剂明显增多。

表4 不同砂比下的平衡参数Table 4 Balance parameters at different sand ratios

综上， 在一定排量下， 适当提高砂比会使低密度支撑剂在复杂裂缝中铺置效果更好， 但过高的排量会导致缝口悬浮态支撑剂扰动效应增强， 动能降低， 不利于裂缝深处的铺砂。

3.5 压裂液黏度的影响

对方案表1 中序号5、 6、 7 进行模拟计算， 研究2、 5、 10 mPa·s 的压裂液对低密度支撑剂铺置的影响， 根据平衡参数表5 可知， 黏度越大， 裂缝中hc越大， 压裂液携砂能力更强。 而αc越大， 说明支撑剂在高黏压裂液中运移性更好， 沉降慢， 使得平衡时间增加。 图10、 图11 为35 s 时不同黏度复杂裂缝中砂堤分布， 可以看出随着压裂液黏度的增加， 悬浮态支撑剂占主体， 黏度5 mPa·s 时，支撑剂几乎全部悬浮在主缝和分支缝中， 以一种稳定的形式向缝网深处推进。 但随着压裂液黏度增加， 低密度支撑剂在缝口处基本无沉降砂堤， 若长时间注砂可能导致缝口无砂， 使得铺砂不均匀。 因此， 采用高黏度压裂液时， 可以适当降低排量使支撑剂铺置更加均匀。

表5 不同压裂液黏度下的平衡参数Table 5 Balance parameters at different viscosities of the fracturing fluid

综上， 高黏压裂液可使低密度支撑剂运移性更好， 在现场施工采用高黏度压裂液对远端缝网铺砂后， 可适当减小泵注排量， 使缝口支撑剂铺置更均匀。

3.6 裂缝夹角的影响

对方案表1 中序号5、 8 进行模拟计算， 研究不同裂缝夹角对低密度支撑剂铺置的影响， 根据表6 可知， 平衡时裂缝夹角60°比90°主缝中hc稍小，而二级缝hc稍大，说明在相同情况下分支缝角度越小， 对支撑剂运移的阻力越小， 在复杂裂缝中铺置更均匀， 而平衡参数差别不大， 说明裂缝夹角对低密度支撑剂的运移影响较小。 图12 为40 s 时裂缝砂堤分布， 在泵注排量不大的情况下， 主缝和分支缝中支撑剂分布基本相同， 而分支缝本身具有分流的作用， 会导致支撑剂在铺置过程中不均匀性增加， 特别是在裂缝相交处， 常规支撑剂堆积过快可能会发生砂堵现象， 采用低密度支撑剂可以减缓沉降， 有效避免砂堵的发生。

表6 不同裂缝夹角的平衡参数Table 6 Balance parameters at different fracture angles

4 结 论

（1） 对比物理模型和数值模型， 两者模拟的铺砂剖面和平衡参数符合度高， 表明欧拉—欧拉两相流数值模型能够较为准确地描述低密度支撑剂在裂缝中的铺置运移情况。

（2） 低密度支撑剂体系相比于常规支撑剂体系在复杂裂缝中整体铺置效果更好， 分布更均匀。

（3） 适当的提高泵注排量， 增加压裂液黏度，选择合适的砂比能使低密度支撑剂在复杂裂缝中运移性更好， 而裂缝夹角对低密度支撑剂铺置影响较小。