第44卷第4期2021年4月
V ol.44，No.4
April2021核技术
NUCLEAR TECHNIQUES
仿真研究
朱光昱1,2闵金坤1,3张丽1元一单1
1（中国核电工程有限公司核电安全研究中心北京100840）
3（清华大学工程物理系北京100084）
摘要基于COMSOL Multiphysics软件的非等温流动模块建立了一种适用于熔融池的仿真计算模型。采用相 网络谣言变材料模型模拟实验中的非共晶摩尔比例20%NaNO
3-80%KNO
3
熔融物，通过将实验台架顶端封闭腔的辐射换
热转化为表面对外界环境辐射边界条件，对COPRA（COrium Pool Research Apparatus）实验的湍流场和传热相变过程进行了仿真模拟。仿真结果表明：熔融池内除了存在明显的自然对流外还存在大量的涡流，且在稳态下熔融池内会出现明显的热分层结构。受自然对流和重力影响，冷却壁面附近熔融物结壳的厚度沿冷却壁面随熔融池圆弧壁面径向角度增加而减小。
关键词熔融池，自然对流，热分层，COMSOL
中图分类号TL334
DOI:10.11889/j.0253-3219.2021.hjs.44.040601
Simulation study on COPRA corium pool based on COMSOL
ZHU Guangyu1,2MIN Jinkun1,3ZHANG Li1YUAN Yidan1
1(Center for Nuclear Power Safety,China Nuclear Power Engineering Co.,Ltd.,Beijing100840,China)
2(Nuclear and Radiation Safety Center,MEE,Beijing100082,China)
3(Department of Engineering Physics,Tsinghua University,Beijing100084,China)
Abstract[Background]COMSOL is considered to be an ideal approach to study the complicated physical
processes in the corium pool because of its good simulation ability in multiphysics field.[Purpose]This study aims to build a computational model by COMSOL to explore the turbulent flow field and phase change in the corium pool.
[Methods]Based on non-isothermal flow calculation module,the non-eutectic binary mixture of20%NaNO
3
-
80%KNO
3
used in COrium Pool Research Apparatus(COPRA)experiment was simulated by using phase-change material model.The radiation heat transfer in the closed cavity at the top of corium pool was translated into the radiation heat transfer between the upper surface and the environment available as boundary condition in COMSOL.
[Results]Simulation results show that in addition to the evident natural convection main flow,a mass of vortexes are existed in the corium pool,which lead to the thermal stratification of the corium pool in steady-state.Under the influence of natural convection main flow and gravity,the crust thickness along the cooling wall surface decreases with the vessel polar angle of corium pool increases.[Conclusions]The current computational model is helpful for severe accident mitigation system design.
Key words Corium pool,Natural convection,Thermal stratification,COMSOL
国家重点研发计划(No.2018YFB1900100)资助
第一作者：朱光昱，男，1989年出生，2015年毕业于哈尔滨工程大学，现为硕士研究生，研究领域为反应堆严重事故
收稿日期：2020-10-20，修回日期：2020-12-16
菊粉酶Supported by National Key R&D Program of China(No.2018YFB1900100)
First author:ZHU Guangyu,male,born in1989,graduated from Harbin Engineering University in2015,master student,focusing on reactor severe accident
Received date:2020-10-20,revised date:2020-12-16
朱光昱等：基于COMSOL的COPRA熔融池实验仿真研究
熔融物堆内滞留技术（In-vessel Retention，IVR）
是目前广泛采用的核反应堆严重事故堆芯熔融物滞
留方案，其有效性取决于反应堆压力容器下封头壁
面是否可以有效导出熔融物的衰变热。为了研究反
应堆严重事故后压力容器内熔融池与冷却壁面的换
热特性，国内外设计了许多压力容器下封头二维切
片实验台架，例如COPO［1］（COrium Pool avec
formation de glace aux parois）、BALI［2］（BAin
Liquide）、SIMECO［3］（Simulation of In-vessel MElt
COolability）以及由中国核电工程有限公司和西安
交通大学针对华龙一号堆型搭建的COPRA实验台
架［4−5］（COrium Pool Research Apparatus）。其中，
COPRA实验采用了穆尔比为20%NaNO
3-80%KNO
3
的二元混合物作为实验工质，该二元混合物与反应
堆内真实熔融物UO
2-ZrO
2
的相图存在相似性［4］，可
以较为准确地模拟下封头熔融池内的熔融物相变行为和高瑞利数自然对流工况，为IVR相关系统设计提供了必要的技术储备。由于熔融池内的物理变化过程十分复杂，限制了COPRA实验中测量手段，需要建立仿真计算模型来弥补实验的不足。对此，本文基于COMSOL Multiphysics多物理场计算软件建立了一种熔融池计算模型，通过耦合湍流场和传热相变过程，实现了对COPRA实验熔融池内的流动传热、相变过程的仿真模拟。
1仿真模型
1.1数学物理模型
采用COMSOL软件中的非等温传热模块耦合实验中的流体传热、相变和湍流流动。穆尔比为
20%NaNO
3-80%KNO
3
的二元混合物具有较大固液
相线温差，适宜采用流体传热模块中的相变材料模型进行模拟，其原理如图1所示，温度低于固相线区域采用固相物性计算，温度高于液相线区域采用液相物性计算，在固液相线之间采用固液两相按体积平均的物性参数进行计算。图1中θ为相体积分数；ΔT为固液相线温差；T
pc
为COMSOL中的相变温度，这里设置为固相线温度和液相线温度的均值。参考Gaus-Liu［6］和Janz［7］的研究结果，20%NaNO
3
-
80%KNO
3
的二元混合物的固相线温度为497.15K，液相线温度为557.15K，由固相线到液相线的相变潜热为161960J·kg−1，其他主要物性参数如表1所示。
湍流流动模块采用低雷诺数k-ε模型进行计算，参考相变材料模型设置思路，该模块中的粘度项设置为固液两项的体积平均值，可通过设置较大的粘度来模拟固相不存在流动性的情况。通过预计算对比，当固相的粘度系数高于103Pa·s后对计算的影响很小。
求解过程采用分离式求解器，时间步进选择向后差分，计算容差因子为0.5，计算加速性和稳定性设置为Anderson加速度。
1.2计算域设置
图2为参考COPRA实验台架的切片结构建立的二维计算模型。台架整体呈半径为2.2m的1/4圆形，加注熔盐后台架可分为底部的熔融池区域和顶部的辐射换热区域，熔盐加注的最大高度为1.9m。其中1/4圆形半径上的边界（c、d和f）为绝热边界，圆弧边界（a和b）外侧为水冷流道，其平均温度约维持在323.15K，可以认为是等温边界条件。
顶部的辐射换热区域可认为是4个表面组成的封闭腔辐射传热问题，其辐射换热网格图如图3所示，其中绝热壁面c、d可认为是重辐射面，则熔融池
上表面（e）向水冷壁面辐射热流Φ
eb
可通过式（1）
计算：
图1相变模型原理图
Fig.1Schematic diagram of phase change model
表120%NaNO
3-80%KNO
3
混合物物性
Table1Property of20%NaNO
3-80%KNO
3
温度Temperature/K ＜497.15 557.15
573.15
623.15密度
Density/kg·m−3
2200
1914
1902
1866
比热容
Specific heat/J·kg−1·K−1sfv
674.25+1.13T
1328
1332
1346
导热系数
Heat conductivity coefficient/W·m−1·K−1
0.5
采油速度0.6
0.439
0.422
动力粘度
Viscosity/Pa·s
—
0.00377
0.00333
0.00252
核技术2021,44:040601
Φeb =
E be -E bb
∑R t
(1)
图3和式（1）中：E bb 、E bc 、E bd 、E be 分别为各表面上的黑体辐射热流密度；J b 、J c 、J d 、J e 分别为各表面上的有效辐射热流密度；ε1和ε2为熔融池上表面和台架内
不锈钢壁面的发射率，参考张卢腾［4］
的成果分别取0.44和0.8；A b 、A c 、A d 、A e 分别为各表面的面积，对于二维模型，可通过将切片台架的厚度进行归一化计算各面的面积；X e ，b 、X e ，d 、X e ，c 、X b ，d 、X d ，c 、X b ，c 为表面e 和b 、e 和d 、e 和c 、b 和d 、d 和c 、b 和c 间的角系数，可将
封闭腔简化为长方形结构通过查表［8］
获得；R t 为图3中所示的各辐射热阻，总等效热阻需根据各热阻串并联关系计算，具体计算方法在文献［8］中有详细介绍。
通过引入黑体辐射公式可将式（1）化为式（2）的形式，此时可将总等效热阻的倒数作为一个等效表面发射率考虑，则辐射换热区域可简化为熔融池上表面对外界环境的辐射换热过程，因此仅对由a 、e 和f 三个边界围成的熔融池区域的建立计算域即
可。其中，将a 设置为边界恒温条件，e 设置为COMSOL 软件中的表面对环境辐射换热边界条件，f 设置为绝热边界条件。
Φeb =1
∑R t
σ(T 4e -T 4b )(2)
式中：σ为黑体辐射常数；T e 为熔融池上表面的温度；
T b 为水冷壁面的温度，这里设置为“表面对环境辐射换热”边界条件中的外界环境温度。
1.3网格敏感性分析
根据网格敏感性分析，当网格数大于10924后，网格变化对计算得到的温度场分布和熔融物结壳厚度变化率降低至1%以内，所以本文仅选取流场中的最大计算流速作为分析依据。图4为加热功率为10kW 时的网格敏感性分析结果，随着网格加密，计算得到的流场中最大流速逐渐增加，在网格数加密至15462后网格宽度对计算得到的最大速度变化率约为2%，因此最终计算采用的网格数为15462。计算采用的初始时间步长为1ms ，由于COMSOL 软件会自动根据局部误差对计算时间步长进行调整，因此未对时间步长的敏感性进行分析。
2计算结果分析
COPRA 台架通过在熔融池中的电加热器模拟
熔融物衰变热，对于实验中的T1工况［4］
，加热功率按照A 阶段（功率15kW ，持续时间23700s ）、B 阶段（功率10kW ，持续时间40578s ）、C 阶段（功率15kW ，持续至稳态）变化。计算过程中，参照文献［5］中的实验台架设计信息，将不同的加热功率折合成体积热源并采用阶梯函数控制，从而实现对实验过程的模拟。图5为A 阶段熔融池的计算温度稳态分布，与实验结果类似，在熔融池内的温度总体呈现随高度增加逐渐升高的趋势，
且不同高度出现了明
建湖实验小学图2计算域和网格划分示意
Fig.2Schematic diagram of calculation domain and mesh
generation
图3辐射传热等效网络图
Fig.3
Equivalent network diagram of radiation heat
transfer
图4不同网格数下的最大计算速度
Fig.4
Maximum simulation velocity under different mesh
number
朱光昱等：基于COMSOL 的COPRA 熔融池实验仿真研究
显的热分层现象，在底部结壳的区域的温度梯度较大，中部和顶部的温度梯度较小。
图6为各阶段稳态下的熔融池温度分布的实验值和计算值。其中，底部区域的计算偏差相对较大，这可能是由于当前计算模型没有引入在熔融物凝固
时会在硬壳和下封头间形成气体间隙导致的［4］
。对于熔融池高度0.5m 以上区域，计算得到的温度分布与实验值十分接近，最大计算误差不超过1.5%。
图7（a ）为B 阶段的流场，熔融池内形成了由冷却壁面向下流动，在中部受热后在浮力作用下缓慢
向熔融池顶部和冷却壁面流动的自然对流趋势。受
自然对流和顶部辐射换热的影响，在熔融池内存在大量的涡流，这些涡流搅混的区域温度分布趋于均匀，因此易于在熔融池内形成热分层结构。图7（b ）为B 阶段的固态体积分数，在冷却壁面附近，受高速向下主流和重力的共同作用，凝固形成的糊状物会更容易向熔融池底部转移并最终结壳，从而使熔融物结壳的厚度呈现沿冷却壁面由熔融池底部到顶部逐渐降低的趋势。
基因敲除小鼠3结语
采用COMSOL Multiphysics 软件的非等温流动模块、相变材料模型建立的熔融池计算模型对COPRA 实验进行了仿真模拟，结果表明：
1）稳态熔融池内自然对流形式为：由冷却壁面向下流动，在中部受热后在浮力作用下缓慢向熔融池顶部和冷却壁面流动。
2）在自然对流和顶部辐射换热的共同影响下，熔融池内会形成了大量的涡流。
3）在冷却壁面向下主流和重力作用下，熔融物凝固硬壳呈现沿冷却壁面由熔融池底部到顶部逐渐
变薄的趋势。参考文献
1
Helle M,Kymäläinen O,Tuomisto H.Experimental COPO II data on natural convection in homogenous and stratified pools[C].The 9th International Topical Meeting on Nuclear Reactor Thermal-Hydraulics (NURETH-9),San Francisco,USA,1999.2
Bonnet J M,Seiler J M.Thermal hydraulic phenomena in corium pools:the BALI experiment[C].The 7th International
Conference
on
Nuclear
Engineering
(ICONE-7),Tokyo,Japan,
1999.
图5A 阶段熔融池的热分层结构
Fig.5
Thermal stratification of corium pool in stage
A
图6各阶段熔融池内温度稳态分布
Fig.6Steady temperature profile in different stages
of corium
pool
图7B 阶段熔融池内速度场(a)和固相体积分数(b)
Fig.7Velocity field (a)and soild volume fraction (b)of corium pool in stage B
核技术2021,44:040601
3Sehgal B R,Bui V A,Dinh T N,et al.SIMECO experiments on in-vessel melt pool formation and heat transfer with and without a metallic layer[C].Proceedings of In-Vessel Core Debris Retention and Coolability Workshop,Garching,Germany,1999.
4张卢腾.反应堆严重事故下高瑞利数熔融池自然对流传热特性研究[D].西安：西安交通大学,2017.
ZHANG Luteng.Research on the natural convection heat transfer characteristics in corium pool with high rayleigh number during reactor severe accident[D].Xi'an:Xi'an Jiaotong University,2017.
5Zhang L T,Zhang Y P,Zhao B,et al.COPRA:a large scale experiment on natural convection heat transfer in corium pools with internal heating[J].Progress in Nuclear
Energy,2016,86:132−140.DOI:10.1016/j.pnucene.
2015.10.006.
6Gaus-Liu X,Mjasoedov A,Cron T,et al.Test and simulation results of LIVE-L4+LIVE-L5L[R].Germany: Karlsruhe Institute of Technology,2011.
7Janz G J,Allen C B,Bansal N P,et al.Physical properties data compilations relevant to energy storage II molten salts:data on single and multi-component salt systems[R].
USA:Rensselaer Polytechnic Institute Cogswell Laboratory,1979.
8杨世铭，陶文铨.传热学[M].北京：高等教育出版社, 2007.
YANG Shiming,TAO Wenquan.Heat transfer[M].
Beijing:Higher Education Press,2007.
豫公网安备41160202000603
站长QQ:729038198
关于我们
投诉建议