核安全是核电厂的生命线，在核电厂的设计、建造和运营中，核安全理念始终贯彻其中，而验证核电厂的工程设计也成为核电厂安全审评的重要环节。核安全分析软件是核安全审评工作的重要辅助手段和工具。使用与设计单位完全不同的安全分析软件，对拟建核电厂设计方案和运行核电厂改进项目进行独立审核计算，以全面评估核电厂的安全性，是核电发达国家核安全审评工作的标准模式之一。福岛核事故后，中国核安全监管机构对核安全审评工作提出了更高的要求，为了进一步深入研究中国新建核电项目和已装料运行核电机组的安全性，为核安全审评工作提供技术支撑，国家核安全局与美国核安全监管当局(NRC)签署了《中美核反应堆安全联合研究协议》，根据协议，美方转让给中方的核安全审评软件包括SNAP2.2、MELCOR2.1、PARCS、TRACE5等，以上核安全审评软件均由环保部核与辐射安全中心负责应用开发，并将其落实到国内核电项目的独立审核计算工作中。

1 分析工具简介

本研究采用的分析软件是MELCOR 2.1和SNAP 2.2。MELCOR程序是NRC委托Sandia国家实验室开发的一体化严重事故分析程序，应用于安全审评工作中。程序的开发始于1982年，初始的考虑是整合并替代源项分析程序包，消除模块间的各种繁杂的手动耦合，并考虑各种反馈信息(包括温度、释放率、衰变热、热源的重新定位等)影响的综合程序。MELCOR程序通过了大量试验数据的验证和严重事故分析程序的相互验证，对事故进程的模拟具有很高的可信度，在许多国家得到了广泛的应用。MELCOR在程序开发过程中与应用反馈紧密结合，目前仍在不断完善中。当前NRC已不再进一步开发用于安全壳分析的CONTAIN程序，而使用MELCOR程序开展安全壳相关分析模型和事故的验证和评价工作^[1-3]。

SNAP程序是在美国NRC支持下开发的一款图形可视化的核能系统工程分析包，方便工程技术人员创建反应堆系统相关的热工水力和物理模型，该程序具有建模错误检测、输入卡自动生成、多计算任务管理、耦合计算中的数据传递、计算结果曲线生成和事故进程可视化动画演示等功能，当前SNAP可以支持包括MELCOR、TRACE、CONTAIN、COBRA、FRAPCON-3、PARCS、RELAP5等分析程序的应用，为每款程序提供独立的计算插件。

2 分析对象

本文针对AP1000进行冷段双端剪切断裂事故，即大破口失水事故(large break loss of coolant accident, LB LOCA)下的安全壳响应分析，AP1000安全壳为双层结构，内层为圆柱形钢制容器，外层为钢筋混凝土屏蔽构筑物^[4]。图 1为AP1000核电厂安全壳示意图。

3 现象与传热模型 3.1 安全壳内热工水力现象

在一、二回路高能管道破裂的过程中，高温高压的水或蒸汽从破口喷放到安全壳内，由于安全壳内压力较低，高温高压水在较低压力下发生闪蒸，安全壳压力升高^[5]。在喷放的强烈扰动下，喷放介质与安全壳内的不可凝气体相互作用，在安全壳内发生相互夹带的流动，在流动过程中，通过对流换热和蒸汽冷凝将热量传递给安全壳内的非能动吸热构件，这些非能动吸热构件包括安全壳建/构筑物、支撑件、未包裹保温层的系统和部件以及其他设备等，表 1给出了安全壳内非能动热阱涉及的相关物项。

发生假想事故后，非能动安全壳冷却系统(passive containment cooling system，PCS)由安全壳高压力信号触发投入，在钢壳外壁面形成喷淋液膜，通过液膜蒸发排出热量，混凝土壳体与钢壳之间有一环形空气通道，空气在该通道中借助于导流叶片的导向作用，在密度差的驱动下进入环腔下降段，然后折流进入环腔上升段，最终通过高位排气孔排出，从而形成自然循环。在AP1000核电厂安全壳冷却过程中，对流传热、蒸汽冷凝及液膜蒸发是安全壳事故降温降压的重要机制。表 2给出了设计基准事故(design basis accident, DBA)下AP1000安全壳中的重要现象。

3.2 传热模型

从事故后AP1000核电厂安全壳的降温降压特点来看，不论是非能动热阱还是非能动安全壳冷却系统，均通过壁面的对流传热和蒸汽冷凝实现其安全功能，这与采用喷淋系统降压的M310核电厂有很大区别^[6]。AP1000安全壳内热工水力瞬态过程的传热是典型的含有不可凝气体的传热问题。工程实践证明，如果水蒸气中含有1%的空气，就会使冷凝传热系数降低约60%^[7]。当前在包括AP1000在内的非能动安全壳冷却系统的设计中，含有不可凝气体的传热问题是影响安全壳分析的关键问题。国内外一些学者对非能动安全壳内热工水力过程开展了大量的研究工作^[8-12]。在现有的核电厂安全壳分析程序中，所采用的冷凝模型主要包括Uchida关系式、Gido-Koestel关系式、Tagami关系式和基于传热传质相似原理的关系式^[8]。

在MELCOR程序中，针对含不可凝气体的冷凝传热问题，传热过程的计算分为传质计算和对流传热计算2部分。传质计算模型建立在传热传质相似原理关系式基础上，由此确定传质系数，由冷凝质量流率方程确定壁面冷凝质量流率。对于气体对流传热，根据流动工况的不同，可选择自然对流换热、强迫对流换热和混合对流换热三类实验关系式^[13]。

根据传质传热相似原理，

$ {h_D} = \frac{{Sh \cdot D}}{{{L_{\rm{c}}}}}, $

(1)

$ Sh = CN{u^a}S{c^b}P{r^d}, $

(2)

$ Sc = \mu /\left( {\rho D} \right). $

(3)

其中：h_D为传质系数，m/s；D为扩散率，m²/s；L_c为特征长度或表面尺寸，m；μ为壁面和相邻气空间平均温度下的气体动态黏性系数，kg/(ms)；ρ为气体的密度，kg/m³；C、a、b、d为常数。

$ {{\dot m}_{\rm{c}}} = {h_D}{\rho _{\rm{v}}}\ln (\Delta {P_{{\rm{srf}}}}/\Delta {P_{{\rm{atm}}}}). $

(4)

其中：$ {\dot m_{\rm{c}}}$为壁面冷凝的质量流率，kg/(m²·s)；ρ_v为饱和温度T_sat(P_tot)下的蒸汽密度，kg/m³；ΔP_srf=P_tot-P_srf, Pa；ΔP_atm=P_tot-P_stm, Pa；P_tot为相邻空间的总压；P_srf为壁面温度下蒸汽的饱和压力；P_stm为蒸汽分压。

自然对流换热关系式，

$ N{u_{{\rm{natural}}}} = CR{a^m} + D, $

(5)

强迫对流换热关系式，

$ N{u_{{\rm{forced}}}} = CR{e^m}P{r^n} + D, $

(6)

混合对流换热关系式，

$ \begin{array}{l} \;N{u_{{\rm{mix}}}} = \left[ {(R{e^2}/Gr - 1/9} \right]\\ \left[ {N{u_{{\rm{forced}}}} - N{u_{{\rm{natural}}}}} \right] + N{u_{{\rm{natural}}}}. \end{array} $

(7)

其中：C、m、n、D为常数，下标mix表示混合对流，forced表示强迫对流，natural表示自然对流。

另一方面，MELCOR程序中提供了液膜跟踪模型，通过将一系列相连接的模拟壁面传热的热构件组成一个相互连通的网络，从而使液膜能在不同的热构件间流动和分配。液膜跟踪模型可用于模拟AP1000非能动安全壳冷却系统钢壳内外壁面的液膜流动特性。液膜跟踪模型中的液膜厚度由式(8)和(9)迭代确定。

$ {\delta _{\rm{f}}} = \left\{ {\begin{array}{*{20}{l}} {{C_{{\rm{f}}, l}} \cdot {{[{{({\mu _{\rm{f}}}/{\rho _{\rm{f}}})}^2}/\left( {g \cdot \sin \theta } \right)]}^{1/3}} \cdot Re_{\rm{f}}^{e{\rm{f}}, l}, }&{R{e_{\rm{f}}} < R{e_{{\rm{LAM}}}};}\\ {{C_{{\rm{f, }}t}} \cdot {{[{{({\mu _{\rm{f}}}/{\rho _{\rm{f}}})}^2}/\left( {g\cdot\sin \theta } \right)]}^{1/3}} \cdot Re_{\rm{f}}^{e{\rm{f}}, l}, }&{R{e_{\rm{f}}} > R{e_{{\rm{TURB}}}};}\\ {插值计算, }&{R{e_{{\rm{LAM}}}} \le R{e_{\rm{f}}} \le R{e_{{\rm{TURB}}}}.} \end{array}} \right. $

(8)

$ {\delta _{\rm{f}}} = \left[ {{m_{{\rm{f}}, 0}} + ({{\dot m}_{{\rm{in}}}} + {{\dot m}_{\rm{c}}} - {{\dot m}_{{\rm{out}}}}) \cdot \Delta t} \right]/({\rho _{\rm{f}}}{A_{{\rm{sur}}}}), $

(9)

$ R{e_{\rm{f}}} = 2({{\dot m}_{{\rm{in}}}} + {{\dot m}_{{\rm{out}}}})/(w{\mu _{\rm{f}}}). $

(10)

其中：C_{f, x}为常数；x为层流或者湍流；μ_f为液膜的黏性系数，kg/(ms)；ρ_f为液膜的密度，kg/m³；θ为表面倾角；m_{f, 0}为初始液膜质量，kg；A_sur为热构件表面积，m²；$ {\dot m_{{\rm{in}}}}$为从上游热构件流入的液膜质量流率，kg/s；$ {\dot m_{{\rm{out}}}}$为从本热构件流出的液膜质量流率，kg/s；w为热构件表面宽度，m。

4 AP1000安全壳分析模型

AP1000安全壳的模拟需要考虑一回路系统相关设备所占据的空间、隔间、隔间之间的门廊、安全壳内的热构件、安全壳消氢系统、安全壳释放类别、堆腔、堆腔注水系统、非能动安全壳冷却系统、安全壳泄漏等^[14]。

从整体上讲，AP1000安全壳以操作平台为界分为操作平台以下的隔间和操作平台以上的大空间。在分析模型的节点划分上，充分考虑了操作平台以下隔间布置的复杂性和非能动安全壳冷却系统的特性，以反映在事故情况下安全壳内的传质传热过程。对于操作平台以下的隔室空间划分，与安全壳内隔间布置保持一致，对于操作平台以上的大空间，径向划分为4环，轴向划分为8层^[15]。图 2为AP1000安全壳节点划分示意图。

PCS是AP1000核电厂重要的专设安全设施，其最主要的功能是在安全壳内发生冷却剂丧失或主蒸汽管道破裂事故时，导出安全壳内热量，降低安全壳的温度和压力，降低安全壳大气和外部环境的压差，限制事故后放射性物质的释放。图 3为采用液膜跟踪模型模化PCS的节点划分。

5 事故瞬态 5.1 分析假设

本研究选取冷段双端剪切断裂事故下的质能释放数据开展瞬态分析，研究AP1000核电厂安全壳瞬态下的性能。分析模型中假设PCS由安全壳压力高2信号触发，并延迟350 s开启PCS给水阀门，以考虑液膜形成过程的延迟，液膜覆盖率假设为90%。

5.2 瞬态分析结果

在反应堆冷却剂系统发生LOCA事故后，由于安全壳背压相对较低，一回路的高温高压冷却剂向安全壳泄放，并发生闪蒸，导致安全壳的温度和压力急剧上升。

在喷放阶段，安全壳内的非能动热阱由于总面积较大，其表面的蒸汽冷凝传热起到了决定作用(最高热量吸收能力可达到600 MW)，其次是钢制安全壳内壁面的蒸汽冷凝传热，而对流传热相对于冷凝传热是很小的。随着非能动热阱的吸热，热量从热阱表面向中心传递，热阱各区温度逐渐上升，升温速率取决于热阱材料的几何尺寸和物理特性。非能动热阱的吸热能力和喷放阶段的质能释放相比，热阱的吸热能力远远小于质能释放。因此，在事故初期的喷放阶段，主要依赖安全壳的大空间来抑制安全壳压力不超过设计压力。

喷放结束后，堆芯进入再灌水、再淹没和长期冷却阶段，随着一回路冷却剂的泄放以及非能动热阱表面温度的升高，通过非能动热阱的冷凝换热减弱。而随着钢制安全壳壁面温度的升高，PCS的自然循环能力逐渐增强。LOCA事故后约200 s，通过钢制安全壳外的空气自然循环冷却已成为带走壳内热量的主要形式。在各种载热形式的共同作用下，安全壳内的大气温度逐渐降低。

在PCS投入以后，上部贮水箱内的水沿着钢制安全壳外壁面向下流动，形成厚度约0.25 mm左右的液膜，在此过程中伴随着液膜的蒸发，形成的蒸汽随着空气的自然循环流动，从而使得PCS出风口蒸汽份额上升，同时，使得液膜的厚度沿流动方向逐渐减薄直至蒸干，位于钢制安全壳下部的区域在事故初期处于蒸干状态，直至事故后约2 000 s才形成液膜覆盖。事故后期，在PCS的作用下，安全壳压力和温度逐渐降低。

在喷放阶段安全壳压力达到第一峰值压力0.408 MPa，事故后期达到最高峰值压力0.482 MPa，低于安全壳的设计压力0.507 MPa，安全壳大空间峰值温度达到190.2 ℃。在此过程中，安全壳内的蒸汽含量迅速上升，在喷放阶段达到60%，事故过程最高达到74.8%。图 4-10给出了事故瞬态过程中主要参数的计算结果。

5.3 液膜覆盖率的影响

考虑到钢制安全壳外壁面环境条件对液膜覆盖影响的不确定性，选取了不同的液膜覆盖率(0.80、0.85、0.90、0.95、1.0)展开分析，以研究不同液膜覆盖率对安全壳峰值压力和PCS载热能力的影响。表 3给出了不同液膜覆盖率下的安全壳峰值压力。图 11-13给出了不同液膜覆盖率条件下主要参数的变化曲线。结果表明：安全壳峰值压力对液膜覆盖率不敏感，最大差异仅为1.9 kPa，约0.4%；液膜覆盖率越高，PCS载热能力越强，且事故后期安全壳压力下降较快。液膜覆盖率越高，液膜厚度越薄，由于液膜的蒸发使得钢制安全壳外壁面的下部存在的干涸期越长。

5.4 形成有效液膜覆盖时间影响

为了评价形成有效液膜覆盖时间对安全壳峰值压力的影响，假设了不同的PCS给水启动延迟时间(250、300、350、400、450 s)展开分析。表 4为不同PCS给水启动时间下的安全壳峰值压力。图 14-16给出了不同PCS给水启动时间下主要参数的变化曲线。结果表明：形成有效液膜覆盖的时间对安全壳峰值压力有较大的影响，液膜有效覆盖形成的时间越晚，安全壳峰值压力越高，最大差异为7.4 kPa，约1.5%。但过早形成液膜覆盖也会导致事故初期部分冷却液膜未完全蒸发而损失部分带热能力。

6 结论

本文针对AP1000核电厂，采用NRC软件MELCOR和SNAP建立了安全壳及非能动安全壳冷却系统模型，并使用液膜跟踪模型模拟非能动安全壳冷却系统的特性，分析了冷段双端剪切断裂事故下的安全壳热工水力瞬态过程，给出了事故后各阶段的计算结果。与NRC的计算结果^[1]相比，主要计算结果的变化趋势是一致的。由于无法获取NRC的建模数据，二者在非能动热阱、质能释放数据、模型节点划分等方面存在一定的差异，因此二者的计算结果存在一定的合理偏差，但从总体上看，分析结果基本一致，偏差在可接受的范围之内。

此外，还对液膜覆盖率以及形成有效液膜覆盖时间对安全壳峰值压力的影响进行了敏感性分析，结果表明：安全壳峰值压力对液膜覆盖率不敏感，但对事故后期的安全壳降压速率有影响。有效液膜覆盖形成时间对安全壳峰值压力影响较大，时间越晚，安全壳峰值压力越高，但过早形成液膜覆盖也会导致事故初期部分冷却液膜未完全蒸发而损失部分带热能力。

MELCOR程序能很好地模拟非能动安全壳冷却系统的热工水力现象，该研究成果可为后续应用程序分析核电厂安全壳系统响应特性提供参考和借鉴。

环保部核与辐射安全中心通过国外合作，建立了一套与设计方完全不同的热工程序系统，能够在一定范围内开展核电厂设计基准事故独立审核计算工作，深入研究和确认国内新建核电机组和已投运核电机组的安全性能。