CN110472270A - 一种钠冷堆非能动余热排出***的建模优化方法 - Google Patents

Abstract

本发明涉及一种钠冷堆非能动余热排出***的建模优化方法，属于钠冷堆余热排出***的技术领域，初步设定热力参数、空气热交换器结构参数、拔风烟囱结构参数等；预设空气侧出口温度，设定传热量等相关计算判断空气侧流动平衡；预设钠侧进口温度与钠侧出口温度设定传热量等相关计算判断钠侧流动平衡；空气侧流动平衡且钠侧流动平衡后经相关计算判断***能量是否平衡，并输出相关设计参数。本发明模型便于电脑语言编程并查找编写错误，可准确高效的完成设备及工艺流程设计工作，并可根据设计需要添加或删减相关模块。本发明的建模优化方法为建立准确的非能动余热排出***模型创造了前提，并显著提高参数估计的精度和速度。

Description

一种钠冷堆非能动余热排出***的建模优化方法

技术领域

本发明涉及钠冷堆余热排出***的技术领域，具体是涉及一种钠冷堆非能动余热排出***的建模优化方法。

背景技术

日本福岛核泄漏事故后，核电***的安全性和可靠性越来越受到人们的重视。非能动安全作为***核技术的关键考虑因素，可有效降低堆芯融毁的概率。钠冷堆相对于传统的压水堆可将铀资源的利用率最大化，同时热响应时间长，冷却剂难发生沸腾现象，反应堆有着更好的安全特性。

目前，对非能动余热排出***进行高效设计并优化尚缺乏完善的解决方法。

发明内容

为了解决非能动余热排出***的高效设计和技术优化的技术问题，本发明提供一种钠冷堆非能动余热排出***的建模优化方法。

为了实现本发明的目的，本发明采用了以下技术方案：一种钠冷堆非能动余热排出***的建模优化方法，钠冷堆非能动余热排出***包括水平高度依次升高的堆芯、事故热交换器、空气热交换器；所述空气热交换器设于拔风烟囱的内部；所述事故热交换器的壳程为堆芯内待冷却的热介质，管程为冷却介质；所述事故热交换器的冷却介质与空气热交换器的加热介质为形成中间自然循环回路的液钠；所述拔风烟囱设有入风口，所述入风口设于空气热交换器的底部；所述空气热交换器的冷却介质为空气，方法包括以下步骤：

S1、特定物性参数为特定条件下的固有参数；初步设定热力参数、空气热交换器结构参数、拔风烟囱结构参数，及事故热交换器液钠出口与空气热交换器液钠出口之间的高度差H₂；

所述热力参数为空气热交换器中的各参数，分别包括空气侧进口温度T_e、设定传热量Q、钠侧进口温度T_s,i与钠侧出口温度T_s,o；

所述空气热交换器为液钠走管程，空气走壳程的翅片管式换热器；翅片管式换热器的换热管为翅片管，所述翅片管包括基管及设置在基管上的翅片；空气热交换器结构参数包括基管结构参数和翅片结构参数；所述基管结构参数包括基管总管数n_s、基管外直径d_o、基管内直径d_i、基管的总管长L_s、基管的排列方式与排列的间距，及根据以上数据得到的沿空气流通方向基管的管排数N_a、单位长度基管外表面积A_o*，基管内流通截面积S_s，空气热交换器翅片管间最小流通面积S_a,min；所述翅片结构参数包括翅片的翅高H_f、节距y，及根据以上数据得到的单位长度翅片表面积A_f与翅化比β；

所述液钠的中间自然循环回路中，钠侧弯头总数为N_s，

所述拔风烟囱结构参数包括空气换热器的壳程空气出口与拔风烟囱出口高度差H₁、烟道内直径D_i、烟道流通截面积S_c、烟道弯头数N_c；

所述特定物性参数包括钠侧比热cp_s与空气侧比热cp_a；不同温度对应下的空气密度ρ_a；不同温度对应下的液钠密度ρ_s；热耗散q、直管阻力系数f、局部阻力系数ξ、进出口阻力系数ε、污垢热阻r、钠侧沿管程的雷诺数Re、导热系数λ、对数平均温差ΔT_m、粘度μ、贝克来数Pe_s与对数平均温差修正系数φ；

S2、预设空气热交换器的空气侧出口温度T_a,o1，根据设定传热量Q确定空气质量流量m_a；计算得到空气侧浮升力F_a和空气侧总阻力降ΔP_a后；计算空气侧浮升力F_a与空气侧总阻力降ΔP_a的相对误差当相对误差大于1％时，则返回调整预设空气侧出口温度T_a,o1；当相对误差小于等于1％时，则判断空气侧流动平衡；

预设钠侧进口温度T_s,i与钠侧出口温度T_s,o，由设定传热量Q确定液钠的质量流量m_s，计算得到钠侧提升力F_s与钠侧总阻力降ΔP_s；计算钠侧提升力F_s与钠侧总阻力降ΔP_s的相对误差当相对误差大于1％时，则返回调整钠侧进口温度T_s,i与钠侧出口温度T_s,o；当相对误差小于等于1％时，则判断钠侧流动平衡；

S3、待S2中的空气侧流动平衡且钠侧流动平衡后，计算得到空气侧传热系数h_a和钠侧传热系数h_s；进而再求得空气热交换器总传热系数K、对数平均温差ΔT_m与总传热量Q_s；再计算总传热量Q_s与设定传热量Q的相对误差当相对误差大于1％时，则判断***能量不平衡，并返回步骤S2重新计算；当相对误差小于等于1％时，则判断***能量平衡，并输出空气质量流量m_a、钠侧质量流量m_s、空气侧进口温度T_e、空气侧出口温度T_a,o1、钠侧进口温度T_s,i、钠侧出口温度T_s,o、空气热交换器总传热系数K、总传热量Q_s、总传热量Q_s与设定传热量Q的相对误差

钠侧比热cp_s是指管路翅片管中液钠的比热容；与空气侧比热cp_a是指管路翅片管中空气的比热容。

空气热换热器壳程空气的入口与拔风烟囱的入风口可以是在同一水平高度，也可以是空气热换热器壳程空气的入口的水平高度高于拔风烟囱的入风口。

进一步地，步骤S2所述空气侧浮升力F_a和空气侧总阻力降ΔP_a的计算步骤如下：

空气热交换器中换热的能量平衡方程如下式(1)，钠侧能量损失量与空气侧能量增加量均分别等于设定传热量Q：

m_scp_s(T_s,i-T_s,o)＝m_acp_a(T_a,o1-T_e)＝Q (1)

其中，T_s,i、T_s,o分别为空气热交换器的钠侧进口温度、钠侧出口温度；T_e、T_a,o1分别为空气热交换器的空气侧进口温度、空气侧出口温度；m_s、m_a分别为空气热交换器的钠侧质量流量、空气质量流量；cp_s、cp_a分别为空气热交换器中的钠侧比热、空气侧比热；

拔风烟囱内的空气侧浮升力F_a的计算公式如下式(2)：

其中，H₁为空气换热器的壳程空气出口与拔风烟囱出口高度差；ρ_e为空气热交换器的空气侧进口温度T_e下的空气密度；ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，式(2)中的ρ_a,om与式(3)中的T_a,om对应，并满足公式T_a,om由式(3)求得，其中，T_a,o1为空气热交换器的空气侧出口温度，q_c为拔风烟囱与外部环境之间的热耗散，m_a为空气热交换器的空气质量流量，cp_a为空气侧比热；g为重力加速度；dh为高度h的微分；

空气侧总阻力降ΔP_a的计算公式如下式(4)：

其中，f_a为翅片管外阻力系数，N_a为沿空气流通方向基管的管排数，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积，ρ_a,m为空气热交换器中的空气平均密度；f_c为烟囱内沿程阻力系数，S_c为烟道流通截面积，ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，D_i为拔风烟囱烟道内直径；N_c为拔风烟囱烟道弯头数，ξ_c为拔风烟囱局部阻力系数，ε_c为拔风烟囱进出口阻力系数；

步骤S2所述的钠侧提升力F_s与钠侧总阻力降ΔP_s计算步骤如下：忽略事故热交换器液钠出口与空气热交换器液钠入口之间管路的热量损耗，事故热交换器液钠出口温度等于空气热交换器液钠入口温度；事故热交换器液钠出口与经过换热后的空气热交换器液钠出口之间存在温度差，并存在密度差，两者之间存在高度差H₂；事故热交换器液钠出口与空气热交换器液钠出口之间钠侧提升力F_s的计算公式如下式(5)：

其中，H₂为事故热交换器液钠出口与空气热交换器液钠出口之间的高度差；ρ_s,o、ρ_s,i分别为空气热交换器的出口液钠密度、进口液钠密度；g为重力加速度；dh为高度h的微分；

钠侧总阻力降ΔP_s的计算公式如下式(6)：

其中，式(6)中的f_s为钠侧沿程阻力系数，钠侧沿程阻力系数f_s与雷诺数Re_s的关系如式(7)；m_s为钠侧质量流量，S_s为基管内流通截面积；ρ_s,m为自然循环回路中液钠的平均密度；L_sur为自然循环回路中除基管外的剩余段管长；L_s为基管的总管长；d_i为基管内直径；N_s为钠侧弯头总数；ξ_s为钠侧局部阻力系数。

进一步地，空气热交换器中的翅片管为三角形错列排布；所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

步骤S3所述空气对流传热系数h_a的计算公式如下式(8-1)：

其中：λ_a为空气导热系数；d_o为基管外直径；m_a为空气质量流量；S_a,min为空气热交换器翅片管间最小流通面积；μ_a为空气粘度；cp_a为空气侧比热；A_o为单位长度翅片管总外表面积；为单位长度基管外表面积；

液钠冷却过程，于均匀热流边界条件下，液钠流动的对流传热系数h_s的计算公式如下式(9)：

其中：λ_s为液钠的导热系数；d_i为基管内直径；Pe_s为贝克来数；

根据式(8-1)和式(9)，空气热交换器总传热系数K的计算公式如下式(10)：

其中：β为翅化比，h_a为空气侧对流传热系数；r_a为空气侧污垢热阻；A_o为单位长度翅片管总外表面积；λ_w为管材导热系数；d_o为基管外直径；d_i基管内直径；r_s为钠侧污垢热阻；η为翅片壁面总效率；h_s为钠侧对流传热系数；管材导热系数是指翅片管选用材料的传导效率，基管和翅片为相同材质，具有相同的传导效率。

***总传热量Q_s的计算公式如下式(11)：

Q_s＝φKn_sL_sA_oΔT_m (11)

其中，φ为对数平均温差修正系数；K为空气热交换器总传热系数；n_s为基管总管数；L_s为基管的总管长；A_o为单位长度翅片管总外表面积；ΔT_m为对数平均温差。

进一步地，所述翅片管为正方形直列排布；所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

步骤S3所述空气对流传热系数h_a的计算公式如下式(8-2)：

其中：λ_a为空气导热系数；d_o为基管外直径；m_a为空气质量流量；S_a,min为空气热交换器翅片管间最小流通面积；μ_a为空气粘度；cp_a为空气侧比热；A_f单位长度翅片表面积；为单位长度基管外表面积；

液钠冷却过程，于均匀热流边界条件下，液钠流动的对流传热系数h_s的计算公式如下式(9)：

其中：λ_s为液钠的导热系数；d_i为基管内直径；Pe_s为贝克来数；

根据式(8-2)和式(9)，空气热交换器总传热系数K的计算公式如下式(10)：

其中：β为翅化比，h_a为空气侧对流传热系数；r_a为空气侧污垢热阻；A_o为单位长度翅片管总外表面积；λ_w为管材导热系数；d_o为基管外直径；d_i基管内直径；r_s为钠侧污垢热阻；η为翅片壁面总效率；h_s为钠侧对流传热系数；

***总传热量Q_s的计算公式如下式(11)：

Q_s＝φKn_sL_sA_oΔT_m (11)

其中，φ为对数平均温差修正系数；K为空气热交换器总传热系数；n_s为基管总管数；L_s为基管的总管长；A_o为单位长度翅片管总外表面积；ΔT_m为对数平均温差。

进一步地，当翅片管为三角形错列排布时，公式(4)中的翅片管外阻力系数f_a为：

当翅片管为正方形直列排布时，公式(4)中的翅片管外阻力系数f_a为：

其中，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积；d_o为基管外直径；μ_a为空气粘度；P_t为管束横向管间距；P_l为管束纵向管间距；y为翅片管节距；H_f为翅高。其中，以翅片管为正方形直列排布的情况考虑；空气是自低处向高处流动，流通方向为竖直方向，其管束纵向管间距即是为管束沿空气流通方向上的相邻两管束自间的距离；液钠的流通在空气热交换器的管程内，其是水平流动的，管束横向管间距是指同一水平高度上的相邻两管束之间的距离。

进一步的，构建物性参数计算模块、空气侧流动平衡计算模块、液态钠侧流动平衡计算模块与***能量平衡计算模块；所述特定物性参数存于物性参数计算模块内；步骤S2从预设空气侧出口温度至判断空气侧流动平衡的过程属于空气侧流动平衡计算模块；步骤2从预设钠侧进口温度与钠侧出口温度至判断钠侧流动平衡的过程属于液态钠侧流动平衡计算模块；步骤S3的过程属于***能量平衡计算模块；所述物性参数计算模块分别与空气侧流动平衡计算模块、液态钠侧流动平衡计算模块连接；所述空气侧流动平衡计算模块判断空气侧流动平衡后进入***能量平衡计算模块；液态钠侧流动平衡计算模块判断钠侧流动平衡后进入***能量平衡计算模块。

本发明的有益效果在于：本发明钠冷堆非能动余热排出***主要是由空气热交换器、事故热交换器、拔风烟囱组成，拔风烟囱底部设有风门(风门为拔风烟囱的入风口)。***包括三个自然循环部分，热钠池中事故热交换器和堆芯构成了一回路自然循环，事故热交换器和空气热交换器之间的液钠循环构成了中间回路自然循环，位于拔风烟囱中的空气热交换器和外部环境构成了空气冷却自然循环。事故发生时，一回路自然循环将堆芯热量转移至事故热交换器中，热量通过事故热交换器经中间回路传递至空气热交换器，最终依靠拔风烟囱中空气冷却循环将热量转移至大气中。本发明主要基于钠冷堆余热排出***的工艺过程进行模块化分析，根据设计传热量与传热元件的物理特性参量、自然动力源高度和内直径等参数，即可获得***各状态点热力参数。计算结果为钠冷反应堆余热排出***的安全设计提供理论依据和实用价值。

数学模型的模块化分析方法依据实现功能的差异将设计模型分为不同模块，便于电脑语言编程并查找编写错误，可准确高效的完成设备及工艺流程设计工作，并可根据设计需要添加或删减相关模块。选择合理的模块化结构是建立准确的非能动余热排出***模型的重要前提；而恰当地确定模型参数的范围，可以显著提高参数估计的精度和速度。

并且，本发明在构建的排出***模型的基础上还验证了传热元件不同的物理特性参量、自然动力源高度和内直径参数对***模型数据选取的影响；本发明为非能动自然循环***，总驱动压等于公式(2)求得的空气侧浮升力F_a；结果表明，1)、***的总传热系数和总驱动压随着环境温度升高而逐渐减小，且环境温度对驱动压力的影响更为明显。2)、随着拔风烟囱高度增加，***总传热系数和总驱动压均增大，且增大趋势不断变缓，存在设计最优值；3)、翅片管翅高减小，***总传热系数及单位压降传热系数大幅增加，对***的传热性能影响明显；4)、拔风烟囱内直径增加，***总传热系数增大，总驱动压减小，且变化趋势不断变缓，内直径增加对***传热和流动有积极影响。5)、翅片厚度增加，***总传热系数及单位压降传热系数大幅增加，对***的传热性能影响明显。以上趋势的总结为之后实际建筑的参数设置有着指导意义。

附图说明

图1为本发明钠冷堆余热排出***模型的结构示意图；

图2为钠冷堆余热排出***模型计算流程图；

图3(a)-3(c)为不同环境温度下的影响实验线形图；

图4(a)-4(d)为不同烟囱高度及翅片高度的影响实验线形图；

图5(a)-5(d)为不同烟囱内直径及翅片厚度的影响实验线形图。

具体实施方式

以下结合附图和具体实施方式对本发明技术方案作进一步详细描述。

图1为钠冷堆余热排出***模型示意图，所述钠冷堆非能动余热排出***包括空气热交换器、拔风烟囱、事故热交换器，及堆芯与事故热交换器的作用形成的对液钠升温提升的热钠池。

表1各公式中字母含义

表2空气物性与空气侧温度T_a(K)的关系

表3液钠物性与钠侧温度T_s(K)的关系

实施例1

非能动余热排出***初始设定传热量与传热元件的物理特性参量、自然动力源烟囱的高度和内直径参数如下表4所示；

表4

翅片管结构

设备结构

如图2钠冷堆余热排出***模型计算流程图所示的，经以下步骤进行计算：一种钠冷堆非能动余热排出***的建模优化方法，钠冷堆非能动余热排出***包括水平高度依次升高的堆芯、事故热交换器、空气热交换器；所述空气热交换器设于拔风烟囱的内部；所述事故热交换器的壳程为堆芯内待冷却的热介质，管程为冷却介质；所述事故热交换器的冷却介质与空气热交换器的加热介质为形成中间自然循环回路的液钠；所述拔风烟囱设有入风口，所述入风口设于空气热交换器的底部；所述空气热交换器的冷却介质为空气，方法包括以下步骤：

S1、特定物性参数为特定条件下的固有参数；初步设定热力参数、空气热交换器结构参数、拔风烟囱结构参数，及事故热交换器液钠出口与空气热交换器液钠出口之间的高度差H₂；

所述热力参数为空气热交换器中的各参数，分别包括空气侧进口温度T_e、设定传热量Q、钠侧进口温度T_s,i与钠侧出口温度T_s,o；

所述空气热交换器为液钠走管程，空气走壳程的翅片管式换热器；翅片管式换热器的换热管为翅片管，所述翅片管包括基管及设置在基管上的翅片；空气热交换器结构参数包括基管结构参数和翅片结构参数；所述基管结构参数包括基管总管数n_s、基管外直径d_o、基管内直径d_i、基管的总管长L_s、基管的排列方式与排列的间距，及根据以上数据得到的沿空气流通方向基管的管排数N_a、单位长度基管外表面积A_o ^*，基管内流通截面积S_s，空气热交换器翅片管间最小流通面积S_a,min；所述翅片结构参数包括翅片的翅高H_f、节距y，及根据以上数据得到的单位长度翅片表面积A_f与翅化比β；

所述液钠的中间自然循环回路中，钠侧弯头总数为N_s，

所述拔风烟囱结构参数包括空气换热器的壳程空气出口与拔风烟囱出口高度差H₁、烟道内直径D_i、烟道流通截面积S_c、烟道弯头数N_c；

所述特定物性参数包括钠侧比热cp_s与空气侧比热cp_a；不同温度对应下的空气密度ρ_a；不同温度对应下的液钠密度ρ_s；热耗散q、直管阻力系数f、局部阻力系数ξ、进出口阻力系数ε、污垢热阻r、钠侧沿管程的雷诺数Re、导热系数λ、对数平均温差ΔT_m、粘度μ、贝克来数Pe_s与对数平均温差修正系数φ；

S2、预设空气热交换器的空气侧出口温度T_a,o1，根据设定传热量Q确定空气质量流量m_a；计算得到空气侧浮升力F_a和空气侧总阻力降ΔP_a后；计算空气侧浮升力F_a与空气侧总阻力降ΔP_a的相对误差当相对误差大于1％时，则返回调整预设空气侧出口温度T_a,o1；当相对误差小于等于1％时，则判断空气侧流动平衡；

预设钠侧进口温度T_s,i与钠侧出口温度T_s,o，由设定传热量Q确定液钠的质量流量m_s，计算得到钠侧提升力F_s与钠侧总阻力降ΔP_s；计算钠侧提升力F_s与钠侧总阻力降ΔP_s的相对误差当相对误差大于1％时，则返回调整钠侧进口温度T_s,i与钠侧出口温度T_s,o；当相对误差小于等于1％时，则判断钠侧流动平衡；

S3、待S2中的空气侧流动平衡且钠侧流动平衡后，计算得到空气侧传热系数h_a和钠侧传热系数h_s；进而再求得空气热交换器总传热系数K、对数平均温差ΔT_m与总传热量Q_s；再计算总传热量Q_s与设定传热量Q的相对误差当相对误差大于1％时，则判断***能量不平衡，并返回步骤S2重新计算；当相对误差小于等于1％时，则判断***能量平衡，并输出空气质量流量m_a、钠侧质量流量m_s、空气侧进口温度T_e、空气侧出口温度T_a,o1、钠侧进口温度T_s,i、钠侧出口温度T_s,o、空气热交换器总传热系数K、总传热量Q_s、总传热量Q_s与设定传热量Q的相对误差

构建物性参数计算模块、空气侧流动平衡计算模块、液态钠侧流动平衡计算模块与***能量平衡计算模块；所述特定物性参数存于物性参数计算模块内；步骤S2从预设空气侧出口温度至判断空气侧流动平衡的过程属于空气侧流动平衡计算模块；步骤2从预设钠侧进口温度与钠侧出口温度至判断钠侧流动平衡的过程属于液态钠侧流动平衡计算模块；步骤S3的过程属于***能量平衡计算模块。

步骤S2所述空气侧浮升力F_a和空气侧总阻力降ΔP_a的计算步骤如下：

空气热交换器中换热的能量平衡方程如下式(1)，钠侧能量损失量与空气侧能量增加量均分别等于设定传热量Q：

m_scp_s(T_s,i-T_s,o)＝m_acp_a(T_a,o1-T_e)＝Q (1)

其中，T_s,i、T_s,o分别为空气热交换器的钠侧进口温度、钠侧出口温度；T_e、T_a,o1分别为空气热交换器的空气侧进口温度、空气侧出口温度；m_s、m_a分别为空气热交换器的钠侧质量流量、空气质量流量；cp_s、cp_a分别为空气热交换器中的钠侧比热、空气侧比热；

拔风烟囱内的空气侧浮升力F_a的计算公式如下式(2)：

其中，H₁为空气换热器的壳程空气出口与拔风烟囱出口高度差；ρ_e为空气热交换器的空气侧进口温度T_e下的空气密度；ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，式(2)中的ρ_a,om与式(3)中的T_a,om对应，并满足公式T_a,om由式(3)求得，其中，T_a,o1为空气热交换器的空气侧出口温度，q_c为拔风烟囱与外部环境之间的热耗散，m_a为空气热交换器的空气质量流量，cp_a为空气侧比热；g为重力加速度；dh为高度h的微分；

空气侧总阻力降ΔP_a的计算公式如下式(4)：

其中，f_a为翅片管外阻力系数，N_a为沿空气流通方向基管的管排数，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积，ρ_a,m为空气热交换器中的空气平均密度；f_c为烟囱内沿程阻力系数，S_c为烟道流通截面积，ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，D_i为拔风烟囱烟道内直径；N_c为拔风烟囱烟道弯头数，ξ_c为拔风烟囱局部阻力系数，ε_c为拔风烟囱进出口阻力系数；

当翅片管为三角形错列排布时，公式(4)中的翅片管外阻力系数f_a为：

当翅片管为正方形直列排布时，公式(4)中的翅片管外阻力系数f_a为：

其中，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积；d_o为基管外直径；μ_a为空气粘度；P_t为管束横向管间距；P_l为管束纵向管间距；y为翅片管节距；H_f为翅高。

步骤S2所述的钠侧提升力F_s与钠侧总阻力降ΔP_s计算步骤如下：忽略事故热交换器液钠出口与空气热交换器液钠入口之间管路的热量损耗，事故热交换器液钠出口温度等于空气热交换器液钠入口温度；事故热交换器液钠出口与经过换热后的空气热交换器液钠出口之间存在温度差，并存在密度差，两者之间存在高度差H₂；事故热交换器液钠出口与空气热交换器液钠出口之间钠侧提升力F_s的计算公式如下式(5)：

其中，H₂为事故热交换器液钠出口与空气热交换器液钠出口之间的高度差；ρ_s,o、ρ_s,i分别为空气热交换器的出口液钠密度、进口液钠密度；g为重力加速度；dh为高度h的微分；

钠侧总阻力降ΔP_s的计算公式如下式(6)：

其中，式(6)中的f_s为钠侧沿程阻力系数，钠侧沿程阻力系数f_s与雷诺数Re_s的关系如式(7)；m_s为钠侧质量流量，S_s为基管内流通截面积；ρ_s,m为自然循环回路中液钠的平均密度；L_sur为自然循环回路中除基管外的剩余段管长；L_s为基管的总管长；d_i为基管内直径；N_s为钠侧弯头总数；ξ_s为钠侧局部阻力系数。

进一步地，所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

空气热交换器中的翅片管为三角形错列排布；步骤S3所述空气对流传热系数h_a的计算公式如下式(8-1)：

所述翅片管为正方形直列排布；步骤S3所述空气对流传热系数h_a的计算公式如下式(8-2)：

其中：λ_a为空气导热系数；d_o为基管外直径；m_a为空气质量流量；S_a,min为空气热交换器翅片管间最小流通面积；μ_a为空气粘度；cp_a为空气侧比热；A_f单位长度翅片表面积；A_o为单位长度翅片管总外表面积；为单位长度基管外表面积；

液钠冷却过程，于均匀热流边界条件下，液钠流动的对流传热系数h_s的计算公式如下式(9)：

其中：λ_s为液钠的导热系数；d_i为基管内直径；Pe_s为贝克来数；

根据式(8-1)和式(9)，或是式(8-2)和式(9)，空气热交换器总传热系数K的计算公式如下式(10)：

其中：β为翅化比，h_a为空气侧对流传热系数；r_a为空气侧污垢热阻；A_o为单位长度翅片管总外表面积；λ_w为管材导热系数；d_o为基管外直径；d_i基管内直径；r_s为钠侧污垢热阻；η为翅片壁面总效率；h_s为钠侧对流传热系数；

***总传热量的计算公式如下式(11)：

Q_s＝φKn_sL_sA_oΔT_m (11)

其中，φ为对数平均温差修正系数；K为空气热交换器总传热系数；n_s为基管总管数；L_s为基管的总管长；A_o为单位长度翅片管总外表面积；ΔT_m为对数平均温差。

对本发明所构建***的合理性与有效性进行验证，以下分别对不同的空气侧进口温度、换热器的翅片参数和烟囱参数进行测试；表4是对应的当翅片管为三角形错列排布时，余热排出***中的设计参数、换热管结构设计参数以及设备机构参数，其本验证实验为单一因素考察实验，除变量外的其他参数均同表4；表4中拔风烟囱高度近似等于空气换热器的壳程空气出口与拔风烟囱出口高度差H₁；事故热交换器与空气热交换器高度差近似等于事故热交换器液钠出口与空气热交换器液钠出口之间的高度差H₂。

环境温度的影响：环境温度即是空气侧进口温度；预设不同空气侧进口温度，并分别进行***热力参数计算结果，不同环境温度下***热力参数计算结果如下表5和图3(a)-3(c)所示。本发明为非能动自然循环***，总驱动压等于公式(2)求得的空气侧浮升力F_a；结果表明，***的总传热系数和总驱动压随着环境温度升高逐渐减小，且环境温度的变化对驱动压力的影响更为明显。如图3所示的，***的总驱动压和总传热系数随着环境温度升高逐渐减小，且环境温度对驱动压力的影响更为明显。

表5不同环境温度下***热力参数计算结果

拔风烟囱高度及翅片高度的影响：于相同翅片高度，不同烟囱高度下计算***热力参数，计算结果如表6-1、表6-2、表6-3与图4(a)-4(d)所示，由表中数据可知，翅片高度固定时，随着拔风烟囱高度增加，***总驱动压和总传热系数均增大，且增大趋势不断变缓，则可预料其总传热系数必然会存在设计最优值；另外，随着翅片高度减小，***总传热系数及单位压降传热系数大幅增加，对***的传热性能也有较高的影响。结果表明，拔风烟囱高度增加，***总传热系数和总驱动压均增大，且增大趋势不断变缓，存在设计最优值；翅片高度减小，***总传热系数及单位压降传热系数大幅增加，对***的传热性能影响明显。

表6-1翅片高度9mm时不同烟囱高度***热力参数

表6-2翅片高度11mm时不同烟囱高度***热力参数

表6-3翅片高度13mm时不同烟囱高度***热力参数

烟囱内直径及翅片厚度的影响：于相同翅片厚度时，不同烟囱内直径下***热力参数的测试，测得结果如下表7-1、表7-2、表7-3与图5(a)-5(d)所示。由表中数据可以明确得出，随着拔风烟囱内直径增加，***总传热系数增大，总驱动压减小，且变化趋势不断变缓，内直径增加对***传热和流动有积极影响。随着翅片厚度增加，***总传热系数及单位压降传热系数大幅增加，对***的传热性能影响明显。

表7-1翅片厚度0.4mm时不同烟囱内直径***热力参数

表7-2翅片厚度0.6mm时不同烟囱内直径***热力参数

表7-3翅片厚度0.8mm时不同烟囱内直径***热力参数

由上测试结果可以看出，相关参数的变化具有规律可寻，本发明的模拟***也具有较高的实用性，可以为实际的核电站施工建筑的各种结构参数进行指导，可以提高工程效率，减少风险，并且本发明的模拟***具有操作简便、计算速度快的优点。

以上实施方式仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施方式对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施方式所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施方式技术方案的精神和范围。

Claims (6)

1.一种钠冷堆非能动余热排出***的建模优化方法，其特征在于，钠冷堆非能动余热排出***包括水平高度依次升高的堆芯、事故热交换器、空气热交换器；所述空气热交换器设于拔风烟囱的内部；所述事故热交换器的壳程为堆芯内待冷却的热介质，管程为冷却介质；所述事故热交换器的冷却介质与空气热交换器的加热介质为形成中间自然循环回路的液钠；所述拔风烟囱设有入风口，所述入风口设于空气热交换器的底部；所述空气热交换器的冷却介质为空气；该方法包括以下步骤：

S1、特定物性参数为特定条件下的固有参数；初步设定热力参数、空气热交换器结构参数、拔风烟囱结构参数，及事故热交换器液钠出口与空气热交换器液钠出口之间的高度差H₂；

所述热力参数为空气热交换器中的各参数，分别包括空气侧进口温度T_e、设定传热量Q、钠侧进口温度T_s,i与钠侧出口温度T_s,o；

所述空气热交换器为液钠走管程，空气走壳程的翅片管式换热器；翅片管式换热器的换热管为翅片管，所述翅片管包括基管及设置在基管上的翅片；空气热交换器结构参数包括基管结构参数和翅片结构参数；所述基管结构参数包括基管总管数n_s、基管外直径d_o、基管内直径d_i、基管的总管长L_s、基管的排列方式与排列的间距，及根据以上数据得到的沿空气流通方向基管的管排数N_a、单位长度基管外表面积A_o ^*，基管内流通截面积S_s，空气热交换器翅片管间最小流通面积S_a,min；所述翅片结构参数包括翅片的翅高H_f、节距y，及根据以上数据得到的单位长度翅片表面积A_f与翅化比β；

所述液钠的中间自然循环回路中，钠侧弯头总数为N_s，

所述拔风烟囱结构参数包括空气换热器的壳程空气出口与拔风烟囱出口高度差H₁、烟道内直径D_i、烟道流通截面积S_c、烟道弯头数N_c；

所述特定物性参数包括钠侧比热cp_s与空气侧比热cp_a；不同温度对应下的空气密度ρ_a；不同温度对应下的液钠密度ρ_s；热耗散q、直管阻力系数f、局部阻力系数ξ、进出口阻力系数ε、污垢热阻r、钠侧沿管程的雷诺数Re、导热系数λ、对数平均温差ΔT_m、粘度μ、贝克来数Pe_s与对数平均温差修正系数φ；

S2、预设空气热交换器的空气侧出口温度T_a,o1，根据设定传热量Q确定空气质量流量m_a；计算得到空气侧浮升力F_a和空气侧总阻力降ΔP_a后；计算空气侧浮升力F_a与空气侧总阻力降ΔP_a的相对误差当相对误差大于1％时，则返回调整预设空气侧出口温度T_a,o1；当相对误差小于等于1％时，则判断空气侧流动平衡；

预设钠侧进口温度T_s,i与钠侧出口温度T_s,o，由设定传热量Q确定液钠的质量流量m_s，计算得到钠侧提升力F_s与钠侧总阻力降ΔP_s；计算钠侧提升力F_s与钠侧总阻力降ΔP_s的相对误差当相对误差大于1％时，则返回调整钠侧进口温度T_s,i与钠侧出口温度T_s,o；当相对误差小于等于1％时，则判断钠侧流动平衡；

S3、待S2中的空气侧流动平衡且钠侧流动平衡后，计算得到空气侧传热系数h_a和钠侧传热系数h_s；进而再求得空气热交换器总传热系数K、对数平均温差ΔT_m与总传热量Q_s；再计算总传热量Q_s与设定传热量Q的相对误差当相对误差大于1％时，则判断***能量不平衡，并返回步骤S2重新计算；当相对误差小于等于1％时，则判断***能量平衡，并输出空气质量流量m_a、钠侧质量流量m_s、空气侧进口温度T_e、空气侧出口温度T_a,o1、钠侧进口温度T_s,i、钠侧出口温度T_s,o、空气热交换器总传热系数K、总传热量Q_s、总传热量Q_s与设定传热量Q的相对误差

2.根据权利要求1所述的钠冷堆非能动余热排出***的建模优化方法，其特征在于：步骤S2所述空气侧浮升力F_a和空气侧总阻力降ΔP_a的计算步骤如下：

空气热交换器中换热的能量平衡方程如下式(1)，钠侧能量损失量与空气侧能量增加量均分别等于设定传热量Q：

m_scp_s(T_s,i-T_s,o)＝m_acp_a(T_a,o1-T_e)＝Q (1)

其中，T_s,i、T_s,o分别为空气热交换器的钠侧进口温度、钠侧出口温度；T_e、T_a,o1分别为空气热交换器的空气侧进口温度、空气侧出口温度；m_s、m_a分别为空气热交换器的钠侧质量流量、空气质量流量；cp_s、cp_a分别为空气热交换器中的钠侧比热、空气侧比热；

拔风烟囱内的空气侧浮升力F_a的计算公式如下式(2)：

其中，H₁为空气换热器的壳程空气出口与拔风烟囱出口高度差；ρ_e为空气热交换器的空气侧进口温度T_e下的空气密度；ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，式(2)中的ρ_a,om与式(3)中的T_a,om对应，并满足公式T_a,om由式(3)求得，其中，T_a,o1为空气热交换器的空气侧出口温度，q_c为拔风烟囱与外部环境之间的热耗散，m_a为空气热交换器的空气质量流量，cp_a为空气侧比热；g为重力加速度；dh为高度h的微分；

空气侧总阻力降ΔP_a的计算公式如下式(4)：

其中，f_a为翅片管外阻力系数，N_a为沿空气流通方向基管的管排数，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积，ρ_a,m为空气热交换器中的空气平均密度；f_c为烟囱内沿程阻力系数，S_c为烟道流通截面积，ρ_a,om为对环境散热后拔风烟囱内平均温度T_a,om下的空气密度，D_i为拔风烟囱烟道内直径；N_c为拔风烟囱烟道弯头数，ξ_c为拔风烟囱局部阻力系数，ε_c为拔风烟囱进出口阻力系数；

步骤S2所述的钠侧提升力F_s与钠侧总阻力降ΔP_s计算步骤如下：事故热交换器液钠出口与空气热交换器液钠出口之间钠侧提升力F_s的计算公式如下式(5)：

其中，H₂为事故热交换器液钠出口与空气热交换器液钠出口之间的高度差；ρ_s,o、ρ_s,i分别为空气热交换器的出口液钠密度、进口液钠密度；g为重力加速度；dh为高度h的微分；

钠侧总阻力降ΔP_s的计算公式如下式(6)：

其中，式(6)中的f_s为钠侧沿程阻力系数，钠侧沿程阻力系数f_s与雷诺数Re_s的关系如式(7)；m_s为钠侧质量流量，S_s为基管内流通截面积；ρ_s,m为自然循环回路中液钠的平均密度；L_sur为自然循环回路中除基管外的剩余段管长；L_s为基管的总管长；d_i为基管内直径；N_s为钠侧弯头总数；ξ_s为钠侧局部阻力系数。

3.根据权利要求1所述的钠冷堆非能动余热排出***的建模优化方法，其特征在于：空气热交换器中的翅片管为三角形错列排布；所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

步骤S3所述空气对流传热系数h_a的计算公式如下式(8-1)：

其中：λ_a为空气导热系数；d_o为基管外直径；m_a为空气质量流量；S_a,min为空气热交换器翅片管间最小流通面积；μ_a为空气粘度；cp_a为空气侧比热；A_o为单位长度翅片管总外表面积；为单位长度基管外表面积；

液钠冷却过程，于均匀热流边界条件下，液钠流动的对流传热系数h_s的计算公式如下式(9)：

其中：λ_s为液钠的导热系数；d_i为基管内直径；Pe_s为贝克来数；

根据式(8-1)和式(9)，空气热交换器总传热系数K的计算公式如下式(10)：

其中：β为翅化比，h_a为空气侧对流传热系数；r_a为空气侧污垢热阻；A_o为单位长度翅片管总外表面积；λ_w为管材导热系数；d_o为基管外直径；d_i基管内直径；r_s为钠侧污垢热阻；η为翅片壁面总效率；h_s为钠侧对流传热系数；

***总传热量Q_s的计算公式如下式(11)：

Q_s＝φKn_sL_sA_oΔT_m (11)

其中，φ为对数平均温差修正系数；K为空气热交换器总传热系数；n_s为基管总管数；L_s为基管的总管长；A_o为单位长度翅片管总外表面积；ΔT_m为对数平均温差。

4.根据权利要求1所述的钠冷堆非能动余热排出***的建模优化方法，其特征在于：所述翅片管为正方形直列排布；所述壳程空气流动方向与管程液钠流动方向垂直；空气为竖向流通，液钠为水平方向流动；

步骤S3所述空气对流传热系数h_a的计算公式如下式(8-2)：

其中：λ_a为空气导热系数；d_o为基管外直径；m_a为空气质量流量；S_a,min为空气热交换器翅片管间最小流通面积；μ_a为空气粘度；cp_a为空气侧比热；A_f单位长度翅片表面积；为单位长度基管外表面积；

液钠冷却过程，于均匀热流边界条件下，液钠流动的对流传热系数h_s的计算公式如下式(9)：

其中：λ_s为液钠的导热系数；d_i为基管内直径；Pe_s为贝克来数；

根据式(8-2)和式(9)，空气热交换器总传热系数K的计算公式如下式(10)：

其中：β为翅化比，h_a为空气侧对流传热系数；r_a为空气侧污垢热阻；A_o为单位长度翅片管总外表面积；λ_w为管材导热系数；d_o为基管外直径；d_i基管内直径；r_s为钠侧污垢热阻；η为翅片壁面总效率；h_s为钠侧对流传热系数；

***总传热量Q_s的计算公式如下式(11)：

Q_s＝φKn_sL_sA_oΔT_m (11)

其中，φ为对数平均温差修正系数；K为空气热交换器总传热系数；n_s为基管总管数；L_s为基管的总管长；A_o为单位长度翅片管总外表面积；ΔT_m为对数平均温差。

5.根据权利要求2所述钠冷堆非能动余热排出***的建模优化方法，其特征在于：

当翅片管为三角形错列排布时，公式(4)中的翅片管外阻力系数f_a为：

当翅片管为正方形直列排布时，公式(4)中的翅片管外阻力系数f_a为：

其中，m_a为空气质量流量，S_a,min为空气热交换器翅片管间最小流通面积；d_o为基管外直径；μ_a为空气粘度；P_t为管束横向管间距；P_l为管束纵向管间距；y为翅片管节距；H_f为翅高。

6.根据权利要求1所述钠冷堆非能动余热排出***的建模优化方法，其特征在于：构建物性参数计算模块、空气侧流动平衡计算模块、液态钠侧流动平衡计算模块与***能量平衡计算模块；所述特定物性参数存于物性参数计算模块内；步骤S2从预设空气侧出口温度至判断空气侧流动平衡的过程属于空气侧流动平衡计算模块；步骤2从预设钠侧进口温度与钠侧出口温度至判断钠侧流动平衡的过程属于液态钠侧流动平衡计算模块；步骤S3的过程属于***能量平衡计算模块；所述物性参数计算模块分别与空气侧流动平衡计算模块、液态钠侧流动平衡计算模块连接；所述空气侧流动平衡计算模块判断空气侧流动平衡后进入***能量平衡计算模块；液态钠侧流动平衡计算模块判断钠侧流动平衡后进入***能量平衡计算模块。