【论文】邓佳佳等：不同雷诺数下的蒸发甲烷液滴阻力特性

0　引言

LNG  是高效、清洁的能源，需求量逐年上升 ^[1-3] 。运输或存储天然气最为安全、经济的方式是将其深冷液化成 LNG 。大型 LNG  储罐投产前的准备工作中，   约有一半的低温液舱事故发生在液舱的预冷过程 ^[4-5] 。液舱预冷是将制冷剂雾化后喷入液舱，利用制冷剂液滴相变制冷 ^[6] 。中小型液舱多用液氮预冷，大型液舱由于预冷所需冷量较大，一般采用 LNG  预冷 ^[7-8] 。为保证液舱及设备安全，预冷前需用制冷剂同种蒸气将舱内的气体置换。

通过分析液滴尾部流场及圆球表面的压力变化并选取阻力系数（ C _d ）作为比较依据，可以更加详细地了解流动情况 ^[9] 。故液滴在气体中的运动需要先了解液滴的阻力特性。液滴阻力主要来自于液滴前缘区和尾部区之间的压力差，由于这种前后压差的变化，会使 LNG  液滴发生变形破碎 ^[10-11]   及流动分离现象。郑德乾等 ^[12]   研究得出流动分离点位置会随雷诺数（ Re ）的变化发生改变，导致屋盖表面风压分布受雷诺数效应影响显著。 Achenbach ^[13]   研究得出了边界层从层流向湍流过渡和边界层分离位置的信息及雷诺数对局部表面摩擦和静压分布的影响。吴以坚等 ^[14]   对小球在圆柱管道的下落过程进行分析，得出小球受到的黏滞阻力是小球半径与圆柱管道半径的比例函数，黏滞阻力、压强力、黏性力都随着球心到轴线的距离先减小后增大，且具有不同的极小值点。学者还对特殊环境下球体的阻力进行了理论研究。 Xiong  等 ^[15]   对超临界水（ Supercritical Water ， SCW ）   中随温度和压力变化而变化的球体的 C _d  和努塞尔数（ Nu ）进行研究，发现其不能很好地拟合传统流体的关联式。 Jin  等 ^[16]   进一步研究 Re ，发现当 10  ＜ Re ≤ 200  时固体球的 C _d  和 Nu 的相关性还与强度有关。 Wu 等 ^[17]   研究结果表明颗粒表面与流入颗粒之间的焓差较大时， C _d  显著增大，且主要由摩擦阻力系数引起。

对于阻力特性的研究，学者们还做了大量的实验研究。张亮等 ^[18]   利用数值模拟与实验相结合的方法对波纹管中流体的传热和阻力特性进行研究，得出流体处于湍流状态时，传热强化幅度随 Re  的增大而降低；平均摩擦系数随着 Re  的增大而减小。麦浚宇 ^[19]   对单个细小颗粒在静止液体中的沉降状态进行实验研究得出， 0.1  ＜   Re  ＜   10  时， C _d  下降幅度大， 10  ＜   Re  ＜   140  时， C _d  下降趋于平缓。并对 C _d — Re  的实验数据进行了拟合。

液滴蒸发也会对阻力特性产生一定的影响。马力等 ^[20]   设计液滴蒸发实验，得出随着液滴直径的增大，内部环流加强，液滴尾涡变大，传质阻力减小，   液滴的蒸发常数变大。 Gutierrez  等 ^[21]   研究发现蒸发液滴表面的温度跃变主要与界面热传导有关而非蒸发，气液界面附近的蒸汽密度是由界面处的传质阻力，即蒸发阻力确定。刘财喜 ^[22]   通过分析不同粒径液滴对阻力和传热效率的影响，得出阻力随液滴质量浓度的增大而增大；阻力和传热效率随着颗粒动量增大而减小。 Wang  等 ^[23]   研究了不同相对方向和粒子距离的相同球形粒子之间的流体与传热特性，得出对于所有构型粒子表面的速度和温度边界层厚度，   进一步表现为 C _d  及 Re  的减小。对于干燥塔内脱硫废水的处理邓佳佳和金宏伟学者通过数值模拟塔内液滴运动与蒸发，发现喷雾液滴速度越大，对烟气的流动阻力越大，烟气下降速度就越快，使烟气流动最终趋于均匀化。且当废水处理量增大时，干燥塔阻力升高，液滴蒸发所需的时间和距离增加 ^[24-25] 。

综上所述，对于蒸发液滴的阻力特性的研究仍然存在一些盲点，特别是液滴蒸发对液滴阻力特性的影响研究仍不够充分。因此以液舱预冷的液滴蒸发过程为研究对象，基于液滴界面能量守恒，建立 LNG  饱和液滴在其蒸气中强制对流蒸发模型，从压力、蒸气蒸发速度、分离点研究蒸发液滴的阻力特性，   并通过蒸发液滴与低温球（不蒸发液滴）的阻力特性对比来分析蒸发对液滴阻力特性的影响，并建立总阻力系数（ C _dt ）与 Re  的无量纲关联式，为 LNG  的实际运输投产提供强有力的数据支撑，减少 LNG  运输前期因实施预冷作业时液滴与周围蒸气巨大温差而可能引起爆炸、冰堵等潜在风险的发生，对实际工程应用提供一定的理论指导价值，以期提高 LNG  运输前的液舱预冷的安全性和高效性。

1 　模型建立

LNG  液舱预冷过程中，喷淋的 LNG  液滴是饱和液滴，因此该液滴吸收的热量完全用于液相蒸发。 LNG  单液滴在其蒸气中蒸发过程较为复杂，为了简化模型，做出如下假设 ^[26-27]   ：

1 ）蒸发过程中液滴保持球形，且气相流场结构为轴对称，采用二维轴对称坐标系。

2 ） LNG  组分较为复杂，但是组分含量最大的为甲烷，因此采用甲烷液滴和气相甲烷作为模型工质。

3 ）气液界面及液滴内部温度均为饱和温度，气液界面相变符合准稳态假设且相变只发生在气液界面处。

4 ）因为 LNG  液滴为饱和液滴且处于稳态蒸发阶段 ^[28] ，液滴内部运动对传热传质的影响很小，暂时可忽略液滴内部的运动带来的影响。

5 ）忽略辐射热带来的影响。

基于以上假设，采用稳态二维轴对称坐标系，   其模型方程包含连续性方程、动量方程、能量方程和湍流模型 ^[29] 。

连续性方程为：

动量守恒方程为：

能量守恒方程为：

由于液滴蒸发时液滴界面边界层的影响非常大，   为准确模拟边界层，采用 SST k - ω  湍流模型 ^[30] 。该模型采用近壁区间与远场区间分开处理的方式，能够准确模拟近壁区间的层流边界层内的流动。

其中湍动能（ k ）方程为：

湍动耗散率（ ε ）方程为：

气液界面处气相传递给液相的热量等于液滴相变所需要的相变潜热，即气液界面处的边界条件表达式为：

选择理想气体模型来计算气相介质密度 ^[31]   ：

以 LNG  液舱喷雾预冷中单个 LNG  液滴在强制对流条件下的换热蒸发建立二维轴对称物理模型。采用 Gambit  进行几何建模和网格划分，整个网格化模型以 x  轴为对称轴，液滴和低温球的几何中心位于坐标轴原点处。 x  轴正方向施加重力加速度为 9.81 m/s ² 。液滴气液界面温度为饱和温度 110 K ，压力设置为 0.1 MPa 。

具体物理模型及网格如图 1  所示。图 1  中黄色椭圆区域为 LNG  储罐，罐内为喷管喷出的 LNG  液滴示意图，左下方为喷嘴处局部放大图，中间为喷出的单液滴局部放大图，右下方为单液滴网格划分图，其中 v ₁  表示来流气体速度， m/s  ；   φ  表示来流夹角，（ ° ）；液滴粒径为 d ， mm 。计算域上边界设为速度进口，下边界设为压力出口，右边界设为壁面，   液滴表面设为速度进口，恒温低温球表面设为壁面，   x  轴为对称边界。液滴和低温球表面温度设为饱和温度 110 K 。选择基于压力的 SIMPLE  算法，采用二阶迎风格式离散。根据运动液滴蒸发的流场特性不同，   分别采用不同的网格。表 1  为所用介质的物理参数表 ^[29]   ；计算工况统计如表 2  所示，主要用来分析不同液滴粒径对 Re  的影响。

2 　模型验证

2.1 　网格无关性验证

计算过程中的网格数量对模拟结果的准确性有着较大的影响。因此划分了 5  组不同数量的网格，即网格 1 、 2 、 3 、 4 、 5 ，网格数量分别为 9 877 、  21 753 、 41 828 、 83 898  和 160 152 。图 2  为液滴在不同网格数量下其表面压力（ p ）随来流夹角（ φ ）的变化曲线图，由图 2  可知，网格 1 、 2 计算结果与网格 3 、 4 、 5  计算结果存在着一些差距，为保证计算精度，   因此选择网格 4 用于后续模拟计算。

用幂函数的形式对液滴的 C _dt  进行拟合，并就麦浚宇 ^[19]   的研究依据 R _e  区分的液滴 C _dt  曲线与岳湘安 ^[32]   研究文章中 Kurten  和 Rumpt  学者的计算模型进行对比，进而对本模型的准确性进行验证，验证结果如图 3  所示，其中实心球代表数值模拟数据，实线代表拟合曲线。由图 3  可得，液滴的数值模拟结果与参考模型的结果基本一致，拟合得到的液滴优度系数（ R ² ）为 0.944 ，且最小偏差为 0.6% ，说明模型计算结果与参考模型结果吻合良好。

3 　结果与讨论

低温球与液滴相比，低温球与周围环境对流换热时，没有蒸气蒸发，不存在 “ 吹拂效应 ” 。首先分析不同 Re  条件下，液滴的压力特性 ，进而分析流场特性即表面蒸气蒸发速度以及分离点；再通过对比分析得到液滴分离对阻力 的影响，进而分析蒸发液滴的阻力特性。

对液滴蒸发及低温球压力特性及换热过程周围流场特性进行对比，分析其表面蒸气蒸发速度及分离特性。

图 4  表示了不同工况下液滴表面压力（ p ）随来流夹角（ φ ）的变化情况。表 3  表示低温球表面压力变化情况。由图 4  可知，当 d  相同时， Re  越大，   p  随 φ  的增大变化越大，但总体压力都呈现先快速下降后缓慢上升趋势。由表 3  可知，当 d  相同、 Re  不同时， φ =0° 时压力取得最大值 p _max ，且随着 Re  的增大，压力到达最小值 p _min  时低温球的来流夹角越来越小，差值 Δ p ₁ （ p _max  与 p _min  的差值）也越来   越大。

结合图 4  可知， d  为 1.0 mm 、 Re  较小时，液滴的表面压力总体变化趋势不大，都随 φ  的增大呈现先下降后缓慢上升的趋势； d =0.1 mm 、 Re =224.3  时，   随 φ  的增大， Δ p ₁  最大为 60.3 Pa 。由图 4-b  可知，   0.1 mm  液滴的 Δ p ₁  比 1.0 mm  液滴的 Δ p ₁  大，这是由于 0.1 mm  的液滴速度（ v ₁ ）较大，所以液滴表面受到的 p  越大。

液滴在强制对流条件下与周围蒸气进行对流换热过程中，液滴表面有蒸气蒸发，产生垂直于液滴表面的蒸气蒸发速度（ v ₂ ），蒸气蒸发产生的吹拂效应对液滴的流场产生变化，因此需要分析液滴表面蒸气蒸发速度的变化规律。

液滴表面的蒸气蒸发速度随来流夹角的变化曲线如图 5  所示。可以看出， d  相同时， v ₂  的变化趋势基本一致，大体呈现两段式， v ₂  随着 φ  的增大先减小后增大，且总体随着 Re  的增大而增大。由图 5  可知，在驻点 φ =0° 处， v ₂  最大，再随 φ  的增大而下降直至达到最小值 v _2min  ；随后再随 φ  的进一步增大逐渐增大，但在 φ =180.0° 时的 v ₂  小于 φ =0° 时的 v ₂ 。 d  相同、 Re  不同时，蒸气蒸发速度达到最小值 v _2min  时的来流夹角 φ _{v 2min}   不同。从图 5  可知， Re =22.4  时， v ₂  随着 φ  的增大呈单调下降趋势，并未出现随 φ  的进一步增大逐渐增大的趋势。

图 6  为不同 Re  下液滴表面蒸气蒸发最小速度值（ v _2min ）及其对应的来流夹角（ φ _{v 2min} ）的曲线图，其中上下角标 1.0  和 0.1  分别代表 d  为 1.0 mm  和 0.1 mm  的液滴。图 7  是 d =1.0 mm 、 Re  分别为 22.4  和 224.3  时液滴和低温球的矢量图，其中图 7-a 、 b  为不同 Re  下液滴表面蒸发速度，图 7-c 、 d  为不同 Re  下液滴外部流场矢量图，图 7-e 、 f  为不同 Re  下低温球外部流场矢量图。由图 6  可知， d  相同时， v _2min 随 Re  的增大总体呈上升趋势； Re  相同时， d =0.1 mm  液滴的 v _2min 比 d =1.0 mm  液滴的 v _2min 大。但不管 d  是否相同， v _2min 对应的 φ _{v 2min}  随 Re  的增大呈下降趋势，   且 Re  ≤ 44.9  时的 φ _{v 2min} 几乎为 180.0 ° ， 而 Re  相同时， d =1.0 mm  和 d =0.1 mm  的 φ _{v 2min}  几乎一致。进一步结合图 7  流场矢量图可知，当 d =1.0 mm ， Re =22.4  时，液滴表面几乎没有发生边界层分离，且 φ _{v 2min}  为 180.0° ，而当 Re =224.3  时，出现了较为明显的边界层分离， φ _{v 2min} =134.6° 明显变小，分析认为液滴表面 v _2min 出现角度可能与液滴表面的流场分离结构有关。

图7 　粒径为 1.0 mm  的液滴表面蒸发速度和液滴、低温球外部流场矢量图

注：液滴表面蒸发速度矢量表示速度分布特征，仅代表速度相对大小和方向。

图 8  是液滴与低温球在不同 Re  下对应的分离角（ φ _s ）即分离点位置 ^[33]   的来流夹角、 φ _{p min}   以及 φ _{v 2min}   的变化图，其中分离点位置为液滴、低温球表面切向速度梯度为 0  对应的点，分别记为 φ _{s - d}  和 φ _{s - b} ，液滴、低温球达到 p _min  时对应的来流夹角分别记为 φ _{p min-d} 、 φ _{p min-b} 。图 9  是不同 d 、 Re  的工况，液滴与部分低温球表面速度梯度（ Lv ）随 φ  变化图。由图 8  可知，   φ _{v 2min} 、 φ _s 、 φ _{p min}   都随 Re  的增大而减小，   即 Re  越大分离点越靠前。当 Re  ≤ 22.4  时， φ _{v 2min} 、 φ _{s - d} 、 φ _{s - b}  均   为 180.0° ； 22.4  ＜   Re  ≤ 44.9  时， φ _{v 2min} 、 φ _{s - b}  为 180.0° ，   φ _{s - d}  ＜ 180.0° ；当 Re ＞ 44.9 时， φ _{s - d}  ＜ φ _{s - b}  ＜ φ _{v 2min}   ＜ 180.0° 。由于液滴表面蒸发的原因，使得 φ _{s - b}  与 φ _{s - b}  产生了差异。结合图 7  以 d =1.0 mm  的液滴和低温球为例得出，   当 Re =22.4  时， φ _s =180.0° ；当 Re =224.3  时， φ _s  随 Re  的增大而减小，且液滴表面有蒸气喷出，流场发生分离，使得 φ _{s - d}  ＜   φ _{s - b} ，即液滴表面的流场分离提前，   但是分离角位置与蒸发最低速度对应位置并不一致，   而是存在一个角度差。

由图 8  可知，当 Re  相同时， φ _{p min}   ＜   φ _s ，结合图 9  中表面 Lv  随 φ  先增大后减小以及图 4  可知，压力在达到 p _min  之前，压力逐渐下降，压力梯度（ Gp ）小于 0 ，即正压区，流体在压力作用下速度增加。达到 p _min  之后，压力逐渐上升， Gp  ＞   0 ，即逆压区，在压力作用下流速降低，最终导致分离，所以   φ _{p min}   ＜   φ _s 。由于液滴表面蒸发的原因，使 φ _{p min-d} 、 φ _{p min-b}   产生了差异，大小具体来说分成 3  段。在 Re  ≤ 224.3  时，   φ _{p min-d}   ＞   φ _{p min-b} 。

绕流物体的阻力由两部分组成： ① 由于流体的黏性引起的切向摩擦阻力； ② 由于边界层分离，物体前后形成压强差而产生的压差阻力。阻力系数定义如下所示 ^[33]   ：

液滴表面压差阻力系数随来流夹角变化如图 11  所示。 C _d  具体变化情况参见表 4 ，其中增长率（ ε ₀ ）   的计算方法为 ^[31]   ：

由图 11 可知，当 Re  相同， d  不同时， C _d  基本一致，   只有在 φ 较小时有略微的差别。且 C _d  都呈现随 φ  的增大，先快速下降直至达到最小值 C _dmin  随后再缓慢上升，且同一工况下，在驻点 φ =0° 处的 C _d  最大。当 Re 不同， d  相同时， C _d  随 Re  的增大而减小，到达某一点之后，随 Re  的增大而增大。结合图 4  和图 9 可知，   正压区， C _d  快速下降，但随着 φ  的增大，到达逆压区，   使 C _d  略有上升。在 C _d  达到 C _dmin  之前 Re 越大， C _d  越小，即 C _d  与 Re  成反比；在 C _d  达到 C _dmin  之后 Re 越大， C _d  越大，即 C _d  与 Re  成正比。当 Re =22.4  较小时，液滴表面受到的压力较小，所以 C _d  较大，由于没有发生边界层分离，只有摩擦阻力，所以达到最小值的来流夹角越大，即 Re  越小，分离角 φ _s  越大。当 Re =224.3  时，边界层发生分离，边界层内的液滴还受到压强差的作用，结合表 4  可知， Re 越大， Δ p ₁  越大，使液滴表面压差阻力系数的 ε ₀  越小。

图 12  是不同 Re  相同 d  的液滴、低温球压差阻力系数对比图，分别表示为 C _d-d 、 C _d-b 。由图 12  可知，   液滴与低温球的 C _d  都随着 φ  的增大呈现先快速下降后缓慢上升的趋势。由于液滴表面有蒸气蒸发，存在压差阻力，使边界层变厚，受到液滴吹拂效应的影响，   使相同工况下 C _d-d  ＞   C _d-b 。

采用幂函数的形式对液滴和低温球的 C _dt  进行拟合如图 13  所示分别记为 C _dt-d 、 C _dt-b ，拟合得到的优度系数（ R ² ）分别为 0.944  和 0.994 。由图 13  可知液滴和低温球的 C _dt  都随着 Re  的增加而减小，直至两者几乎重叠，且 C _dt-d  ＜   C _dt-b 。由液滴与低温球的分离特性并结合图 9  可知，相同 Re  下液滴表面速度梯度（ Lv _-d ）小于低温球表面速度梯度（ Lv _-b ），即液滴的摩擦阻力小于低温球的摩擦阻力，使 φ _{s - d}  ＜   φ _{s - b} ，所以同种工况下，液滴比低温球提前分离，使 C _dt-d  ＜   C _dt-b 。

4 　结论

p  随 Re  增大而增大， Re  相同时随 φ  的变化趋势基本一致，均呈现两段式，即 p  随 φ  的增大先减小后增大，存在最小值 p _min ，当 Re =224.3  时，   d =0.1 mm  时， p _min-d  为－   19.4 Pa ， p _min-b  为－   22.0 Pa ，   p _min-d  ＞   p _min-b ，   且 0.1 mm  的 Δ p  比 1.0 mm  的 Δ p  大，   液滴 Δ p ₁  最大为 60.3 Pa ，低温球 Δ p ₁  最大为 62.7 Pa 。 p _min  随 Re  的增大而减小，由于达到 p _min  之后，存在逆压区，在此作用下使流速（ v ₁ ）降低，最终导致分离，分离造成 v ₂  的变化， v ₂  随着 Re  的变化趋势与 p  一致，存在 v _2min ， v _2min  随 Re  的增大而增大，且   v _2min-0.1  ＞   v _2min-1.0 。

φ _s 、 φ _{v 2min} 、 φ _{p min} 随 Re 的增大而减小，但三者略有差异。当 Re ≤22.4 时， φ _{v 2min} 、 φ _s  均为 180.0° ； 22.4 ＜ Re ≤ 44.9  时，   φ _{v 2min} 、 φ _{s - b}  为 180.0° ， φ _{s - d}  ＜   180.0° ；当 Re  ＞   44.9  时，   φ _{p min}   ＜   φ _s  ＜   φ _{v 2min} ，由于分离作用且液滴表面受蒸气蒸发影响，使 φ _{s - d}  ＜   φ _{s - b} ，即 Re  越大，分离点越靠前。

对液滴和低温球在强制对流蒸发下的压差阻力系数进行分析，得出 C _d-d  和 C _d-b  都随 φ  的增大先减小后增大，且同一工况下 C _d-d  ＞   C _d-b ，当 Re =224.3 、 d =1.0 mm  时， C _d-d  随 φ  的增大下降最快，   减小了 141.7% ，随后再缓慢上升。由于绕流物体受到阻力包括压差阻力和黏性阻力，而 Lv _-d  ＜   Lv _-b ，即液滴的黏性阻力小于低温球的黏性阻力，使 C _dt-d  ＜   C _dt-b 。并对其分别建立无量纲关联式为 C _dt-d =9.8 Re ^{－   0.51} 、   C _dt-b =31.4 Re ^{－   0.62} ， R ²   分别为 0.944  和 0.994 。