不同喷射条件下闪蒸喷雾液滴尺寸特性数值分析

柯炳正， 高璞珍， 王博， 李茹， 卢川， 田瑞峰

(1.哈尔滨工程大学 核安全与仿真技术国防重点学科实验在室，黑龙江 哈尔滨 150001；2.中国核动力研究设计院 核反应堆系统设计技术重点实验室，四川 成都 610213)

闪蒸喷雾行为存在于很多工业领域中，如当核电站发生一回路冷却剂丧失事故(loss of coolant accident, LOCA)时，大量高温高压冷却剂从承压管路的破口处以闪蒸喷雾的形式快速向外界喷出，造成安全壳内压力温度急剧升高，安全壳的完整性受到严重威胁[1]。闪蒸喷雾属于非平衡热力学过程，其相变情况十分复杂，很多学者对闪蒸喷雾现象进行了数值模拟研究。Marsh等[2]使用6方程模型模拟喷嘴闪蒸流动，考虑了由成核和相变引起的相间质量、动量以及能量的传递。Janet等[3]通过使用CFX 14.5中的5方程模型研究了闪蒸喷嘴流中各种成核模型的性能。Mimouni等[4]使用NEPTUNE CFD中的6方程模型模拟了闪蒸空化流动，通过使用确定的传热系数确保蒸汽温度接近饱和温度。季璨[5]利用离散相的DPM模型和用户自定义函数对闪蒸罐内的外部闪蒸喷雾现象进行研究。聂永广等[6]结合混合物模型与用户自定义函数，计算出沥青喷雾造粒喷嘴中的戊烷闪蒸情况。在闪蒸喷雾的研究领域中，文献[6-11]对燃料在内燃机中的闪蒸喷雾以及对真空中制冷剂闪蒸关注较多，而在液态水闪蒸喷雾的研究领域中，虽然有学者进行了一定的研究[5]，但对于闪蒸喷雾液滴的外部尺寸特性研究偏少。本文利用计算流体力学软件FLUENT对不同喷射条件下闪蒸喷雾过程进行计算模拟研究，对闪蒸喷雾过程中液滴尺寸的变化特性进行分析，得到不同初始压力，喷口直径和喷射距离时，闪蒸喷雾液滴尺寸的变化规律。

1 闪蒸喷雾数值方法及验证

1.1 闪蒸喷雾数值方法介绍

在本文建立了一个150×150×1 200的计算区域并划分结构化网格，如图1所示。为了使计算更准确，对计算区域进行了网格无关性验证。在计算区域内计算了不同网格数量下，流动空气在200 mm处的截面平均流速，空气的入口流速设为1 m/s。其计算结果如图2所示。

图1 计算区域网格图Fig.1 Grid of calculation domain

图2 网格无关性验证结果Fig.2 Grid independence verification results

从验证结果中可以看出，当网格数量达到150万以上时，网格数量对计算结果的影响很小。本研究中计算区域的整体网格数量为165万。

研究利用离散相模型中的过热喷射器模型模拟闪蒸喷射，连续相入口和出口条件均为压力边界，离散相边界条件为逃逸边界条件。初始喷射环境温度设置为23 ℃，初始环境压力为大气压，初始水蒸气湿度设置为0.5%。离散相考虑了液滴颗粒受到加热冷却及相变的影响，同时考虑了压力变化导致的饱和温度变化的影响，液滴曳力的模型为动力阻力模型，液滴二次破碎模型采用了开尔文亥姆霍兹-瑞利泰勒破碎(Kelvin-Helmholtz and Rayleigh-Taylor hybrid wave, KHRT)模型[12]。液滴温度设为初始计算压力下的饱和温度，液滴尺寸分布为Rosin-Rammler分布。压力速度耦合算法、梯度插值分别采用压力耦合半隐式算法以及基于节点的格林高斯插值方法，其他插值采用二阶迎风格式。

1.2 液滴特性研究

本研究中的液滴尺寸特性以在计算域不同截面处计算得到的索特平均直径(Sauter mean diameter, SMD)进行分析。索特平均直径d32是一种描述喷雾液滴尺寸的参数，在流体力学领域中，该参数被广泛应用于颗粒粒径分布的研究中[13]。其定义是与所研究的颗粒群的总体积与总表面积之比相同的一群同等大小的颗粒群中，单个颗粒的直径大小。 索特平均直径d32为：

(1)

式中：dv和ds分别代表体积定义的直径与表面积定义的颗粒直径；Vp为液滴总体积；Ap为液滴总表面积。

在闪蒸喷射的研究中，选择合适的喷雾液滴破裂模型是成功模拟喷雾过程，得到合理的液滴尺寸分布结果的关键。本文闪蒸喷雾过程属于高速条件下的射流喷雾过程，在高速条件下，射流液滴的破裂是空气动力、表面张力和液滴黏滞力相互作用的结果。这些力的相互作用可以以韦伯数We和奥内佐格数Oh来表示。

韦伯数We以及奥内佐格数Oh的分别为：

(2)

(3)

式中：ρ为液滴密度；u为特征流速；d为特征长度；σ为表面张力；μ为动力粘度。

从式(2)、(3)中可以看出，韦伯数We表示的是空气动力与表面张力的关系，奥内佐格数Oh则表示液滴黏滞力和表面张力的关系。在闪蒸喷射中，由于射流流体速度快，射流的韦伯数We很大，射流破裂形式主要为雾化形式[14]。

目前在射流喷雾雾化模型的研究中，比较常用的有泰勒类比破碎(Taylor analogy breakup, TAB)模型[15]和KHRT模型[12]。

1.2.1 泰勒类比破碎模型

TAB模型将液滴的行为与振荡弹簧质量系统相结合，将液体受到的表面张力、液体阻力、液体黏性力分别与振荡弹簧质量系统中弹簧回复力、外力以及弹簧阻尼力进行类比。液滴在气动力的作用下产生形变，当形变达到临界时，液滴发生破裂产生子液滴。液滴形变为：

(4)

式中：y表示液滴发生形变的无量纲量；ρl和ρg是离散相和连续相密度；u是液滴的相对速度，r是未产生形变的液滴半径；σ是液滴表面张力；μl是液滴粘度；CF、CK和Cd为无量纲常数。

当y≤1时，液滴不断振荡变形；当y>1时，液滴发生破裂，产生子液滴。子液滴的速度为：

(5)

式中CV为常数。

1.2.2 开尔文亥姆霍兹-瑞利泰勒破碎模型

KHRT模型将开尔文-亥姆霍兹(Kelvin-Helmholtz)不稳定波理论与瑞利泰勒(Rayleigh-Taylor)不稳定波理论相结合，共同预测射流中液滴的破碎情况。该模型将射流破碎分为一次破碎和二次破碎这2个过程，一次破碎过程中，将射流整体当作是直径与喷口直径相当的大液滴，利用KH不稳定波理论对该液滴的破裂进行模拟计算。

破裂产生的子液滴半径为：

rc=B0ΛKH

(6)

式中ΛKH是对应于具有最大生长速率的KH波的波长；B0是等于0.61的常数。

在二次破碎过程中，模型将RT模型与KH模型结合使用来预测液滴的破碎。

RT模型中增长最快的波的频率与波数为：

(7)

(8)

式中：gt是行进方向上的加速度；ρl和ρg分别为液相和气相的密度。

由RT模型计算得到的新的子液滴的半径为：

(9)

式中CRT为常数，CRT=1。

在KHRT模型中引入了液体核心破碎长度的概念，当液滴破裂发生在液体核心破碎长度内，破裂只考虑KH模型的影响。当破裂发生在液体核心破碎长度外，液滴的破裂则考虑KH模型RT模型的共同作用。

TAB模型虽然在喷雾雷诺数较低时计算结果较为理想，而当雷诺数较大时，计算结果偏差较大。而KHRT模型在雷诺数较大的喷雾情况中计算结果较为理想[16]。在本次研究中，由于液体的温度和压力高于大气压，当液体喷放至大气环境中时，液体与环境的相对流动速度大，雷诺数很大。因此TAB模型并不适用，最终本研究的液滴破碎模型选择的是KHRT模型。

1.3 数值计算方法验证

为了验证本研究中采用的数值计算方法是否合理，利用文献[17]中测量得到的液滴索特平均直径的实验数据和数值计算的结果进行了对比分析。

图3展示了喷嘴直径d为0.75 mm和1mm时，液滴的索特平均直径随喷雾距离变化的实验数据与计算结果对比情况。

图3 实验与计算结果对比Fig.3 Comparison of experiment and calculation results

通过对实验与计算结果的分析可知，在250 mm到750 mm的喷雾距离的范围内，在2种工况下计算结果与实验数据在数值上存在一定偏差，在d为0.75 mm的条件时计算结果与实验结果的最大偏差为13.98%；d为1 mm时计算结果与实验数据的最大偏差为24.90%。

虽然计算结果与实验数据依然有一定数值上的差距，但是2者之间偏差在合理区间之内，因此在本次研究中利用该数值计算方法是合理的。

2 数值计算结果与分析

本文研究了喷雾初始压强范围在0.2～2.0 MPa的计算工况下，液滴索特平均直径变化规律；喷雾喷口直径范围在2～4 mm的计算工况下，液滴索特平均直径的变化规律；喷雾距离范围200～1 000 mm的计算工况下，液滴索特平均直径的变化规律。

2.1 不同喷雾初始压力的影响

图4为喷雾喷口直径为1 mm，喷雾初始压强在0.2～2.0 MPa，喷雾温度为120 ℃～212 ℃，环境压强为0.1 MPa时，喷雾距离分别为200、400、600、800和1 000 mm时，液滴SMD的变化规律。

图4 不同初始压力的SMD结果Fig.4 SMD results of different initial pressures

随着初始喷雾压力的提高，喷雾液滴的索特平均直径不断降低。当喷雾距离为600、800和1 000 mm时，在初始喷雾压力1.0 MPa之前，喷雾液滴的SMD下降较快，当初始喷雾压强在1.0～2.0 MPa时，索特平均直径的变化趋势趋于平缓。而当喷雾距离为200～400 mm时，索特平均直径在0.2～2.0 MPa的压强区间内一直呈现平缓的下降趋势。在过热喷放过程中，喷放流体与环境中的气体存在较大的速度差，环境气体对喷放流体存在气动力的作用。在该气动力的作用下，喷放流体产生流动不稳定性，破裂的不同尺寸的子液滴。同时由于液相流体的温度高于在环境下的饱和温度，喷放流体以及子液滴中产生了大量的汽化核心，部分液体发生闪蒸相变，成为蒸汽。随着闪蒸相变的不断发生，液相流体尺寸逐渐减小。而在闪蒸相变和外部气动力的共同影响下，喷放流体以及子液滴继续发生破碎，产生了更加细小的液滴。

当喷雾初始压力提高时，喷口处两相流体质量流量升高，汽化核心数量提高，使得更多的液滴发生闪蒸相变。同时大量闪蒸相变产生的蒸汽和液滴混合，使得射流更加剧烈，液滴的韦伯数增大，气动力对液滴破裂的影响增大，射流中的液滴发生破裂的几率更大，液滴的平均粒径随着喷雾初始压力的提高不断减小。当压力持续增大时，在喷口附近两相流以临界流的形式流出，由于在临界流中，流体以临界流速流动，所以当压力继续增大时，流体的流量不变，外部气动力对射流的影响不变，因此只有液滴的自身参数对破裂产生影响，液滴的平均粒径减小的趋势趋于平稳。

2.2 不同喷雾喷口直径的影响

图5展示了喷雾压强1 MPa，喷雾温度为179 ℃，喷雾喷口直径在2.0～4.0 mm，环境压强为0.1 MPa的范围时，喷雾距离分别为400、500、600、700、800和900 mm时，液滴索特平均直径的变化规律。

图5 不同喷口直径的SMD结果Fig.5 SMD results of different injector diameters

当喷口直径增大时，喷雾液滴的索特平均直径呈现一种先增大后降低的趋势。而且当喷雾距离越大时，变化趋势越明显，喷雾液滴索特平均直径越大。当喷雾距离为400 mm和500 mm时，液滴索特平均直径基本无变化。由于喷雾压强为1 MPa，在喷口处的流动属于两相临界流。喷口直径的增大导致喷射面积增大，射流流体的整体流量增大。整体流量的增大导致流体的相对体积增大，射流内部与环境的换热发生恶化，射流中汽化核心密度相对减少，闪蒸汽化的影响减弱，使得SMD有增大的趋势。但是由于射流处于临界流动状态，流体仍以临界流速流动，随着喷口直径的增大，流体的韦伯数增大，气动力对液滴破裂的影响提高，由气动力作用产生的液滴破裂增多，从而导致液滴的SMD变小。在喷射直径增大的初期过程中，射流流体的流量增大，此时闪蒸汽化对液滴尺寸变化的影响起到了主要作用，使得液滴SMD相对增大。而当喷射直径持续增大时，射流流体的韦伯数持续增大，此时气动力对液滴尺寸变化的影响增大，成为主要的影响因素，导致大量液滴破裂，SMD随之减小。

2.3 不同喷雾距离的影响

图6(a)展示了喷雾喷口直径为1 mm，喷雾温度为饱和温度，喷雾距离为200～1 000 mm，环境压强为0.1 MPa的范围时，喷雾初始压强分别为0.4、0.8、1.2、1.6和2.0 MPa时，液滴SMD的变化规律。

图6(b)展示了喷雾初始压力为1 MPa，喷雾温度为179 ℃，喷雾距离在200 mm至1 000 mm的范围时，喷雾喷口直径分别为2、2.5、3、3.5和4 mm时，液滴SMD的变化规律。

图6 不同喷雾距离的SMD结果Fig.6 SMD results of different injector downstream distances

随着喷雾距离的增加，喷雾液滴的SMD不断增大，在600 mm之前，SMD变化趋势较为平缓，而且液滴SMD较小，在60 μm以下。当喷雾距离达到600 mm时，液滴SMD的增加趋势明显，达到了100～200 μm。在喷雾距离为600 mm之前，喷口喷放的两相流体大部分处于过热状态，当射流喷放至环境中时，在喷放流体及子液滴中发生闪蒸现象，产生大量的蒸汽，使得液滴的尺寸处于较小的范围。同时由于射流流速较快，韦伯数较大，气动力对液滴破裂的影响也很大，在两者的共同作用下，流体内的液滴SMD在较小的数值区间内变化。当喷雾距离达到600～800 mm时，射流整体的压力和温度大幅下降，闪蒸相变对液滴尺寸变化的影响减少，同时射流流速降低，液滴的韦伯数减小，气动力对液滴破裂的影响也减小。由于冷凝作用，在之前相变产生的部分蒸汽冷凝成液滴聚集在一起，从而导致液滴SMD不断增大，且变化趋势明显。

3 结论

1)初始喷雾压力越高，喷射流体受到的气动力和闪蒸相变的影响越大，液滴的索特平均直径越小；当初始喷雾压力增大至一定程度后，气动力对流体的影响不再变化，液滴直径减小趋势变缓。

2)随着喷雾喷口直径的增大，闪蒸汽化的影响减弱，但气动力对液滴破裂的影响提高，喷雾液滴的索特平均直径先增大后降低。且当喷雾距离越大时，变增大趋势越明显。

3)当喷雾距离增加时，闪蒸相变作用和气动力的影响减小，喷雾液滴的索特平均直径整体变化趋势呈现先平稳波动后剧烈上升的情况。