西南石油大学学报(自然科学版)  2016, Vol. 38 Issue (3): 135-143
钻井环空中赫-巴流体内气泡的上升速度    [PDF全文]
尹浚羽1,2 , 周英操1, 张辉2, 蒋宏伟1, 朱磊1,2    
1. 中国石油集团钻井工程技术研究院, 北京 昌平 102206;
2. 中国石油大学(北京)石油工程学院, 北京 昌平 102249
摘要: 气泡在钻井环空中的运动规律是井筒多相流的一个研究重点,气泡上升速度则是其中的一个关键参数。井筒内各流型内气泡/气体上升速度的准确性很大程度上影响了井筒多相流描述的精确性。为此,按照含气率对环空流型进行了划分,对各流型气泡的运动行为进行了分析。系统地总结分析了单气泡在无限流域中的滑脱速度及气泡纵横比的计算方法。利用赫-巴流变模式对气泡周围表观黏度进行了修正,通过实例对比分析单气泡在赫-巴流体内滑脱速度几个计算公式,认为Rodrigue公式相对于Harmathy公式来说是极大的提高,能适应不同形状的气泡,可作为泡状流流域气泡滑脱速度计算的首选公式。通过对单气泡滑脱速度进行含气率、尺寸和井斜修正,得出了环空中各流型的气体滑脱速度计算方法,过渡流型的气泡滑脱速度则采用线性处理,并对环空气泡滑脱速度随含气率变化关系进行了实例分析。最终通过漂移流模型给出了计算环空气泡上升速度的方法。
关键词: 钻井环空流型     赫-巴流体     气泡滑脱速度     漂移流模型     气泡上升速度    
Bubble Rise Velocity in Herchel-Bulkley Fluid of Drilling Annulus
YIN Junyu1,2 , ZHOU Yingcao1, ZHANG Hui2, JIANG Hongwei1, ZHU Lei1,2    
1. CNPC Drilling Research Institute, Changping, Beijing 102206, China;
2. School of Petroleum Engineering, China University of Petroleum(Beijing), Changping, Beijing 102249, China
Abstract: Bubble motion laws in drilling annulus is an important part of borehole multiphase flow research. Among which bubble rise velocity is a key parameter. The accuracy of borehole multiphase flow description highly depends on the accuracy of bubble/gas rise velocity of each flow pattern. In the view of this, according to the void fraction, the drilling annulus flow is divided into four patterns and then the analysis of bubble or gas motion characteristic of every pattern is carried out. The calculation methods of single bubble slip velocity and aspect ratio in Infinite media are summarized. Apparent viscosity round bubbles has been corrected through the Herchel-Bulkley reheology model. Through the contrastive analysis of case study among several formulas of single bubble's slip velocity in Herchel-Bulkley fluid, the Rodrigue's formula is considered as an improvement compared to the Harmathy's formula and can be the optimal formula of bubble rise velocity in bubbly flow pattern. The modification of calculation formulas of bubble slip velocities with regard to void fraction,bubble size and the angle of borehole has been adopted to correct the slip velocities in four flow patterns. For the bubble slip velocity in the transforming pattern,the liner method is proposed. The case study is applied in the relationship between bubble rise velocity and void fraction. Finally, bubble rise velocity can be acquired through the method of drift model.
Key words: flow pattern in driling annulus     Herchel-Bulkley fluid     bubble slip velocity     drift model     bubble rise velocity    
引 言

气泡上升运动规律一直是热点研究问题。初期的研究重点是单个气泡的上升过程。单气泡在上升过程中受到浮力、黏滞力、表面张力、惯性力、气泡尺寸、气泡变形率、不稳定扰动量等因素影响。前人通过大量的实验研究和理论分析发展了众多单气泡滑脱速度计算公式。目前钻井环空中泡状流流域的气泡滑脱速度大多采用Harmathy 无限流域中单气泡上升速度公式[1]。该公式没有考虑气泡尺寸等因素对上升速度的影响,仅对表面张力占主导地位的流域有较高的精度。若将该公式用于整个泡状流域的气泡滑脱速度计算会造成较大误差。因此,寻找考虑因素较完整、适用范围较广的单气泡上升速度公式是本文的一个重点。

除了以上因素以外,井筒中气泡上升过程会受到含气率的影响。气泡群在上升过程中可能相互碰撞、挤压或融合为一个气泡从而改变气泡尺寸形状,前行气泡形成的尾流会使得后续气泡的速度增大。对于大尺寸气泡(Taylor 气泡),井筒尺寸也是一个重要的影响因素。这两个方面的影响,本文会给出相应的影响因子对单气泡上升速度公式进行修正。值得一提的是大部分钻井液是非牛顿流体,流变模式考虑为赫-巴流变模式。该模式的表观黏度是随着剪切速率变化的。而大部分气泡上升速度计算方法认为气泡周围黏度是常数,故需要进行黏度修正。井筒气液混合流动过程中,针对同一截面液流速度剖面的非线性和气泡分布不均等造成的气泡上升速度变化,本文会采用漂移流模型进行讨论。

综上所述,气泡上升速度是一个多因素影响的函数。本文给出气泡上升速度计算的一般思路是:首先计算出单个气泡在无限流域中的滑脱上升速度,再通过与含气率、流道尺寸和倾角相关的修正因子改进,最后通过漂移流模型得出气泡上升速度。

1 环空气泡形状和流型

图 1 展示了钻井环空可能出现的气泡形状和流型[2-3]。同一口气井可能自下而上在不同井段分别经历这4 个过程。下面将根据Julia J E 等[2-3] 及Waltrich P J 等[4-7] 对环空流型的研究介绍各流型特性。

图1 气体形状和流型 Fig. 1 Flow regimes definition
1.1 泡状流

在体积含气率${\alpha}$ 较小的情况下(一般${\alpha}$<15%或者${\alpha}$<20%),液体为连续相,气体以气泡形式分散在液相中。在泡状流型下,气泡有3 个基本形状:球形,椭球和球冠状。Clift R[5] 利用气泡的Reynolds(Re)数、Eötvös(Eo)数和Morton(M)数绘制了气泡的形状版图以方便使用。

$\left\{ \begin{array}{l} Eo = {\rm{g}}\Delta \rho d_{\rm{e}}^2/\sigma \\ Re = {\rho _{\rm{l}}}{d_{\rm{e}}}{U_\infty }/{\mu _a}\\ M = {\rm g}{\mu _{\rm{a}}}^4\Delta \rho /(\rho _{\rm{l}}^2{\sigma ^3}) \end{array} \right.$ (1)

式中:Eo Eötvös 数,无因次;

M—Morton数,无因次;

g—重力加速度,g = 9.8 m/s2

$\Delta \rho$—气液密度差,kg/m3

$d_{\rm{e}}$—气泡等效直径(当量直径),m;

$\sigma$—气液表面张力,N/m;

${\rho _{\rm{l}}}$—液体密度,kg/m3

${U_\infty }$—气泡在无限流域中稳定滑脱速度,m/s;

${\mu _{\rm{a}}}$—气泡界面附近的黏度,Pa·s。

当气泡$Re< 1$时,气泡近似为球形;随着$Re$和$Eo$($1<Eo< 40$)增大,气泡横纵比开始增大,气泡形状为椭球形和类椭球形;当气泡$Eo>40$时,气泡变形继续增强,气泡形状为球冠形或类球冠形。

1.2 冠状-断塞流

气泡在环空上升过程中,液柱压力逐渐减小,气泡体积膨胀增大,使得气泡含气率增大(20%<${\alpha}$<52%)。气泡之间的相互作用也逐渐加强,小气泡开始融合成大气泡。由于管壁的限制,气泡横向膨胀受到限制,纵向则变长。当体积足够大时,气泡便成为了断塞状(图 1b)。值得一提的是,环空断塞气泡形状又与圆形截面流道不同:圆形截面流道形成的断塞气泡顶部为球冠状,底部类似为圆柱;而环空中的断塞气泡横向左侧膨胀受到内管壁阻碍,右侧受到外管壁限制,气泡只能沿着环空流道周向膨胀变长,使得气泡截面变为月牙状或者圆弧环状(图 1b)。

1.3 搅拌流

气泡继续上升膨胀,气泡含气率进一步加大(52%<${\alpha}$<90%),气泡之间的相互作用也更加剧烈,小气泡时而融合形成断塞,断塞时而被湍流震荡冲散为小气泡(图 1c)。在搅拌流域,气液交互也十分混乱,液体时而形成段桥阻碍气泡的运动,又被气泡抬升驱散。

1.4 环雾流

气泡含气率达到一定程度(${\alpha}$>90%),气体流速会非常大,气体和液体运动变为分层流动状态(气泡处于环空中心,液体则在内外管壁上形成液膜),部分液体则以液珠的形式分散在气相中(图 1d)。

2 单气泡在赫-巴流体无限流域中的滑脱速度 2.1 气泡表观黏度

钻井液的常用体系绝大部分都是非牛顿流体,而应用范围较广的钻井液流变模式是赫-巴流变模式

$\tau = {\tau _0} + k{\gamma ^n}$ (2)

式中:${\tau _0}$—屈服应力,Pa;

k—黏度系数,Pa$\cdots^n$;

$\gamma$—剪切速率,s$^{-1}$;

n—流性指数,无因次。

参照Mendelson H D 等对气泡周围幂律流体表观黏度的定义[8],赫-巴流体中气泡表面的平均剪切速率等于${\mu _{\rm{a}}}/a$(a—气泡的横向最大尺寸,m),气泡附近流体表观黏度可写为

${\mu _{\rm{a}}} = {\tau _0}{\left (\dfrac{{{U_\infty }}}{a} \right)^{ - 1}} + k{\left( {\dfrac{{{U_\infty }}}{a}} \right)^{n - 1}}$ (3)
2.2 球形气泡滑脱速度和阻力系数

当气泡尺寸较小和气泡雷诺数Re < 1 时,表面张力对气泡形状影响占主导地位,气泡形状成球形,而黏性力是气泡最主要的阻力。Clift R 等[5] 总结了两个适用于球形气泡的两种计算方法。

第一种是Stokes 方法,该方法认为小尺寸气泡的气液界面表现出固壁特性。基于这个前提,球形气泡的滑脱速度与球形固体的滑落速度类似

${U_\infty } = {\rm g}d_{\rm{b}}^2\Delta \rho /\left( {18{\mu _a}} \right)$ (4)

式中:$d_{\rm{b}}$ 气泡直径,m。

第二种方法是Hadamard-Rybczynski 方法,该方法考虑了气泡上升过程中,黏性切力引起的气泡内部气体循环对气泡滑脱速度的影响。计算公式为

${U_\infty } = {\rm g}d_{\rm{b}}^2\Delta \rho /\left( {12{\mu _{\rm{a}}}} \right)$ (5)

对于稳定上升的气泡来说,气泡处于浮力和上升阻力受力平衡状态,气泡阻力系数$C_{\rm{D}}$表示为

${C_{\rm{D}}} = \dfrac{{4{\rm{g}}d_{\rm{e}}^3\Delta \rho }}{{3{\rho _{\rm{l}}}{a^2}U_\infty ^2}}$ (6)

对于球形气泡来说,$d_{\rm{b}}=d_{\rm{e}}=a$,将式(6)代入式(4),得到Stokes 阻力系数

${C_{\rm{D}}} = \dfrac{{Re }}{{24}}$ (7)

将式(6)代入式(5),得到Hadamard-Rybczynski阻力系数

${C_{\rm{D}}} = \dfrac{{Re }}{{16}}$ (8)
2.3 椭球形气泡滑脱速度、纵横比和阻力系数

(1)椭球形气泡滑脱速度

椭球形气泡一般处于$Re >0.1$和$Eo< 40$的区域。该区域气泡的上升阻力中黏性力比例减弱,表面张力占主导地位,惯性力比例开始增加。忽略黏性阻力的作用,Medelsonb H D[8] 创新性地提出了预测气泡上升速度的波动理论,他认为气泡界面扰动的运动特性和理想流体中波的运动特性十分相似,从而导出了气泡的上升速度公式

${U_\infty } = \sqrt {\dfrac{{2\sigma }}{{{\rho _{\rm{l}}}{d_{\rm{e}}}}} + \dfrac{{{\rm{g}}{d_{\rm{e}}}}}{2}}$ (9)

文献[9][10]认为该公式特别适合椭球形气泡滑脱速度计算。Lahrer I H 认为上升过程中气泡势能会转换为动能,随后又在波动中耗散掉一部分,这刚好跟Medelson 波动理论契合,并对Medelson气泡滑脱速度公式进行了改进[11]

${U_\infty } = \sqrt {\dfrac{{3\sigma }}{{{\rho _{\rm{l}}}{d_{\rm{e}}}}} + \dfrac{{\Delta \rho {\rm{g}}{d_{\rm{e}}}}}{{2\rho }}}$ (10)

式(10)根号下第一项表征表面张力对气泡滑脱速度的作用,第二项则表征惯性力对气泡滑脱速度的作用。式(10)相对于式(9)在第一项有所强化。

除此之外,Clift R 等[5] 在前人研究成果基础上,给出了适用于$M< 10^{-3}$,$Eo< 40$,$Re>0.1$条件的椭球形气泡滑脱速度计算公式

${U_\infty } = \dfrac{{{\mu _a}}}{{{\rho _{\rm{l}}}{d_{\rm{e}}}}}{M^{ - 0.149}}(J - 0.857)$ (11)

其中

$\left\{ \begin{array}{l} J = 0.94{H^{0.757}}{\kern 10pt}{\left(2 < H \leqslant 59.3 \right )}\\ J = 3.42{H^{0.441}}{\kern 10pt} \left ( H > 59.3 \right )\\ H = \dfrac{4}{3}Eo{M^{ - 0.149({\mu _{\rm{a}}}/{\mu _{\rm{w}}})}}\\ {\rm{ }}{\mu _{\rm{w}}} = 0.0009 \end{array} \right.$ (12)

(2)纵横比和阻力系数

变形气泡的等效直径定义为

${d_{\rm{e}}} = \sqrt[3]{{{a^2}b}}$ (13)

式中:

a—气泡的横向最大尺寸,m;

b—气泡的纵向最大尺寸,m。

$E = \dfrac{b}{a}$ (14)

式中:E—气泡纵横比,无因次。

则有

$\left\{ \begin{array}{l} a = {E^{ - 1/3}}{d_{\rm{e}}}\\ b = {E^{2/3}}{d_{\rm{e}}} \end{array} \right.$ (15)

式(15)代入式(7)得到变形气泡的阻力系数

${C_{\rm{D}}} = \dfrac{{4{E^{2/3}}{\rm{g}}{d_{\rm{e}}}\Delta \rho }}{{3{\rho _{\rm{l}}}U_\infty ^2}}$ (16)

Clift R 等[5] 总结分析了大量的前人实验数据并得出了气泡纵横比版图,该版图表明,气泡的纵横比不但随着Eo 增大而逐渐增大,还跟Morton 数密切相关。大Morton 数下,黏度大,纵横比小,气泡更不易变形。Bozzano G[12] 提出了类似的气泡变形率(a=de)2 计算方式(式(17)),从而得到Bozzano 和Dente 气泡纵横比版图(图 2)。

$E = {\left[ {\dfrac{{10(1 + 1.3M^{1/6}) + Eo}}{{10(1 + 1.3M^{1/6}) + 3.1Eo}}} \right]^{3/2}}$ (17)
图2 Bozzano 和Dente 气泡纵横比版图 Fig. 2 Bozzano and Dente bubble aspect ratio E Vs. Eo

Tadaki和Maeda则认为气泡变形率是与$Re$和$Eo$有关的$T_{\rm{a}}=ReMo^{0.23}$的函数,Vakrushev和Efremov扩展了这一方法并得到了适用于$M< 10^{-3}$的气泡横纵比的计算公式

$E = \left\{ \begin{array}{l} 1{\kern 10pt}{\rm{ }}({T_{\rm{a}}} < 1)\\ {\left[ {0.81 + 0.206\tanh \left( {1.6 - 2\lg {T_{\rm{a}}}} \right)} \right]^3}\\{\kern 21pt}(1 \leqslant {T_{\rm{a}}} < 39.8)\\ 0.24 {\kern 10pt} {\rm{ }}({T_{\rm{a}}} \geqslant 39.8) \end{array} \right.$ (18)

Vakrushev 和Efremov 气泡纵横比版图见图 3

图3 Vakrushev 和Efremov 气泡纵横比版图 Fig. 3 Vakrushevand Efremov bubble aspect ratio E Vs. Re
2.4 球冠形气泡滑脱速度

球冠形气泡一般处于$Eo >40$的区域。在该区域,气泡所受的上升阻力中,惯性力占主导地位,黏性力十分微小,表面张力也较小。忽虑黏性力和表面张力的作用,球冠形气泡滑脱速度可以写为[5]

${U_\infty } \!=\! 0.707\sqrt {{\rm{g}}{d_{\rm{e}}}\Delta \rho /{\rho _{\rm{l}}}} {\kern 9pt}\left( {Re\! > \! 150{,}Eo \!\geqslant \!40} \right)$ (19)
2.5 单气泡滑脱速度的几个通用公式

(1)Karamanev 公式

对于小气泡,气泡的滑脱速度遵从Stokes 公式。但是随着气泡形状的变化,气泡的阻力系数也随之改变。Karamanev 阻力系数公式[13] 是利用Reynolds 数对Stokes 公式进行修正的。

纯净液体中的气泡在上升过程中气泡内部还可能出现内部循环,气泡的阻力系数可利用Hadamard-Rybczynski 公式进行修正。而井筒内的钻井液包含了很多表面活性剂和其他杂质,这大大地阻碍了气泡的内部循环,钻井液中的气泡表面更加倾向为“固壁”。阻力系数的修正也应针对于Stokes 公式

${C_{\rm{D}}} = \left\{ \begin{array}{l} \dfrac{{24}}{{R_{\rm{et}}}}\left (1 + 0.173{R_{\rm{et}}}^{0.657} \right) + \\{\kern 21pt}\dfrac{{0.413}}{{1 + 163000{R_{\rm{et}}}^{ - 1.09}}}{\kern 10pt}(Re \leqslant 135)\\ 0.95{\kern 10pt}(Re > 135) \end{array} \right.$ (20)

其中

${R_{\rm{e}}}_{\rm{t}} = \dfrac{{{\rho _{\rm{l}}}a{U_\infty }}}{{{\mu _{\rm{a}}}}} = \dfrac{{{\rho _{\rm{l}}}{E^{ - 1/3}}{d_{\rm{e}}}{U_\infty }}}{{{\mu _{\rm{a}}}}}$

将式(20)带入式(16)可得到气泡的滑脱速度

${U_\infty } = \sqrt {\dfrac{{4{E^{2/3}}\Delta \rho {\rm{g}}{d_{\rm{e}}}}}{{3{\rho _{\rm{l}}}{C_{\rm{D}}}}}}$ (21)

(2)Jamaialahamdi 公式

Jamaialahamdi 等综合考虑Stokes 公式和Medelson 公式,得出了无限流域中气泡滑脱速度的计算公式[14]

$\left\{ \begin{array}{l} {U_\infty } = \sqrt {\dfrac{{4{E^{2/3}}\Delta \rho {\rm{g}}{d_{\rm{e}}}}}{{3{\rho _{\rm{l}}}{C_{\rm{D}}}}}} \\[6pt] {U_{{\rm{b1}}}} = {\rm{g}}d_{\rm{b}}^2\Delta \rho /\left (18{\mu _{\rm{a}}}\right )\\[6pt] {U_{{\rm{b2}}}} = \sqrt {\dfrac{{2\sigma }}{{{\rho _{\rm{l}}}{d_{\rm{e}}}}} + \dfrac{{\Delta \rho g{d_{\rm{e}}}}}{{2{\rho _{\rm{l}}}}}} \end{array} \right.$ (22)

气泡尺寸较小时,Stokes 公式占主导地位;当气泡处于大尺寸时,Medelson 公式占主导地位。

(3)Bozzano 和Dente 公式

Bozzano G 和Dente M 认为气泡阻力系数是摩擦因子和气泡变形率的函数[12]

${C_{\rm{D}}} = f{\left( {\dfrac{a}{{{d_{\rm{e}}}}}} \right)^2}$ (23)

式中:f—摩擦因子,无因次。

摩擦因子f是$Eo$和M的函数

$f = \dfrac{{48}}{{{\mathop{ Re}\nolimits} }}\left( {\dfrac{{1 + 12{M^{1/3}}}}{{1 + 36{M^{1/3}}}}} \right) + \\{\kern 42pt} 0.9\dfrac{{E{o^{3/2}}}}{{1.4\left( {1 + 30{M^{1/6}}} \right) + E{o^{3/2}}}}$ (24)

气泡的变形率用$(a/ d_{\rm{e}})^2$来标定

${\left( {\dfrac{a}{{{d_{\rm e}}}}} \right)^2} = \dfrac{{10\left( {1 + 1.3{M^{1/6}}} \right) + 3.1Eo}}{{10\left( {1 + 1.3{M^{1/6}}} \right) + Eo}}$ (25)

利用受力平衡可得出气泡的稳定滑脱速度

${U_\infty } = \sqrt {\dfrac{{4\Delta \rho {\rm{g}}{d_{\rm{e}}}}}{{3{\rho _{\rm{l}}}{C_{\rm{D}}}}}}$ (26)

(4)Rodrigue 公式

Rodrigue[15] 通过无量纲分析方法提出了一个很好的单气泡上升速度计算公式

$\left\{ \begin{array}{l} U \!=\! {\dfrac{{F{{\left( {1 \!+\! 1.31 \!\times \!{{10}^{ - 5}}{M^{11/20}}{F^{73/33}}} \right)}^{21/176}}}}{{12{{\left( {1 + 0.02{F^{10/11}}} \right)}^{10/11}}}}} \\[9pt] U \!=\! {U_\infty }{\left( {\dfrac{{\rho _{\rm{l}}^2d_{\rm{e}}^2}}{{\sigma {\mu _{\rm{a}}}}}} \right)^{1/3}}\\[9pt] F \!=\! {\rm{g}} {\dfrac{{\Delta {\rho ^5}d_{\rm{e}}^8}}{{\sigma \mu _{\rm{a}}^4}}} \end{array} \right.$ (27)

从式(27)可见,该方法考虑到了气泡尺寸、表面张力、黏性和密度等因素,能够在适应非常宽的气泡变形区域和不同的流体介质。

3 不同含气率区域内气泡/气体滑脱速度和井斜校正 3.1 泡状流流域气泡滑脱速度

气泡之间的相互作用会影响气泡的滑脱速度,在泡状流流域(${\alpha}$<20%),含气率对气泡滑脱速度的影响可用含气率修正因子来表征。根据Richardson和Zaki 的启发,气泡的修正因子可写为[14]

${F_\alpha } = {(1 - \alpha )^{n'}}$ (28)

式中:F${\alpha}$ 含气率修正因子,无因次;

$n'$ 系数,无因次。

$n'$ 在1.5~2.0 取值[16],本文取为1.75。从而泡状流流域的气泡滑脱速度经含气率修正因子修正后可表示为

${U_{\rm{S}}} = {U_\infty }{(1 - \alpha )^{1.75}}{\kern 6pt}(\alpha < 20\% )$ (29)

式中:${U_{\rm{S}}}$—滑脱速度,m/s。

3.2 环空Taylor 气泡滑脱速度

气泡尺寸和流道尺寸的比例关系也会对气泡的滑脱速度产生影响,气泡尺寸越大,管壁对气泡的上升阻力越大。类似于含气率,流道尺寸对气泡滑脱速度的影响可以流道尺寸修正系数来表征。Krishna R 给出了Collins 圆形流道尺寸修正系数[17]

$\left\{ \begin{array}{l} {F_{\rm{p}}} = 1 {\kern 66pt}\left(\dfrac{{{d_{\rm{e}}}}}{{{D_{\rm{T}}}}} < 0.125 \right)\\ {F_{\rm{p}}} = 1.13{{\rm{e}}^{ - {d_{\rm{e}}}/{D_{\rm{T}}}}} {\kern 26pt} \left( 0.125 < \dfrac{{{d_{\rm{e}}}}}{{{D_{\rm{T}}}}} < 0.6 \right)\\ {F_{\rm{p}}} = 0.496\sqrt {{D_{\rm{T}}}/{d_{\rm{e}}}} {\kern 16pt} \left(\dfrac{{{d_{\rm{e}}}}}{{{D_{\rm{T}}}}} > 0.6\right) \end{array} \right.$ (30)

式中:${F_{\rm{p}}}$—流道尺寸修正系数,无因次;

${D_{\rm{T}}}$—圆管内径,m。

当$d_{\rm{e}}$趋近于${D_{\rm{T}}}$时,气泡变为Taylor气泡,此时${F_{\rm{p}}} \thickapprox 0.496$。Taylor气泡的前端为球冠状,该结构是影响气泡滑脱速度最主要的因素。因此,圆管内其滑脱速度可以利用流道尺寸修正系数0.496乘以球冠状气泡进行修正得到

${U_{\rm{S}}} = 0.35\sqrt {{\rm g}{D_{\rm{T}}}\Delta \rho /{\rho _{\rm{l}}}}$ (31)

Hasan A R 等[18] 发现环空中的Taylor 气泡滑脱速度也和环空几何有关,并扩展得出了环空中Taylor 气泡的滑脱速度计算公式(40% < ${\alpha}$ < 52%)

${U_{\rm{S}}} = \left(0.345 + 0.1\dfrac{{{D_{\rm{i}}}}}{{{D_{\rm{o}}}}}\right)\sqrt {\dfrac{{\Delta \rho {\rm g}{D_{\rm{o}}}}}{{{\rho _{\rm{l}}}}}}$ (32)

式中:${D_{\rm{i}}}$—环空的内径,m;

${D_{\rm{o}}}$—环空的外径,m。

断塞流流域的气泡不仅仅是Taylor 气泡,还有一些体积较小的气泡。因此式(32)仅可作为以Taylor气泡为主的断塞流流域($40\% < \alpha < 52\%$)气泡滑脱速度。

3.3 搅拌流流域的气泡

搅拌流流域的气泡运动则较为不稳定,湍流波动使得气泡在大量破碎的同时又大量融合。对于该流域气泡的滑脱速度,Zuber 认为它和Harmathy 气泡上升速度公式相同[1, 16]

${U_{\rm{S}}} = 1.53{\left( {\dfrac{{\Delta \rho \sigma {\rm{g}}}}{{\rho _{\rm{l}}^2}}} \right)^{1/4}} {\kern 10pt}(60\% < \alpha < 90\%)$ (33)
3.4 环雾流流域气体滑脱速度

当环空流域为环雾流流域时,气体为连续相,当含气率趋进于100% 时,滑脱速度接近于0。

$\mathop {\lim ({U_{\rm{S}}})}\limits_{\alpha \to 100\% } = 0$ (34)
3.5 过渡含气率区域的气泡滑脱速度

过渡含气率区域为式(31)~ 式(34)没有覆盖所有的含气率区域。针对这一问题,参考Shi H 等[19]的处理方法,对处在过渡区域内的气泡滑脱速度看作是含气率的线性函数,函数的前后两端点分别服从前一含气率区域的气泡滑脱速度公式和后一含气率区域的滑脱速度公式。

3.6 气泡/气体滑脱速度的井斜校正

钻井技术的飞速发展,使得井眼形状多样,钻井环空不再局限于垂直状态,定向井、水平井等迅速发展。这些井型的井斜角必然会对气泡的滑脱速度造成影响。基于实验数据的观察,Hasan A R等[18] 得到了气泡/气体滑脱速度的井斜修正因子,修正后气泡/气体的滑脱速度可表示为

${U_{\rm{S}}}(\theta ) = {U_{\rm{S}}}{(\cos \theta )^{0.5}}{(1 + \sin \theta )^2}$ (35)

式中:θ 井斜角,(◦)。

4 漂移流模型

两相流漂移流模型认为气体和液体之间的滑脱是两个机理联合作用的结果。第一个机理是流体速度剖面非均匀性和流道剖面上气体的非均匀分布使得流道中心处的气体/液体流量最高,混合物的局部速度也最大,对流道剖面积分后得出气体的平均速度会比液体大。第二个机理是由于浮力等作用使得气体超越液体流动。

Zuber N 和Findlay J A[20] 给出了描述气泡/气体上升速度的漂移流模型

${U_{\rm{g}}} = {C_{\rm{o}}}{U_{\rm{m}}} + {U_{\rm{S}}}(\theta )$ (36)

式中:${U_{\rm{g}}}$—气泡/气体上升速度,m/s;

${C_{\rm{o}}}$—速度分布因子,无因次;

${U_{\rm{m}}}$—气液混合物流速,m/s。

气液混合物流速计算公式为

${U_{\rm{m}}} = \dfrac{{{Q_{\rm{l}}} + {Q_{\rm{g}}}}}{A}$ (37)

式中:${Q_{\rm{l}}}$—液体流量,m3/s;

${Q_{\rm{g}}}$—气体流量,m3/s;

A—横截面积,m2

Ishii M 给出了圆管的速度分布系数[21]

${C_{\rm{o}}} = 1.2 - 0.2\sqrt {{\rho _{\rm{g}}}/{\rho _{\rm{l}}}}$ (38)

式中:ρg 气体密度,kg/m3

环空中的速度分布系数则要考虑环空尺寸。对式(38)进行扩展,高永海[6] 给出了环空中各流域气泡的速度分布系数

泡状流域

${C_{\rm{o}}} = 1.2 + 0.371\left( {{D_{\rm{i}}}/{D_{\rm{o}}}} \right)$ (39)

段塞流域

${C_{\rm{o}}} = 1.182 + 0.9\left( {{D_{\rm{i}}}/{D_{\rm{o}}}} \right)$ (40)
5 实例分析

根据实际情况,计算钻井过程中的滑脱速度。计算使用的参数见表 1

表1 基本参数列表 Table 1 Basic parameter list
5.1 泡状流流域中赫-巴流体内单气泡滑脱速度的几个通用公式比较

钻井井筒泡状流流域的气泡滑脱速度多采用Harmathy 公式[1, 22-23],在基本参数不变的情况下,该公式计算出的滑脱速度为一常数,如图 4 虚线所示。与其他气泡计算公式相比较,Harmathy 气泡计算公式表征的仅仅是中等尺寸气泡的滑脱速度。因此,利用Harmathy 公式计算整个泡状流流域气泡的滑脱速度必然会产生很大的误差。

图4 泡状流流域中单气泡滑脱速度 Fig. 4 Single bubble slip velocity of bubbly flow pattern

图 4 展示了4个通用公式计算的滑脱速度随气泡当量直径的变化关系。在气泡当量直径很小($d_{\rm{e}}$< 1~2 mm)的情况下,气泡的滑脱速度接近于0。主要原因在于赫-巴钻井液特殊的流变模式,赫-巴流体的流变曲线具有较大的屈服值,使得小尺寸气泡周围的表观黏度很高($\mu_{\rm{a}}>100$~Pa·s),此时气泡所受到的黏滞阻力会很大,气泡很难向上运动。随着气泡膨胀或气泡之间的融合,气泡当量直径会逐渐增大,赫-巴流体的剪切稀释性显著增强,气泡周围的表观黏度迅速减小($\mu_{\rm{a}}< 1$~Pa·s),气泡在浮力的作用下克服黏滞力向上运动,速度增大幅度也较强。对于小气泡(0~7 mm)来说,气泡形状近似为圆球体,黏滞力是最主要的滑脱阻力,气泡滑脱速度遵从Stokes定律(图 4 中Stokes 曲线)。在这一区域与Stokes公式符合最好的通用公式是Rodrigue公式,其次是Bozzano & Dente公式和Jamialahmadi公式,最差的是Karamaev公式。随着气泡的当量直径继续增大,气泡不再为圆球形,气泡形变增大,形状近似为椭球体,此时,气泡不但受到黏滞力的作用,还受到表面张力和惯性力的影响,气泡滑脱速度的增大幅度减弱。在这区域,Medelson & Lahrer方法计算出的椭球形气泡滑脱速度较之4个通用公式都偏大,原因在于Medelson & Lahrer方法没有考虑黏滞阻力对气泡滑脱速度的影响。这一区域,根据实例计算出的气泡周围的表观黏度量级大于等于10$^{-1}$,M的量级则大于等于10$^{-2}$,黏滞力对气泡滑脱运动的阻力不能忽略, Medelson & Lahrer方法计算值偏高。而这种情况不满足Clift关于椭球形气泡滑脱速度计算公式的应用条件($M< 10^{-3}$),比较符合的是Rodrigue公式和Bozzano & Dente公式。当气泡体积增大到一定程度,气泡变形率趋于恒定,气泡形状类似于球冠状,黏滞力和表面张力对气泡滑脱运动影响变小,惯性力成为主要的影响因素。气泡滑脱速度接近球冠形气泡的滑脱速度公式计算值(如图 4 中球冠形气泡滑脱速度公式)。在这一区域符合最好的通用公式是Jamialahmadi公式,其次是Rodrigue公式和Bozzano & Dente公式。而在各区域表现较差的是Karamaev公式。

综合考虑,本文建议泡状流流域赫-巴流体内单气泡滑脱速度计算使用Rodrigue 公式。该公式在各区域的表现都较好,特别是在小气泡区域内计算值和Stokes 公式计算值重合度极高。

5.2 钻井环空中赫-巴流体内气泡/气体滑脱速度和上升速度

泡状流流域的气泡运动较为稳定,气泡之间的融合率也很低。气泡尺寸的增大主要是由于气泡的上升过程中周围压力降低引起的膨胀作用。因此,在泡状流流域可以认为气泡含气率和气泡平均体积成正比关系,即气泡的等效直径与含气率的1/3 次幂成正线性关系。假设气泡在形成后经过一个不稳定过程变为5 mm 左右的气泡,此时的含气率为1%。而当含气率1%~20% 时,气泡滑脱速度可使用Rodrigue 公式,再经含气率校正得到。再分别计算出其他不同含气率区域环空气泡/气体滑脱速度最终得到了气泡含气率与环空气泡滑脱速度的关系曲线(图 5)。

图 5 中可以看出整个泡状流流域的气泡滑脱速度均小于0.17 m/s,如果采用Harmathy 公式则计算值在0.25 m/s 左右。当含气率大于20%时,Taylor气泡开始形成,直到含气率达到40% 左右时,Taylor气泡成为环空中气泡的主要组成形式。此时气泡滑脱速度计算采用环空Taylor 气泡计算公式较为合理。在含气率大于20% 和40% 的区间,滑脱速度则认为随含气率线性变化。而一些文献视含气率一旦达到泡状流向段塞流转换边界(${\alpha}$=20%),滑脱速度就可采用环空Taylor 气泡公式计算的做法欠考究。当含气率达到52% 时,湍流扰动突然加强,Taylor 气泡的破碎概率增大,同时伴随着小气泡融合成Taylor 气泡过程。当含气率达到60% 左右时,这一种扰动现象成为主流,环空流域从段塞流真正过渡到搅拌流。对于过渡区域的处理方法与泡状流向段塞流转换的处理类似。当含气率达到90% 时,环空中气体为连续相,部分液体则靠近环空外壁形成液膜,部分液体成为液滴分散在气体中。当含气率达到趋近于1 时,环空中为气体单相流动,气体滑脱速度为0。这一区域内,气泡的滑脱速度是随着含气率的升高逐步减小的(图 5 环雾流区域)。

图5 气体含气率与环空气泡滑脱速度的关系 Fig. 5 Void fraction vs. bubble slip velocity in annulus
6 结语

(1)根据含气率把环空流型进行了划分:含气率在0~20% 为泡状流,20%~52% 为断塞流,52%~90% 为搅拌流,90%~100% 为环雾流。并对各流型的气泡/气体运动行为进行了分析。

(2)对单气泡在赫-巴流体无限流域中的滑脱速度进行了广泛深入的研究。根据赫-巴流变模式对气泡周围的表观黏度进行了修正,分别给出了适合球形气泡、椭球形气泡和球冠形气泡滑脱速度的计算公式以及几个通用公式。通过实例对比分析发现泡状流流域赫-巴流体内单气泡滑脱速度使用Rodrigue 公式计算较好。该公式在各区域的表现都较优,特别是在小气泡区域内计算值和Stokes 公式计算值重合度极高。该公式相对于Harmathy 公式来说有了长足的进步,能作为精确描述井筒泡状流的一个有力工具。

(3)给出了不同含气率区域内气泡/气体滑脱速度的计算公式,对于过度区域的气泡/气体滑脱速度采用线性处理。通过环空气泡滑脱速度随含气率变化关系的实例分析发现这种处理方式能改善一些文献视含气率一旦达到泡状流向段塞流转换边界就采用Taylor 气泡计算公式计算滑脱速度等处理方法的缺陷。在得到环空中气泡滑脱速度以后,如果存在井斜,首先利用(35)式进行井斜校正,最后通漂移流模型便可求出气泡的上升速度。

(4)钻井环空中流型的转换规则还需进一步研究。

参考文献
[1] HARMATHY T Z. Velocity of large drops and bubbles in media of infinite or restricted extent[J]. AIChE Journal, 1960, 6 (2) : 281 –288. DOI:10.1002/(ISSN)1547-5905
[2] JULIA J E, OZAR B, JEONG J J, et al. Flow regime development analysis in adiabatic upward two-phase flow in a vertical annulus[J]. International Journal of Heat and Fluid Flow, 2011, 32 (1) : 164 –175. DOI:10.1016/j.ijheatfluidflow.2010.09.003
[3] JULIA J E, HIBIKI T. Flow regime transition criteria for two-phase flow in a vertical annulus[J]. International Journal of Heat and Fluid Flow, 2011, 32 (5) : 993 –1004. DOI:10.1016/j.ijheatfluidflow.2011.06.001
[4] WALTRICH P J, FALCONE G, BARBOSA J R. Axial development of annular, churn and slug flows in a long vertical tube[J]. International Journal of Multiphase Flow, 2013, 57 : 38 –48. DOI:10.1016/j.ijmultiphaseflow.2013.06.008
[5] CLIFT R, GRACE J R, WEBER M E. Bubbles, drops, and particles[M]. New York: ACademic Press,Inc, 1978 .
[6] 高永海. 深水钻探井筒多相流流动与井控研究[D]. 青岛:中国石油大学(华东), 2007.
[7] DRANCHUK P M, PURVIS R A, ROBINSON D B. Computer calculation of natural gas compressibility factors using the Standing and Katz correlation[C]//Annual Technical Meeting, Petroleum Society of Canada, 1973.
[8] MENDELSON H D. The prediction of bubble terminal velocities from wave theory[J]. AIChE Journal, 1967, 13 (2) : 250 –253. DOI:10.1002/(ISSN)1547-5905
[9] Hassan N M S, Khan M M K, RASUL M G. Bubble rise velocity and drag co-efficient at high reynolds number in power-law fluids[C]. International Conference on Fluid Mechanics and Aerodynamics, Athens, Greece, 2007.
[10] NEDELTCHEV S, SCHUMPE A. New approaches for theoretical estimation of mass transfer parameters in both gas-liquid and slurry bubble columns[M]. InTech Open Access Publisher, 2011 .
[11] LEHRER I H. A rational terminal velocity equation for bubbles and drops at intermediate and high Reynolds numbers[J]. Journal of Chemical Engineering of Japan, 1976, 9 (3) : 237 –240. DOI:10.1252/jcej.9.237
[12] BOZZANO G, DENTE M. Shape and terminal velocity of single bubble motion:A novel approach[J]. Computers & Chemical Engineering, 2001, 25 (4) : 571 –576.
[13] KARAMANEV D G, NIKOLOV L N. Free rising spheres do not obey Newton's law for free settling[J]. AIChE Journal, 1992, 38 (11) : 1843 –1846. DOI:10.1002/(ISSN)1547-5905
[14] SIMONNET M, GENTRIC C, OLMOS E, et al. Experimental determination of the drag coefficient in a swarm of bubbles[J]. Chemical Engineering Science, 2007, 62 (3) : 858 –866. DOI:10.1016/j.ces.2006.10.012
[15] RODRIGUE D. A general correlation for the rise velocity of single gas bubbles[J]. The Canadian Journal of Chemical Engineering, 2004, 82 (2) : 382 –386.
[16] WALLIS G B. One-dimensional two-phase flow[M]. America: McGraw-Hill, 1969 .
[17] KRISHNA R, URSEANU M I, VAN BATEN J M, et al. Rise velocity of a swarm of large gas bubbles in liquids[J]. Chemical Engineering Science, 1999, 54 (2) : 171 –183. DOI:10.1016/S0009-2509(98)00245-0
[18] HASAN A R, KABIR C S. Two-phase flow in vertical and inclined annuli[J]. International Journal of Multiphase Flow, 1992, 18 (2) : 279 –293. DOI:10.1016/0301-9322(92)90089-Y
[19] SHI H, HOLMES J A, DURLOFSKY L J, et al. Drift-flux modeling of two-phase flow in wellbores[J]. SPE Journal, 2005, 10 (1) : 24 –33. DOI:10.2118/84228-PA
[20] ZUBER N, FINDLAY J A. Average volumetric concentration in two-phase flow systems[J]. Journal of Heat Transfer, 1965, 87 (4) : 453 –468. DOI:10.1115/1.3689137
[21] SUN X, KURAN S, ISHII M. Cap bubbly-to-slug flow regime transition in a vertical annulus[J]. Experiments in Fluids, 2004, 37 (3) : 458 –464.
[22] KULKARNI A A, JOSHI J B. Bubble formation and bubble rise velocity in gas-liquid systems:A review[J]. Industrial & Engineering Chemistry Research, 2005, 44 (16) : 5873 –5931.
[23] WANG Z, PEDEN J M, LEMANCZYK R Z. Gas kick simulation study for horizontal wells[C]. IADC/SPE 27498, 1994.