1 数值模型

汽轮机低压级级内的自发凝结流动通过商用软件 Ansys CFX 进行模拟，汽液两相均采用欧拉方法进行描述。由于自发凝结生成的一次液滴粒径较小，液滴能很好地随蒸汽运动。因此，计算过程中忽略蒸汽与液滴之间的速度滑移。对于自发凝结流动，该假设与真实情况差异较小，可大幅减少计算量。在模拟过程中，液滴的形成采用 Kantrowitz 提出的非等温修正成核模型计算，表达式如下：

J = (qc / (1+ε)) √(2σ/(πm³)) (ρg²/ρl) exp( - (4πσ r*²) / (3KTg) )(1)

其中，

r* = 2σ / (ρl R Tg ln S)(2)

通过引入表面张力修正系数来修正液滴表面张力σ：

σ = a × σ0(3)

当a=1.0时，即表示采用平面水表面张力计算液滴表面张力，表达式如下：

σ0 = 0.2358 (1 - Tl/Tc)1.256 [1 - 0.625 (1 - Tl/Tc)](4)

可见，通过改变a可以改变液滴表面张力的大小。

在液滴形成后，其尺寸的变化采用 Gyarmathy 提出的液滴增长率模型来描述：

dr/dt = (kg / (r ρl (1+3.18Kn))) ((Tl - Tg) / (hg - hl))(5)

湿蒸汽流动采用 CFX 求解器中高阶格式求解汽液两相流动控制方程组，湍流模型采用 SST 湍流模型，非平衡凝结流动求解过程中水蒸气热力学性质选用 IAPWS-IF 97 标准数据。Ansys CFX 中非平衡凝结流体动力学模型可参考文献。

2 计算结果与分析

2.1 表面张力修正系数与蒸汽膨胀速率的相关性

汽轮机低压级中蒸汽凝结总是伴随着蒸汽的膨胀做功而产生，因此凝结的发生与蒸汽膨胀速率密切相关。笔者基于 Dykas 等的 2 个喷管实验数据，研究蒸汽膨胀速率变化时，表面张力修正系数的取值对湿蒸汽流动数值模拟精度的影响，探讨表面张力修正系数的最佳取值与蒸汽膨胀速率的关系。

D1、D2 喷管为半弧形非对称拉法尔喷管，喉部高度分别为 5.75 cm 和 2.484 cm，壁面曲率半径分别为 70 cm 和 52.5 cm，坐标原点位于喷管进口位置，喉部位于 x=20 cm 处，如图 1 所示。显然，D2 喷管内的蒸汽膨胀速率大于 D1 喷管。

图 1 D1、D2 喷管型线

数值计算时，两喷管的进、出口条件与 Dykas 的实验工况完全一致，进口总压为 98 kPa、总温为 378.15 K，出口静压为 35 kPa，壁面均采用无滑移绝热壁面边界条件，计算网格采用结构化网格。为更好地捕捉到壁面边界层的流动及喷管喉部附近蒸汽的凝结现象，对喷管壁面及喉部的网格进行加密。为保证模拟精度，2 个喷管壁面最大 y+ 均小于 5，网格质量在 0.85 以上。D1、D2 喷管网格数量分别为 42 万和 38 万。其中，D1 喷管的网格示意图见图 2。

图 2 D1喷管网格示意图

图 3 给出了不同表面张力修正系数下，2 个喷管平板壁面上的蒸汽静压分布，图中实验数据来自文献。从图 3 可以看出，蒸汽静压先沿 x 轴方向快速下降，并分别在距喉部 6 cm 和 4 cm 处出现压力突跳，这是由于蒸汽在喷管下游发生凝结释放大量潜热导致压力升高。随后，蒸汽继续膨胀，压力不断下降。D2 喷管的蒸汽压力在 x=35 cm 左右再次发生突跳。这是由于对于 D2 喷管，其蒸汽膨胀速率高，蒸汽出口压力远低于设计压力使得激波向喷管内移动，但激波与凝结冲波相距甚远，不存在干涉作用。

图 3 喷管平板壁面静压分布

对比图 3(a) 和图 3(b) 可以发现，当 a=1.00，即直接采用平面水表面张力计算成核速率时，2 个喷管的模拟结果与实验数据均有较大偏差。当 a 从 1.15 增大至 1.30 时，模拟预测所得凝结发生位置改变并向喷管出口移动，凝结冲波强度也不断减弱。由式(1)~式(3)可知，随着 a 的增大，成核速率公式中的液滴表面张力增大，从而导致液滴临界半径增加，因此需要较高的过冷度来形成凝结核心，成核位置则不断向下游移动，且凝结冲波也减弱。而当 a=1.20 时，模拟结果较好地预测了 2 个喷管内的蒸汽静压分布。这表明通过修正成核模型中液滴表面张力的大小可以提高蒸汽凝结流动预测的准确性。

表 1 给出了不同表面张力修正系数下 2 个喷管内蒸汽凝结位置的实验数据和模拟结果。由表 1 可知，a=1.00 时，D1 喷管、D2 喷管内蒸汽凝结位置相对误差的绝对值分别达到 19.23% 和 15.63%，与实验数据相差较大。当 a=1.20 时，D1 喷管、D2 喷管内蒸汽凝结位置相对误差的绝对值分别减小至 0.12% 和 0.96%。表明调整 a 的值能够使模拟结果准确地捕捉到蒸汽凝结位置，这与图 3 的分析相对应。同时，由蒸汽膨胀速率公式 $\dot{p}(x)=-\frac{v}{p}\frac{dp}{dx}$ 可知，D1 喷管和 D2 喷管在凝结区前的蒸汽膨胀速率分别为 2000 s-1 和 4000 s-1，二者差异较大，但在 a=1.20 时，两喷管内蒸汽凝结位置相对误差的绝对值均小于 1%，在可接受范围内。因此，可认为在一定的误差范围内，蒸汽膨胀速率不是影响表面张力修正系数最佳取值的主要因素。

表 1 蒸汽凝结位置