裂缝扩展与拉格朗日支撑剂运移的双向耦合数值模拟方法

本发明涉及油气田开发工程水力压裂，特别涉及一种裂缝扩展与拉格朗日支撑剂运移的双向耦合数值模拟方法。

背景技术：

1、水力压裂是开发页岩油气藏、致密砂岩油气藏等非常规油气资源的核心主体技术。水力压裂的关键在于形成具有一定裂缝导流能力的支撑裂缝，从而为油气从储层流入井筒提供“高速通道”。在水力压裂作业中，滑溜水压裂液体系因其低储层损伤、低作业成本、强造缝能力已被广泛用于非常规油气储层改造过程。然而，滑溜水压裂液因黏度低导致其携砂性能差，常常造成水力裂缝支撑不足，制约了水力压裂效果，成为非常规油气开发的一个瓶颈问题。因此，研究滑溜水压裂过程中支撑剂在裂缝中的运移行为和铺置特征可为水力压裂施工参数优化提供重要理论支撑与技术指导。

2、目前，绝大部分支撑剂运移数值模拟研究假设裂缝壁面固定(静态裂缝)，与实际水力压裂存在的裂缝扩展过程有着较大差异。水力压裂涉及流体速度场、流体压力场以及固体应力场等多场之间的耦合。离散元(discrete element method)、有限元(finiteelement method)、位移不连续法(displacement discontinuity method)以及上述三种方法的变体和组合是模拟裂缝(变形)扩展最常见的数值方法。相较于离散元和有限元方法，位移不连续法具有网格数量较少、数据准备简单、求解速度较快的显著优点。更重要地，位移不连续法以裂缝开度为直接求解未知量，为裂缝中支撑剂运移求解提供了便利。此外，多相质点网格法已被成功应用于工程尺度(百米量级长度裂缝)的支撑剂运移数值模拟，具有较大的应用前景。因此，位移不连续法与多相质点网格法的结合是开展裂缝扩展和支撑剂运移研究的一个重要研究方向。

3、然而，现有基于位移不连续法的耦合模型大都采用浓度输运模型刻画缝中支撑剂运移过程。浓度输运模型是一种欧拉-欧拉方法，其会破坏颗粒相的“离散特性”。其次，浓度输运模型适用于高黏压裂液携砂情况，而不适合于低黏压裂液(滑溜水)工况。再次，浓度输运模型需要使用拟合实验结果的经验公式，相关经验公式的适用性与可靠性仍有待证实。

4、另外，虽然现有技术基于多相质点网格法实现了支撑剂运移与pkn型裂缝扩展的单向耦合(one-way coupling)。但是，该单向耦合模型实际上并未严格实现流体流动场和固体应力场之间的耦合，即双向耦合(two-way coupling)，存在一定的局限性。此外，该单向耦合模型难以实现在具有地应力差异储层中的水力裂缝扩展及支撑剂运移数值模拟。

技术实现思路

1、针对上述问题，本发明旨在提供一种裂缝扩展与拉格朗日支撑剂运移的双向耦合数值模拟方法，基于位移不连续法和多相质点网格法实现裂缝扩展与支撑剂运移的高效双向耦合，能够实现工程尺度的裂缝扩展与拉格朗日支撑剂运移模拟，为水力压裂泵注程序优化提供数值工具。

2、本发明的技术方案如下：

3、一种裂缝扩展与拉格朗日支撑剂运移的双向耦合数值模拟方法，包括以下步骤：

4、s1：定义数值模拟区域，输入基本参数；

5、s2：建立压裂液流动控制方程与平面三维位移不连续法的水力裂缝变形控制方程，并对所述压裂液流动控制方程与所述水力裂缝变形控制方程进行耦合求解，获得流场参数与裂缝宽度；

6、s3：判断裂缝是否扩展：若扩展，则在现有缝尖单元相邻位置添加新的位移不连续单元；若不扩展，则保持裂缝长度不变；

7、s4：在预定时间节点向数值模拟区域放置支撑剂颗粒包，建立支撑剂运移控制方程；

8、s5：计算所述支撑剂颗粒包所受颗粒应力项、拖曳力项、压力梯度项、净浮力项以及碰撞阻尼项；

9、s6：求解所述支撑剂运移控制方程，更新所述支撑剂颗粒包的动力学属性；

10、s7：计算流体-颗粒相互作用项，并根据计算结果更新步骤s2的流场参数；

11、s8：重复步骤s2-s7，直至达到预设模拟时间。

12、作为优选，步骤s2中，所述压裂液流动控制方程包括：

13、质量守恒方程：

14、

15、式中：w为裂缝宽度，m；θf为流体体积分数，m3/m3；t为时间，s；为散度算子；uf为压裂液流动速度矢量，m/s；cl为压裂液滤失系数，m/s0.5；τ(x)为某一位置x压裂液开始滤失的时间，s；

16、动量守恒方程：

17、

18、式中：μf为压裂液粘度，pa·s；为梯度算子；p为压裂液压力，pa；ffp为流体-颗粒相互作用项，pa/m。

19、作为优选，步骤s2中，所述平面三维位移不连续法的水力裂缝变形控制方程为：

20、

21、式中：为平面三维位移不连续单元i所受法向应力，pa；n为目标系统中平面三维位移不连续单元总数；为平面三维位移不连续单元j的法向位移不连续量对平面三维位移不连续单元i的法向应力影响系数，pa/m；为平面三维位移不连续单元j的法向位移不连续量，m；为原位地应力。

22、作为优选，步骤s2中，采用牛顿-拉夫森迭代方法进行耦合求解。

23、作为优选，步骤s3中，判断裂缝是否扩展时，根据缝尖单元ⅰ型应力强度因子进行判断；若缝尖单元ⅰ型应力强度因子大于等于岩石断裂韧性，则裂缝扩展；若缝尖单元ⅰ型应力强度因子小于岩石断裂韧性，则裂缝不扩展。

24、作为优选，所述缝尖单元ⅰ型应力强度因子通过下式进行计算：

25、

26、式中：kⅰ为ⅰ型应力强度因子，pa·m1/2；e为杨氏模量，pa；v为泊松比；d为裂缝尖端单元长度，m；dn为裂缝尖端单元的法向位移不连续量，m。

27、作为优选，步骤s4中，所述支撑剂运移控制方程为：

28、

29、式中：ap为支撑剂颗粒包加速度矢量，m/s2；up为支撑剂颗粒包速度矢量，m/s；dp为拖曳力系数，1/s；ρp为支撑剂颗粒密度，kg/m3；ρf为压裂液密度，kg/m3；g为重力加速度，m/s2；θp为颗粒体积分数，m3/m3；τp为颗粒应力，pa；<up>为某一支撑剂颗粒包周围所有支撑剂颗粒包的平均速度矢量，m/s；τd为碰撞松弛时间，s。

30、作为优选，步骤s6中，基于向前欧拉法、蛙跳法或预测-修正法更新所述支撑剂颗粒包的动力学属性。

31、作为优选，步骤s7中，所述流体-颗粒相互作用项通过下式进行计算：

32、

33、式中：f为颗粒相概率分布函数；vp为支撑剂颗粒包的体积，m3。

34、本发明的有益效果是：

35、首先，通过所述压裂液流动控制方程，能够使本发明实现水力压裂裂缝变形(扩展)与支撑剂运移的双向耦合；其次，本发明保留了颗粒的离散特性(拉格朗日特征)，能够灵活方便地实现较为复杂的支撑剂颗粒注入(变密度、变粒径支撑剂颗粒注入)；最后，通过颗粒包与颗粒应力的引入，能够使本发明实现实际工程尺度的裂缝(变形)扩展与支撑剂运移数值模拟。