2. 中国科学院核能安全技术研究所 中子输运理论与辐射安全重点实验室 合肥 230031
2. Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China
根据国际原子能机构(International Atomic Energy Agency,IAEA)的定义，非能动系统是不需要任何输入即可运行的系统[1]。非能动设计已广泛用于先进反应堆的设计中，非能动余热排出系统受到了第四代核能系统国际论坛(Generation IX International Forum,GIF)推荐的包括铅冷快堆在内的6种第四代新堆型的青睐[2]，如美国的先进嬗变快堆(Advanced Breeder Reactor,ABR)、铅冷快堆(Small,Sealed,Transportable,Autonomous Reactor,SSTAR)、先进液态金属堆(Advanced Liquid Metal Reactor,ALMR)，欧盟的实验加速器驱动系统(The Xperimental Accelerator Driven System,XADS)、高科技应用多功能混合动力研究堆(Multi-purpose Ybrid Research Reactor for High-tech Applications,MYRRHA)，俄罗斯的模块式铅铋快堆(Svintsovo Vismutny Bystryi Reactor,SVBR)等均采用了非能动余热排出系统[3, 4]。虽然非能动设计可以避免能动失效和人因失效，但非能动系统功能受环境因素和材料参数影响大，存在失效的可能。为了更客观地评价非能动系统，需要对非能动系统的可靠性进行定量化分析[5]。
非能动系统的失效分为两类：一类是硬件失效；另一类是设计功能的不确定性，即功能失效[5, 6]。近年来，国内外对非能动安全系统可靠性评价进行了大量研究。欧洲核能机构(European Nuclear Energy Agency,ENEA)提出非能动安全系统可靠性评估方法(Reliability Evaluation of Passive Safety System,REPAS)，并应用于两相流自然循环系统[7]；法国核能委员会开发出非能动安全功能可靠性方法(Reliability Methods for Passive Safety Functions,RMPS)[8]；国内的清华大学分别采用蒙特卡罗模拟方法[9]和响应面法[10]估计高温气冷堆HTR-10的余热排除系统失效概率；中国科学院核能安全技术研究所的FDS团队采用人工神经网络方法估计中国铅基研究实验堆(China lead-based research reactor,CLEAR-I)容器空气冷却系统(Reactor Vessel Air Cooling System,RVACS)的可靠性[11]。
由于人工神经网络法存在过拟合和欠拟合问题，计算结果存在不确定性，需要建立其他模型进行验证。对于非能动系统的可靠性分析，目前常用的分析方法有一次二阶矩法、响应面法和蒙特卡罗模拟[12]。因一次二阶矩法不能求解隐式的性能函数，而非能动系统的性能函数通常是模拟程序，所以一次二阶矩法不能直接应用于非能动系统的可靠性分析中。响应面法通过近似构造一个具有明确表达形式的多项式来表达这种隐式性能函数[13]，然后用插值方法来确定响应面性能函数的系数。由于非能动系统的模拟程序计算耗时一般较长，此方法当随机变量个数较大时，运用模拟程序次数多，计算代价大，且该方法的误差无法估计[10]。蒙特卡罗模拟需要进行多次抽样，且每次抽样后都需要运用模拟程序，同样存在计算代价很大的问题。
本文提出一种响应面拟合方法用于分析中国铅基研究实验堆CLEAR-I的RVACS的可靠性，并验证文献[11]的结论。响应面拟合方法采用构造响应面性能函数表达反应堆非能动安全系统的数学模型，通过增加失效面附近的抽样频率提高失效面的拟合精度，利用最小二乘法和bootstrap方法求解响应面性能函数，使用得到的响应面模型代替Fluent模型估计RVACS系统的非能动失效概率。
1 分析方法影响系统可靠性的不确定性输入参数与输出的可靠性参数之间的作用关系用数学模型进行描述，其输入向量为由n个不确定输入变量组成的X，输出量为Y，则其函数映射关系为：
| Y=y(X)，X=X1,X2,…,Xn | (1) |
假设系统失效的失效准则是Y>A或Y<A，A为定义的可靠性参数阈值。为方便讨论，取前者，定义性能函数为：
| g(x)=A-y(X) | (2) |
根据Weierstress多项式最佳逼近定理，任何类型的复杂函数都可以用多项式逼近。在实际应用中，多采用二次多项式的响应面函数来近似表达隐式性能函数：
| \[\begin{align} & G=g(X)=g({{X}_{1}},{{X}_{2}},\ldots ,{{X}_{n}}) \\ & =a+\sum\nolimits_{i=1}^{n}{{{b}_{i}}{{X}_{i}}}+{{\sum\nolimits_{i=1}^{n}{{{c}_{i}}{{X}_{i}}}}^{2}}+ \\ & _{{}}^{{}}\sum\nolimits_{i=1}^{n}{{}}\sum\nolimits_{j>i}^{n}{{{d}_{ij}}{{X}_{i}}}{{X}_{j}} \end{align}\] | (3) |
其中需要求解的系数个数为(n2+3n)/2+1。当不确定性输入变量个数多时，求解系数太多，因而需要的数据样本太多，需要对响应面形式进行简化，通常舍去交叉项，得二次多项式响应面函数：
\[\begin{align}
& G=g(X)=g({{X}_{1}},{{X}_{2}},\ldots ,{{X}_{n}}) \\
& =a+\sum\nolimits_{i=1}^{n}{{{b}_{i}}{{X}_{i}}}+{{\sum\nolimits_{i=1}^{n}{{{c}_{i}}{{X}_{i}}}}^{2}}+ \\
& _{{}}^{{}}\sum\nolimits_{i=1}^{n}{{}}\sum\nolimits_{j>i}^{n}{{{d}_{ij}}{{X}_{i}}}{{X}_{j}}
\end{align}\]
(4)
若用响应面法直接求解，需要运用模拟程序的k(2n+1)次，当迭代次数k较大时，计算代价很大。为减少计算代价，本文采用多维拟合方法求解响应面方法构造出的性能函数（式(4)），用响应面性能函数代替模拟程序。响应面拟合方法流程图如图1所示。
|
图1 响应面拟合方法的流程图 Fig.1 Flowchart of response-surface-fitting method. |
具体实施步骤如下：
1) 获取样本点：为使拟合的响应面在失效面局部更加准确，本文中的样本点包括：第一部分样本点采用蒙特卡罗抽样方法在不确定性输入参数的分布区间中抽取得到输入向量，带入系统模拟程序计算得到相应的输出量。由于物理状态的连续性，失效样本点附近出现其它失效点的概率较高。第二部分样本点是通过在失效样本点附近进行蒙特卡罗抽样，并由系统模拟程序计算得到相应的输出量。第二部分样本点中有较多的失效样本点。
2) 估计响应面性能函数的系数：为减小响应面性能函数的系数估计的方差，本文采用bootstrap方法从步骤1)获得的样本点中产生自举样本作为求解响应面性能函数（式(4)）的输入参数，利用MATLAB中自带的最小二乘法多维曲面拟合回归函数LSQCURVEFIT计算出每组响应面性能函数的系数，并以不包含在自举样本中的第二部分样本点作为测试样本点，计算测试样本点的均方根误差(Root Mean Square Error,RMSE)。抽样次数为t (t >1000)次，得到t组数据。删除测试样本点均方根误差最大的i (i<t)组，以剩下t-i组响应面性能函数的系数的均值为响应面性能函数的系数，确定响应面性能函数。
3) 非能动系统可靠性分析：从不确定性输入量的分布中进行蒙特卡罗随机抽样，得q组输入向量，代入步骤2)得到的响应面性能函数中计算，得到q组输出量。对q组输出量进行统计和计算可得到非能动系统的失效概率。
响应面拟合方法中的样本点个数m即是该方法需要运用模拟程序的次数。m一般小于k(2n+1)，响应面拟合方法需要运用模拟程序的次数比响应面方法少。所以相较于响应面方法，响应面拟合方法可在消耗较小的计算代价前提下完成非能动系统的可靠性分析。
2 RVACS系统可靠性分析 2.1 RVACS系统模型中国铅基研究实验堆CLEAR-I是FDS团队设计的一种铅铋冷却研究堆，其RVACS系统采用非能动自然循环设计。详细系统模型介绍参见文献[14, 15, 16]。
RVACS系统包含4组并联排热管道，设计分析表明三组排热管正常运行即可保证余热可靠排出[16, 17]。本文考虑在运行过程中可能发生的更加严峻的一种情况：一组维修，一组失效，分析两组排热运行时，输入参数不确定性对系统余热排出可靠性的影响。
RVACS系统的成功准则为：主容器壁面温度和安全容器壁面平均温度均不超过设计值450 ºC和400ºC[17]。但考虑物理模型近似带来的误差，定义失效准则：主容器壁面平均温度超过435 ºC，安全容器壁面平均温度超过385 ºC。
文献[18]确定了RVACS系统的6个不确定性输入变量：堆芯余热功率(t=0 h)、外部空气温度、等效摩擦压降系数、主容器壁面发射率、安全容器壁面发射率、圆柱形隔热层壁面发射率。各参数具体分布类型及分布区间见表1。采用Fluent软件建立RVACS系统物理模型。
|
表1 不确定参数具体分布类型及区间 Table 1 Distribution pattern and distribution interval of uncertainty parameters. |
近似RVACS系统的Fluent模型的响应面性能函数为：
\[G=a+\sum\nolimits_{i=1}^{6}{{{b}_{i}}{{X}_{i}}}+{{\sum\nolimits_{i=1}^{6}{{{c}_{i}}{{X}_{i}}}}^{2}}\]
(5)
用于拟合响应面性能函数的样本点包括：第一部分采用蒙特卡罗抽样方法在6个不确定性输入参数的分布区间中抽取得到100组输入向量，带入RVACS系统的Fluent模型计算得到100组输出量。第二部分样本点是以其中的失效样本点（主容器壁面平均温度超过435 ºC或安全容器壁面温度超过385 ºC）为圆心，各区间波动范围的10%为半径共抽取50组输入向量，带入Fluent模型中计算出50组输出量。总共得到150组样本点。
根据§1步骤2)，利用150组样本点求解响应面性能函数的系数，确定响应面性能函数。
用响应面性能函数代替RVACS系统的Fluent模型进行可靠性分析。从不确定性输入量的分布中(表1)进行蒙特卡罗随机抽样，得10000组输入向量，代入响应面性能函数，得到10000组输出量。对10000组输出量进行统计可得到RVACS系统的失效概率。计算结果如图2、表2所示。
|
图2 10000组样本对应的反应堆主容器(a)和安全容器(b)壁面平均温度 Fig.2 Mean temperature of main vessel (a) and safety vessel (b) of 10000 samples. |
|
|
表2 主容器与安全容器平均温度的累积概率 Table 2 Cumulative probability of mean temperature of main vessel and safety vessel. |
分析结果得出：1) 只有两组余热排出系统成功开启时，99.38%概率下RVACS能够成功将余热导出；2) 共有62组样本对应的工况造成RVACS系统失效，其中只有两组样本对应的工况造成安全容器失效，这两组样本对应的工况也造成主容器失效。主容器温度先于安全容器达到极限温度；3) 主容器壁面平均温度的置信水平1-α为90%的上下限置信区间为[397,429]，不确定性分布区间较大。
3 结语本文提出一种用于非能动系统可靠性分析的响应面拟合方法，相较于响应面方法，该方法可在消耗较小的计算代价前提下完成非能动系统的可靠性分析。将该方法应用于中国铅基研究实验堆反应堆容器空气冷却系统的可靠性分析，分析结果表明，非能动系统RVACS具有较大的冗余性，在所有能动余热排除系统不可用的情况下，RVACS四组并联排热管中的两组也能够可靠地导出反应堆余热。
致谢 感谢中广核研究院有限公司的夏少雄同志对本工作的支持，同时非常感谢核能安全技术研究所FDS团队其他成员对本工作的大力支持。
| 1 | Atomique E. Safety related terms for advanced nuclear plants[R]. International Atomic Energy Agency, 1991( 1) |
| 2 | 郑明光, 叶成, 韩旭. 新能源中的核电发展[J]. 核技术, 2010, 33(2):81-86 ZHENG Mingguang, YE Cheng, HAN Xu. The development of nuclear power as an alternative energy[J]. Nuclear Techniques, 2010, 33(2):81-86(1) |
| 3 | Maes D. Mechanical design of the small-scale experimental ADS:MYRRHA[J]. Energy Conversion and Management, 2006, 47(17):2710-2723(1) |
| 4 | Zrodnikov A V, Toshinsky G I, Komlev O G, et al. Nuclear power development in market conditions with use of multi-purpose modular fast reactors SVBR-75/100[J]. Nuclear Engineering and Design, 2006, 236(14-16):1490-1502(1) |
| 5 | Mathews T S, Ramakrishnan M, Parthasarathy U, et al. Functional reliability analysis of safety grade decay heat removal system of Indian 500 MWe PFBR[J]. Nuclear Engineering and Design, 2008, 238(9):2369-2376(1) |
| 6 | Burgazzi L. Reliability evaluation of passive systems through functional reliability assessment[J]. Nuclear Technology, 2003, 144(2):145-151(1) |
| 7 | Ricotti M E, Bianchi F, Burgazzi L, et al. The REPAS study:reliability evaluation of passive safety systems[C]. 10th International Conference on Nuclear Engineering, Arlington, 2002(1) |
| 8 | Marques M, Pignatel J F, Saignes P, et al. Methodology for the reliability evaluation of a passive system and its integration into a probabilistic safety assessment[J]. Nuclear Engineering and Design, 2005, 235(24):2612-2631(1) |
| 9 | 谢国锋, 童节娟, 何旭洪, 等. 用Monte Carlo方法计算HTR-10余热排出系统物理过程的失效概率[J]. 核动力工程, 2008, 29(2):85-87 XIE Guofeng, TONG Jiejuan, HE Xuhong, et al. Calculation of physical failure probability of HTR-10 residual heat removal system by Monte Carlo method[J]. Nuclear Power Engineering, 2008, 29(2):85-87(1) |
| 10 | 谢国锋, 何旭洪, 童节娟, 等. 响应面方法计算HTR-10余热排出系统物理过程的失效概率[J]. 物理学报, 2007, 56(6):3192-3197 XIE Guofeng, HE Xuhong, TONG Jiejuan, et al. Calculating physical failure probability of HTR-10's residual heat removal system by response surface method[J]. Acta Physica Sinica, 2007, 56(6):3192-3197(2) |
| 11 | 夏少雄, 王家群, 潘晓磊, 等. 中国铅基研究实验堆非能动余热排出系统可靠性分析[J]. 核技术, 2015, 38(2):020605. DOI:10.11889/j.0253-3219.2015.hjs.38.020605 XIA Shaoxiong, WANG Jiaqun, PAN Xiaolei, et al. Reliability analysis of passive decay heat removal system of China lead-based research reactor[J]. Nuclear Techniques, 2015, 38(2):020605. DOI:10.11889/j.0253-3219.2015.hjs.38.020605(2) |
| 12 | 赵国藩. 工程结构可靠性理论与应用[M]. 大连:大连理工大学出版社, 1996 ZHAO Guofan. Reliability theory and application for engineering structure[M]. Dalian:Dalian University of Technology Press, 1996(1) |
| 13 | 王永菲, 王成国. 响应面法的理论与应用[J]. 中央民族大学学报(自然科学版), 2005, 14(3):236-240 WANG Yongfei, WANG Chengguo. The application of response surface methodology[J]. Journal of Minzu University of China(Natural Sciences Edition), 2005, 14(3):236-240(1) |
| 14 | 吴宜灿, 柏云清, 宋勇. 中国铅基研究反应堆概念设计研究[J]. 核科学与工程, 2014, 34(2):201-208 WU Yican, BAI Yunqing, SONG Yong. Conceptual design of China lead-based research reactor CLEAR-I[J]. Nuclear Science and Engineering, 2014, 34(2):201-208(1) |
| 15 | Wu Y C, Bai Y Q, Song Y, et al. Development strategy and conceptual design of China lead-based research reactor[J]. Annals Nuclear Energy, 2016, 87:511-516(1) |
| 16 | 盛美玲. 铅合金冷却反应堆非能动余热排出系统优化设计研究[D]. 合肥:中国科学技术大学, 2013 SHENG Meiling. Optimal design and study of passive decay heat removal system for lead alloy cooled reactor[D]. Hefei:University of Science and Technology of China, 2013(1) |
| 17 | 盛美玲, 金鸣, 柏云清, 等. 中国铅合金冷却研究堆事故余热排出系统概念设计与分析[J]. 核科学与工程, 2014, 34(1):91-96 SHENG Meiling, JIN Ming, BAI Yunqing, et al. Preliminary design and analysis of emergency decay heat removal system for China lead alloy cooled research reactor[J]. Nuclear Science and Engineering, 2014, 34(1):91-96(2) |
| 18 | 王强龙. 加速器驱动次临界堆始发事件及非能动安全特性初步研究[D]. 合肥:中国科学技术大学, 2013 WANG Qianglong. Preliminary reseach on initial events and passive safety characteristics of accelerator friven system reactor[D]. Hefei:University of Science and Technology of China, 2013(1) |