单一环空窄间隙条件下固井液当量循环密度及其数值模拟计算

李 旭， 任胜利， 刘文成， 赵丹汇， 廖茂林， 林黎明

1中海油田服务股份有限公司油田化学事业部 2北京科技大学机械工程学院3中国科学院力学研究所流固耦合系统力学重点实验室

0 引言

对固井水泥浆体系当量循环密度(以下简称“ECD”)的准确计算是保障固井作业顺利实施的关键，因为合理的ECD将确保环空压力系统的稳定[1]。窄间隙是固井作业中遇到的一种特殊而复杂的工况，由窄间隙所引起的ECD的显著变化会大大增加压漏地层或引起井壁坍塌的风险，因此在固井作业中应该给予特别关注。

在固井水泥浆体系ECD的计算中，摩擦压降的计算精度直接决定了ECD的准确性。在通用摩擦压降计算模型[2]的基础上，国外的学者们对窄间隙环空的摩擦压降的计算提出了一系列的修正模型[3-5]，并针对固井水泥浆的不同流态[3]、不同流变模型[4]以及不同作业条件[5]分别进行了讨论，并依据现场实测数据开展了半经验模型的验证工作[6]。虽然大多数学者习惯于以钻井液为对象开展窄间隙条件下的ECD修正研究，但基本的修正思路和修正方法是可以被固井水泥浆的ECD计算所借鉴的。国内学者们对窄间隙环空的研究侧重于理论与实验的结合，明确了其具有摩阻大、井底压力高的特点[7]，同时也归纳出了一套提高窄间隙水泥浆顶替效率的有效措施[8]。除了通过理论模型来计算ECD之外，一些学者还通过有限元软件[9-10]来近似模拟固井过程，并基于此计算摩擦压降，进而推算固井水泥浆的ECD[11-12]，该方法有助于提升对固井作业的感性认识，计算流程也较为直接。此外，在油田实际进行固井设计时，常常采用专业的钻完井软件来进行固井水泥浆的ECD计算，主流的软件包括[13]：Schlumberger公司的Drill Bench软件，Holiburton公司的Landmark软件，以及PVI公司的PVI CEMPRO Plus软件，但不同软件对窄间隙条件下ECD的计算结果并不一致，而且有时差距还较大，这也给实际的固井设计带来了诸多不确定因素，急需进一步明确。

1 MR理论

以赫巴(Herschel-Bulkley)模型为例，其固井水泥浆具有如下的流变模型：

(1)

若用于宾汉(Bingham)塑性流体模型，则n=1，且K=μp，其中μp为塑性黏度,Pa·s；若应用于幂律(Power-Law)流体模型，则τy=0。

当前计算环空ECD包括由静水压力引起的ESD和摩擦压降两个部分，即：

(2)

式中:Δpf—摩擦压降,Pa;g—重力加速度，9.81 m/s2；H—井底垂直深度，m;ECD—当量循环密度，g/cm3;ESD—当量静态密度,g/cm3。

因此，环空ECD计算中最为重要的是得到摩擦压降Δpf，或对应的Fanning摩阻因子f(无量纲)，定义如下：

(3)

式中:τw—壁面剪切应力，Pa；ρ—流体密度，kg/m3;Dw和De—分别为环空外径和内径,m;L—井段长度,m;V—剖面平均流速,m/s。

V定义为：

(4)

式中：Q—体积排量，m3/s。

当前环空内流动假设为层流，从而有如下的摩阻系数和MR雷诺数的关系式：

(5)

对应的湍流临界雷诺数为：

Recr=4 150-1 150nL

(6)

因此当计算出该摩阻系数后，即可根据Fanning摩阻因子定义式(3)求出摩擦压降，从而最终求出井底ECD。

2 基于Fluent的ECD计算

为了证明MR法计算ECD的准确性，在缺乏实验条件的前提下，本文考虑采用直接数值模拟的方式来进行验证。

2.1 软件设置

本文采用ANSYS 17.0的Workbench平台进行有限元结构模型的搭建。首先，使用Geometry模块建构井筒环空中水泥浆的3D几何模型；然后，采用ICEM CFD模块来实现几何模型的网格划分和边界层加密；最后，调用Fluent模块(有限体积方法)完成直接数值模拟。在调用Fluent模块对所建立的有限元结构模型进行数值模拟之前，还需要进行网格尺寸与计算区域的无关性验证。其中，网格尺寸的调整参数包括横截面边界层的网格尺寸δ(mm)和环空轴向上的网格尺寸Δy(mm)。关于计算区域的调整参数为环空长度与其水力当量直径的比值，L/Dh。在进行无关性检验时，主要的验证指标是出口界面的最大轴向流速vz_max(m/s)。通过比较发现，当采用L/Dh=30，δ=0.01 mm，Δy=5 mm时，得到的vz_max能较好地满足计算精度的要求，同时模型总体的网格数量可控。

2.2 摩擦压降的计算

当模型轴向尺寸为环空水力当量直径的30倍时，即L=30Dh时，流体在出口处的流速就已经可以到达稳定状态；而当模型轴向尺寸进一步增加至60Dh甚至120Dh时，并不会对出口处的流速造成影响。这是因为在长度为30Dh的环空中，流体的流动已经可以实现充分发展，而模型轴向尺寸的进一步增加只是加长了环空中部流体稳定流动区域的范围，因而对出口处的流体流速没有影响。为了验证这一分析，对轴向尺寸为环空水力当量直径30倍、60倍和120倍的有限元模型进行模拟，得到的流体在轴向上的流速分布如图1所示，图中的4个子图展示了重点区域的流速变化。

由图1可知，流体在进入环空入口后的大约15Dh的长度范围内有一个流动发展过程，流速随着所处环空高度的增加，先是显著加快，然后逐渐趋于稳定。而在距离环空出口大约10Dh的长度范围内，流体流速呈现缓慢下降的趋势，并在出口附近有较为明显的流速放缓现象。除了在入口和出口的一定范围内流速会出现变化外，在剩余的环空中部区域，流速则保持稳定，并且环空模型轴向尺寸越大，流速稳定区域也相应地变大。基于以上认识，在分析长达数千米的油气井筒环空中的流体流动时，则可以通过采用一个轴向尺寸极大收缩的有限元模型，例如L=30Dh，来代替实际的井筒环空进行直接数值模拟，这样就能在保持计算精度的情况下避免了建立全尺寸的环空结构模型，从而能极大的提高数值模拟的效率。

图1 L=30Dh、60Dh、120Dh的环空中流体在x=0，z=0处的轴向流速的分布比较

对流体流速分布规律的分析是为了准确地计算摩擦压降，因为根据前述的MR理论，流速将决定MR雷诺数，从而改变摩阻系数，并最终影响ECD的计算结果。因此，为了进一步明确环空中流体流速变化的原因，对流体在轴向尺寸为L=120Dh的环空中不同轴向位置处的轴向、径向和切向流速在横截面单侧上的分布进行了分析，结果如图2所示。

图2 L=120Dh的环空中流体在x=0，y=(a)15Dh、(b)60Dh、(c)110Dh、(d)120Dh处各向流速在单侧横截面分布

由图2(a)可知，当流体从底部的环空入口流动至15Dh的环空高度时，轴向流速剖面稳定，径向和切向流速的变化也已经极为微弱，因此，流体的流速也基本上达到稳定值，这段环空区域就对应了图1中所提到的流速发展段。而如图2(b)、图2(c)所示，当流体进一步向上流动至60Dh、甚至110Dh的环空高度时，各个方向上的流速都保持稳定，这就对应了图1中所提到的流速稳定段。最后，如图2(d)所示，当流体到达环空出口处时，其径向和切向流速都有明显的波动，而这也是出口流速突然下降的原因，这对应了图1中所提到的流速减缓段。

综上所述，在计算大轴向尺度的环空摩擦压降时，应排除流速发展段和流速减缓段的影响，而以流速稳定段的平均摩擦压降为基准进行计算。例如，当环空结构模型的轴向尺寸取为L=30Dh时，只需要计算入口下游15Dh处的横截面平均压力pinlet+15Dh和出口上游10Dh处的横截面平均压力poutlet-10Dh，就可计算出在长度为(L-25Dh)的流速稳定段的平均摩擦压降：

(7)

然后，基于Δpmiddle，在确保环空径向尺寸不变和流体条件相同的情况下，就可以进一步推算出任意轴向尺寸环空的整体摩擦压降：

Δpf=Δpmiddle·L

(8)

最后，将式(8)代入式(2)就可以实现对井底ECD的计算。

2.3 ECD计算

对井深1 000 m的单一尺寸环空中窄间隙条件下的固井水泥浆进行ECD计算，为了分析流变模型的影响，选取的三种水泥浆体系(如表1所示的A、B和C)的流变特性分别采用宾汉、幂律和赫巴模型进行表示。然后通过Fluent软件对轴向尺寸L=30Dh的环空结构模型进行数值模拟，模拟时保证流动全部处于层流状态。相同水泥浆体系在不同流变模型下计算的ECD结果如图3所示。由图3可知，幂律模型与赫巴模型的ECD计算结果较为接近，而宾汉模型的结果与它们之间的差别较为明显，特别是当环空的窄间隙特性显著时，宾汉模型的ECD计算结果明显高于其他两种模型，但随着环空间隙尺寸的增大和水泥浆排量的增加，水泥浆的流变模型对井底ECD逐渐减小，三种流变模型的计算结果也逐步趋于一致。由此可知，在环空窄间隙条件下，水泥浆的流度模型对井底ECD的计算结果十分敏感，固井设计时应该给予特别关注。

表1 三种典型流体的密度和每分钟对应转速下的流变仪读数

图3 三种水泥浆体系(A，B，C)在不同流变模型下计算的井底ECD的比较

3 不同商业软件计算ECD对比

为了对不同软件计算的ECD结果进行比较，使用三款钻完井软件(Drill Bench、Landmark、PVI)和MR法分别对相同窄间隙算例进行计算和对比。由于室内试验测试较为复杂，业界对层流计算采用直接数值模拟的方式是认可的，鉴于此，以DNS法计算的ECD为标准值，将Drill Bench、Landmark、PVI以及MR法计算的ECD的对应误差按照宾汉模型、幂律模型和赫巴模型进行分类，并分别绘制于图4(a)、4(b)、4(c)中。基于图4所展示的误差分析可知，采用MR法迭代计算的ECD与通过DNS法模拟得到的ECD基本一致，其计算误差控制在3%以内，这首先证明了MR理论计算方法的准确性。而三款钻完井软件针对不同流变模型的适应性则存在明显差异，具体而言：

图4 采用不同方法计算的ECD与DNS法计算的ECD之间的误差比较

(1)对于宾汉模型，PVI的计算误差最小，基本与DNS法的计算精度一致。而Landmark的计算结果都明显小于DNS法的计算结果，误差水平处于-40%～30%之间。Drill Bench的计算误差波动较大，特别是对水泥浆体系B，其计算的ECD接近Landmark的计算结果，即误差水平接近-40%；但Drill Bench在分析水泥浆系统A和C时，计算误差又较小，与DNS法的计算精度相近。

(2)对于幂律模型，计算误差最小的软件变为Landmark，而PVI的计算结果普遍偏高20%左右。Drill Bench的计算误差还是呈现较大的波动性，它对水泥浆体系B的计算误差明显偏大，误差水平处于20%～40%之间；而对水泥浆体系C的计算结果又总是小于DNS法的计算结果。

(3)对于赫巴模型，由于Drill Bench和Landmark软件的使用条件限制，笔者只采用了PVI软件来进行ECD计算，其计算结果均略小于MR法，且与DNS法的计算误差保持在3%以内。

4 结论与建议

(1)采用Fluent 17.0商用流体计算软件对单一环空窄间隙条件下的层流流动进行模拟时，当环空轴向尺寸与其水力当量直径的比值为L/Dh=30、横截面边界层的网格尺寸为δ=0.01 mm、环空轴向上的网格尺寸为Δy=5 mm时，数值模拟的结果就能较好地满足计算精度的要求，同时流动计算区域内的网格数量适中，计算费用合理。

(2)通过不同计算软件的对比分析，DNS法的计算结果与MR法的计算结果基本一致，平均误差小于3%，因此MR法是一种计算窄间隙条件下固井ECD的有效方法。

(3)窄间隙环空的井底ECD对水泥浆的流变模型十分敏感，而不同流变模型下ECD的显著变化会大大增加压漏地层或引起井壁坍塌的风险，因此在窄间隙固井作业中应该特别关注水泥浆流变性对井底ECD的影响。