气流中液化天然气液滴破碎数值模拟研究

樊玉光，宋光辉，袁淑霞，刘家豪，雷 瑶

（西安石油大学 机械工程学院，陕西 西安 710065）

在液化天然气（LNG）生产过程中，存在LNG 液滴在气相作用下的两相流动问题。在不同冷却级间需进行气液分离，而液滴的大小影响分离效率。液滴直径越大，惯性力越大，越易分离。当液滴破碎时，直径变小，惯性力变小，不易分离。因此LNG 液滴的变形破碎将会对气液分离器的分离效率产生重要影响。许多学者通过数值模拟法研究液滴的变形破碎过程，其中用于处理相界面的有格子玻尔兹曼（LBM）法［1-2］、水平集（Level Set）法［3-4］和流体体积（VOF）法［5-10］等。LBM法相比VOF 法不够成熟，且会受到壁面边界处理问题的限制，涉及到相变流动时，VOF 法捕捉界面的动态行为较准确［11］。Level Set 法质量守恒性较差，VOF 法质量守恒性较好，在追踪LNG 液滴相界面的变形破碎等问题上具有良好的优越性。张文英等［6］采用VOF 法对不同We下的水滴进行了模拟，分析典型的袋状破碎形态。Kékesi 等［7-8］采用VOF 法研究了初始球形液滴在稳态流动和剪切流动下的变形和剪切破碎状态。Minakov 等［9］研究了We范围为7 ～212 的水煤浆液滴的二次破碎过程。梁伟等［12］研究了常温常压下We范围为0.5～200的海水液滴碰撞过程。廖达雄等［13］的研究表明高速气流下，气流速度决定氮液滴发生碰撞时的破碎形态。LNG 生产过程中，关于LNG 液滴变形破碎的研究较少，且LNG 液滴的物性与水滴、油滴各有不同，如LNG 液滴的密度和表面张力均小于水滴的密度和表面张力，相界面的变化及破碎形态可能会有所不同。

本工作研究了LNG 生产过程中液滴在气流作用下变形破碎的过程及影响其形态的因素。采用VOF 法与湍流模型realizablek-ε相结合，对气相作用下LNG 液滴变形破碎的过程进行数值模拟，研究液滴的初始直径、气液初始相对速度及We对液滴形态变化及过程的影响，分析不同We下液滴的形态特征及破碎模式，进一步研究液滴在不同流场中的变形破碎过程，探索LNG 生产过程中液滴发生剪切破碎的临界We。

1 数值计算模型及验证

1.1 模型假设

在VOF 模型求解过程中，追踪LNG 液滴的相界面是通过求解相连续方程来完成的，通过求出体积分量中急剧变化的点来确定LNG 气液两相分界面的位置。由于计算过程的复杂性，对模型做如下假设：1）LNG 主要成分为甲烷，采用甲烷工质进行模拟，连续相和离散相分别为甲烷的气相和液相；2）该流场模型为二维流场，液滴的初始形状为圆形；3）气液两相均为等温且为不可压缩流体。气液相界面的表面张力模型采用连续表面张力模型。考虑流场中的湍流，选取基于RANS 的realizablek-ε湍流模型，液滴破碎过程中的主要无量纲参数有We，Oh，Re等。

1.2 计算模型和求解设置

模型计算区域见图1，采用笛卡尔坐标系下的二维结构化网格进行划分。气相以一定初速度从左边垂直入口进入，给予液滴一定的初速度；左边界为速度入口，右边界为压力出口，上下边界为wall。选择基于压力的求解器，压力速度耦合方程选用PISO 算法，压力计算采用PRESTO 格式。

图1 计算区域Fig.1 Computational domain.

1.3 网格独立性验证

为了确保计算结果的准确性，选择四种不同的网格分辨率（D0/h，D0为液滴的初始直径，h为全局网格尺寸下的宽度，D0/h=8，16，21，24）对VOF 模型进行独立性验证。图2a 为液滴在不同时刻、不同D0/h下的形态变化；图2b 为液滴中心坐标点的压力随时间的变化曲线。由图2 可知，在D0/h≥16 时，液滴的形态及其中心处的速度均无明显改变。因此，通过多个不同D0/h的计算及液滴形态对比，选取D0/h=16 的网格来进行计算。

1.4 模拟有效性验证

采用前述数值模拟方法对文献［14］中的数据进行数值模拟，并与已有实验结果及文献［14］中数据的数值模拟结果进行对比，结果见图3。数值模拟及实验中采用的液滴数据均为：液滴直径1.2 mm，We=78，Re=2 850。由图3 可知，模拟结果与实验结果吻合良好，从而证明了该数值模拟方法的有效性。

2 数值计算结果分析

在不同初始直径，不同的初始相对速度下，对LNG 液滴变形和破碎过程进行数值模拟。天然气在-161.5 ～-82.6 ℃下为气液共存状态，生产工艺中压力范围为0 ～5 MPa。以-100 ℃、0.1 MPa状态下液滴为例展开研究［15］，连续相密度（ρg）和液相密度（ρl）分别为1.141 5，302.5 kg/m3，黏度分别为6.77×10-6，4.2×10-5kg/（m·s），表面张力系数为2.278 mN/m。根据已有研究及气液分离器设计规范［16］，选取液滴直径分别为50，100，200，400 μm；为出现多种不同的破碎模式，初始的气液相对速度（Ur）为20，30，40 m/s。液滴在低Oh（Oh＜0.1）时，不同We下的破碎模式见表1［17-18］。为覆盖所需的破碎模式，选择对We范围为0 ～360 内的液滴形态和特征进行分析。图4 为不同粒径下LNG 液滴的破碎模式。由图4 可知，各种破碎模式与表1 较为吻合。

图2 不同D0/h 液滴形态变化（a）及液滴中心处坐标速度曲线（b）Fig.2 Droplet morphology changes with different grid resolutions(D0/h)(a) and velocity curve at droplet center coordinate(b).

图3 数值模拟结果及实验结果对比Fig.3 Comparison of numerical simulation results and experimental results.

表1 液滴在不同We 下的破碎机制（Oh＜0.1）［17-18］Table 1 Breakup mechanisms of droplet at different We(Oh＜0.1)［17-18］

图4 液滴破碎模式分布Fig.4 Droplet breakup pattern distribution.

2.1 不同We 的破碎机制

We和Oh是描述液滴破碎过程的重要参数，We为气体惯性力和液滴表面张力的比值，惯性力作用使液滴破碎，而表面张力作用使液滴维持自身形态。气体惯性力和液滴表面张力的相对关系决定液滴是否发生破碎及其破碎模式。当Oh＜0.1 时，随着We的增大，液滴会表现出不同的破碎模式；随着Ur增加，液滴的临界破碎时间也会变化。

2.1.1D0对液滴形态的影响

由数值模拟结果可知，Ur一定时，随着液滴直径的增大We增大，液滴破碎模式会发生不同的转变。以Ur=30 m/s 为例，对D0分别为100，200，400 μm 的液滴进行分析。图5 给出了三种不同D0及Wer（气液相对速度的韦伯数）在Ur=30 m/s 和D0=400 μm 不同Ur时部分时刻的液滴形态。由图5a 可知，200，400 μm 的液滴破碎为剪切破碎，100 μm 液滴为多模式破碎。液滴直径越小，变形程度越大。同种破碎模式下，直径较小的液滴变形时间会缩短。液滴在受到气相作用时，形变主要受气体惯性力、表面张力和阻力的影响。液滴直径越小所受阻力越小，反之阻力越大。

图5 不同D0（Ur=30 m/s）（a）及不同Ur（D0=400 μm）（b）的液滴形态Fig.5 Droplet morphologies with different D0(Ur=30 m/s)(a) and different Ur(D0=400 μm)(b).

为研究液滴运动特性，对其时间特性进行分析。液滴的变形破碎时间分为变形时间和总破碎时间。Pilch 和Erdman 对其进行了定义，并拟合了不同范围的总破碎时间和We的相关性［18］。变形时间定义为当一个完整的液滴不再存在，即从初始阶段到液滴首次破碎的时间。总破碎时间为液滴及其部分不再进一步破碎的时间。图6 为液滴的无量纲变形时间（T）与Wer关系曲线。其中t为液滴临界破碎时间。由图6 可知，同一Wer下，T随D0的增大而减小，D0较小的液滴T区间较广。随着Wer的增大，液滴的T先减小后趋向于平缓。

2.1.2Ur对液滴形态的影响

由图5b 可知，D0一定时，随着初始Ur的增加Wer增大。同一时刻液滴的变形程度逐渐增大，变形破碎时间也会缩短，破碎模式也发生了转变。随着Ur的增大，液滴所受的气体惯性力增大，变形破碎程度也会发生变化。为得到剪切破碎的临界条件，只需增大Ur即可。在D0=400 μm 时，分别增大Ur进行数值模拟，根据模拟结果可知，液滴发生剪切破碎的临界We为80 左右。

图6 T 随Wer 的变化曲线Fig.6 Dimensionless deformation time(T) change curves along with Wer.

液滴在变形破碎过程中，随着We的增大，阻力也会发生变化［18-19］。引入阻力系数（CD），CD=4D0ρla/3ρgUr，a为液滴加速度。图7a 为液滴速度和加速度随时间的变化曲线。由图7a 可知，液滴速度随着时间的延长一直增大，但增加幅度有所不同，液滴加速度随着时间的延长先增大后减小。液滴处于变形阶段时，液滴速度缓慢增加；在0.15 ms 时，液滴发生了剪切破碎，随着子液滴的增多，母液滴的质量随之减小，液滴的速度快速增大；在0.36 ms 液滴加速度达到最大值，此时液滴速度变化值最大。随后液滴加速度下降，液滴速度增量减小，液滴速度缓慢增加。图7b 为液滴阻力系数随T的变化曲线，并给出了变形破碎过程的部分液滴形态。由图7b 可知，液滴的阻力系数随着T的增加呈现先增大后减小的趋势。T约为0.7 时，液滴上下方刚好剥离出子液滴，发生剪切破碎，在这之前液滴处于变形阶段，与液滴速度缓慢增加阶段相对应。随着剥离出子液滴的增多，液滴的阻力系数快速增大。T约为1.9时，液滴阻力系数达到最大值，随后液滴阻力系数减小。

2.2 不同流场下液滴的破碎情况分析

为研究湍流涡对液滴破碎的影响，进一步探索了带挡板的流场下液滴破碎情况。折流板分离元件是气液分离中常用的结构，而挡板结构产生的湍流进一步加快了液滴的破碎。在流场中，挡板的边界类型设置为wall，其余求解设置和边界条件同图1保持一致。图8 为液滴在不同时刻的体积分数云图和速度矢量图。

图7 液滴速度和加速度随时间变化曲线（a）和液滴阻力系数随T 变化曲线（b）Fig.7 Droplet velocity and acceleration change curves with time(a) and drag coefficient change curves with T(b).

图8 液滴在不同时刻的体积分数和速度矢量云图Fig.8 Contour of volume fraction and velocity vector of droplet at different time.

由图8 可知，在挡板的作用下，气流的流向发生了改变，流场中产生不同的涡团。进入流场的气流已不再均匀，两侧的气流先向中间流去，相邻挡板间及液滴背风面形成了涡团。作用在液滴上的气体惯性力增大，气体惯性力的增大加速了液滴的变形破碎过程，液滴的形态变化有很大差异。同一We下，图5 中0.05 ms 时刻所对应的流场中液滴由圆变扁，而图8 流场中的液滴发生了剪切变形，液滴的上下两端即将剥离。图5 中0.10 ms 时液滴仍在变形中，直至0.15 ms 液滴发生了剪切剥离，但图8 液滴在0.10 ms 时剪切剥离程度已经显著，剥离出的子液滴也有跟随气流流向漩涡的趋势。由模拟结果可知，在0.06 ms 液滴就已经发生了破碎；0.15 ms 时液滴剥离出更多子液滴，部分子液滴由于漩涡的作用已经附着在挡板壁面上；0.23 ms 时液滴又剥离出了较大的子液滴，且之前附着在挡板上的子液滴仍受到漩涡影响运动。由此可知，同样Wer下，液滴在有挡板流场的破碎时间比无挡板流场的破碎时间更短，且变形破碎程度更显著。

3 结论

1）液滴的D0和初始Ur是影响LNG 液滴变形及破碎的重要因素。随着We的增大，液滴呈现出不同的破碎模式。

2）同一Wer下，液滴的T随着D0的增大而减小。随着Wer的增大，液滴的T先减小后趋于平缓。D0=400 μm 时，液滴发生剪切破碎的临界We在80左右。

3）LNG 液滴速度随着时间的延长一直增大，但增加幅度有所不同，加速度随着时间的延长先增后减小，在0.36 ms 液滴加速度达到最大值。阻力系数随T的增加呈现先增大后减小的趋势。T约为1.9 时，液滴阻力系数达到最大值。这是由于液滴经历了从变形到首次破碎再到二次破碎的过程。

4）分析了带挡板流场下的LNG 液滴破碎情况，由于挡板的作用，作用在液滴上的惯性力增大，液滴破碎时间更短，变形破碎程度更显著。