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 蒸汽凝結位置