引言在工业应用中许多地方都会出现水蒸气在壁面上冷凝的现象[1-4]。冷凝的发生影响设备的安全性和系统的效率，需要对水蒸气在壁面的冷凝过程进行研究。研究的方法主要包括实验研究、理论研究和数值模拟研究等3种。近年来，随着计算流体力学CFD(computational fluid dynamics)的快速发展，促进了数值模拟方法的广泛应用[5-7]。数值模拟方法不仅可以提供整个流场内的详细信息和参数的可视化输出，同时可以不受模型尺寸的限制，降低研究成本，缩短研究周期。利用数值模拟方法实现冷凝的关键是在控制方程中添加源项[8-10]，源项的设置体现了冷凝的传热和传质过程。源项的选取不仅对数值模拟的准确性至关重要，同时也对方程求解过程中的稳定性有重要影响。本文对数值模拟中的控制方程进行阐述，介绍源项的作用，着重介绍应用范围较广的几种源项模型，对各个模型之间的异同点和适用场景进行介绍，指出水蒸气壁面冷凝数值模拟中源项模型存在的问题以及对今后的发展方向进行展望。1控制方程数值模拟是求解控制方程（以下简称方程）的过程。式（1）代表方程的通用形式，方程从左至右，第一项代表非稳态项，第二项表示对流项，第三项表示扩散项，第四项代表源项。前3项的表达式固定，第四项是为了满足不同计算模型之间的差别可以自定义的源项。冷凝现象涉及气相和液相之间的传热传质，必须在质量方程中添加质量源项Sm、能量方程中添加能量源项Se、动量方程中添加动量源项Su、组分方程中添加组分源项Si，具体形式见式（2）~式（5）[11]。通用形式：∂(ρϕ)∂t+∇(ρuϕ)=∇(F∇ϕ)+S (1)质量守恒方程：∂(ρ)∂t+∇(ρu)=Sm (2)能量守恒方程：∂(ρcT)∂t+∇(ρucT)=∇(k∇T)+S (3)动量守恒方程：∂(ρu)∂t+∇(ρuu)=∇(ρu∇u)-∇p+ρg+Su (4)组分守恒方程：∂(ρYi)∂t+∇(ρuYi)=-∇Jt+Si (5)式中：ρ——密度，kg/m3；ϕ——通用变量；u——速度张量，m/s；F——广义扩散系数；Sm——质量源相，kg/(m3∙s)；Se——能量源项，W/m3；Su——动量源项，n/m3；Si——组分源项，kg/(m3∙s)；T——温度，K；g——重力加速度，9.8 m/s。冷凝按照气相的组分不同可以分为纯工质（纯水蒸气）的冷凝和含不凝性气体的混合工质冷凝。在纯工质的情况下不涉及气相中的组分输运，所以不用求式（5）。基于冷凝的原理分析，数值模拟中源项的设置是准确模拟冷凝的关键，而上述4种源项中质量源项和组分源项是冷凝产生的源头，能量源项和动量源项则是由于相变传质和组分扩散所引起的能量和动量的变化，所以能量和动量源项都可以从质量源项或组分源项计算得到。基于此，本文涉及的冷凝模型主要侧重在质量源项的计算。2源项模型2.1相变系数模型1953年，Schrage在Knudsen方程的基础之上[12-13]，根据气体分子动力学理论发展出一种计算质量转移量的模型。该模型认为，蒸发和冷凝现象是一些分子从气液界面逃逸或吸收的过程，当分子从界面逃逸和吸收的分子数相等时，则处于气液平衡状态，且假设相界面是平衡的。质量转移量简化后的方程式如下：m=2σ2-σM2πRpvTv-p1T1 (6)式中：m——气液界面的质量转移量，kg/(m2∙s)；M——水蒸气的摩尔质量；R——通用气体常数，8.314 J/(mol∙K)；p——压力，Pa；T——温度，K；Tv——气相；Tl——液相；σ——冷凝（蒸发）系数。σ代表冷凝的难易程度。如果取值过大，导致数值计算难以收敛；如果取值过小，导致气液交界面温度与饱和温度偏差过大，不能反映实际情况。适宜的值需要进行相关的实验或计算得到。Marek[14]等通过分析已有的实验数据，建议对于液相液面不断变化的情况σ=0.1~1，对于液面稳定的情况σ0.1。Komnos[15]认为冷凝（蒸发）系数与pv和p0下的蒸汽密度有关，表达式为：σp=0.05ρv,p0ρv,pv (7)式中：ρv,p0——p0压力下水蒸气的密度，kg/m3；ρv,pv——pv压力下水蒸气的密度，kg/m3。Finkelstein[16]等对蒸汽在液膜上直接接触冷凝的实验研究中，分析冷凝系数随着温度和压力的变化情况，Marek[14]等拟合实验数据得出冷凝系数与压力的关系为：σp=5.9083×pvp0-1.3686 (8)-式中：p0——参考压力，p0=1×105 Pa。式（7）和式（8）计算的冷凝（蒸发）系数有差异，随着p的增大，两者计算值的差别越来越大。Tanasawa[17]做了进一步的简化，假设气液界面达到饱和状态，计算公式为：m=2σ2-σM2πRρvhfgTsat-TTsat32 (9)式中：hfg——汽化潜热，J/kg；Tsat——饱和温度，K。Lee[18]在前人的基础之上提出了一种更加简洁的模型，该模型没有相界面的限制，将流动看成是直径不变的分散流动（例如细小的液滴或气泡流），非饱和相的温度和饱和温度的差值是相变的动力，冷凝（蒸发）速率也与差值成正比。蒸发过程，TlTsatSl=-Sv=rα1ρ1Tsat-TlTsat (10)冷凝过程，TvTsat-S1=Sv=rαvρvTv-TsatTsat (11)式中：Sl、Sv——相变过程的质量转移量，kg/(m3∙s)。α——各相的体积分数；ρ——各相的密度；T——流体温度；Tsat——饱和温度，K；r——控制相变发生强弱的相变系数1/s。相变系数r的计算公式为：r=AirkM2πRTsathfgα1ρvρ1-ρv (12)式中：Ai——界面浓度，即单位体积内的相界面面积；rk——调节系数；M——水蒸气的摩尔质量，kg/mol；R——通用气体常数，8.314 J/(mol∙K)；hfg——汽化潜热，J/kg；Tsat——流体饱和温度，K。为了简化计算的复杂度，在使用Lee模型时，相变系数r取为一个定值。r的取值与气体成分、温度、压强、模型的结构、网格的大小、时间步长等都有关系[19]，取值的大小直接影响计算结果的准确性。已有研究都是与实验数据进行比较从而确定r值的大小。r的取值非常广，从0.000 1到107都有。Riva[20]等研究环状管内冷凝，r随时间的变化在7.5×105~7.5×107之间变化；Lee[21]等研究环状管中冷凝时，r取10 000；Chen[22]等在研究矩形通道冷凝时取100；另有许多研究中，r取0.1[18,23-25]。Lee模型虽然是在纯蒸汽的基础之上推导得到的，但是也有研究[26-28]在蒸汽中含有不凝性气体的冷凝过程中使用该模型，通过与实验数据进行比较来确定相变系数。邱国栋[29]等对Lee模型中相变系数的取值合理性做了研究。提出采用2个指标来评价相变系数的合理性：潜热换热量占总换热量的份额、流体温度和饱和温度的差值。通过稳态模型推导得出潜热份额和饱和温差的理论表达式。推荐在常见工况下相变系数应该在104以上。Lee模型相较于其他模型其计算公式简洁，模型的适用性广，与各类数值模拟软件结合较好，所以在模拟冷凝时多采用Lee模型。2.2能量平衡模型为了避免相变系数的局限性，发展了一种解决相变问题的纯理论公式[30-31]。该公式运用Fourier定律求解气液交界面处的传热量，并依据潜热确定相应的传质大小。界面处的传热量和质量转移量为：qI⃗=(-λl∂T∂nl)-(-λv∂T∂nv)n⃗ (13)Sl=-Sv=qI⃗n⃗AI,ChfgVC (14)式中：qI⃗——界面指向液相的单位法向量；AI,C——网格单元的面积，m2；Vc——单元格的体积, m3；λ——导热系数，W/(m∙K)；qI⃗——界面处的能量跃迁值，W/m2。许多研究都是在这两个公式的基础之上编写各自的程序代码来模拟相变[32-33]。该模型的关键是准确计算相界面的热流密度，但是在Fluent的UDF中要准确求解相界面的热流密度是很困难的，因此也就限制了该模型在Fluent中的应用。对该模型进一步的简化推导[34-36]，得出如下公式：Sl=-Sv=αvλv,sat+αlλl,sat∇αl∇Thfg (15)式中：λv,sat——饱和气相导热系数，W/(m∙K)；λl,sat——饱和液相的导热系数，W/(m∙K)；α——体积分数。该模型在Fluent的UDF中容易实现。但是实际情况下冷凝发生时，气相是饱和相，液相是不饱和相，而式（15）中气相和液相的导热系数都是按照饱和状态计算的，这就导致其计算精度低于式（14）。为克服式（15）的缺点，Sun[37]等对该模型进行改进，他们将饱和相的热导率按λv,sat=0计算，具体表达式如下：Sl=-Sv=2λl,unsat∇αl,unsat∇Thfg (16)式中：λl,unsat——非饱和相液相的热导率，W/(m∙K)。为了论证该模型的准确性，Sun将该模型运用于一维Stefan问题，发现相界面的变化情况和理论解吻合；又将该模型运用于二维膜态沸腾，验证了其计算结果中气相体积的正确性。但是没有将该模型运用于冷凝的验证。式（15）和式（16）与实际冷凝（蒸发）过程的区别在于，模型中的传质通量发生在相界面附近的局部区域，而实际过程是发生在整个相界面上。为了进一步提高计算精度，Sun[38]等基于Fluent中的VOF模型建立了相变源项计算关系式，提出适用于气液两相中一相饱和、另一相不饱和情况下的新模型。具体计算公式如下：质量源项：Sl=-Sv=QCI,outhfgVCI+QCIs,outhfgVCIs (17)能量源项：Sh=QCI,outVCI网格CI1030Tsat-1030TCIs网格CIs (18)式中：CI——存在相界面的网格；CIs——紧靠CI全部是饱和相的网格；VCI——CI网格的体积，m3；VCIs——CIs网格的体积，m3；QCI,out——网格CI传出的热量，J；QCIs,out——网格CIs传出的热量，J；T——温度，K。Sun[38]等为了验证该模型的精确性，分别使用该模型模拟一维Stefan沸腾的相界面位置、二维膜态沸腾的中气泡的变化情况以及水中单气泡的冷凝过程，模拟结果与前人的模拟结果相吻合。式（15）、（16）、（17）和（18）不涉及相变系数选取的不确定性，但是其在推导的过程中假设气液两相中一相为饱和相、另一相为非饱和相，所以对于过热状态的冷凝还需要在使用时通过试验验证。且该模型需要有初始的相界面，因此不能用于一开始没有冷凝液的模拟。2.3壁面反应模型当水蒸气中有不凝性气体存在时，发生相变的传热传质规律和纯蒸气的过程有非常大的区别。当不凝性气体中的蒸汽在冷却壁面凝结，形成液膜，在液膜附近不凝性气体的浓度会增加，形成不凝性气体层，此时混合气体与壁面交换的热量不仅需要穿过液膜层还需要穿过不凝性气体层；由于壁面不凝性气体的浓度增加，则蒸汽冷凝必须以扩散的方式穿过不凝性气体层到达液膜冷凝[39]。Comini[40]等根据纳维斯托克斯方程并且考虑气体扩散的作用，忽略液膜的厚度，将水蒸气的总通量表示为扩散和对流通量的总和：mv=mdv+mav=-ρD∂ω∂nW+ρωWνn (19)式中：m——单位面积的质量流量，kg/(m2∙s)；ρ——密度，kg/m3；D——扩散系数，m2/s；ω——壁面水蒸气的质量分数；νn——壁面附近混合气体的速度，方向与壁面垂直，具体表达式如下：νn=-D1-ωW∂ω∂nW (20)如果研究的模型中νn比较小，对水蒸气的总通量影响很小，则式可以简化为：mv=mdv=-ρD∂ω∂nW (21)但是若νn的影响不可以忽略，表达式简化后计算公式为：mv=-ρD1-ωW∂ω∂nW (22)上述关系式在推导的过程中忽略了液膜的存在，李成[41]等基于式（21）和式（22）对湿空气流过竖直冷却壁面的冷凝过程进行数值模拟，验证分析了忽略凝结液的合理性，同时探讨两者的适用范围。结果表明，式（21）只适用于壁面温度较低的工况，如制冷空调的蒸发冷凝设备等；当壁面温度比较高时，传热传质规律会发生变化，此时按照式（21）忽略壁面的法向速度νn，则与实际情况出现大的偏差，就必须用式（22）综合考虑扩散和对流对总的传质过程的影响。2.4实验关联式模型许多研究总结归纳了许多含有不凝性气体蒸气冷凝的实验关联式[42-47]，但是这些关联式大部分是对实验数据的简单拟合，只适用于特定的场合并没有发展出通用的形式。2015年Dehbi[45]基于传热传质的类比，对上述实验的原始数据进行分析，提出一个通用的冷凝计算关系式：mv=0.185φD23ρw+ρ∞ρw-ρ∞μ13ln1-ωv,W1-ωv,∞ (23)式中：φ——抽吸作用调节系数；D——质量扩散系数，m2/s；ρw——壁面上混合气体；ρ∞——主流区混合气体的密度，kg/(m∙s)；μ——动力黏度；ωv,W——水蒸气浓度；ωv,∞——主流区水蒸气浓度。该模型适用于压力小于2 000 Pa，水蒸气质量分数在0.1~0.95的湍流和自由对流的情况。陈富强[49]等运用该模型研究湍流和层流情况下含不凝性气体水蒸气在竖直壁面的凝结。研究表明，数值计算结果获得的温度场和不凝性气体的浓度场符合实际情况，体现出不凝性气体在冷却壁面的积聚现象，与理论分析结果一致。还比较该模型与壁面反应模型两者计算结果的差异性，发现在物理量的变化趋势上两者表现出一致性，但是在局部细节方面，壁面反应模型更加的准确，原因是实验关联式模型受到实验条件的限制，在测量数据上会出现误差，从整体上考虑凝结问题，而忽视了一些微观上面的变化情况，导致在细节上没有纯理论推导的壁面反应模型精确。3结语运用数值模拟方法研究水蒸气在壁面冷凝，影响其最终计算结果的因素有很多，其中源项的计算对其结果的影响非常大。目前运用较多的几种源项模型都有其局限性，并不能适用于所有的情况，应针对特定的情况进行合理的假设。相变系数模型是适用范围比较广的模型，但是由于涉及相变系数的不确定性，需要有实验验证来确定相变系数的合理性。虽然相变系数模型是在纯蒸汽的基础之上推导得到的，但是也有研究将其应用在水蒸气中含不凝性气体的冷凝。能量平衡模型虽然不涉及相变系数，但是在理论推导的过程中假设了一相是饱和相、另一相是非饱和相，如果涉及过饱和的情况需要验证其正确性，必须事先存在初始的气液界面，因此不能用于一开始没有冷凝液的数值模拟，且收敛性较差。壁面反应模型考虑到气相为混合相时水蒸气在混合气体中扩散过程的影响，所以适用于含不凝性气体的混合气体中的冷凝，但是忽略了壁面上液膜的厚度，实际上就是忽略了热量在传递过程中液膜的导热热阻。因此当液膜的厚度比较厚，不凝性气体的含量较低，液膜的热阻相对于不凝性气体层的热阻不能忽略时，壁面反应模型就不再适用。实验关联式模型是通过试验数据总结归纳得来，所以适用范围有限，如果模拟的工况或物理模型与实验中的差别很大，不能够保证结果的正确性，同时实验关联式模型在局部细节方面精度不是很好。基于以上所述，对于今后水蒸气壁面冷凝数值模拟中冷凝源项的研究可以从以下几个方面入手：(1) 对现有的模型进行更加细致的研究，找出引起误差的参数，对参数进行合理的优化，提高源项计算的准确性和稳定性。(2) 将理论研究与数值模拟软件结合起来，针对不同软件的特点，对源项模型进行改进，使两者可以更好地结合使用。(3) 针对数值模拟结果，需要开发更多的方法和参数来验证其正确性。(4) 继续深入研究冷凝的传热传质理论机理，综合考虑多种因素的共同作用，为冷凝源项的开发提供理论支持。