水平通气横流中气泡的运动特性数值模拟

钱凯凯， 王治云

（上海理工大学 能源与动力工程学院，上海 杨浦 200093）

气液两相流动及其涉及到的气泡运动规律的研究广泛应用在化工、船舶、医药技术、石油运输等各种工程技术领域[1]。气泡两相流的特征是气泡分散在连续的液相中，气泡的最大尺寸远小于容器或管道的直径。在文献中，竖直流动比水平流动受到更多的关注，但水平流动往往在工业和工程应用中很常见，比如化工生产中在反应池中水平通入气体来去除液相中的杂质[1]，使用鼓泡器通入气体来增大管道内气液的反应面积与效率[2]。气体通入后引起的气泡的运动特性和参数会影响反应的速率和效率以及安全。

目前研究气泡在液相中的运动主要有两类方法，实验研究与数值模拟研究。其他的方法如理论分析、经验拟合等相对较少。其中实验研究大多是利用光学、声学来测量气泡的各种参数，如Lee等[3]通过使用高速摄像机与双针探头来获取气泡的图像和测量气泡的大小；Ran Kong等[4]通过高速摄像和图像处理技术来获取气泡的界面和局部速度；Nikitopoulos等[5]利用相位多普勒和视频成像技术获得了气泡速度。数值模拟研究的方法则主要通过求解气泡在两相流中的运动方程来获取所需的气泡参数，如Yeoh等[6]通过MUSIG和DQMOM耦合的双流体模型来预测水平管内的气泡的粒径分布；Tsui等[7]将VOF（volume of fluid，流体体积法）和Level set进行耦合，对静水中上升的气泡进行了三维数值模拟。

为研究气泡流动特性，学者们通常会在开放流场（不受边界影响）的情况下进行单气泡的动力学研究，这种方式比较基础但能够更简单与准确地对气泡形状、速度、分布、大小等进行分析，代表性综述有[8-9]。气泡的运动特性主要的影响因素有表面张力、惯性、黏性、密度等。本文研究的气泡运动则是连续的气体通入有着横流的水平管道内，研究各种因素对气泡运动特性的影响。这种气泡的运动状态十分复杂，目前相关的文献研究较少且需要更多的研究。

前文提到VOF方法具有求解简便、精度较高[7]、追踪界面锐利性好、计算成本低等优点，已经被广泛地运用于气泡两相流的模拟与分析研究中[8]。刘娟等[10]使用VOF方法模拟了离心式喷嘴内部的流动过程；Valencia等[11]使用VOF方法研究了通气速度等对气泡速度和大小的影响。

本文的研究来自于化工领域的反应器中，在具有流动的液相的管内通入气体，通入的气体及生成气泡的运动特性会影响反应器的反应速率与效率以及反应过程中的安全。本文基于VOF方法以及界面重构法（PLIC），通过数值模拟的方法，对水平圆管内在恒定液相速度和气相速度下气泡的运动状态进行追踪。然后通过自编后处理程序获得结果数据并处理，得到整个过程中气泡粒径的大小、数量与分布，分析气泡运动中的粒径和数量变化，从而为相关的气液两相流中气泡运动特性的研究提供参考。

1 仿真模拟和计算方法

1.1 数值计算模型

如图1所示为计算模型示意图，计算区域中采用液相、气相的速度入口边界条件，出口采用压力出口边界条件，其他壁面采用无滑移边界条件。主要模型参数：长1.5 m，直径0.1 m的水平管道，中间的进气管长0.05 m，直径为0.005 m。气相和液相分别从中间的进气管和从进气管的两侧进入水平圆管。计算过程为，开始阶段通入液相并充满管道，等管道内液相流动达到动态平衡后打开进气管，开始通入气相。

图1 模型示意图

1.2 计算方法

本项目模拟的流动是气液两相流动，需采用多相流计算模型。VOF方法是一种前端捕捉的界面追踪方法，由Hirt和Nichols在1981年提出，提出了相体积分数的概念，其值是控制体积单元中流体所占体积与该控制单元体积之比（0～1）[12]。考虑到气泡数量相对较少，为了获得气泡的直径等参数，本文采用VOF这种模型，可以捕捉和追踪气液两相分界面。速度与压力的耦合求解采用的是PISO算法[13]。

1.3 控制方程

1）连续性方程

2）动量方程

VOF模型中的表面张力模型是Brackbill[14]等提出的连续表面力模型（Continuum Surface Force，CSF）。在VOF模型中需要考虑表面张力时，需要在动量方程中增加源项。该模型把表面张力定义为一种作用于相界面区域网格单元内流体的体积力Fvol，并将其引入到动量方程（2）中，如式（3）所示。

其中，p为压强，ρ为流体密度，∇为哈密顿算子，v为速度矢量，μ为粘性系数，Fvol为表面张力源项，g为重力加速度，σij为i、j相流体表面应力，α i、α j分别为i、j相流体体积分数，ρ i、ρ j分别为i、j相流体密度，ki、kj分别为i、j相流体表面曲率。

1.4 网格无关性验证和模型验证

选择三组典型工况，重点考察网格数量对脱离气泡平均直径的影响，综合对比来选择较为合理的网格数量。图2是不同工况下脱离气泡直径与网格数量的关系，从图2可以看出，当达到60万个时，直径不再有较大的变化，所以选用60万个规模的网格可以获得不错的准确性和效率。本文采用区域结构化划分网格，生成网格都是四边形网格，总共59.9万个且网格质量优良。

图2 脱离气泡直径和网格数量的关系

关于模型的验证和分析，马斗[15]重点考察了气泡上升时的速度和直径变化，与众多的文献结果进行对比，符合程度比较好，目前有很多文献[7,12-13,15-16]都验证了该模型对相关问题求解的可行性与准确性。

2 物性参数对气泡个数及粒径的影响

2.1 气泡运动中气泡个数和粒径计算方法

本文采用用户自定义方法（UDF）内嵌入模拟软件（fluent）中提取出所需数据。通过此方法可以提取每个气泡的气泡面积、几何中心、当量直径等气泡运动过程中的特性参数。其中气泡直径大小采用索特等效直径，如式（4）所示。

其中，Deq为气泡当量直径（索特等效直径），A为气泡截面面积。

在统计管道内的气泡运动时，本文将只占有4个以下的网格数目的特殊气泡、管道内靠近上壁面2 mm内的气泡和距离管道进气孔0.025 m内的气泡（即圆柱射流）剔除。本文的瞬时相分布图，黑色部分是气相，白色部分是液相，时间上会选取气泡平均脱离粒径、频率、位置保持稳定的时间点。此时相分布图可以比较准确地体现主流区的气泡分布特性，具有一定的代表性，更易在对比中体现它的一些规律。

2.2 液相密度对气泡运动的影响

本节模拟计算主要研究液相密度的改变对管道内气泡个数和气泡直径的影响。图3是进口气液两相速度分别为3 m/s、0.8 m/s，在1.08 s时管道内的瞬时相分布图。其中气相密度为1.225 kg/m3，表面张力为0.072 N/m。从图3中可以看出，不同的液相密度下都表现出相似的特征，长勺状气泡从射流上脱离，然后几个小气泡迅速从气泡上脱落，大气泡旋转、变形、上浮到上壁面。同时随着液相密度的增加，管道内破碎的小气泡数量增加。提取该时刻下气泡从进气口到靠近上壁面之间的大气泡数据即主流区直径大于0.03 m的气泡，如图4所示，不同液相密度下的大气泡个数都是3个，计算出大气泡的平均的直径，分别是0.036 4 m、0.036 9 m、0.036 6 m、0.038 m，考虑第三张图含有刚脱落且将要分裂的“不稳定双气泡”，脱离后的大气泡的平均直径基本相同，都是0.036 6 m左右，液相密度对脱离后的大气泡的直径影响不明显[9]。

图3 1.08 s时刻气泡运动相分布图

图4 1.08 s时不同液相密度下气泡直径散点图

从前面的相分布图3可以看出，不同条件的管道内总气泡个数相差较大，图5所示为不同液相密度在0～1.2 s时间间隔内对管道内气泡个数的影响。从图3可以明显看出，在0～1.2 s间隔内，随时间的增加，气泡个数不断增加；分别经过一定时间后，气泡个数都开始呈周期性地线性增长，气泡个数的平均增长频率基本都是0.02 hz。这一线性增长的现象会一直持续直到在气泡开始流出管道，接下来气泡数量将维持在一个常数附近小范围波动，液相密度越小，周期性波动大。同时也可以看到气泡个数前期的增长率都较大，分析是前期小气泡不够密集，不容易聚合形成更大的气泡，这也是气泡个数的增长率接下来趋于相同的原因，所以当气泡密集程度上升时，气泡个数的增长率出现降低。而液相密度越大，气泡个数的增长越早降低，然后趋于稳定。对气泡最终的数量统计时，进行时均化，液相密度增加2 500 kg/m3内，每增加500 kg/m3，小气泡的数量增加50个，总共增加了80%，基本呈线性增加。

图5 不同液相密度下气泡个数与时间的关系

2.3 液相粘度对气泡的影响

本文设置了液相流体粘性力分别为μ＝0.000 5 Pa·s、μ＝0.001 Pa·s、μ＝0.005 Pa·s时的管道中气泡运动情况。其中μ＝0.001 Pa·s为标准情况下，水的粘度（粘性系数）。图6所示为在0.94 s时刻，进口气液两相速度分别为4 m/s、0.8 m/s时以及进口气液两相速度分别为6 m/s、1 m/s时管道内瞬时相分布情况。

图6 0.94 s时刻气泡运动相分布图

采用控制变量法，仅更改液相粘度，保持其他的物性参数不变，进行数值模拟计算。其中液相密度为998.2 kg/m3，气相密度为1.225 kg/m3，表面张力为0.072 N/m。通过用户自定义方法计算得0.94 s时管道内所有气泡的直径，相对应的气泡个数和直径的对比图7。

观察图6和图7并综合比较分析，在标准情况的液相粘度下，两个不同的速度条件下，气泡的破碎程度都较高，大气泡的平均粒径相对较小，气泡的数量最多。除此之外，从相分布图上可以明显看到，随着液相粘度增加，气泡脱离时的气相射流的长度减少，气泡脱离周期变小，气泡平均脱离粒径减小，但其趋势都越来越小。分析是液相粘度的增加会加快气泡的脱离，气泡的平均脱离直径减小，但是也会抑制气泡脱离后的破裂，所以最后液相粘性对气泡直径影响的综合作用结果是比较小的[15-16]。而总体上，减小液相粘度时气泡的直径和数量变化更明显；通过统计管道气泡运动稳定后的气泡直径、数量并时均化，减小液相粘度到标准情况一半时，大气泡的平均直径都增大超过了10%，气泡数量都减少了至少20%。

图7 0.94 s时不同液相粘度下气泡个数和直径的对比

2.4 表面张力对气泡个数影响

气泡在运动过程中会受到多方面的因素影响，其中气泡表面所受到的表面张力对气泡运动过程中维持气泡形状具有十分重要的作用。国内外学者[17]经常采用无量纲常数奥托斯数E0（Eotvos number）来描述气泡运动过程中浮力和表面张力相对强弱的关系。E0表示气泡所受浮力和气泡表面张力的比值。其中气泡粒径大小de在多气泡运动过程中无法准确代表管道内多气泡在运动过程中所受浮力，本文基于奥托斯数的计算公式进行简单的近似修正，如式（6）所示。

式中，n表示为多气泡运动过程中的气泡个数，E0'表示为修正后的奥斯托数。

为了更深入地了解表面张力在多气泡运动中对气泡形状变化的影响大小，本文采用控制变量方法，将表面张力分别设置为σ＝0.0072 N/m、σ＝0.072 N/m、σ＝0.72 N/m等三个不同的数值，并进行相关数值模拟。图8为2 s时，进口气液两相进口速度分别为1、0.8 m/s，以及进口气液两相进口速度分别为2、0.6 m/s时气泡运动情况。其中液相密度为998.2 kg/m3，气相密度为1.225 kg/m3，液相粘度为0.001 Pa‧s。

图8 2 s时气泡运动的相分布图

本文选择的气泡直径d设置为多个气泡的平均直径，计算区域为进气管进气开始到第一个气泡开始靠近管道上壁面截止。为了方便对比，表1为气相自进管道后1 s和2 s时气液两相流体的有关物性参数。

表1 考察表面张力影响所用物性参数

表中de1为1 s时气泡平均粒径；de2为2 s时气泡平均粒径；为1 s时的奥托斯数；为2 s时的奥托斯数。结合表1与图8分析，表面张力的改变使气泡脱离时的形态和大小发生改变。随着表面张力的增加，奥托斯数变小，小气泡的个数减少，气泡的平均直径增加，气泡脱离时的射流长度减小；气泡脱离时受到的表面张力变大，气泡的平均脱离直径变小。当表面张力达到0.72 N/m时，气泡脱离由表面张力主导，呈椭球状。为了获得更大的气泡或者减少小气泡的数量，可以考虑加入无机盐来增大表面张力；反之，可以添加一些表面活性剂来减小表面张力。

图9 为不同条件下的管道中气泡个数随时间的变化，由图中可以观察到，气泡个数随时间的推移，当气泡运动达到稳定时，气泡个数呈现出基本线性增长；当气泡流到出口处后，分散在管道中的小气泡伴随着液相流出管道，气泡个数呈现出低频波动的稳定不变的状态。

党的十九大报告中，习总书记在“加快生态文明体制改革，建设美丽中国”章节中，明确提出“推进资源全面节约和循环利用”、“加强固体废弃物和垃圾处置”。党中央国务院在《关于加快推进生态文明建设的意见》中进一步部署了推进固体废弃物综合利用，发展循环经济的任务。

当表面张力相同，进口气相流速和进口液相速度不同时，由图9可以观察到气泡个数有差距，但其总体趋势基本相同。如表面张力为0.072 N/m，当进口气液两相流体速度不同时，气泡个数随时间的变化曲线趋势仍然相同，并且气泡个数稳定时的气泡个数相差较小。而当表面张力不同时，气泡个数却有很大的差距，即表面张力越小，将加剧气泡破碎效果，产生管道内气泡个数越多的效果在水平管道内气泡运动过程中，表面张力和奥斯托数呈负相关；表面张力越小，奥斯托数越大，并且气泡越容易变形。

图9 不同表面张力下气泡个数对比

对表面张力和气泡数量的关系进行统计，随着表面张力从0.007 2 N/m增长到0.72 N/m的过程中，两个工况下，气泡数量不断减小，大气泡的直径不断增加，但都呈现出逐渐趋于平缓的状态。

3 结论

本文以标准状况下的水作为基础液相介质，选取几个不同的气相速度和液相速度，作为典型工况条件，分别改变液相密度、液相粘度、表面张力进行数值模拟，来观察、统计和分析单一物性参数的改变对水平管道内气泡破碎情况、粒径分布、气泡数量以及气泡的脱离周期、直径、数量等的影响。虽然针对水平圆管内液相物性的影响的研究较浅，但是物性参数对气泡的影响在一定范围内仍有很多共同的规律可循，本文的研究结果可为相关的研究提供参考依据。相关物性参数范围：液相密度（998.2～2 500 kg/m3），液相粘度（0.000 5～0.05 Pa·s），表面张力（0.007 2～0.72 N/m）。

1）在气泡运动中，增大液相密度到2 500 kg/m3的过程中，会加剧气泡的破碎，气泡数量不断增加，且基本呈现出线性增长的趋势。同时这个范围内的液相密度的改变对大气泡直径的影响不明显。

2）减小或增大液相粘度时，气泡特性受到的影响是复杂且不明确的，液相粘度对气泡直径和数量的影响会存在极大值。增大液相粘度，会减小气泡脱离的平均直径，也会抑制气泡脱离后的破裂，综合作用下改变液相粘度后对气泡的影响较小。减小液相粘度至0.000 5 Pa·s时，变化较明显，气泡个数减少了20%，大气泡的直径减小了10%；而增大至0.005 Pa·s才达到了相同程度的影响。

3）在水平管道内气泡运动过程中，表面张力和奥斯托数呈负相关；表面张力越小，奥斯托数越大，并且气泡越容易变形，并加剧气泡运动过程中的破碎效果；但是表面张力的增大也会减少气泡的脱离直径，对气泡脱离形态产生影响。表面张力从0.007 2 N/m增长到0.72 N/m的过程中，两个工况下，气泡数量不断减小，大气泡的直径不断增加，但都呈现出逐渐趋于平缓的状态。当表面张力大于0.072 N/m后，其对气泡数量和直径的影响趋于平缓。