高级检索

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

沸腾传热数值模拟方法及多管耦合应用特性研究

张振国 谭思超 李小畅 刘思超 冯艺 黄玉健 田瑞峰

张振国, 谭思超, 李小畅, 刘思超, 冯艺, 黄玉健, 田瑞峰. 沸腾传热数值模拟方法及多管耦合应用特性研究[J]. 核动力工程, 2024, 45(2): 130-138. doi: 10.13832/j.jnpe.2024.02.0130
引用本文: 张振国, 谭思超, 李小畅, 刘思超, 冯艺, 黄玉健, 田瑞峰. 沸腾传热数值模拟方法及多管耦合应用特性研究[J]. 核动力工程, 2024, 45(2): 130-138. doi: 10.13832/j.jnpe.2024.02.0130
Zhang Zhenguo, Tan Sichao, Li Xiaochang, Liu Sichao, Feng Yi, Huang Yujian, Tian Ruifeng. Numerical Simulation Method of Boiling Heat Transfer and Its Application Characteristics under Multi-tube Coupled Heat Transfer[J]. Nuclear Power Engineering, 2024, 45(2): 130-138. doi: 10.13832/j.jnpe.2024.02.0130
Citation: Zhang Zhenguo, Tan Sichao, Li Xiaochang, Liu Sichao, Feng Yi, Huang Yujian, Tian Ruifeng. Numerical Simulation Method of Boiling Heat Transfer and Its Application Characteristics under Multi-tube Coupled Heat Transfer[J]. Nuclear Power Engineering, 2024, 45(2): 130-138. doi: 10.13832/j.jnpe.2024.02.0130

沸腾传热数值模拟方法及多管耦合应用特性研究

doi: 10.13832/j.jnpe.2024.02.0130
基金项目: 国家自然科学基金(U22B2093,11805044)
详细信息
    作者简介:

    张振国(1976—),男,博士研究生,现从事核动力装置与反应堆工程研究,E-mail: zhangzhenguoheu@163.com

    通讯作者:

    李小畅,E-mail: lixiaochang@hrbeu.edu.cn

  • 中图分类号: TL334;O368

Numerical Simulation Method of Boiling Heat Transfer and Its Application Characteristics under Multi-tube Coupled Heat Transfer

  • 摘要: 针对直流蒸汽发生器(OTSG)中全流型沸腾传热及一、二次侧耦合换热等复杂物理现象,计算流体动力学(CFD)数值分析普遍面临计算难度大、计算效率低及不确定性大等问题。基于欧拉两流体多相流模型与临界热流密度(CHF)壁面沸腾模型,建立了管内全流型流动沸腾传热数值分析模型,并验证了模型的有效性。基于所验证的模型,开展了数值模型在多管耦合传热下的应用特性研究,明确了该数值模拟方法在多管耦合下的可靠性,并对温度与相分布计算结果对相间作用力模型的敏感性进行了数值分析。研究结果表明:基于欧拉两流体多相流模型与CHF壁面沸腾模型,能够较准确地预测管内水介质由过冷到过热的全流型流动沸腾传热过程,计算的“干涸”点位置及壁面峰值温度与实验值符合较好,最大误差小于10%;基于欧拉两流体多相流模型与CHF壁面沸腾模型的数值方法对多管耦合工况有较好的适用性,计算的二次侧温度与实验结果吻合良好;两相间曳力对壁面温度及空泡份额的计算结果有较明显的影响,但非曳力对壁面温度的影响较小,因此对于大规模工程应用计算,可在分析中不考虑部分相间非曳力的影响。本文研究结果可为OSTG的三维精细化数值分析的模型选择提供有益参考。

     

  • 直流蒸汽发生器(OTSG)[1]既是一体化反应堆及快堆的一、二次侧热交换器,又是两侧的压力边界,其热工水力特性直接影响快堆电站的安全性与经济性,一直以来是核工程领域关注的重点问题。OTSG中存在两侧耦合传热及单相液到过热蒸汽的全流型沸腾传热过程[2],使得其三维精细化数值分析普遍面临计算难度高、计算量大、稳定性差等问题。因此,基于计算流体动力学(CFD)的三维精细化数值分析方法一直未能广泛应用于实际工程中的全尺寸OTSG。

    近年来,国内外学者对于蒸汽发生器的数值模拟做了很多研究。丛腾龙等[3]基于漂移流单流体多相流模型和多孔介质方法,对蒸汽发生器二次侧流场进行了分析。干依燃等[4]基于一维均相流数学模型和传热分区模型,预测了不同工况下管壁面温度跃升的幅度、位置及各传热区域传热管的长度,得到了OTSG工作恶劣区域,该方法较好地模拟了蒸汽发生器一、二次侧的温度分布,但仅能模拟管长方向的参数变化。张振国等[5]基于跨尺度耦合方法,对钠冷快堆OTSG开展了一维与三维跨尺度耦合的精细化数值方法研究。同时,不少学者[6-12]对管内两相流动换热流型识别及相变传热问题进行了数值研究。

    当前研究主要针对一维模型或单管三维模型,壁面加热方式大多为定热流密度加热,且主要关注流动沸腾不同流型段的数值模拟方法,对于三维精细化沸腾传热数值分析方法在多管耦合传热工况下的适应性研究相对较少,而实际工程中的OSTG均为三维多管耦合传热情形。此外,不同沸腾流型段传热机理和传热效率的不同主要由管内的两相份额及分布特性决定,而两相分布又主要由相间作用力决定,因此相间作用力不仅影响相分布和阻力特性,还影响三维温度场特性。然而,目前鲜有针对相间作用力对传热特性计算结果影响的研究。

    目前,针对管内全流型流动沸腾较为精细的数值模拟技术为欧拉两流体多相流模型结合临界热流密度(CHF)壁面沸腾模型的方法,该方法能够较全面地考虑两相间的相互作用力和两相间的传热传质现象。为了明确欧拉两流体多相流模型结合CHF壁面沸腾模型在多管耦合传热下的可靠性及沸腾传热计算结果对相间作用力模型的敏感性,为今后相关研究的模型选取提供参考,本文首先基于这2个模型建立了管内全流型流动沸腾传热数值分析模型并验证了模型的有效性;之后基于所验证的模型,建立了针对一、二次侧多管耦合传热的三维精细化数值分析模型,对二次侧流体两相沸腾传热及“蒸干”现象进行数值模拟研究;最后对管内流动沸腾开展了壁面温度与空泡份额计算结果对相间作用力模型的敏感性分析,并给出了相间作用力模型的选取建议。

    本文采用欧拉两流体多相流模型耦合CHF壁面沸腾模型的方法来模拟全流型沸腾传热问题。气液两相流动与传热分别由各自的连续方程、动量方程及能量方程来描述,并以源项的形式考虑相间质量传递、相间动量交换和相间传热等物理现象,各方程间存在较强的耦合关系。由于模型中各相均拥有独立的控制方程组,相比于其他多相流模型,能够更精细地考虑两相之间的质量、动量与能量交换,但计算量相对较大。

    欧拉两流体多相流模型控制方程如下。

    连续方程:

    $$ \left\{\begin{gathered}\frac{\partial}{\partial t}\left(\alpha_{\text{l}}\rho_{\text{l}}\right)+\nabla\cdot\left(\alpha_{\text{l}}\rho_{\text{l}}\boldsymbol{U}_{\text{l}}\right)=\dot{m}_{\text{vl}}-\dot{m}_{\text{lv}} \\ \frac{\partial}{\partial t}\left(\alpha\mathrm{_v}\rho_{\mathrm{v}}\right)+\nabla\cdot\left(\alpha_{\mathrm{v}}\rho\mathrm{_v}\boldsymbol{U}_{\mathrm{v}}\right)=\dot{m}_{\text{lv}}-\dot{m}_{\text{vl}} \\ \end{gathered}\right. $$ (1)

    动量方程:

    $$ \begin{split} &\frac{\partial }{{\partial t}}\left( {{\alpha _{\text{l}}}{\rho _{\text{l}}}{{\boldsymbol{U}}_{\text{l}}}} \right) + \nabla \cdot \left( {{\alpha _{\text{l}}}{\rho _{\text{l}}}{{\boldsymbol{U}}_{\text{l}}}{{\boldsymbol{U}}_{\text{l}}}} \right) = - {\alpha _{\text{l}}}\nabla p + {\alpha _{\text{l}}}{\rho _{\text{l}}}{\boldsymbol{g}} + \\ & \nabla \cdot {{ {\textit{τ}}}_{\text{l}}} + {C_{\text{D}}}\left( {{{\boldsymbol{U}}_{\text{v}}} - {{\boldsymbol{U}}_{\text{l}}}} \right) + \left( {{{\dot m}_{{\text{vl}}}}{{\boldsymbol{U}}_{{\text{vl}}}} - {{\dot m}_{{\text{l}} {\text{v}}}}{{\boldsymbol{U}}_{{\text{lv}}}}} \right) + {{\boldsymbol{F}}_{{\text{nd}}}} \end{split}\tag{2-1} $$
    $$ \begin{split} &\frac{\partial }{{\partial t}}\left( {{\alpha _{\text{v}}}{\rho _{\text{v}}}{{\boldsymbol{U}}_{\text{v}}}} \right) + \nabla \cdot \left( {{\alpha _{\text{v}}}{\rho _{\text{v}}}{{\boldsymbol{U}}_{\text{v}}}{{\boldsymbol{U}}_{\text{v}}}} \right) = - {\alpha _{\text{v}}}\nabla p + {\alpha _{\text{v}}}{\rho _{\text{v}}}{\boldsymbol{g}} + \\ &\nabla \cdot {{ {\textit{τ}}}_{\text{v}}} + {C_{\text{D}}}\left( {{{\boldsymbol{U}}_l} - {{\boldsymbol{U}}_{\text{v}}}} \right) + \left( {{{\dot m}_{{\text{l}} {\text{v}}}}{{\boldsymbol{U}}_{{\text{lv}}}} - {{\dot m}_{{\text{vl}}}}{{\boldsymbol{U}}_{{\text{vl}}}}} \right) + {{\boldsymbol{F}}_{{\text{nd}}}} \end{split}\tag{2-2} $$

    能量方程:

    $$ \begin{split} &\frac{\partial }{{\partial t}}\left( {{\alpha _{\text{l}}}{\rho _{\text{l}}}{h_{\text{l}}}} \right) + \nabla \cdot \left( {{\alpha _{\text{l}}}{\rho _{\text{l}}}{h_{\text{l}}}{U_{\text{l}}}} \right) \\ &= \nabla \cdot \left( {{\alpha _{\text{l}}}{k_{\text{l}}}\nabla {T_{\text{l}}}} \right) + \left( {{{\dot m}_{{\text{vl}}}}{h_{{\text{vl}}}} - {{\dot m}_{{\text{lv}}}}{h_{{\text{lv}}}}} \right) \end{split}\tag{3-1} $$
    $$ \begin{split} &\frac{\partial }{{\partial t}}\left( {{\alpha _{\text{v}}}{\rho _{\text{v}}}{h_{\text{v}}}} \right) + \nabla \cdot \left( {{\alpha _{\text{v}}}{\rho _{\text{v}}}{h_{\text{v}}}{U_{\text{v}}}} \right) \\ &= \nabla \cdot \left( {{\alpha _{\text{v}}}{k_{\text{v}}}\nabla {T_{\text{v}}}} \right) + \left( {{{\dot m}_{{\text{lv}}}}{h_{{\text{lv}}}} - {{\dot m}_{{\text{vl}}}}{h_{{\text{vl}}}}} \right) \end{split}\tag{3-2} $$

    式中,$ \alpha $为体积分数;ρ为密度,kg/m3U为流速,m/s;Fnd为除曳力以外相间作用力,N/m3p为压强,Pa;g为重力加速度,m/s2τ为剪切应力张量,Pa;CD为曳力系数;h为比焓,J/kg;T为温度,K;k为导热系数,W/(m·K);${\dot m_{{\text{vl}}}}$为气相向液相质量传递速率,kg/(m³·s);${\dot m_{{\text{lv}}}}$为液相向气相质量传递速率,kg/(m³·s);$ {h_{{\text{lv}}}} $与$ {h_{{\text{vl}}}} $分别为蒸发及冷凝时相界面温度下的相变潜热,J/kg;t为时间,s;下标:l表示液相;v表示气相。

    管内流体从过冷态到过热蒸汽态,将经历一系列不同的流动沸腾流型,不同沸腾流型流体的传热特性不同,本文采用CHF壁面沸腾模型模拟流体在管内不同沸腾流型下的传热现象。CHF壁面沸腾模型基于Kurul和Podowski提出的 RPI壁面沸腾数值模型[13],通过引入壁面热流密度分配函数来考虑高空泡份额的纯气相对流换热及气相温度的变化。其核心机理为:将壁面传递给流体的总壁面热流密度$ q\mathrm{_W} $分为单相液对流传热热流密度$ {q_{\text{C}}} $、蒸发传热热流密度$ {q_{\text{E}}} $、气液置换传热热流密度$ {q_{\text{Q}}} $及单相气对流传热热流密度$ {q_{\text{V}}} $四部分,如式(4)所示。

    $$ {q_{\text{W}}} = \left( {{q_{\text{C}}} + {q_{\text{Q}}} + {q_{\text{E}}}} \right)\left[ {1 - f({\alpha _{\text{l}}})} \right] + f({\alpha _{\text{l}}}){q_{\text{V}}} $$ (4)

    式中,$ {\alpha _{\text{l}}} $为液相体积分数;$ f\left( {{\alpha _{\text{l}}}} \right) $为壁面热流量分配函数,可由式(5)求解:

    $$ f\left(\alpha_{\text{l}}\right)=\left\{\begin{array}{ll}0 & \alpha_{\text{l}} < \alpha_{\text{l},1} \\ \dfrac{1}{2}\left[1-\mathrm{cos}\left(\text{π}\dfrac{\alpha_{\text{l}}-\alpha_{\text{l},1}}{\alpha_{\text{l},2}-\alpha_{\text{l},1}}\right)\right] & \alpha_{\text{l},1}\leqslant\alpha_{\text{l}}\leqslant\alpha_{\text{l}\text{,}2} \\ 1 & \alpha_{\text{l}} > \alpha_{\text{l},2}\end{array}\right. $$ (5)

    式中,$ {\alpha _{{\text{l}},1}} $、$ {\alpha _{{\text{l}},2}} $为参考液相体积分数,分别取0.9、0.95[14]。从式(4)与式(5)可以看出,$ f\left( {{\alpha _{\text{l}}}} \right) $主要用来控制模型中近壁面两相份额对壁面热流密度分配的影响,以表征不同沸腾模式下的不同传热机制:当近壁面空泡份额较低时,$ f\left( {{\alpha _{\text{1}}}} \right) $=0,此时的壁面热流量中不存在气相对流传热;当近壁面空泡份额较高时,$ f\left( {{\alpha _{\text{1}}}} \right) $=1,此时壁面传热主要由单相气的对流换热主导。

    CHF壁面沸腾模型中各热流密度的求解模型如式(6)所示:

    $$ \left\{\begin{gathered}q_{\text{C}}=h_{\text{l}}(T_{\text{w}}-T\mathrm{_l})(1-A_{\text{b}}) \\ q_{\text{Q}}=C_{\text{wt}}\frac{2k_{\text{l}}}{\sqrt{\text{π }\lambda_{\text{l}}t}}(T_{\text{w}}-T_{\text{l}})A_{\text{b}} \\ q_{\text{E}}=V_{\text{d}}N_{\text{w}}\rho_{\mathrm{v}}h_{\text{fv}}f \\ q_{\text{V}}=h_{\text{g}}(T_{\text{w}}-T_{\text{g}}) \\ \end{gathered}\right. $$ (6)

    式中,$ {h_{\text{l}}} $与$ {h_{\text{g}}} $分别为单相液与单相气的对流换热系数,由使用的湍流模型对应的近壁面处理方法获得,W/(m2·K);Tw为壁面温度,℃;Ab为壁面气泡的单位覆盖面积,m2;$ {C_{{\text{wt}}}} $为气泡等待系数;$ {\lambda _{\text{l}}} $为液相热扩散系数;$ {V_{\text{d}}} $为汽泡脱离壁面时的体积,m3Nw为壁面成核面密度,m−2f为气泡脱离频率,Hz;t为气泡脱离周期,s;hfv为汽化潜热,J/kg;AbNwfDw的封闭关系式参见文献[15]。

    相间动量传输的本质是两相间的相互作用力。相间作用力决定了气相的分布,而近壁面区域的空泡份额决定了流体与壁面间的传热机制,是各个沸腾流型存在不同传热特性的根本原因。欧拉两流体多相流模型相比于其他多相流模型能更精细地模拟两相间相互作用力,理论上能更准确地预测多相速度场与相分布场。流动沸腾过程中,气相的分布对沸腾传热特性有显著影响,因此欧拉两流体多相流模型中的相间作用力模型对关键热工水力参数的预测有重要意义。

    相间作用力主要包含曳力和非曳力两种,其中非曳力主要考虑升力、湍流耗散力与壁面润滑力。对于各个相间作用力,本文选用如下相间作用力模型[16]进行求解:

    $$ \left\{\begin{gathered}\boldsymbol{F}_{\text{D}}=\frac{\rho_{\mathrm{p}}f\mathrm{_d}}{6\tau\mathrm{_p}}D_{\text{b}}A_{\mathrm{i}}\left(\boldsymbol{U}_{\text{v}}-\boldsymbol{U}_{\text{l}}\right) \\ \boldsymbol{F}_{\text{L}}=C_{\text{L}}\rho_{\text{l}}\alpha_{\mathrm{v}}(\boldsymbol{U}_{\text{v}}-\boldsymbol{U}_{\text{l}})\times(\nabla\times\boldsymbol{U}_{\text{l}}) \\ \boldsymbol{F}_{\text{TD}}=-C_{\text{TD}}\rho_{\text{l}}k_{\text{l}}\times\nabla\alpha_{\text{v}} \\ \boldsymbol{F}_{\text{WL}}=C_{\text{wl}}\alpha_{\text{v}}\rho_{\text{l}}\left|\boldsymbol{U}_{\text{v}}-\boldsymbol{U}_{\text{l}}\right|^2\boldsymbol{n}_{\text{w}} \\ \end{gathered}\right. $$ (7)
    $$ \begin{gathered}\text{ }C_{\text{wl}}=\text{max}\left(0,\frac{C_{\text{W1}}}{D_{\mathrm{b}}}+\frac{C_{\text{W2}}}{y_{\text{w}}}\right) \\ f_{\text{d}}=\frac{C_DRe}{24};\tau_{\text{p}}=\frac{\rho_{\mathrm{v}}D_{\mathrm{b}}^2}{18\mu\mathrm{_l}}\end{gathered} $$

    式中,$ {{{\boldsymbol{F}}}_{\text{D}}} $、$ {{{\boldsymbol{F}}}_{\text{L}}} $、$ {{{\boldsymbol{F}}}_{{\text{TD}}}} $及$ {{{\boldsymbol{F}}}_{{\text{WL}}}} $分别为曳力、升力、湍流耗散力及壁面润滑力,N/m3;$ \boldsymbol{n}_{\text{w}} $为壁面法向单位向量;Ai为相界面积浓度,m2/m3CL为升力系数;CDCL可选用多个模型得到[17],本文后续将对比不同系数模型对计算结果的影响;CTD为湍流耗散力系数,该参数为经验参数,取1或其他经验证的值[18]Cw为壁面润滑力系数,CW1CW2的推荐值分别为−0.01与0.05[19] ;$ {f_{\text{d}}} $为阻力函数;$ {\tau _{\text{p}}} $为气泡“弛豫时间”,s;Re为以两相速度差及气泡直径计算的雷诺数;Db为气泡直径,m;μl为流体的动力粘性系数。

    基于欧拉两流体多相流模型和CHF壁面沸腾模型的流动沸腾数值模拟方法包含一系列封闭子模型和模型参数,对工况依赖性较强,本节首先基于Becker的经典圆管流动沸腾实验数据[20]对数值模型进行验证,以明确该数值方法的有效性。表1为Becker实验的单管几何参数与实验工况,实验中单管轴向高度与径向直径的尺寸比约为700∶1,是典型的大尺寸跨度窄通道,实验条件为高温高压,其几何特征与工况均与OTSG的传热管一致,具有较好的验证价值。

    表  1  圆管几何模型与实验工况
    Table  1.  Geometric Model of Circular Tube and Experimental Conditions
    参数 圆管轴
    向高度/m
    管内径/
    mm
    入口过
    冷度/K
    系统压
    力/MPa
    质量流速/
    (kg·m−2·s−1)
    壁面热流密度/
    (kW·m−2)
    数值 7 10 10 7.02 14.95 797
    下载: 导出CSV 
    | 显示表格

    本文数值模型验证是基于ANSYS FLUENT 求解器进行的,该求解器提供了欧拉两流体多相流模型框架,在加热壁面边界条件的处理中嵌入了CHF壁面沸腾模型,并提供了部分相间作用力模型的模型系数,未提供的模型系数利用用户自定义函数(UDF)加入到求解器中。对于高压条件下流动沸腾数值分析的最佳湍流模型目前尚无定论,但对于直管内的流动沸腾,推荐使用更为简单的标准k-ε模型结合标准壁面函数的方法[21]

    单管几何模型采用全结构化网格进行建模,网格模型如图1所示,数值计算采用SIMPLE算法,动量及湍流相关量采用二阶格式。为了保证$ q\mathrm{_Q} $的计算准确性,根据文献[22]推荐的方法,应确保壁面第一层网格无量纲高度y+落在(200, 250)的区间来确定壁面第一层网格高度[22],因此本计算中网格单元尺寸相对较大。对于全流型流动沸腾传热,两相介质的物性均随压强和温度变化而发生变化。由于管内流体的流动压降相对于系统压力非常小,因此仅考虑物性参数随温度变化。本文中水及水蒸气物性参数均基于国际通用的IAPWS-IF97公式得到。

    图  1  单管几何模型建模及网格划分示意图
    Figure  1.  Single Tube Geometry Modeling and Meshing

    为了明确本文采用的方法预测管内全流型流动沸腾关键热工参数的可靠性,图2给出了壁面温度及质量含气率随轴向高度变化的计算结果与实验值的对比。从图2可以看出,沿着轴向流动方向,质量含气率持续升高,在出口处接近1,变成饱和蒸汽;壁面温度先缓慢上升,在轴向高度约4 m处出现“跃升”现象,随后壁面温度出现一定程度的回落。在轴向高度4 m前,管内主要经历了单相液对流换热与核态沸腾换热过程,由于入口过冷度只有10 K,这一阶段主要为两相段,传热效率较高,因此壁面温度上升非常缓慢。随着壁面液相被“蒸干”,加热壁面逐渐暴露在气相中,传热出现恶化,此时壁面温度“跃升”。随着近壁面完全被气相覆盖,形成了稳定的蒸汽对流换热,对流换热系数增大,壁面温度出现一定程度的回落。从计算值与实验值的对比结果来看,本文采用的欧拉两流体多相流模型耦合CHF壁面沸腾模型的方法很好地预测了全流型沸腾传热各个阶段的传热现象和传热“特征点”,计算的壁面温度和质量含气率与实验数据的最大误差小于10%。

    图  2  壁面温度及质量含气率沿轴向高度变化
    Figure  2.  Variation of Wall Temperature and Mass Quality along Axial Height

    图3进一步给出了圆管壁面的温度分布情况,可以明显看到液相“蒸干”点处温度发生了突变,证明了温度的“跃升”。

    图  3  壁面温度分布
    Figure  3.  Distribution of the Wall Temperature

    以上分析验证了欧拉两流体多相流模型耦合CHF壁面沸腾模型的方法对全流型流动沸腾传热数值分析的有效性。

    图4为圆管内壁处第一层网格节点的y+沿轴向高度的分布。从图中可以看出,整个圆管轴向范围内,y +的取值在80~240之间,符合采用壁面函数法第一层网格高度的要求。此外,结合图2可以看出,管内核态沸腾两相区约在轴向高度1~4 m之间,该区域气液置换传热比较显著,因此该区域y+应落在(200, 250)的区间以保证$ q\mathrm{_Q} $的求解准确性。从图4可以看出,在1~4 m的轴向高度之间,y+大约在180~240之间,能够满足气液置换换热区液相温度求解的要求[21]

    图  4  壁面第一层网格y+沿轴向高度的变化
    Figure  4.  Variation of y+ of the First Layer Mesh Near Wall

    实际OSTG中的换热管束均为多管耦合换热情况,即一、二次侧逆向流动的耦合传热,与单管定热流密度加热的热工水力特性存在明显不同。本文以印度钠冷快中子原型堆(PFBR)模型蒸汽发生器为研究对象,建立一次侧流体域(钠)、传热管固体域、二次侧流体域(水)的三维精细化网格模型。建模时忽略过热段,只考虑预热段和蒸发段,即仅保留一次侧钠与二次侧流体耦合换热的主体部分,简化掉流体的入口和出口腔室、流量分配罩、管板以及支撑板。最终建立的几何模型与网格模型如图5所示,其中网格模型采用全结构化网格建模,多管几何模型结构及运行参数如表2所示。

    图  5  多管耦合的几何模型与网格模型
    Figure  5.  Geometric Model and Mesh Model of Multi-tube Coupling
    表  2  几何模型尺寸与设计参数
    Table  2.  Geometric Model Size and Design Parameters
    参数 设计值
    管外径/mm 17.2
    壁厚/mm 2.3
    管节距/mm 32
    轴向高度/m 14
    外壳直径/mm 164
    钠侧系统压力/MPa 0.1
    水侧系统压力/MPa 17.2
    钠侧质量流量/(kg·s−1) 25.6
    水侧质量流量/(kg·s−1) 2.45
    下载: 导出CSV 
    | 显示表格

    图6给出了蒸汽发生器中心管二次侧流体不同位置的截面平均温度的CFD模拟结果与PFBR模型实验数据[23-24]的对比。可以看出,二次侧流体温度的计算值与实验结果较为吻合,二次侧过冷水进入管内被加热,温度逐渐升高,传热模式逐渐从单相水对流传热转变为过冷沸腾传热,当二次侧流域液体温度达到饱和温度时,传热模式由过冷沸腾转变为饱和沸腾,在轴向高度约9 m处,液膜逐渐蒸发,在轴向高度约11.5 m处,液膜被完全“蒸干”,二次侧换热进入膜态沸腾与蒸汽对流阶段。

    图  6  多管耦合情况下二次侧流体温度
    Figure  6.  Secondary Side Fluid Temperature under Multi-Tube Coupling

    图7进一步给出了二次侧流体与壁面温度沿轴向高度的变化。可以明显看出,在约9 m处壁面温度逐渐“上扬”,此时壁面附近液膜逐渐被蒸发。在多管耦合情况下,由于一、二次侧之间为逆向流动的耦合传热,一侧温度始终受到另一侧温度的限制,因此没有出现单管定热流密度计算时的壁面温度大幅“跃升”现象。

    图  7  多管耦合下二次侧流体与壁面温度变化
    Figure  7.  Variation of Secondary Side Fluid and Wall Temperatures under Multi-tube Coupling

    以上对比分析表明,本文模型得到的数值结果具有较高可信度,数值精度在可接受范围内,本文使用的模型和方法可以较为准确地预测OTSG多管耦合下的一次侧、传热管及二次侧间的耦合传热特性。

    为了明确相间作用力对二次侧管内流动沸腾和相变换热的影响,本文进一步开展了计算结果对相间作用力模型的敏感性分析。仍以Becker经典圆管实验为研究对象,重点分析曳力及三个非曳力(升力、湍流耗散力与壁面润滑力)模型对计算结果的影响。

    图8~图11分别给出了两相间曳力、升力、壁面润滑力及湍流耗散力模型或模型参数对壁面温度及气泡径向分布的影响,其中空泡份额的径向分布结果取自轴向位置3 m处,半径0 mm表示管轴线处,半径5 mm表示管壁附近。

    图  8  曳力系数模型对计算结果的影响
    Figure  8.  Influence of Drag Coefficient Model on Calculation Results
    图  9  升力模型对计算结果的影响
    Figure  9.  Influence of Lift Force Model on Calculation Results
    图  10  壁面润滑力模型对计算结果的影响
    Figure  10.  Influence of Wall Lubrication Force Model on Calculation Results 
    图  11  湍流耗散力系数模型对计算结果的影响
    Figure  11.  Influence of Turbulent Dispersion Force Model on Calculation Results

    图8给出了不同曳力系数模型的计算结果可以看出不同曳力系数模型对壁面温度计算结果有一定影响,Symmetric模型相比于其他模型的误差更大。同时,Symmetric模型计算的空泡份额径向分布与其他曳力模型的差异较大,这主要是由于Symmetric模型较为简单,模型中的物性参数直接选取了两相体积份额的加权平均值。可见两相间曳力影响了空泡份额,而空泡份额将影响壁面与流体之间的换热机制,当壁面附近的空泡份额高时,传热能力相对较低,此时壁面温度更高。由于现有的曳力系数模型是由针对不同气泡形状和气泡数量、采用理论分析或实验研究的方法得到的,这些曳力系数模型对流动沸腾工况和几何结构均较为敏感,因此对于实际工程问题,应进行模型敏感性分析以确定最优曳力系数模型。

    图9为升力模型的敏感性计算结果,从图中可以看出,不同升力模型计算的壁面温度基本相同,但计算的空泡份额径向分布存在较明显的差异。图10给出了壁面润滑力模型对壁面温度和空泡份额计算结果的影响,该计算结果与升力模型的分析结果相似,壁面润滑力对壁面温度的影响较小,对空泡份额有一定的影响。尤其当不考虑壁面润滑力时,计算的壁面温度与采用不同壁面润滑力模型系数的计算结果非常接近,但空泡份额存在一定差异。图11为湍流耗散力系数模型对计算结果的影响,结果表明不同湍流耗散力系数对壁面温度和空泡份额均有一定影响,CTD取0.5时的误差相对较大。由于非曳力封闭模型的推荐模型系数大多为经验值,有一定局限性,因此对于具体的工况需进行模型验证,以确定最佳模型系数。

    以上分析结果表明,曳力对流动沸腾相变传热计算结果的影响相对较大,非曳力的影响相对较小,尤其是两相间的升力和壁面润滑力对数值结果的影响不明显。这说明,对于管内流动沸腾,两相间的主导作用力是曳力。因此对于包含大规模换热管束的实际工程问题,在计算中可以不考虑升力和壁面润滑力,以获得更好的求解稳定性和计算效率。

    本文基于欧拉两流体多相流模型与CHF壁面沸腾模型,建立了管内全流型沸腾传热的数值分析模型,验证了其在多管耦合换热下的有效性,并对多管耦合情况下两相间作用力对沸腾传热特性及空泡份额的影响开展了数值研究,得到如下主要结论:

    (1)基于欧拉两流体多相流模型与CHF壁面沸腾模型的管内流动沸腾数值求解方法,很好地预测到了全流型沸腾传热各个阶段的传热现象和传热“特征点”,计算的壁面温度、质量含气率及“蒸干”点的位置均与实验结果符合良好,最大误差小于10%。

    (2)欧拉两流体多相流模型结合CHF壁面沸腾模型的数值方法对多管耦合传热管束具有较好的适应性,较准确地预测了一、二次侧的耦合传热特性,计算的二次侧温度与实验结果吻合良好:多管耦合传热工况下,由于一、二次侧温度的相互限制,二次侧“蒸干”点处的壁面温度“跃升”幅度远小于定热流密度工况。

    (3)采用相间力模型模拟管内全流型沸腾传热时,两相间曳力对流动沸腾与相变传热计算结果的影响相对较大,非曳力的影响相对较小。对于实际工程问题,部分非曳力模型可以不予考虑,以兼顾数值分析的效率与稳定性。

  • 图  1  单管几何模型建模及网格划分示意图

    Figure  1.  Single Tube Geometry Modeling and Meshing

    图  2  壁面温度及质量含气率沿轴向高度变化

    Figure  2.  Variation of Wall Temperature and Mass Quality along Axial Height

    图  3  壁面温度分布

    Figure  3.  Distribution of the Wall Temperature

    图  4  壁面第一层网格y+沿轴向高度的变化

    Figure  4.  Variation of y+ of the First Layer Mesh Near Wall

    图  5  多管耦合的几何模型与网格模型

    Figure  5.  Geometric Model and Mesh Model of Multi-tube Coupling

    图  6  多管耦合情况下二次侧流体温度

    Figure  6.  Secondary Side Fluid Temperature under Multi-Tube Coupling

    图  7  多管耦合下二次侧流体与壁面温度变化

    Figure  7.  Variation of Secondary Side Fluid and Wall Temperatures under Multi-tube Coupling

    图  8  曳力系数模型对计算结果的影响

    Figure  8.  Influence of Drag Coefficient Model on Calculation Results

    图  9  升力模型对计算结果的影响

    Figure  9.  Influence of Lift Force Model on Calculation Results

    图  10  壁面润滑力模型对计算结果的影响

    Figure  10.  Influence of Wall Lubrication Force Model on Calculation Results 

    图  11  湍流耗散力系数模型对计算结果的影响

    Figure  11.  Influence of Turbulent Dispersion Force Model on Calculation Results

    表  1  圆管几何模型与实验工况

    Table  1.   Geometric Model of Circular Tube and Experimental Conditions

    参数 圆管轴
    向高度/m
    管内径/
    mm
    入口过
    冷度/K
    系统压
    力/MPa
    质量流速/
    (kg·m−2·s−1)
    壁面热流密度/
    (kW·m−2)
    数值 7 10 10 7.02 14.95 797
    下载: 导出CSV

    表  2  几何模型尺寸与设计参数

    Table  2.   Geometric Model Size and Design Parameters

    参数 设计值
    管外径/mm 17.2
    壁厚/mm 2.3
    管节距/mm 32
    轴向高度/m 14
    外壳直径/mm 164
    钠侧系统压力/MPa 0.1
    水侧系统压力/MPa 17.2
    钠侧质量流量/(kg·s−1) 25.6
    水侧质量流量/(kg·s−1) 2.45
    下载: 导出CSV
  • [1] 孙中宁,阎昌琪. 核动力设备[M]. 哈尔滨: 哈尔滨工程大学出版社,2004: 20-23.
    [2] 王冲. 中国实验快堆蒸汽发生器稳态特性分析及热工水力设计[D]. 大连: 大连理工大学,2013.
    [3] CONG T L, TIAN W X, QIU S Z, et al. Study on secondary side flow of steam generator with coupled heat transfer from primary to secondary side[J]. Applied Thermal Engineering, 2013, 61(2): 519-530. doi: 10.1016/j.applthermaleng.2013.08.024
    [4] 干依燃,孙宝芝,齐洪亮,等. 基于传热分区的直流蒸汽发生器换热性能仿真[J]. 化工学报,2015, 66(Z1): 123-129.
    [5] 张振国,李小畅,谭思超,等. 钠冷快堆直流蒸汽发生器跨尺度耦合分析[J]. 原子能科学技术,2022, 56(11): 2354-2364.
    [6] TALEBI S, KAZEMINEJAD H, DAVILU H. Prediction of dryout and post-dryout wall temperatures using film thickness model[J]. Nuclear Engineering and Design, 2012, 244: 73-82. doi: 10.1016/j.nucengdes.2011.12.032
    [7] LI H P, ANGLART H. Prediction of dryout and post-dryout heat transfer using a two-phase CFD model[J]. International Journal of Heat and Mass Transfer, 2016, 99: 839-850. doi: 10.1016/j.ijheatmasstransfer.2016.04.021
    [8] ZHANG H B, HEWITT G F. Phenomenological modelling of CHF in annular flow in annuli using new models of droplet deposition and entrainment[J]. Nuclear Engineering and Design, 2016, 305: 284-292. doi: 10.1016/j.nucengdes.2016.05.029
    [9] ZHANG R, CONG T L, TIAN W X, et al. Prediction of CHF in vertical heated tubes based on CFD methodology[J]. Progress in Nuclear Energy, 2015, 78: 196-200. doi: 10.1016/j.pnucene.2014.10.001
    [10] MALI C R, VINOD V, PATWARDHAN A W. Comparison of phase interaction models for high pressure subcooled boiling flow in long vertical tubes[J]. Nuclear Engineering and Design, 2017, 324: 337-359. doi: 10.1016/j.nucengdes.2017.09.010
    [11] MALI C R, VINOD V, PATWARDHA A W. New methodology for modeling pressure drop and thermal hydraulic characteristics in long vertical boiler tubes at high pressure[J]. Progress in Nuclear Energy, 2019, 113: 215-229. doi: 10.1016/j.pnucene.2019.01.026
    [12] SHI J X, SUN B Z, ZHANG G L, et al. Prediction of dryout and post-dryout wall temperature at different operating parameters for once-through steam generators[J]. International Journal of Heat and Mass Transfer, 2016, 103: 66-76. doi: 10.1016/j.ijheatmasstransfer.2016.07.027
    [13] KURUL N, PODOWSKI M Z. On the modeling of multidimensional effects in boiling channels[C]//27th National Heat Transfer Conference. Minneapolis, Minnesota, USA: ANS, 1991: 28-31.
    [14] IOILEV A, SAMIGULIN M, USTINENKO V, et al. Advances in the modeling of cladding heat transfer and critical heat flux in boiling water reactor fuel assemblies[C]//Proceedings of 12th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-12). US: Pittsburgh, 2007: 1-16.
    [15] 李小畅,郜冶. 压水堆子通道欠热沸腾数值验证及交混翼研究[J]. 原子能科学技术,2013, 47(12): 2208-2215.
    [16] 李小畅. 棒束通道单相及两相流动与传热数值研究[D]. 哈尔滨: 哈尔滨工程大学,2015.
    [17] 霍飞鹏,陈纲,丁炜堃,等. 用CFD方法评价低压过冷沸腾相间传递模型[J]. 高技术通讯,2015, 25(7): 734-745.
    [18] DE BERTODANO L. Turbulent bubbly two-phase flow in a triangular duct[D]. Troy: Rensselaer Polytechnic Institute, 1991.
    [19] ANTAL S P, LAHEY R T JR, FLAHERTY J E. Analysis of phase distribution in fully developed laminar bubbly two-phase flow[J]. International Journal of Multiphase Flow, 1991, 17(5): 635-652. doi: 10.1016/0301-9322(91)90029-3
    [20] BECKER K M, LING C H, HEDBERG S, et al. An experimental investigation of post dryout heat transfer: KTH-NEL-33[R]. Sweden: Dept. of Reactor Technology, 1983.
    [21] MURALLIDHARAN J S, PRASAD B V S S S, PATNAIK B S V, et al. CFD investigation and assessment of wall heat flux partitioning model for the prediction of high pressure subcooled flow boiling[J]. International Journal of Heat and Mass Transfer, 2016, 103: 211-230. doi: 10.1016/j.ijheatmasstransfer.2016.06.050
    [22] EGOROV Y, MENTER F. Experimental Implementation of the RPI Wall Boiling Model in CFX-5.6: TR-04-10[R]. Pennsylvania: ANSYS GmbH, 2004.
    [23] NANDAKUMAR R, SELVARAJ P, ATHMALINGAM S, et al. Thermal simulation of sodium heated once through steam generator for a fast reactor[J]. International Journal of Advances in Engineering Sciences and Applied Mathematics, 2012, 4(3): 127-137. doi: 10.1007/s12572-012-0063-1
    [24] VINOD V, SIVAKUMAR L S, KUMAR V A S, et al. Experimental evaluation of the heat transfer performance of sodium heated once through steam generator[J]. Nuclear Engineering and Design, 2014, 273: 412-420. doi: 10.1016/j.nucengdes.2014.03.034
  • 加载中
图(11) / 表(2)
计量
  • 文章访问数:  498
  • HTML全文浏览量:  59
  • PDF下载量:  51
  • 被引次数: 0
出版历程
  • 收稿日期:  2023-05-03
  • 修回日期:  2023-06-26
  • 刊出日期:  2024-04-12

目录

/

返回文章
返回