Vol. 55,No. 4
Apr 2021
第55卷第4期
2021年4月
原子能科学技术
AtomicEnergyScienceandTechnology
典型压水堆堆芯物理-热工耦合 稳态计算软件的开发与验证
李治刚丄2，安萍* *2潘俊杰*2卢川*2芦餠丄2"，杨洪润12
收稿日期2020-05-06；修回日期.020-06-26
作者简介：李治刚(1989-)，男，四川成都人，工程师，硕士 ，从事反应堆堆芯计算软件研发
* 通信作者：芦 餠,E-mail : ***************
"1中国核动力研究设计院，四川成都610213；
2.核反应堆系统设计技术重点实验室，四川成都610213)
摘要：为能更加准确地模拟典型压水堆中强烈的物理-热工耦合现象，研制了压水堆堆芯物理-热工耦合 计算软件ARMcc 。其中物理计算模块基于四阶节块展开法(NEM)和格林函数节块法(NGFM )，热工计
算模块基于一维的单相单通道换热模型和一维圆柱导热计算模型，在程序中采用有限体积法和有限差 分法求解一维圆柱导热模型$基于典型压水堆基准题NEACRP-L-335对程序的稳态耦合计算能力进
行了验证，程序计算的堆芯关键参数如临界硼浓度、堆芯多普勒温度等参数与参考结果符合良好，临界
硼浓度与参考结果的相对偏差均小于0. 5% $另外研究4种计算模式对模拟堆芯物理-热工耦合过程的 影响，选择PARCS 程序计算结果为对比，发现NGFM+DIF 模式能更加准确地模拟堆芯燃料多普勒温
度和堆芯功率分布；NGFM+VOL 模式能更加准确地模拟临界硼浓度；NEM+VOL 模式能更加准确地
模拟堆芯燃料最高温度$
关键词：压水堆堆芯；物理-热工耦合；稳态；圆柱导热模型
中图分类号:TL33
文献标志码:A 文章编号:1000-6931(2021)04-0685-08
doi ：10. 7538/yzk. 2020. youxian. 0296
Development and Verification of Typical PWR
Core Physical and Thermal-hydraulic Steady Coupling Code
LI Zhigang 12 , AN Ping 12 , PAN Junjie 12 , LU Chuan *'2 , LU Wei *'2'* , YANG Hongrun 12
(1. Nuclear Poxver Institute of China , Chengdu 610213 , China ；
2. Science and Technology on Reactor System Design Technology Laboratory , Chengdu 610213 , China )
Abstract : In order to more accurately simulate the strong neutronics physical and the?
mal-hydrauliccouplingphenomenoninatypicalPWR ARMcc asoftwareforthephys-
icalandthermal-hydrauliccouplingcalculation ofPWR core was developed Inthe
ARMccprogram thephysicalcalculation module is based on the fourth-order nodal expansion method (NEM ) and Nodal Green's function method (NGFM ) , the thermal-
hydrauliccalculation moduleis based on one-dimensional single-phase single-channel
heattransfer modeland one-dimensional cylinder heat conduction calculation model
686原子能科学技术第55卷
The finite volume method and finite difference method were used to solve heat conduction model in ARMcc program.Based on the typical PWR benchmark NEACRP-L-335, the ability of steady-state coupling calculation of the program was verified.The key parameters of the program,such as critical boron concentration and core Doppler temperature,are in good agreement with reference results.The relative deviation between the critical boron concentration and the reference results is less than0.5%.In addition, theinfluencesofthefinitevolumemethodandthefinitedi f erencemethodontheresults ofthecouplingprogram werestudied.ThePARCSprogram wasselectedasthecompar-ison program.It is found that NGFM+DIF mode can more accurately simulate the core fuel Doppler temperature and core power distribution,NGFM+VOL mode can more accurately simulate the critical boron concentration,and NEM+VOL mode can more a//uratelysimulatethe/orefuelmaximumtemperature.
Key words:PWR core；physical and thermal-hydraulic coupling；steady；cylinder heat conducion model
在压水堆中燃料温度、慢化剂密度、冷却剂温度等热工参数会影响中子的截面，进而影响堆芯中子通量和功率分布，而功率又反过来影响堆芯燃料温度、冷却剂温度等热工参数，即压水堆中子物理与热工之间存在强烈的反馈现象因此目前在压水堆的设计和安全分析中均须考虑物理-热工耦合计算，以便更加真实地模拟堆芯稳态和瞬态过程$自20世纪80年代以来，国内外针对压水堆物理-热工耦合现象开展了大量研究，提出了较多有效的耦合计算方法，如Picard迭代法*2、JFNK迭代法并研发了大量三维物理-热工耦合计算软件，如CRONOS/FLICA4、PARCS/ TRACE⑷、CSSS5+JDYN3D/ATHLET C6］、NAL-SANMT/CORBA-JV［7］等$
目前典型的堆芯物理-热工耦合计算程序常采用外耦合方式实现，且常采用1种中子物理或热工水力计算方法求解守恒方程$为研究不同中子物理、热工水力计算方法对压水堆物理-热工耦合计算的影响，本文通过内耦合方式实现物理-热工耦合计算，研制典型压水堆堆芯物理-热工耦合计算软件ARMcc，其中采用节块展开法（NEM）™和格林函数节块法（NGFM）9求解中子扩散方程，分别采用单通道计算模型［10］和一维圆柱导热计算［11］（采用有限差分法或有限体积法求解）冷却剂温度和燃料温度，并采用典型压水堆耦合基准题NEACRP1213弹棒初始参数对软件稳态计算功能进行验证$
1理论模型
1.1稳态中子扩散方程求解方法
多稳态中子扩散方程「1415+可表达为：—▽•Dg P"+/d"/—,/s/f"/=
///
G
；g=1，…,G（1）g'-l
式中：g为裂变谱;Dg为第g中子扩散系数;/g，为第g'裂变宏观截面；"g为第g中子注量率；/为第g，平均裂变中子数；/g o g为第'到第g的散射截面；/R g为第g的移除截面；为特征值或有效增殖因数$
目前式（1）的典型求解方法包括有限差分法、节块法、有限元方法等，其中节块法由于具有较高的效率和精度，在商业软件中被广泛应用$节块法的重要特点是对式（1）进行横向积分，将三维问题的求解变成联立求解3个一维问题$在网格8内，对给定的坐标方法"（交替等于X,yY），对式（1）沿与6方向垂直的另两个方向#和W进行积分，得到3个一维方程（为便于描述，全文略去节块编号8）:
g—1,…,G；6—D,y,#/F（2）式中，"/”、Q g”和L g”分别为横向积分通量、源项（包括裂变源项和散射源项）和横向泄漏项$典型的节块法有NEM、解析节块法（ANM）
第4期李治刚等：典型压水堆堆芯物理-热工耦合稳态计算软件的开发与验证687
和NGFM等，不同节块法的差异在于%O 和=gu等近似处理方式。本文将采用四阶NEM和第二类边界条件NGFM进行多稳态中子扩散方程的求解，具体的理论详见文献［8-9］。
=2热工水力计算方法
1)单通道计算模型
冷却剂单相单通道的守恒方程为：
&6=q(3)
d z
式中:&为冷却剂密度，kg/m31为冷却剂焓, J/kg；6为冷却剂流动速度，m/s；q为体积热流密度，W/m3$
对式(3)在网格m进行z方向的积分得到：
(1%)out—(1%)in=-(4)
△z
式中：下标out、in分别表示网格m的入口和出口界面;1为网格m的平均体积热流密度$根据入口边界条件，采用式(4)可依次计算得到每个网格出口焓l out，则网格m的平均焓为：
1=(10ut+l n)/2(5)
2)一维圆柱导热模型
一维圆柱导热稳态模型的守恒方程为：
--"0)-)+q=o⑹
式中：为导热系数，W/(cm・H］1)0为温度，°c。
典型的压水堆燃料元件几何形式如图1所示$在径向上燃料芯体划分为&个网格，气隙为1个网格，包壳划分为2个网格$
图1典型压水堆燃料元件几何形式Fig1GeometryoffuelelementintypicalPWR
(1)有限差分法
采用有限差分法(DIF)对式(6)在网格"处进行离散，得到网格"的温度计算关系式：
(—-—120—1+
(—+1/2+--1/2)0"——+1/20汁1)=1(7)式中：为材料的导热系数-为第"个网格的半径；△厂为第"个网格的间距$关于燃料中心、芯体外表面、包壳内外表面的处理方式参见文献［16］$
(2)有限体积法
采用有限体积法对式(6)进行离散，表达形式如下：
丄(———*+
(――1/2——1/2+—+1/2-汁1/2)・
0——+1/2—+1/2T l+i)=1(2+1/2——2T/2)(8)式中，—+1/2和——1/2分别为第"个网格左右边界的半径$
式(7)、(8)均为三对角矩阵，采用高斯-赛德尔迭代计算得到燃料元件的径向温度分布，采用罗兰公式加权计算燃料有效温度$
2验证与分析
基于上述中子物理和热工水力计算理论，采用C/C++语言开发了反应堆堆芯多物理耦合计算软件(ARMcc)。本文分别采用三维LWR基准题「17+和NEACRP-L-335基准题［1-13,8］对程序中子扩散方程求解和物理-热工耦合计算进行验证$
2.1中子扩散方程求解验证
三维LWR基准题的几何布置及材料截面参数见文献［17］$采用NEM和NGFM对三维LWR基准题进行模拟，堆芯有效增殖因数化ff结果及偏差列于表1，三维LWR基准题相对功率分布及偏差如图2所示$
表1Uff结果及偏差
Table1Result and deviation of k ff
计算方法；eff偏差/pcm
NEM 1.02908 5.4
NGFM 1.029052
参考值 1.02903
由表1和图2可知，本软件采用的NEM 和NGFM方法模拟三维LWR
基准题均具有
688原子能科学技术第55卷
较高的精度。与参考结果相比，NEM和
NGFM方法计算的；eff偏差均在10pcm以内；
堆芯径向相对功率与参考结果相比，采用NEM
方法的最大偏差为0.553%，采用NGFM方法的
最大偏差为103%$
相对功率参考值
NEM方法相对偏差/%
NGFM方法相对偏差/%
0.729 -0.261 0.110
1.397
-0.322
—0.490
1.368
-0.241
—0.270
1.178
-0.136
—0.130
0.476
0.042
0.210
0.597
0.553
0.270
1.281
-0.312
—0.540
1.432
-0.335
—0.360
1.311
-
0.214
—0.110
0.972
-0.093
—0.020
1.291
-0.271
—0.360
1.181
-0.127
—0.150
1.193
-0.201
—0.370
1.072
-0.187
—0.200
0.700
0.186
0.610
0.923
0.087
0.110
0.611
0.376
1.030
1.089
0.000
-0.090
0.610
-0.016
0.030
0.976
0.041
—0.160
0.953
0.010
—0.050
0.866
0.289
0.420
1.000
0.100
0.070
0.711
0.352
0.300
0.959
0.146
—0.060
0.757
0.251
0.120
0.777
0.219
—
0.060
图2三维LWR基准题相对功率分布及偏差Fig.2Relative power distribution and
deviation of3D LWR benchmark
2.2物理-热工耦合计算验证
NEACRP弹棒基准题是1991年由Fine-main等「*2「*3+建立的，包含压水堆和沸水堆两种堆型，主要用于轻水堆堆芯三维物理-热工耦合软件的验证$其中压水堆基准题参考了典型压水堆的几何尺寸和运行状态，堆芯布置了157个燃料组件和1圈径向反射层，轴向分为18层，如图3所示「*8+。根据初始功率水平和控
reactor technology反射层
控制棒组件B工况弹出控制棒
CR A工况弹出控制棒C工况弹出控制棒
图3PWR NEACRP弹棒基准题堆芯模型
Fig.3Core model of PWR NEACRP
rod ejection benchmark 制棒弹出位置的不同，共包含6种工况:A1，1/4堆芯、中心控制棒在热态零功率（HZP）弹出；A2,1/4堆芯、中心控制棒在热态满功率（HFP）弹出；\1,1/4堆芯、外围1组控制棒在热态零功率弹出；\2,1/4堆芯、外围1组控制棒在热态满功率弹出；C1,1/2堆芯、外围1组控制棒在热态零功率弹出；C2,1/2堆芯、外围1组控制棒在热态满功率弹出$
关于NEACRP基准题几何尺寸划分及材料截面参数可详见文献[12]，基准题的参考结果由PANTHER程序采用精细时空网格计算，在1993年发布初版.1997年发布修订版，计算结果包括稳态和瞬态关键结果参数$本文将采用1997年修订版的稳态计算结果对本软件物理-热工稳态耦合计算模块进行验证，而对于PANTHER未提供的结果参数（如堆芯组件径向功率分布等）将采用PARCS软件计算结果作为参考$2种物理计算方法和2种热工计算方法组合后形成4种计算模式（NEM i
VOL对应节块展开法和有限体积法，NEM+ DIF对应节块展开法和有限差分法，NGFM+ VOL对应格林函数法和有限体积法，NGFM+ DIF对应格林函数法和有限差分法）$
1）临界硼浓度
采用ARMcc针对NEACRP基准题计算的临界硼浓度对比结果列于表2$与PAN-THERC1997）相比，4种
模式计算的6种工况的临界硼浓度与参考结果符合较好，最大相对偏差在0.5%以内，其中在A2和C2工况中，NEM方法计算临界硼浓度的精度高于NGFM 方法，而在剩余4个工况中，NGFM方法计算临界硼浓度的精度高于NEM方法$
2）堆芯组件径向相对功率最大值
径向相对功率最大值F y,max列于表3$与PARCS结果相比，ARMcc计算的堆芯组件径向相对功率最大值与参考结果符合较好!相对偏差在0.5%以内$采用NEM方法计算的最大相对偏差一0.39%出现在A2算例中，而采用NGFM方法计算的最大相对偏差0.27%出现在C2算例中$采用有限体积法或有限差分法对F y,max的影响较小$
B1算例堆芯组件径向功率分布如图4所示，其中在同一种物理方法中采用VOL和
DIF
第4期 李治刚等：典型压水堆堆芯物理-热工耦合稳态计算软件的开发与验证
689
不同工况下的临界硼浓度/ppm
表2不同工况下的临界硼浓度
Table2 C%iticalbo%onconcent%ationofdife%entcases
注：括号内为相对偏差，％
方法
A1
A2Bl B2C1C2PANTHER(1997)
561. 21 156. 61248.
1183.81128.29
1156.63PARCS 561. 261 154. 361248.211182.87
1128.641154.36
NEM +VOL
562. 82(0. 29)1 158. 37(0. 12)1 252. 88(0. 39)1 186. 97(0. 23)1133. 20.42)1 157. 66(0. 11)NEM+DIF
562. 82(0. 29)
1 157. 94(0. 12)1 252. 88(0. 39)1 186. 54(0. 23)1133. 20.42)
1 157. 87(0. 11)NGFM+VOL 560. 54( — 0. 12)1154.15 ( — .21 )
1247.45 — . 4)1182.63— .1 )1127.72— . 5 )1154.14— .22)NGFM+DIF
561.61 ( . 7)
1 153. 73( — 0 - 25)
1247.45 — . 4)
1182.21 — .13)
1127.72— . 5 )
1153.72— .25 )
表3不同工况下的堆芯组件径向相对功率最大值
Table 3 Maximum radial relative power of core assembly of different cases
不同工况下的F y ,—ax
方法
A1
A2B1B2C1C2
PARCS
1.9113 1.2 15
1.2765 1.17)9 1.4447
1.2)15NEM+VOL 1.9 59— .28) 1.1968— .39) 1.275 — .12) 1.1678—).26)
1.446 8(0. 15) 1.1971 —).37)
NEM+DIF 1.9 59— .28)
1.1968— .39) 1.275
— .12)
1. 167 8( — 0. 24) 1.446 8(0. 15)
1. 196 9( — 0. 38)NGFM+VOL 1. 913 3(0. 10 ) 1.2 460.26) 1. 277 4(0. 07) 1. 171 1(0. 0 2) 1.445 9(0. 08) 1.2)480.27)NGFM+DIF
1. 913 0 (0. 0 9)
1.2 460.26)
1. 277 4(0. 07)
1.17 7— . 2)
1.445 9(0. 08)
1.2)460.26)
注：括号内为相对偏差，％
的相对功率分布基本相同。从图4可知，相对 功率V0. 9的组件功率最大偏差情况为:NEM
方法为2.09% ,NGFM 方法为一0. 14%，相对功 率＞0. 9的组件功率最大相对偏差情况为:NEM
0.4080.7240.976Zl.086 1.261 1.1270.4890.5590.4000.7110.961 1.071 1.249 1.1200.4890.5670.408
0.7240.976 1.086 1.261 1.1260.4890.558
0.8240.885 1.121 1.250 1.1830.9070.5030.8110.871 1.108 1.238 1.1790.9100.5140.825
0.885 1.121 1.250
1.184
0.9070.503
0.549 1.060 1.277 1.199 1.0700.541 1.049 1.270 1.197 1.0820.548
1.059 1.277
1.198 1.070PARCS 1.204 1.231 1.2320.792NEM+VOL 1.197 1.228 1.2370.807NGFM+VOL
1.205
1.230 1.231
0.791
1.1410.8971.1460.9141.141
0.897
图4 B1算例堆芯组件径向功率分布(1/8堆芯)
Fig.4 Radicalpowerdistribution
of Bl case (1/8 core )
方法为一* 48% ,NGFM 方法为一0. 06%。
3） 堆芯轴向相对功率最大值
堆芯轴向相对功率峰值F
y —
ax 列于表4。
与 PARCS 结果相比 ,4 种计算模式得到的堆芯 轴向相对功率最大值与参考结果符合较好，最
大相对偏差不超过0.1%。采用NEM + VOL 模式计算时在满功率工况A2、\2和C2中相对 偏 差达 到 最大 , 为 0.1% $
4） 堆芯燃料最高温度
表 5 列出不 同 工 况 下 的 堆 芯 燃 料 最 高 温 度。与PARCS 结果相比，因零功率算例（A1、
B1和C1）堆芯功率较小，燃料最高温度与初
始温度相同，相对偏差为0. 0% ；在满功率算 例中，不同计算模式得到的燃料最高温度存
在明显差异：1）采用相同物理计算方法时，采
用有限体积法得到的燃料最高温度低于采用有 限差分法；2）采用相同热工计算方法时，NEM
方法计算的结果低于NGFM 方法；3） NEM +
VOL 模式相对偏差最小（\2算例中为0. 77%），
版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。
发表评论