引用本文
王武杰, 崔国民, 魏耀奇, 潘杰. 倾斜气井临界携液流速预测新模型[J]. 石油勘探与开发, 2021,48(5): 1053-1060.
WANG Wujie, CUI Guomin, WEI Yaoqi, PAN Jie. A new model for predicting the critical liquid-carrying velocity in inclined gas wells[J]. Petroleum Exploration & Development, 2021,48(5): 1053-1060.  
Doi: 10.11698/PED.2021.02.17
Permissions
Copyright©2021, 《石油勘探与开发》编辑部
《石油勘探与开发》编辑部 所有
倾斜气井临界携液流速预测新模型
王武杰1, 崔国民1, 魏耀奇2, 潘杰3
1.上海理工大学能源与动力工程学院,上海市动力工程多相流动与传热重点实验室,上海 200093
2.国家石油天然气管材工程技术研究中心,陕西宝鸡 721008
3.西安石油大学石油工程学院,西安 710065
联系作者简介：崔国民（1969-）,男,吉林双辽人,博士,上海理工大学教授,主要从事系统优化、仿真和控制技术研究。地址：上海市杨浦区,上海理工大学能源与动力工程学院,邮政编码：200093。E-mail: cgm@usst.edu.cn

第一作者简介：王武杰（1988-）,男,陕西宝鸡人,硕士,上海理工大学博士研究生,主要从事倾斜管道内气液两相流动特性研究。地址：上海市杨浦区,上海理工大学能源与动力工程学院,邮政编码：200093。E-mail: wwj_usst@163.com

收稿日期: 2020-11-11 修回日期: 2021-08-17
基金项目: 国家自然科学基金资助项目（21978171）
摘要

基于倾斜气井气液两相分层流假设,考虑润湿性和表面张力对液膜沿井筒内壁周向分布的影响,通过考察气液相界面形状变化对单位管长气液两相系统势能、动能和表面自由能的影响,利用能量最小原理建立了临界条件下的气液相分布计算模型;考虑倾角的影响,建立了相界面摩擦因子计算模型并最终闭合控制方程。研究表明：低持液率流动条件下,相界面形状更容易受到壁面润湿性和表面张力的影响呈现弯曲,管径越小、持液率越小、井斜角越小、气体流速越大、气相密度越大时,相界面弯曲越明显;气井临界携液流速随井斜角的增加呈现非线性先增后减的趋势;最大临界携液流速对应的井斜角随井筒直径的增加而增大,同时也受气液两相物性的影响。该模型预测临界携液流速与临界压力梯度的平均相对误差分别为1.19%和3.02%,现场倾斜气井积液误判率2.38%,可对倾斜气井积液进行有效判断。图10表3参17

关键词: 倾斜气井; 气液相分布; 相界面摩擦因子; 临界携液流速; 井底积液
中图分类号:TE355 文献标志码:A 文章编号:1000-0747(2021)05-1053-08
A new model for predicting the critical liquid-carrying velocity in inclined gas wells
WANG Wujie1, CUI Guomin1, WEI Yaoqi2, PAN Jie3
1. Shanghai Key Laboratory of Multiphase Flow and Heat Transfer in Power Engineering, School of Energy and Power Engineering, University of Shanghai for Science and Technology, Shanghai 200093, China
2. Chinese National Engineering Research Center for Petroleum and Natural Gas Tubular Goods, Baoji 721008, China
3. College of Petroleum Engineering, Xi’an Shiyou University, Xi’an 710065, China
Abstract

Based on the assumption of gas-liquid stratified flow pattern in inclined gas wells, considering the influence of wettability and surface tension on the circumferential distribution of liquid film along the wellbore wall, the influence of the change of the gas-liquid interface configuration on the potential energy, kinetic energy and surface free energy of the two-phase system per unit length of the tube is investigated, and a new model for calculating the gas-liquid distribution at critical conditions is developed by using the principle of minimum energy. Considering the influence of the inclination angle, the calculation model of interfacial friction factor on the phase interface is established, and finally closed the governing equations. The interface shape is more vulnerable to wettability and surface tension at a low liquid holdup, resulting in a curved interface configuration. The interface is more curved when the smaller is the pipe diameter, or the smaller the liquid holdup, or the smaller the deviation angle, or the greater gas velocity, or the greater the gas density. The critical liquid-carrying velocity increases nonlinearly and then decreases with the increase of deviation angle. The inclination corresponding to the maximum critical liquid-carrying velocity increases with the increase of the diameter of the wellbore, and it is also affected by the fluid properties of the gas phase and liquid phase. The mean relative errors for critical liquid-carrying velocity and critical pressure gradient are 1.19% and 3.02%, respectively, and the misclassification rate is 2.38% in the field trial, implying the new model can provide a valid judgement on the liquid loading in inclined gas wells.

Keyword: inclined gas well; gas-liquid phase distribution; interfacial friction factor; critical liquid-carrying velocity; bottom-hole liquid loading
0 引言

准确预测临界携液流速对诊断气井积液状况、解决含水气井稳产问题、挖掘煤炭地下气化开发潜力至关重要[1, 2, 3], 同时, 其对促进致密油气快速规模发展, 早日实现中国“ 能源独立” 战略具有积极意义[4, 5]。越来越多的实验数据表明, 倾斜气井的携液能力最差, 预测也最为复杂和困难[6, 7]。因此, 有必要对倾斜气井临界携液流速的预测展开深入研究。分层流作为倾斜管道气液两相流动的主要流型被普遍认可。分层流流型下, 井筒内气液两相的周向分布对流动特性有至关重要的影响, 尤其对液相流量较小而气相流量较大的情形更甚[8, 9]。目前对临界携液状态下的气液相分布缺乏量化研究, 且尚无倾斜管道气液两相临界流动条件下的相界面摩擦因子计算模型。基于此, 本文分别考察单位管长气液两相系统重力势能、表面自由能和气相动能的改变量, 并根据能量最小原理确定临界携液状态下气液相分布的结构参数; 同时采用实验数据修正相界面摩擦因子计算模型, 使预测模型封闭; 并将模型应用于倾斜气井临界携液流速的预测, 以评估模型的应用效果。

1 实验系统

本文通过可视化实验, 模拟气井携液的动态过程, 实验装置如图1所示。

图1 实验装置图

分别采用水和甘油、水的混合溶液（二者体积分数均为50%）模拟井下积液, 利用空气模拟天然气。实验条件下空气密度ρ a=1.184 kg/m3, 空气动力黏度μ a=1.849× 10-5 Pa· s; 水密度ρ w=997.05 kg/m3, 水动力黏度μ w=8.900 8× 10-4 Pa· s; 水、甘油混合溶液密度ρ gw=1 139.56 kg/m3, 水、甘油混合溶液动力黏度μ gw= 6.882× 10-3 Pa· s; 空气、水表面张力σ aw=0.071 97 N/m, 空气与水、甘油混合溶液表面张力σ agw=0.067 11 N/m。

分别选用外径50 mm、壁厚5 mm和外径70 mm、壁厚5 mm的有机玻璃管进行实验, 步骤为：①在管内低洼处注入一定体积的液体, 通过风机增压向玻璃管内注入压缩空气, 计量压力和流量; ②采用变频器实现气体流速的连续可调, 直至低洼处的积液完全铺展于上倾段并与管壁相对静止（此时达临界气体流速）。

2 模型建立与求解方法
2.1 受力分析与控制方程

临界气体流速条件下, 倾斜管内流型为波状分层流, 液膜在重力、相界面和壁面切应力的作用下处于平衡状态。此时, 近气相侧的液膜在相界面切应力的作用下向下游运动, 而近壁面侧的液膜在重力作用下回流, 即液膜内部存在循环流动。倾斜井筒微元段内液膜受力如图2所示。

图2 倾斜井筒内液膜流动与受力分析

上述流动过程可简化为稳态一维流动, 并忽略加速压降, 可对气、液相分别建立动量方程：

${{A}_{\text{g}}}{{\left( -\frac{\text{d}p}{\text{d}x} \right)}_{\text{g}}}={{\tau }_{\text{wg}}}{{S}_{\text{g}}}+{{\tau }_{\text{i}}}{{S}_{\text{i}}}+{{A}_{\text{g}}}{{\rho }_{\text{g}}}g\cos \beta $ （1）

${{A}_{\text{l}}}{{\left( -\frac{\text{d}p}{\text{d}x} \right)}_{\text{l}}}=-{{\tau }_{\text{wl}}}{{S}_{\text{l}}}-{{\tau }_{\text{i}}}{{S}_{\text{i}}}+{{A}_{\text{l}}}{{\rho }_{\text{l}}}g\cos \beta $ （2）

流动过程中近似认为气相与液相的平均压力梯度相等, 联立（1）式和（2）式可得控制方程：

${{\tau }_{\text{wg}}}{{S}_{\text{g}}}\frac{{{A}_{\text{l}}}}{{{A}_{\text{g}}}}+{{\tau }_{\text{wl}}}{{S}_{\text{l}}}+{{\tau }_{\text{i}}}{{S}_{\text{i}}}\left( 1+\frac{{{A}_{\text{l}}}}{{{A}_{\text{g}}}} \right)-\left( {{\rho }_{\text{l}}}-{{\rho }_{\text{g}}} \right){{A}_{\text{l}}}g\cos \beta =0$ （3）

上式中, 几何参数与气液相分布有关, 应力参数与持液率、各相速度和摩擦因子有关。通过闭合几何参数和应力参数可使控制方程闭合。

2.2 气液相分布

根据Brauner等[10]的研究, 液膜的周向分布如图3所示。其中, O和O1分别表示井筒流道的截面轴心和虚拟圆圆心。控制方程中的几何参数采用下式计算：

${{A}_{\text{l}}}={{R}^{2}}\left\{ {{\varphi }_{0}}-\frac{1}{2}\sin \left( 2{{\varphi }_{0}} \right)-\frac{{{\sin }^{2}}{{\varphi }_{0}}}{{{\sin }^{2}}{{\varphi }_{\text{PA}}}}\left[ {{\varphi }_{\text{PA}}}-\text{ }\!\!\pi\!\!\text{ }-\frac{1}{2}\sin \left( 2{{\varphi }_{\text{PA}}} \right) \right] \right\}$ （4）

${{A}_{\text{g}}}={{A}_{\text{T}}}-{{A}_{\text{l}}}$ （5）

${{S}_{\text{l}}}=2R{{\varphi }_{0}}$ （6）

${{S}_{\text{g}}}=\text{ }\!\!\pi\!\!\text{ }D-{{S}_{\text{l}}}$ （7）

${{S}_{\text{i}}}=R\frac{2\left( \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{\text{PA}}} \right)\sin {{\varphi }_{0}}}{\sin {{\varphi }_{\text{PA}}}}$ （8）

图3 井筒周向的液膜分布

不同相界面假设下的几何参数的计算式见表1。

表1 不同相界面下的几何参数计算关系式
2.3 能量最小原理

在相界面形状由水平变为弯曲的过程中, 重力、切应力和表面张力将引起系统势能、动能和表面自由能的变化。临界条件下, 单位管长气液两相系统总能的改变量可表示为：

$\Delta {{E}_{\text{tot}}}=\Delta {{E}_{\text{p}}}+\Delta {{E}_{\text{s}}}+\Delta {{E}_{\text{k, g}}}$ （9）

基于能量最小原理, 单位管长系统总能量最小时, 相界面结构最稳定, 即：

$\frac{\text{d}\Delta {{E}_{\text{tot}}}}{\text{d}{{\varphi }_{\text{PA}}}}=0$ （10）

①单位管长重力势能的改变量[10]

$\Delta {{E}_{\text{p}}}=\left( {{A}_{\text{g}}}{{\rho }_{\text{g}}}+{{A}_{\text{l}}}{{\rho }_{\text{l}}} \right)\left[ g\left( {{Y}_{\text{GPC}}}-{{Y}_{\text{GP}}} \right)\cos \beta \right]$ （11）

其中

${{Y}_{\text{GP}}}=-\frac{2}{3}R\frac{{{\sin }^{3}}{{\varphi }_{0\text{P}}}}{{{\varphi }_{0\text{P}}}-\frac{1}{2}\sin \left( 2{{\varphi }_{0\text{P}}} \right)}$ （12）

${{Y}_{\text{GPC}}}=\frac{R}{\text{ }\!\!\pi\!\!\text{ }}\frac{1}{\varepsilon }\frac{{{\sin }^{3}}{{\varphi }_{0}}}{{{\sin }^{2}}{{\varphi }_{\text{PA}}}}\left( \cot {{\varphi }_{\text{PA}}}- \right.\ \left. \cot {{\varphi }_{0}} \right)\left[ \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{\text{PA}}}+\frac{1}{2}\sin \left( 2{{\varphi }_{\text{PA}}} \right) \right]$ （13）

②单位管长表面自由能的改变量

$\Delta {{E}_{\text{s}}}={{\sigma }_{\text{gl}}}\left( \Delta {{A}_{\text{gl}}}-\Delta {{A}_{\text{lw}}}\cos \theta \right)$ （14）

其中

$\Delta {{A}_{\text{lw}}}=2R\left( {{\varphi }_{0}}-{{\varphi }_{0\text{P}}} \right)$ （15）

$\Delta {{A}_{\text{gl}}}=2R\left[ \sin {{\varphi }_{0}}\frac{\left( \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{\text{PA}}} \right)}{\sin {{\varphi }_{\text{PA}}}}-\sin {{\varphi }_{0\text{P}}} \right]$ （16）

③单位管长气相动能改变量

根据七分之一次方定律, 圆管内湍流气相的速度分布表示如下：

$\frac{{{u}_{\text{g}}}}{{{u}_{\text{g, }\max }}}={{\left( 1-\frac{r}{R} \right)}^{\frac{1}{7}}}$ （17）

其中

${{u}_{\text{g}, \max }}=\frac{\overline{{{u}_{\text{g}}}}}{0.817}$ （18）

单位管长气相动能改变量采用下式计算：

$\Delta {{E}_{\text{k, g}}}=\frac{1}{2}{{\rho }_{\text{g}}}\left[ \int_{-{{\varphi }_{0}}}^{{{\varphi }_{0}}}{\int_{0}^{R-{{\delta }_{\text{l}}}}{u_{\text{g}}^{2}r\text{d}r\text{d}\varphi }}-\int_{-{{\varphi }_{0\text{P}}}}^{{{\varphi }_{0\text{P}}}}{\int_{0}^{R-{{\delta }_{\text{l}}}}{u_{\text{g}}^{2}r\text{d}r\text{d}\varphi }} \right]$ （19）

2.4 剪切应力计算

Taitel等[11]给出了相关应力计算表达式：

${{\tau }_{\text{wg}}}=\frac{1}{2}{{f}_{\text{g}}}{{\rho }_{\text{g}}}u_{\text{g}}^{2}$ （20）

${{\tau }_{\text{wl}}}=\frac{1}{2}{{f}_{\text{l}}}{{\rho }_{\text{l}}}u_{\text{l}}^{2}$ （21）

${{\tau }_{\text{i}}}=\frac{1}{2}{{f}_{\text{i}}}{{\rho }_{\text{g}}}{{\left( {{u}_{\text{g}}}-{{u}_{\text{i}}} \right)}^{2}}$ （22）

Banafi等[12]实验研究指出Taitel模型可准确预测气液两相低持液率流动条件下的气壁面摩擦因子, 表达式如下：

${{f}_{\text{g}}}={{C}_{\text{g}}}Re_{\text{g}}^{-m}={{C}_{\text{g}}}{{\left[ \frac{{{\rho }_{\text{g}}}{{D}_{\text{H}, \text{g}}}{{u}_{\text{sg}}}}{{{\mu }_{\text{g}}}\left( 1-\varepsilon \right)} \right]}^{-m}}$ （23）

其中

${{D}_{\text{H}, \text{g}}}=\frac{{{A}_{\text{g}}}}{{{S}_{\text{g}}}+{{S}_{\text{i}}}}$ （24）

层流时, Cg=16, m=1; 湍流时, Cg=0.046, m=0.2。

Biberg[13]指出, 倾斜管道光滑分层流条件下, 平均液壁面切应力与平均相界面切应力满足下式：

$\overline{{{\tau }_{\text{wl}}}}=f\left( {{\varphi }_{0\text{P}}} \right)\overline{{{\tau }_{\text{i}}}}-{{\left. \frac{8{{\mu }_{\text{l}}}{{u}_{\text{l}}}}{{{D}_{\text{l}}}\left( {{\varphi }_{0\text{P}}} \right)} \right|}_{{{u}_{\text{l}}}=0}}=f\left( {{\varphi }_{0\text{P}}} \right)\overline{{{\tau }_{\text{i}}}}$ （25）

上式中$f\left( {{\varphi }_{0\text{P}}} \right)$是难以获得解析解的积分项。Biberg[13]给出以下近似计算关系式：

$f\left( {{\varphi }_{0\text{P}}} \right)\approx \frac{{{a}_{0}}+{{a}_{1}}\left( \frac{{{\varphi }_{0\text{P}}}}{\text{ }\!\!\pi\!\!\text{ }} \right)+...+{{a}_{4}}{{\left( \frac{{{\varphi }_{0\text{P}}}}{\text{ }\!\!\pi\!\!\text{ }} \right)}^{4}}}{{{b}_{0}}+{{b}_{1}}\left( \frac{{{\varphi }_{0\text{P}}}}{\text{ }\!\!\pi\!\!\text{ }} \right)+...+{{b}_{7}}{{\left( \frac{{{\varphi }_{0\text{P}}}}{\text{ }\!\!\pi\!\!\text{ }} \right)}^{7}}}$ （26）

式中, an和bn的取值见表2[13]。该近似值的最大绝对误差不超过2× 10-6。

表2 近似关系式（（26）式）中参数an和bn的取值

准确计算相界面摩擦因子是预测相界面切应力的关键。倾斜管道中, 界面波和持液率对相界面摩擦因子的影响较大。界面波动特性主要取决于机械能从气相向液相的传递速率以及液体中黏性耗散速率, 而这一过程受倾角的影响[14]。随着管道倾角的增加, 重力对界面波生长的约束减弱, 使得倾斜管中的液层表面更易于观察到界面波。本文认为在倾斜管道中, 重力对倾斜管内临界状态下的界面波动作用更加显著。因此, 一方面将界面波等效为管壁的相对粗糙度, 另一方面, 引入倾角对相界面摩擦因子进行修正[14], 具体过程如下：

$F=\frac{{{f}_{\text{i}}}}{\left( 0.05+{{f}_{\text{i}}} \right){{\left( 1-\varepsilon \right)}^{1.5}}}\left( \frac{{{u}_{\text{sg}}}}{\sqrt{gD}} \right){{\left( \frac{{{\sigma }_{\text{gl}}}}{{{\mu }_{\text{l}}}\sqrt{gD}} \right)}^{0.04}}{{\left( \frac{{{\rho }_{\text{l}}}g{{D}^{2}}}{{{\sigma }_{\text{gl}}}} \right)}^{0.22}}$ （27）

$\frac{k}{D}=\frac{0.514\ 5\varepsilon }{s_{\text{i}}^{1.5}}\left\{ \tanh \left[ 0.057\ 6\left( F-33.74 \right) \right]+0.945 \right\}$ （28）

${{f}_{\text{i}}}=\frac{0.074\ 52}{{{\left\{ \lg \left[ \frac{15}{R{{e}_{\text{g}}}}+\frac{1}{3.71}\left( \frac{k}{D} \right) \right] \right\}}^{2}}}$ （29）

其中, 无因次相界面湿周定义为：

${{s}_{\text{i}}}=\frac{{{S}_{\text{i}}}}{\text{ }\!\!\pi\!\!\text{ }D}$ （30）

再考虑倾角的影响对其进行修正：

${{f}_{\text{i}}}=f\left( \alpha \right)\frac{0.074\ 52}{{{\left\{ \lg \left[ \frac{15}{R{{e}_{\text{g}}}}+\frac{1}{3.71}\left( \frac{k}{D} \right) \right] \right\}}^{2}}}$ （31）

$f\left( \alpha \right)=A\alpha +B$ （32）

式中, A和B与流体物性和管道条件有关[15]。

2.5 模型求解方法

控制方程（（3）式）可以视为持液率的函数, 即：$f\left( \varepsilon \right)=0$。零始终是方程$f\left( \varepsilon \right)=0$的一个解。$f{{\left( \varepsilon \right)}^{\frac{1}{3}}}$的函数图像如图4所示。当usg=usg, cr时, $f\left( \varepsilon \right)=0$只有一个非零解, 且该非零解在（0, 0.1）区间求得, 求解步骤如图5所示。

图4 函数$f{{\left( \varepsilon \right)}^{1/3}}$在不同流动条件下的变化曲线

图5 模型求解过程示意图

3 结果分析与实例验证
3.1 气液相分布预测结果与分析

对控制方程无因次化处理后可得：

$\frac{2\varepsilon }{\text{ }\!\!\pi\!\!\text{ }{{\left( 1-\varepsilon \right)}^{3}}}\left[ {{f}_{\text{g}}}\left( \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{0}} \right)+{{f}_{\text{i}}}\frac{\left( \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{\text{PA}}} \right)\sin {{\varphi }_{0}}}{\sin {{\varphi }_{\text{PA}}}} \right]+ \frac{2{{f}_{\text{i}}}}{\text{ }\!\!\pi\!\!\text{ }{{\left( 1-\varepsilon \right)}^{2}}}\left[ f\left( {{\varphi }_{0\text{P}}} \right){{\varphi }_{0}}+\frac{\left( \text{ }\!\!\pi\!\!\text{ }-{{\varphi }_{\text{PA}}} \right)\sin {{\varphi }_{0}}}{\sin {{\varphi }_{\text{PA}}}} \right]=\frac{Bo}{W{{e}_{\text{sg}}}}\varepsilon \cos \beta $ （33）

其中

$Bo=\frac{\left( {{\rho }_{1}}-{{\rho }_{\text{g}}} \right)g{{D}^{2}}}{{{\sigma }_{\text{gl}}}}$ （34）

$W{{e}_{\text{sg}}}=\frac{{{\rho }_{\text{g}}}u_{\text{sg}}^{2}D}{{{\sigma }_{\text{gl}}}}$ （35）

由上式可知, 临界状态下的气液相分布与邦德数、气相表观韦伯数、持液率和井斜角有关。基于能量最小原理, 量化气液两相系统总能量随相界面结构的变化, 并以空气和水为例, 分别考察上述4个因素对管内气液相分布的影响, 计算结果见图6（图中红色数据点为各工况下的能量最小值点）。可以看到：①随着管径的增加, 最稳定的气液相分布（对应能量最小值）所对应的相界面弯曲角从125° （微凸）逐渐趋于180° （水平）; 反之, 管径越小, 在润湿性和表面张力的影响下相界面弯曲角越小, 呈现凸型结构的液相分布（见图6a）; ②随着井斜角的增加, 重力抵消润湿性和表面张力的作用效果更加显著, 使最稳定的气液相分布所对应的相界面弯曲角从120° （微凸）趋近于180° （水平）（见图6b）; ③结合图6c、图6d可知, 气体流速和气相密度增加, 使气相动能也增加, 在润湿性和表面张力的影响下相界面弯曲角减小（变凸）有利于降低系统表面自由能; ④随着持液率的增加, 液相重力势能的影响增大, 有利于稳定原来水平相界面的形状, 同时界面微凸也可减小持液率增加导致的系统表面自由能的增加。但是在液相密度远大于气相密度的低压系统下, 液相重力势能起主导作用（见图6e）, 随着持液率的增加, 最稳定的气液相分布所对应的相界面弯曲角从105° 逐渐趋于165° （近水平）。综上所述, 低持液率流动条件下, 受表面张力和润湿性的影响, 相界面呈现微凸, 符合如图7所示的实验观察结果。

图6 液相分布的模型预测结果

图7 实验观察的临界状态下的气液相分布示意图

3.2 实验结果与分析

图8为临界携液气体流速和临界压力梯度的实验结果。由图可知, 临界携液流速和临界压力梯度随着倾角的增加呈非线性增加。在相同的管道倾角下, 液相密度或者黏度越大, 临界携液流速越大。这与Birvalski等[16]的实验结果一致。同时, 管径增加使临界携液气体流速增加。这是因为在相同持液率下, 界面形状的变化不及Bo和Wesg数对临界携液流速的影响显著, 要使控制方程依然成立, 必须提高气体流速, 以抵消管径增加引起的Bo/Wesg的增量。

图8 临界携液流速和临界压力梯度的实验结果

3.3 倾斜管内相界面摩擦因子计算模型

采用Birvalski等[16]和本文的实验数据, 基于简化后的水平相界面形状封闭几何参数, 对不同倾角下的相界面摩擦因子修正参数进行反算, 并通过数据拟合得到相界面摩擦因子修正参数（$f\left( \alpha \right)$）, 计算结果如表3所示。采用表3中相界面摩擦因子修正参数对相界面摩擦因子计算公式（（31）式）进行修正后, 模型的预测值达到了非常高的精度。皮尔逊相关系数相对较高。

表3 不同实验条件下的相界面摩擦因子修正参数

由表3可知, $f\left( \alpha \right)$的斜率为正值, 即表明相界面摩擦因子随着倾角的增大整体呈现递增的趋势。主要是因为倾角的增加使得约束界面波动的重力分量减小, 导致界面波动更加剧烈和明显。同时, $f\left( \alpha \right)$小于1, 表明在微倾管道内低持液率流动条件下, 相界面摩擦因子要小于相同气体流速下水平管道相界面摩擦因子预测值, 主要原因是当管道接近水平时, 倾角对界面波的影响不如持液率的影响显著。本文基于实验数据, 结合无因次化分析, 给出相界面摩擦因子修正参数的计算式：

$f\left( \alpha \right)=0.092\ 9{{\left( \frac{{{\mu }_{\text{l}}}\sqrt{gD}}{{{\sigma }_{\text{gl}}}} \right)}^{-0.116}}\alpha +0.002\ 9B{{o}^{0.860\ 4}}$ （36）

3.4 模型对比与工程实例验证

将模型预测结果与实验数据进行对比（见图9）。本文模型预测的临界携液流速随管径的增加而增加, 随液相密度和黏度的增加而增加（同管径下, 甘油溶液的临界携液流速大于水的临界携液流速）, 这一规律与符合Birvalski等[16]和Rastogi等[17]的实验结果; 临界压力梯度的预测结果随管径的增加而减小, 符合实验结果, 管径的增加使液相水力直径增加是压力梯度减小的重要原因。对比本文实验数据, 新模型预测临界携液流速和临界压力梯度的平均相对误差分别为1.19%和3.02%。同时, 临界携液流速随倾角的增加先增大后减小, 最大临界携液流速对应的管道倾角随管径的增加而增加, 同时也受流动介质物性的影响。本文模型模拟计算发现, 标况下天然气携矿化水时, 最大临界携液流速对应的管道倾角分别为30° （D=50.67 mm）、37° （D=62.00 mm）、40° （D=75.9 mm）和45° （D=100.53 mm）。

图9 本文模型预测结果与实验结果对比

本文模型应用于高压系统时需要对预测结果进行转换。本文从能量角度出发, 认为高、低压管道系统内, 临界条件下单位体积气液两相间传递的机械能应相等, 即：

$u_{\text{sg, cr, hp}}^{{}}=\sqrt{\frac{u_{\text{sg, cr}}^{\text{2}}{{\rho }_{\text{g}}}}{\rho _{\text{g, hp}}^{{}}}}$ （37）

为便于比较, 将临界携液流速转换为标况下的临界携液流量：

${{Q}_{\text{cr}}}=2.5\times {{10}^{2}}\frac{{{A}_{\text{T}}}pu_{\text{sg, cr, hp}}^{{}}}{ZT}$ （38）

当天然气的实际流速大于倾斜气井中临界携液流速最大值时, 气井才能连续携液。气田现场需结合方位测井资料来确定倾斜气井的最大临界携液流速。选取中国普光气田和延长气田42口气井的生产数据, 并与实际产气量绘制交会图（见图10）, 可以看到本文模型仅误判1口气井, 误判率2.38%, 说明本文模型的预测结果具有较高的可信度, 可对倾斜气井积液进行有效判断。

图10 本文模型预测结果与对比

4 结论

低持液率流动条件下, 相界面形状更容易受到壁面润湿性和表面张力的影响呈现弯曲, 管径越小、持液率越小、井斜角越小、气体流速越大、气相密度越大时, 相界面弯曲越明显。

气井临界携液流速随井斜角的增加呈现非线性先增后减的趋势; 最大临界携液流速对应的井斜角随井筒直径的增加而增大, 同时也受气液两相物性的影响。

本文模型预测临界携液流速与临界压力梯度的平均相对误差分别为1.19%和3.02%, 现场倾斜气井积液误判率2.38%, 可对倾斜气井积液进行有效判断。

符号注释：

an, bn, A, B—常系数, 无因次; AT—油管截面面积, m2; Ag—气相流动通道面积, m2; Δ Agl—单位管长气液相接触面面积的变化, m2/m; Al—液相流动通道面积, m2; Δ Alw—单位管长液壁接触面面积的变化, m2/m; Bo—邦德数, 无因次; Cg, m—与流态有关的系数, 无因次; D—油管内径, m; DH, g—气相水力直径, m; ${{D}_{\text{l}}}\left( {{\varphi }_{0\text{P}}} \right)$—不同液相分布角下的等效层流直径, m; fg—气相摩擦因子, 无因次; fi—相界面摩擦因子, 无因次; fl—液相摩擦因子, 无因次; $f\left( \alpha \right)$—相界面摩擦因子修正参数, 无因次; $f\left( {{\varphi }_{0\text{P}}} \right)$—与液相分布角有关的函数, 无因次; F—中间变量, 无因次; g—重力加速度, m/s2; k—相界面绝对粗糙度, m; n—常系数编号; p—压力, Pa; Qcr—临界携液流量, m3/d; r—距离井筒轴心的径向位置, m; R—油管半径, m; Reg—气相雷诺数, 无因次; R1—虚拟圆半径, m; Sg—气相湿周, m; Si—相界面湿周, m; si—无因次相界面湿周, 无因次; Sl—液相湿周, m; T—热力学温度, K; ug—气相真实流速, m/s; ${{\bar{u}}_{\text{g}}}$—平均气相真实速度, m/s; ug, max—最大气相真实速度, m/s; ui—相界面真实流速, m/s; ul—液相真实流速, m/s; usg—气相表观流速, m/s; Δ usg—气相表观流速迭代步长, m/s; usg, cr—临界携液气相表观流速, m/s; usg, cr, hp—高压系统下的临界携液气相表观流速, m/s; Wesg—气相表观韦伯数; x—流动方向坐标, m; YGPC—弯曲相界面的重心位置, m; YGP—水平相界面的重心位置, m; Z—压缩因子; α —管道倾角, （° ）; β —井斜角, （° ）; Δ Ek, g—单位管长气相的动量改变量, J/m; Δ Ep—单位管长气液两相系统的重力势能改变量, J/m; Δ Es—单位管长气液两相系统的表面自由能改变量, J/m; Δ Etot—单位管长气液两相系统的总能量改变量, J/m; δ l—液膜厚度, m; ${{\delta }_{\text{l}}}\left( \varphi \right)$—液膜厚度分布, m; ε —持液率; ε 0—持液率迭代初值; Δ ε —持液率迭代步长; θ —润湿角, rad; μ a—空气动力黏度, Pa· s; μ g—气相动力黏度, Pa· s; μ gw—水、甘油混合溶液动力黏度, Pa· s; μ l—液相动力黏度, Pa· s; μ w—水动力黏度, Pa· s; ρ a—空气密度, kg/m3; ρ g—气相密度, kg/m3; ρ g, hp—高压系统下的气相密度, kg/m3; ρ gw—水、甘油混合溶液密度, kg/m3; ρ l—液相密度, kg/m3; ρ w—水密度, kg/m3; σ aw—空气、水界面的表面张力, N/m; σ agw—空气与水、甘油混合溶液的表面张力, N/m; σ gl—气液界面的表面张力, N/m; τ i—相界面剪切应力, Pa; ${{\bar{\tau }}_{\text{i}}}$—平均相界面切应力, Pa; τ ig—相界面反作用于气相的相界面切应力, Pa; τ gi—气相作用于相界面的相界面切应力, Pa; τ wg—气壁面剪切应力, Pa; τ wl—液壁面剪切应力, Pa; $\overline{{{\tau }_{\text{wl}}}}$—平均液壁面切应力, Pa; φ —任意液相分布角, rad; φ 0—弯曲相界面的液相分布角, rad; φ 0P—水平相界面的液相分布角, rad; φ PA—相界面弯曲角, rad。下标：g— 气相; l— 液相。

（编辑 唐俊伟）

参考文献
[1] 潘杰, 王武杰, 王亮亮, 等. 考虑液滴夹带的气井连续携液预测模型[J]. 石油学报, 2019, 40(3): 332-336.
PAN Jie, WANG Wujie, WANG Liangliang, et al. A prediction model for continuous liquid-carrying in gas wells considering droplet entrainment[J]. Acta Petrolei Sinica, 2019, 40(3): 332-336. [本文引用:1]
[2] 雷群, 翁定为, 罗健辉, 等. 中国石油油气开采工程技术进展与发展方向[J]. 石油勘探与开发, 2019, 46(1): 139-145.
LEI Qun, WENG Dingwei, LUO Jianhui, et al. Achievements and future work of oil and gas production engineering of CNPC[J]. Petroleum Exploration and Development, 2019, 46(1): 139-145. [本文引用:1]
[3] 邹才能, 陈艳鹏, 孔令峰, 等. 煤炭地下气化及对中国天然气发展的战略意义型[J]. 石油勘探与开发, 2019, 46(2): 195-204.
ZOU Caineng, CHEN Yanpeng, KONG Lingfeng, et al. Underground coal gasification and its strategic to the development of natural gas industry in China[J]. Petroleum Exploration and Development, 2019, 46(2): 195-204. [本文引用:1]
[4] 孙龙德, 邹才能, 贾爱林, 等. 中国致密油气发展特征与方向[J]. 石油勘探与开发, 2019, 46(6): 1015-1026.
SUN Longde, ZOU Caineng, JIA Ailin, et al. Development characteristics and orientation of tight oil and gas in China[J]. Petroleum Exploration and Development, 2019, 46(6): 1015-1026. [本文引用:1]
[5] 邹才能, 潘松圻, 赵群. 论中国“能源独立”战略的内涵、挑战及意义[J]. 石油勘探与开发, 2020, 47(2): 416-426.
ZOU Caineng, PAN Songqi, ZHAO Qun. On the connotation, challenge and significance of China’s “energy independence” strategy[J]. Petroleum Exploration and Development, 2020, 47(2): 416-426. [本文引用:1]
[6] 王琦, 李颖川, 王志彬, 等. 水平气井连续携液实验研究及模型评价[J]. 西南石油大学学报(自然科学版), 2014, 36(3): 139-145.
WANG Qi, LI Yingchuan, WANG Zhibin, et al. Experimental study and model evaluation on continuous liquid removal in horizontal gas well[J]. Journal of Southwest Petroleum University (Natural Science Edition), 2014, 36(3): 139-145. [本文引用:1]
[7] CHEN D, YAO Y, FU G, et al. A new model for predicting liquid loading in deviated gas wells[J]. Journal of Natural Gas Science and Engineering, 2016, 34: 178-184. [本文引用:1]
[8] ROVINSKY J, BRAUNER N, MARON D M. Analytical solution for laminar two-phase flow in a fully eccentric core-annular configuration[J]. International Journal of Multiphase Flow, 1997, 23(3): 523-543. [本文引用:1]
[9] BARAJAS A M, PANTON R L. The effects of contact angle on two-phase flow in capillary tubes[J]. International Journal of Multiphase Flow, 1993, 19(2): 337-346. [本文引用:1]
[10] BRAUNER N, MARON D M, ROVINSKY J. A two-fluid model for stratified flows with curved interfaces[J]. International Journal of Multiphase Flow, 1998, 24(6): 975-1004. [本文引用:2]
[11] TAITEL Y, DUKLER A E. A model for predicting flow regime transitions in horizontal and near horizontal gas-liquid flow[J]. AIChE Journal, 1976, 22: 47-55. [本文引用:1]
[12] BANAFI A, TALAIE M R. A new mechanistic model to predict gas-liquid interface shape of gas-liquid flow through pipes with low liquid loading[J]. AIChE Journal, 2015, 61(3): 1043-1053. [本文引用:1]
[13] BIBERG D. Liquid wall friction in two-phase turbulent gas laminar liquid stratified pipe flow[J]. The Canadian Journal of Chemical Engineering, 1999, 77(6): 1073-1082. [本文引用:3]
[14] GROLMAN E, FORTUIN J M H. Gas-liquid flow in slightly inclined pipes[J]. Chemical Engineering Science, 1997, 52(24): 4461-4471. [本文引用:2]
[15] 王武杰, 崔国民, 陈家星, 等. 微倾斜管道零液流量临界气速实验与建模[J]. 工程热物理学报, 2020, 41(12): 3068-3076.
WANG Wujie, CUI Guomin, CHEN Jiaxing, et al. Experiments and modelling for predicting critical gas velocity with zero net liquid flow in slightly inclined pipes[J]. Journal of Engineering Thermophysics, 2020, 41(12): 3068-3076. [本文引用:1]
[16] BIRVALSKI M, KOREN G B, HENKES R A W M. Experiments and modelling of liquid accumulation in the low elbow of a gas/liquid pipeline[R]. Banff: North American Conference on Multiphase Technology, 2014. [本文引用:3]
[17] RASTOGI A, FAN Y. Experimental and modeling study of onset of liquid accumulation[J]. Journal of Natural Gas Science and Engineering, 2019, 73: 103064. [本文引用:1]