2. 中国核动力研究设计院 成都 610213
2. Nuclear Power Institute of China, Chengdu 610213, China
运行瞬态主要是指核电厂在调试或功率运行期间经常或定期发生的瞬态事件,比如核电厂升温升压、阶跃负荷变化、甩负荷等,它属于典型的Ⅰ类工况,在瞬态过程中不允许触发反应堆保护或专设安全设施。设计必须事先开展运行瞬态分析,以验证核电厂控制系统能够进行有效的调节,并满足运行要求。而对于运行瞬态分析,通常依赖于专门的分析软件,比如用于法国M310核电厂的CATIA2[1]软件、用于美国AP1000核电厂的CENTS[2]软件等。
为了全面实现核电设计自主化,中核集团于2011年启动核电工程软件自主化研发项目[3],其中包括运行瞬态分析软件PANTO (Program for Analysis of Normal Transient and Overpressure) 的研发。该软件针对压水堆核电厂,模拟的范围包括反应堆堆芯、一回路系统、部分二回路系统、控制及保护系统等,它可用于Ⅰ类和Ⅱ类工况瞬态分析、控制系统优化、稳压器及蒸汽发生器安全阀的设计评价等。
本文主要介绍PANTO软件的研发概况,包括物理模型、软件设计、软件编码、软件测试与软件验证等。
1 物理模型作为系统级的设计与分析软件,PANTO涉及的物理模型较多。研发中一方面借鉴了同类软件中成熟可靠的模型,即采用集总参数方法,假设流动和传热为一维过程;另一方面,对于存在两相流的稳压器和蒸汽发生器,分别建立了部件模型。
1.1 一回路冷却剂流量计算模型假设在给定的时间步长内冷却剂密度为常数,并且系统压力保持不变,对一回路建立一维动量
方程:
| $p = \frac{m}{A}\frac{{{\rm{d}}v}}{{{\rm{d}}t}} = \frac{V}{{{A^2}{g_{\rm{c}}}}}\frac{{{\rm{d}}W}}{{{\rm{d}}t}}$ | (1) |
式中:p、m、A、V分别为控制体内的压力、质量、流通面积和自由容积;W为质量流量;v为流体速度;gc为重力加速度。
对于封闭的一回路有:
| ${P_{\rm{H}}} + \oint {{B_{\rm{H}}}} - \oint {{F_{\rm{H}}}} - \oint {{M_{\rm{H}}}} = 0$ | (2) |
式中:PH为泵压头;BH为提升压降;FH为摩擦压降;MH为局部压降。
将一回路划分为j段控制体,总的摩擦压降为:
| $\oint{{{F}_{\text{H}}}}=\sum\limits_{j}{{{{K}_{j}}W_{j}^{2}}/{{{\rho }_{j}}}\;}$ | (3) |
式中:Kj为第j段的阻力系数;Wj为质量流量;ρj为流体密度。
由此,便得到一回路流量计算的控制方程:
| $\frac{\text{d}W}{\text{d}t}=\frac{{{P}_{\text{H}}}+\oint{{{B}_{\text{H}}}}-\sum\limits_{j}{{{{K}_{j}}W_{j}^{2}}/{{{\rho }_{j}}}\;}}{\sum\limits_{j}{\left( {{{V}_{j}}}/{A_{j}^{2}{{g}_{\text{c}}}}\; \right)}}$ | (4) |
除了上述计算模型外,PANTO软件还包括系统流量计算的经验模型,或用户直接输入随时间变化的流量曲线。
1.2 中子动力学模型对于堆芯核功率的计算,采用点堆中子动力学方程,考虑6组缓发中子[4]:
| $\frac{{{\rm{d}}N\left( t \right)}}{{{\rm{d}}t}} = \frac{{\delta K - \beta }}{\Lambda }N\left( t \right) + \sum\limits_{i = 1}^6 {\lambda _i^{}} C_i^{}\left( t \right) + S$ | (5) |
| $\frac{{{\rm{d}}C_i^{}\left( t \right)}}{{{\rm{d}}t}} = \frac{{\beta _i^{}}}{\Lambda }N\left( t \right) - \lambda _i^{}C_i^{}\left( t \right){\kern 1pt} {\kern 1pt} ,{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} i = 1,6$ | (6) |
| $\beta = \sum\limits_{i = 1}^6 {\beta _i^{}} $ | (7) |
式中:N(t) 为中子通量;δK为反应性;β为总的缓发中子份额;βi为第i组缓发中子份额;λi为第i组缓发中子先驱核衰变时间常数;Ci为第i组缓发中子先驱核浓度;S为中子源项;Λ为中子代时间。
1.3 堆芯传热计算模型假设燃料在径向只有一个节点,堆芯换热方程如下[5]:
| ${m_{\rm{f}}}{C_{{\rm{pf}}}}\frac{{{\rm{d}}{T_{\rm{f}}}}}{{{\rm{d}}t}} = kq - {U_{\rm{A}}}({T_{\rm{f}}} - {T_{{\rm{c,w}}}})$ | (8) |
式中:mf为燃料质量;Cpf为燃料的比热容;Tf为燃料的温度;q为总的核功率;k为燃料产生核功率的系数,k=0.974;UA为等效热导;Tc, w为冷却剂温度。
采用有限差分方法,在时间步长t时的燃料温度Tf可以写成Tc, w的函数,假设Tc, w取前一时间步长的值,在第j个轴向节点上有:
| $\begin{align} & {{T}_{\text{f}j}}(t+\Delta t)={{T}_{\text{f}j}}(t){{\text{e}}^{-\Delta t/{{\tau }_{\text{f}j}}}}+ \\ & \quad \quad \quad \quad \quad \left[ \frac{k{{q}_{j}}}{{{U}_{{{\text{A}}_{j}}}}}+{{{\tilde{T}}}_{\text{c},\text{w}j}} \right]\left[ 1-{{\text{e}}^{-\Delta t/{{\tau }_{\text{f}j}}}} \right] \\ \end{align}$ | (9) |
| ${\tilde T_{{\rm{c,w}}j}} = {T_{{\rm{c,w}}j}}(t) + \frac{{{T_{{\rm{c,w}}j}}(t) - {T_{{\rm{c,w}}j}}(t - \Delta t)}}{2}$ | (10) |
式中:
稳压器的模拟采用两区模型[6],并假设气区和液区都是均匀混合的。
1) 压力方程
| $P_{{\rm{pr}}}^{} = P - \frac{k}{{{A_{{\rm{SU}}}}}} \cdot Q_{{\rm{SU}}}^{} \cdot \left| {W_{{\rm{SU}}}^{}} \right|$ | (11) |
式中:Ppr为稳压器压力;P为一回路系统压力;k为波动管摩擦阻力系数;ASU为波动管流通面积;QSU为波动管体积流量,流向稳压器为正方向;WSU为波动管质量流量。
2) 质量守恒方程
每个时间步长结束时稳压器内水和蒸汽的质量分别为:
| $W' = {W^0} + \Delta W$ | (12) |
| $S' = {S^0} + \Delta S$ | (13) |
式中:W0、S0为t-Δt时刻水和蒸汽质量;W'、S'为t时刻的水和蒸汽质量。
液相质量变化为:
| $\Delta W = \left( {\dot S_{\rm{C}}^{} + \dot W_{{\rm{SP,W}}}^{} - \dot W_{{\rm{SU}}}^{} - \dot W_{{\rm{FL}}}^{} - \dot W_{{\rm{SV}}}^{}} \right) \cdot \Delta t$ | (14) |
气相质量变化为:
| $\Delta S = \left( {\dot W_{{\rm{FL}}}^{} - \dot S_{\rm{C}}^{} - \dot S_{{\rm{SV}}}^{} - \dot S_{{\rm{SU}}}^{} + \dot W_{{\rm{SP,S}}}^{}} \right) \cdot \Delta t$ | (15) |
式中:
3) 能量守恒方程
每个时间步长结束时水和蒸汽的比焓分别为:
| $\begin{array}{l} h_{\rm{w}}^\prime {W^\prime } = h_{\rm{w}}^0W + \left( {{{\dot S}_{\rm{C}}}{h_{\rm{F}}} - {{\dot W}_{{\rm{FL}}}}{h_{\rm{G}}} - {{\dot W}_{{\rm{SU}}}}{h_{{\rm{SU}}}}} \right. - \\ \left. {\quad \quad \quad {{\dot W}_{{\rm{SV}}}}{{\bar h}_{\rm{W}}} + {{\dot W}_{{\rm{SP,W}}}}{h_{{\rm{SP}}}} + Q} \right)\Delta t + {{\bar V}_{\rm{W}}}\Delta P \end{array}$ | (16) |
| $\begin{array}{l} h_{\rm{s}}^\prime {S^\prime } = h_{\rm{s}}^0S + \left( {{{\dot W}_{{\rm{FL}}}}{h_{\rm{G}}} - {{\dot S}_{\rm{C}}}{h_{\rm{F}}} - {{\dot S}_{{\rm{SU}}}}{h_{{\rm{SU}}}} - } \right.\\ \left. {\quad \quad \quad {{\dot S}_{{\rm{SV}}}}{{\bar h}_{\rm{S}}} + {{\dot W}_{{\rm{SP,S}}}}{h_{{\rm{SP}}}}} \right)\Delta t + {{\bar V}_{\rm{S}}}\Delta P \end{array}$ | (17) |
式中:hw0、hs0为t-Δt时刻水和蒸汽比焓;hw'、hs'为t时刻的水和蒸汽比焓;hF为水的饱和比焓;hG为蒸汽的饱和比焓;hSU为波动流体的比焓;hSP为喷淋流体的比焓;Q为电加热器的输出功率;VW、VS为水和蒸汽的平均体积;ΔP为压力变化值。
1.5 蒸汽发生器模型蒸汽发生器模型在一次侧包括多个传热管区,二次侧采用点模型,假设里面充满饱和汽液混合物。一次侧每个节点上的换热量为:
| $Q_j^{} = {U_{\rm{A}}} \cdot \left( {T_{{\rm{P}},j}^{} - T_{\rm{S}}^{}} \right)$ | (18) |
式中:Qj为每个节点的热流量;TP, j为一次侧每个节点的温度;TS为二次侧蒸汽饱和温度。
假定二次侧处于饱和状态。若二次侧流体温度已知,则流体压力、比焓与比体积等可通过水的物性计算得到。根据质量守恒方程和能量守恒方程推导二次侧饱和温度变化率,由温度变化率就可以求得每个时刻的二次侧饱和温度。
质量守恒方程:
| $\frac{{{\rm{d}}M}}{{{\rm{d}}t}} = w_{{\rm{fw}}}^{} - w_{\rm{s}}^{}$ | (19) |
能量守恒方程:
| $\frac{\text{d}\left( Mu \right)}{\text{d}t}={{q}_{\text{SG}}}+{{w}_{\text{fw}}}{{h}_{\text{fw}}}-{{w}_{\text{s}}}{{h}_{\text{s}}}$ | (20) |
式中:M为蒸汽发生器流体质量;qSG为进入蒸汽发生器的热流量;wfw、ws为给水流量和蒸汽流量;hfw、hs为给水比焓和蒸汽比焓;u为蒸汽发生器流体平均内能。
联合质量、能量方程可推导得到二次侧饱和温度变化率:
| $\frac{\text{d}T}{\text{d}t}=\frac{{{q}_{\text{SG}}}+{{w}_{\text{fw}}}\left( {{h}_{\text{fw}}}-h_{\text{f}}^{\prime } \right)-{{w}_{\text{s}}}\left( {{h}_{\text{s}}}-h_{\text{f}}^{\prime } \right)}{M\frac{\partial h_{\text{f}}^{\prime }}{\partial T}+V\frac{\partial }{\partial T}\left( \frac{{{h}_{\text{f}}}-{{h}_{\text{g}}}}{{{v}_{\text{f}}}-{{v}_{\text{g}}}}-P \right)+{{C}_{\text{T}}}{{M}_{\text{T}}}}$ | (21) |
式中:V为蒸汽发生器二次侧总体积;hf和hg分别为饱和水与饱和蒸汽比焓;vf与vg分别为饱和水比容与饱和蒸汽比容;P为蒸汽压力。
对于式 (21) 中与物性相关的两个变量:
| $h_{\text{f}}^{\prime }\text{=}\left( {{h}_{\text{f}}}-{{v}_{\text{f}}}\frac{{{h}_{\text{fg}}}}{{{v}_{\text{fg}}}} \right)\text{=}0.992{{T}_{\text{s}}}-30.25$ | (22) |
| $\frac{\partial }{{\partial T}}\left( {\frac{{{h_{\rm{f}}} - {h_{\rm{g}}}}}{{v_{\rm{f}}^{} - {v_{\rm{g}}}}} - P} \right) = 0.033T_{\rm{s}}^{} - 8.96$ | (23) |
可根据水物性参数,建立与饱和温度Ts相关的拟合多项式。
1.6 控制系统模型PANTO软件开发考虑了核电厂一、二回路主要的控制系统,包括[7]:1) 反应堆功率控制系统;2) 反应堆温度控制系统;3) 稳压器压力与液位控制系统;4) 蒸汽发生器液位控制系统;5) 蒸汽排放控制系统;6) 给水流量控制系统。
这些控制系统在运行瞬态中根据设定逻辑调节相关参数使得核电厂趋于稳定运行。目前,PANTO中内置了针对国内M310电厂的控制系统逻辑,瞬态模拟时只需要输入相关的控制参数,比如时间常数、增益系数和目标值等。
2 软件设计PANTO软件总体流程设计如图 1所示。采用模块化设计理念,软件设计上考虑了三个层次。第一层次分为4大模块:计算输入(读入)、数据初始化、核心计算和结果输出;第二层次主要针对核心计算模块进行了细分,包括:一回路热工水力计算、蒸汽发生器相关计算、堆芯中子动力学计算、给水设备热工水力计算、蒸汽排放系统计算、控制保护系统及辅助系统计算;第三层次则是对第二层次各计算模块的再分解,比如对于一回路热工水力计算模块分为:硼浓度计算、冷却剂流量计算、稳压器计算、热段与波动管流量计算、堆芯温度和热流密度计算等。
|
图 1 PANTO软件总体流程图 Figure 1 Flow chart of PANTO. |
另一方面,为了方便用户使用,设计了图形化的输入界面,如图 2所示。采用列表框和标签页相结合的模式,包括8个分页面:1) 名义工况与模型选择;2) 热工水力参数;3) 蒸汽发生器参数;4) 堆芯物理和中子学参数;5) 保护和控制系统参数;6) 蒸汽管线模型;7) 给水设备特性参数;8) 初始工况和瞬态特性数据。
|
图 2 PANTO软件图形化输入界面 Figure 2 Input GUI of PANTO. |
考虑到计算部分和图形化输入界面部分的设计目标不同,PANTO软件计算部分采用了C++语言编写,开发环境为Oracle Solaris Studio。对于输入界面,为了同其它自主化软件保持统一的风格,同时有利于相似模块的直接继承使用,该部分采用Java语言编写,开发环境为NetBeans IDE 7.1.2。PANTO软件代码行总数约为91900行,其中计算部分代码约为52500行,输入界面部分代码约为39400行。
软件编译后运行在中国核动力研究设计院的核反应堆设计计算软件平台之上,操作系统为SUN Solaris 5.10,支撑软件为Java Runtime Environment 6u10或更高版本,硬件环境为SUN SPARC-Enterprise。
4 软件测试根据软件开发规范,必须开展代码级的软件测试,以验证软件实现的功能是否满足需求说明书和设计说明书的要求,代码编写是否满足编码规范。为此,分三个阶段进行了PANTO软件的测试。第一阶段为单元测试,测试各单元代码对相应编码规范的遵循情况;第二阶段为集成测试,验证软件计算结果的正确性和人机界面的正确性;第三阶段为系统测试,对软件进行总体测试。
PANTO软件测试过程中共发现了12个缺陷。其中单元测试过程中发现两个,系统测试过程中发现10个。这些缺陷经过开发人员修复后再进行回归测试,实现遗留缺陷数目为零的目标。最终测试的结果表明,PANTO软件源代码结构清晰,可读性强,编码格式满足规范要求;软件计算结果准确,对异常计算能进行相应错误提示并终止计算过程;软件图形界面友好,具备良好的人机交互能力;软件功能完整,实现了软件需求分析中所列功能点。
5 软件验证以秦山二期核电厂为对象,选取了调试试验中的阶跃负荷增大和额定功率全部甩负荷瞬态进行了验证计算。秦山二期核电厂为两环路设计,额定工况下的主要参数如表 1[8]所示。
| 表 1 秦山二期核电厂主参数 Table 1 Main parameters of Qinshan phase Ⅱ. |
两个瞬态计算所需的边界条件取自调试试验结果,包括一回路流量、汽机第一级压力、汽机进气流量、给水流量和给水温度,如图 3-7所示。对于堆芯中子学参数(如反应性反馈系数、调节棒价值等),由于试验中没有实测结果,验证计算中采用了设计值。
|
图 3 一回路流量 Figure 3 Coolant flow rate of primary loop. |
|
图 4 汽机第一级压力 Figure 4 First stage pressure of turbine. |
|
图 5 汽机进汽流量 Figure 5 Steam flow rate of turbine inlet. |
|
图 6 给水流量 Figure 6 Feed water flowrate. |
|
图 7 给水温度 Figure 7 Feed water temperature. |
阶跃负荷增大瞬态的计算结果对比如图 8-10所示。前200 s为稳态计算,假设瞬态从200 s开始(对应于边界条件的零时刻),瞬态计算360 s。可以看出,由于二回路汽机负荷突然增大,反应堆冷却剂平均温度和稳压器压力略有降低,但同时由于功率控制系统发出调节信号,功率调节棒快速提升,导致堆芯核功率快速升高,反应堆冷却剂平均温度逐渐上升,并达到满功率运行设定值,反应堆最终稳定运行在100%满功率水平。
|
图 8 堆芯核功率 Figure 8 Core nuclear power. |
|
图 9 反应堆冷却剂平均温度 Figure 9 Reactor average coolant temperature. |
|
图 10 稳压器压力 Figure 10 Pressurizer pressure. |
额定功率全部甩负荷瞬态的计算结果对比如图 11-13所示。同样,前200 s为稳态计算,瞬态计算360 s。瞬态开始,由于汽机负荷迅速下降,蒸汽发生器排热量急剧减少,反应堆冷却剂平均温度略有升高,其与设定值间的温差信号将触发蒸汽旁排系
|
图 11 堆芯核功率 Figure 11 Core nuclear power. |
|
图 12 反应堆冷却剂平均温度 Figure 12 Reactor average coolant temperature. |
|
图 13 稳压器压力 Figure 13 Pressurizer pressure. |
统开启,以缓解一回路、二回路间的能量失配。瞬态参数最终将稳定在厂用负荷水平。在320s后,相对于试验结果,PANTO计算得到的堆芯和功率偏大,主要原因是由于在该时刻后功率控制的调节量较小,而计算又采用的是理论的中子学参数,特别是调节棒的价值跟实际值存在一定的差异(比实际值小)。堆芯核功率偏大也导致了反应堆冷却剂平均温度和稳压器压力偏大。
从上述两个瞬态计算结果跟试验数据的对比中可以看出,关键参数的变化趋势符合较好,该软件准确地反映瞬态特征,在参数值上与试验结果的偏差较小,满足工程设计应用的要求。由此验证PANTO软件开发的正确性。
6 结语PANTO软件的物理模型完备,它采用模块化的软件设计理念,应用了面向对象的C++语言和java语言,开发过程规范,并通过软件测试。针对秦山二期核电厂阶跃负荷增大和额定功率全部甩负荷瞬态,完成软件验证分析。计算结果表明,PANTO软件能够较好地模拟瞬态中反应堆堆芯核功率、冷却剂平均温度和稳压器压力的变化,计算精度满足工程应用要求,该软件适用于压水堆核电厂正常运行瞬态分析。
| [1] |
张英, 陈智, 周祖鉴, 等. 用于核电厂数字化仪表控制系统优化设计的CATIA2程序的改进[J].
核动力工程, 2008, 29(1): 19–23.
ZHANG Ying, CHEN Zhi, ZHOU Zujian, et al. Improvement of CATIA2 code used in optimization design for digital instrument and control system of NPP[J]. Nuclear Power Engineering, 2008, 29(1): 19–23. |
| [2] |
刘立欣, 郑利民, 周全福. AP1000核电厂典型的运行瞬态分析[J].
核技术, 2012, 35(11): 869–876.
LIU Lixin, ZHENG Limin, ZHOU Quanfu. Preliminary study on operational transient analysis for AP1000[J]. Nuclear Techniques, 2012, 35(11): 869–876. |
| [3] |
中国核动力院. 中核集团召开核电软件自主化项目启动会[EB/OL]. 2011-03-11[2016-05-13]. http://www.cnnc.com.cn/publish/portal0/tab426/info52484.htm.
Nuclear Power Institute of China. CNNC held the kick-off meeting of self-reliance nuclear power plant design software project[EB/OL]. 2011-03-11[2016-05-13]. http://www.cnnc.com.cn/publish/portal0/tab426/info52484.htm. |
| [4] |
谢仲生.
核反应堆物理分析[M]. 北京: 原子能出版社, 1981: 310-311.
XIE Zhongsheng. Nuclear reactor physics analysis[M]. Beijing: Atomic Energy Press, 1981: 310-311. |
| [5] |
朱继洲.
核反应堆安全分析[M]. 西安: 西安交通大学出版社, 2004: 42-43.
ZHU Jizhou. Nuclear reactor safety analysis[M]. Xi'an: Xi'an Jiaotong University Press, 2004: 42-43. |
| [6] | Abdallah M A, Ahmed H M, Mohammad A R. Pressurizer transients dynamic model[J]. Nuclear Engineering and Design, 1982, 73(3): 447–453. DOI: 10.1016/0029-5493(82)90018-8 |
| [7] |
张森如, 王建渝. 核电厂控制调节系统在运行瞬态分析中的数值计算[J].
核科学与工程, 1988, 8(4): 324–329.
ZHANG Senru, WANG Jianyu. Digital calculation of control and adjustment system in operational transient analysis for nuclear power plant[J]. Chinese Journal of Nuclear Science and Engineering, 1988, 8(4): 324–329. |
| [8] |
李经纬. 秦山核电二期工程反应堆热工水力设计[J].
核动力工程, 1999, 20(4): 308–312.
LI Jingwei. Reactor thermal hydraulic design for Qinshan Phase Ⅱ nuclear power project[J]. Nuclear Power Engineering, 1999, 20(4): 308–312. |