鉴于超声速旋流分离器在脱水、脱重烃方面的优良性能[1-2]，研究人员尝试将其应用到天然气液化领域。刘恒等[3-4]对分离器内液化过程进行了探索性理论研究，指出其效率高，前景好。Wen等[5]从宏观上分析了技术的可行性。喷管作为旋流分离器的核心部件，是实现气体液化的关键部位。为此，杨文等[6-8]研究了甲烷在Laval喷管内的液化过程与参数。这些研究的成核模型均引用了内部一致经典理论(以下简称ICCT)[9-10]，该理论在预测近似理想气体的成核率方面因其准确度相对较高而受到广泛的关注与应用。不过，ICCT应用于真实气体会存在一定偏差。
天然气液化通常是在高压低温条件下进行的，气体状态远远偏离理想气体条件，而ICCT是基于理想气体假设推导而来，如果直接将其应用于天然气液化成核研究，将不可避免地引入一定误差。为此，需要对液化成核模型进行修正。笔者针对当前成核模型的不足进行了修正，一方面考虑真实气体效应的影响，另一方面考虑曲率对表面张力的影响，并对低温下液滴表面张力公式进行筛选，以期获得较为准确的成核模型，为基于高速膨胀的气体液化研究提高理论参考。
描述气体凝结成核的理论主要有：经典成核理论(CNT)、半现象学理论(SPT)以及内部一致经典成核理论(ICCT)。经典成核理论出现最早，是其他成核理论的基础，其理论公式最早由Volmer等[11]提出。CNT假设凝结核服从统计热力学Bolzmann分布，认为成核率是Gibbs自由能的指数函数，经过不断改进和完善，经典成核理论推导的成核率计算模型如下：
式中：J为单位时间单位体积内形成的凝结核数，1/(m3·s1)；K为蒸汽的动力学前因子，1/(m3·s1)；ΔG*为临界吉布斯自由能，J；kB为Boltzmann常数，kB =1.38×10-23J/K；Tv为热力学温度，K。
式中：ρl、ρv为液、气体密度，kg/m3；m为液滴质量，kg；σ为液滴表面张力，N/m。
根据CNT理论，在没有异质颗粒、离子的条件下，不断压缩气体使其饱和度S连续升高(定义S=pv/ps，式中pv为气体分压，ps为气体在T时的饱和蒸汽压)，当气体达到饱和状态(S=1)时，凝结不会马上发生，继续压缩气体将进入过饱和状态，当达到临界过饱和度时，凝结才会发生。根据热力学关系式，系统的表面自由能差为：
式中：pl、pv为液、气相压力，Pa；Vl为液滴体积，m3；A为液滴表面积，m2；n为液滴数目；μv, μl分别为气、液相化学势，J。
自由能差不是一个严格单调的函数，而是在某一半径处存在一个极大值，此时对应的液滴半径就是临界半径rc：
临界半径的求解通常可以选用Kelvin方程：
经典成核理论基于理想气体假设，采用理想气体方程求气相分子体积Vl*：
式(6)第1项远小于积分项，可忽略。将式(7)代入式(6)得化学势差：
整理以上各式可得经典成核率计算公式如下：
式中:Θ=σa0/kBT, a0=(36π)1/3vl2/3; vl为分子体积，m3。
式(8)在低温下预测的成核率较实验值偏小[12]，究其原因为：①CNT是基于理想气体假设推导而来的，而低温下气体严重偏离理想状态，所以直接应用成核公式时，误差不可避免；②在公式中表面张力是以3次指数的形式出现，表面张力细微的偏差都会对计算的成核速率产生很大影响，且对于微米级液滴很难获得准确的表面张力，故低温下表面张力公式选择对计算结果影响也很大；③经典成核理论将凝结成核与液滴生长分开处理，忽略了相变过程的传热传质问题。
基于以上问题，Girshick等[9-10]重新修正了自由能计算公式，基于Katz-Wiedersich动力学理论推导了内部一致经典成核理论(ICCT)，形式如下：
程万[13]研究指出，ICCT计算的成核率随温度变化趋势比CNT预测的更准确，建议将ICCT公式作为数值模拟过程中计算成核率的基准公式。
ICCT公式虽然在预测结果上优于CNT公式，但仍没有完全避免CNT的不足：一方面在计算化学势差时仍基于理想气体假设，对真实气体情况显然存在偏差; 另一方面公式使用的表面张力为平面表面张力，而凝结形成的微小液滴并非平面，忽略液滴曲率对表面张力的影响，将不可避免地对预测结果产生较大影响。
基于以上分析，从以下4方面对ICCT公式进行修正，具体包括：考虑非理想气体效应的化学势差修正，求解气体物性的状态方程筛选，考虑曲率对液滴表面张力的修正以及低温下平面表面张力公式的筛选。
对式(6)化学势差求解时，ICCT公式进一步展开后，首先略去了第一项，引入了舍入误差；而后对积分项求解时采用了理想气体方程，对真实气体存在计算误差。因此，本文针对这两点进行修正，引入逸度来进行化学势差的计算，既未作任何忽略又考虑了真实气体效应，气液达到相平衡，两相的化学势差分别为：
式中:fv、fl分别为气、液相逸度, Pa; fv0、fl0分别为气液在标准态下的逸度, Pa。
根据热力学理论，处于气液平衡的两个标准态有如下关系：
将式(11)、式(12)代入式(6)，可得逸度表示的气液相间化学势差为：
气液逸度的计算都要选用真实气体状态方程。在几种状态方程中，SRK和LKP方程式应用最广。本文通过比较选用SRK方程，原因如下：①SRK为立方型，计算相对简单，而LKP方程是超越方程, 计算难度大；②成核模型中，摩尔体积是非常重要的一项参数，因为它在表面张力系数和开尔文方程式中均有应用。根据文献[14]，当采用Peneloux校正时，SRK在预测摩尔体积方面远远优于LKP(SRK得出的结果与实验数据仅有0.8%的误差，而LKP的误差则达6%)。SRK方程形式如下：
式中：R为通用气体常数，J/(mol·K)；VSRK为摩尔体积，m3/mol；ω为离心因子；Tr为对比温度；脚标c为临界点。
当前，计算微观液滴表面张力的方法有：①用平面表面张力近似微小液滴的表面张力，这种近似简化了计算过程，所以被ICCT模型采用；②考虑液滴曲率的修正，使用与液滴半径相关的表面张力计算公式，此时表面张力是液滴半径和平面表面张力的函数，这种模型更贴合实际，但计算相对复杂，此类模型文献中出现最多的有Benson & Shuttleworth模型(BS模型)[15]、Kashchiev模型[16]和Tolman模型[17]。由于BS模型缺少对多种物质在各温度范围内的实验验证(微观测量困难)，引入BS模型可能造成更大的误差；而Kashchiev模型能有效表示全范围(0 < r < ∞)的表面张力，但公式是在理想气体假设下推导的，对于非理想过程较复杂。Tolman模型形式相对简单，但在液滴半径较小区域误差偏大，且精度与Tolman系数δ有关。由于缺少微观液滴表面张力的实验数据，王美利[18]提出基于Kashchiev模型选择合适的δ以提高Tolman模型精度。
式中：σ∞为平面表面张力，N/m；r为液滴半径, m；δ为Tolman系数；s表示气体处于饱和状态。
由于低温或超低温下平面表面张力实验数据较少，很多学者推导了不同温度范围的理论公式，本文将基于Hacker[19]的低温实验数据对这些公式作对比分析，并将不同公式的计算结果与实验数据绘制在一张图里进行分析(见图 1)。几个理论公式分别为：Luijten-Prast修正式[20]、Schnerr-Dohrmann修正式[21]、新表面张力修正式(LD)[22]、Strey等表面张力公式[23]，相应形式如下：
不同公式的预测结果与实验数据对比见图 1。
从图 1可知，当T>265 K时，Strey模型与实验结果贴合最好，但在T < 265 K的低温范围里，预测值与实验值偏差逐渐拉大，在低温下Strey模型已不能很好地预测表面张力。在整个实验温度范围内，SD和LD模型均能与实验值保持较好的吻合且变化趋势一致，不过在更低的温度下(245 K以下)，2个模型的预测趋势反向，即SD模型出现拐点，低温下预测值随温度升高而升高，这显然和其他模型预测趋势相反。值得注意的是，LP模型在T < 268 K范围内与实验值贴合度较高，在更低温度下(245 K以下)预测值与LD模型基本重合，这进一步印证了LD模型在低温区的准确性。因此，本文选择LD模型。
用式(19)计算表面张力时，除了要选择合适的σ∞，还需要确定合适的Tolman系数δ以提高计算成核率的精度。基于文献[16]，本文选取δ=0.05 r0、0.1 r0、0.16 r0、0.3 r0进行对比分析，并将计算结果与Kashchiev模型对比，如图 2所示。
由图 2可知，当0 < r < r0时，δ=0.3 r0与Kashchiev模型较吻合，r>3.5r0时，δ=0.05 r0与Kashchiev模型较吻合，而根据Kashchiev的研究，δ=0.16 r0时，2个模型的总体预测值较为接近，平均偏差较小。因此，选择Tolman系数δ=0.16 r0，且可根据成核率计算结果适当调整。
基于以上对ICCT模型的修正，笔者用C语言编制了成核率的计算程序，将计算结果分别与低温下水的实验数据以及超低温下烃类的实验数据作对比，更加全面地检验修正结果的准确性和有效性。通常成核率的数量级较大，为便于误差评判，作图时以log10J为纵轴，衡量理论结果和实验结果的数量级差别。
对于不同温度下水的成核实验，wolk等[24]做了很多重复性实验，积累了大量的实验数据，且数据可靠性高，重复性好，可作为成核模型研究的重要资料。
图 3为低温下水的实验与计算对比图。由图 3可知，修正前后理论值均高于实验值，这与Luijten[20]研究结果一致。进一步分析得知，修正后的结果比ICCT更接近实验结果，且T>240 K时模型修正前后预测值相对接近，T < 240 K时模型修正前后结构偏差逐渐扩大，而修正后的模型更接近实验结果。这是因为微低温下气体偏离理想气体性质较小，随着温度的降低气体与理想气体偏差变大，忽略非理想气体效应的影响显然与实际不符。另外，还应注意到T>245 K时，各表面张力公式均有实验数据作支撑，预测结果都很吻合，所以模型修正前后预测值差别并不明显；而T < 245 K时，不同表面张力模型偏差逐渐变大，预测趋势甚至相反，因此，低温时筛选平面表面张力模型显得很有必要。
图 4为修正前后理论计算结果与实验数据的平均量级偏差。由图 4可知，ICCT值与实验值偏差在2.5~4个数量级，而修正后偏差降至1~2个数量级。还应注意到，随着温度的降低，ICCT预测偏差逐渐扩大，而修正值仅略微升高，足见考虑真实气体效应和表面张力修正的有效性。因此，从低温下水的数据分析可知，修正后的成核模型比ICCT模型能获得更好的预测效果。
美国Wyslouzil课题组[25]做了一系列关于低温下烃类的成核实验。Ghosh等利用美国国家重点实验室的小角度X射线散射测量装置，先后研究了C7H16、C8H18、C9H20、C10H22在喷管内的超声速液化过程，数据准确度高，稳定性好。
图 5是烃类理论预测值和实验数据的对比。由图 5可知，与水的情况类似，ICCT模型预测值普遍偏高，高于实验值2个数量级左右；而修正后只有正庚烷预测值高于实验值，且仅高1个数量级左右；对于其余烷烃修正值大都略低于实验值，偏差在1个数量级以内。所以，修正后的模型比ICCT模型能取得更好的预测结果。进一步比较ICCT预测值和修正值可以发现，ICCT预测值普遍高于修正值2~3个数量级，这说明对ICCT进行修正的有效性。因此，考虑真实气体效应和液滴曲率对表面张力影响的模型修正比ICCT模型准确性高，适用性强，可作为天然气高速膨胀液化研究的成核基准公式。
针对内部一致经典成核理论应用中存在的不足，本文在ICCT公式推导过程中引入了真实气体效应的影响(逸度计算化学势差)和液滴曲率对表面张力的影响，对公式应用过程中需要用到的状态方程和低温下平面表面张力公式进行了筛选。修正后的模型采用逸度求解化学势差，选用Tolman模型修正表面张力，且在应用时状态方程推荐SRK方程，计算低温下平面表面张力推荐LD公式，Tolman系数δ=0.16 r0为最优。
通过将ICCT模型、修正模型与低温下水和烷烃的成核实验数据对比发现，低温时ICCT预测值均高于实验值，偏差最大可达3~4个数量级；而修正后可将ICCT预测值偏差降低1~2个数量级，将预测值与实验的偏差控制在2个数量级以内，特别对烷烃类偏差均在1个数量级以下。因此，修正后的模型比ICCT模型能获得较好的预测结果，可作为天然气高速膨胀液化研究的成核率基准公式。
2017, Vol. 46 Issue (3): 36-41,60
2017, Vol. 46 Issue (3): 36-41,60