海洋条件下反应堆热工水力参数的子通道计算

李志威，张小英，陈焕栋，白 宁，历井钢

(1.华南理工大学 电力学院，广东 广州 510640；2.中科华核电技术研究院有限公司，广东 深圳 518026)

海洋条件下反应堆热工水力参数的子通道计算

李志威1，张小英1，陈焕栋1，白 宁2，历井钢2

(1.华南理工大学 电力学院，广东 广州 510640；2.中科华核电技术研究院有限公司，广东 深圳 518026)

针对海洋条件下反应堆的子通道热工水力分析，建立了海洋运动附加力模型和瞬态入口边界，将起伏、摇摆及复合运动的附加力关系式用于子通道模型的轴向和横向动量方程，并应用到COBRA ⅢC程序将其改造为适应海洋条件的反应堆子通道分析程序。作为验证，计算了加热实验通道和“奥陆”堆在起伏运动情况下热通道的临界热流密度比(CHFR)、出口空泡份额和冷却剂流量，并与文献结果对比。还详细计算了“奥陆”堆在起伏、不同摇摆中心和复合运动情况下，热通道的CHFR和不同位置子通道出口的热工水力参数。研究表明：海洋条件下反应堆的子通道热工水力参数随运动呈周期性变化；起伏运动对子通道的压降影响较大，摇摆运动对子通道冷却剂的流量和温度影响较大。

海洋条件；压水堆；运动附加力；子通道程序

核动力舰船在海洋条件下工作，运动条件下流过反应堆堆芯的冷却剂受到各种惯性力作用，使得堆芯子通道的热工水力参数呈现与陆地反应堆不同的特征，针对最小临界热流密度比(CHFR)和临界热流密度(CHF)的评价方法也不同，进而对反应堆的安全设计也提出不同的要求。目前国际通用的子通道程序一般都是基于陆地反应堆而开发的[1-3]，不能适用于海洋条件下的反应堆热工水力分析。针对深水潜探、海上作业平台及舰船核动力装置设计需要，研究海洋条件下反应堆的热工水力特征具有现实而深远的意义。

目前国内外针对海洋条件反应堆的热工水力研究注重分析附加力模型和反应堆的自然循环能力。钱立波等[4]提出了海洋条件下，运动对冷却剂流动的影响可归结为动量方程中附加力的改变，基于非惯性系海洋条件下几种基本运动特征(前进、倾斜、起伏、摇摆)，建立了典型海洋运动条件下用于动量方程中的附加力模型。谭长禄等[5]研究了倾斜、起伏和摇摆海洋条件的附加力模型，基于系统程序RELAP5开发了适用于海洋条件的反应堆热工水力系统分析程序，并对海洋条件下简单两回路系统的自然循环特性进行了计算分析。Kim等[6]针对海洋条件下船体的三维运动研究了运动附加力模型，并用于系统程序RETRAN的改造，计算摇摆、起伏和倾斜条件下船载反应堆的自然循环能力、冷却剂流量和蒸汽发生器水位的波动。将海洋运动附加力用于开发子通道程序的公开文献较少，Ishida等[7]针对海洋条件的起伏运动，修改了COBRA-Ⅳ-Ⅰ程序中轴向动量方程的重力加速度，采用重力加速度关联的CHF关系式，计算了起伏运动条件下，日本“奥陆”堆的CHF点位置及数值。

海洋条件下，舰船运动呈现单一运动或复合运动的多样性，反应堆子通道的热工水力特征变化复杂，目前尚无公开文献全面地研究海洋各种运动条件下反应堆的子通道热工水力参数。为此，本文将综合起伏、摇摆、起伏/摇摆复合运动下的运动附加力模型，用于目前陆地反应堆子通道程序COBRA ⅢC，将其改造为海洋条件下的程序。用改造后的程序，计算“奥陆”堆在地面和3种海洋运动条件下，堆芯子通道的热工水力参数及CHF点参数，并将起伏运动下的结果与文献[7]的结果进行比较。

1 子通道程序的动量方程

COBRA ⅢC程序的子通道模型方程采用均相流模型[3]，本文改造海洋条件的子通道程序时，主要考虑将相应运动下的附加力模型用于子通道模型的轴向和横向动量方程中。原子通道模型的轴向动量方程增加附加力模型后成为：

(1)

式中：mi为通道中冷却剂质量流量；Ai为通道横截面积；ρi为通道中冷却剂密度；ui为通道中冷却剂流速；pi为通道中的压力；fi为单相摩擦因子；φi为两相摩擦修正系数；Ki为格架摩擦系数；Vi为冷却剂比体积；g0为重力加速度；θ为轴向通道方向与重位方向之间的夹角；u*为相邻通道间的平均轴向流速；Wij为相邻通道的横向交混；Fx为海洋运动轴向附加力；下标i表示通道号。

原子通道模型的横向动量方程增加附加力模型后成为：

(2)

2 海洋运动的附加力模型和边界条件

海洋运动具有起伏、摇摆、起伏/摇摆复合的特征，令舰船上的非惯性坐标系为Oxyz，陆地上的惯性坐标系为O′x′y′z′，舰船在静止条件下非惯性坐标系成为惯性坐标系。本节以堆芯的单个冷却剂通道为对象，研究各种运动条件下，通道的轴向及横向动量方程中的附加力模型。

2.1 起伏运动的附加力模型

舰船的起伏运动通常具有简谐运动的特征，如图1a所示。假设起伏运动时船体呈现一随时间呈正弦曲线变化的运动加速度a(t)：

(3)

式中：t为起伏运动对应时刻；a(t)为对应t时刻起伏运动加速度；am为起伏加速度幅值；T为起伏加速度的波动周期。由此，起伏运动加速度对冷却剂轴向运动产生的附加力为：

(4)

图1 海洋条件下的各种运动Fig.1 Motions under ocean condition

2.2 摇摆运动的附加力模型

(5)

根据文献[4]提出的海洋运动附加力模型，可得摇摆运动的总附加力为：

(6)

式中：ω(t)为摇摆运动的角速度；u为冷却剂相对惯性系的流速；r为冷却剂至摇摆中心点的距离。

轴向x的附加力分量为：

(7)

沿通道横向z的附加力分量为：

(8)

摇摆中心位置不同，计算点的轴向和横向坐标不同，相应的轴向附加力也不同，摇摆中心距计算点越远，轴向和横向附加力就越大。

2.3 起伏/摇摆复合运动的附加力模型

由于舰船的起伏和摇摆运动均为简谐运动，起伏/摇摆的复合运动可考虑为两种运动的叠加。因此，得到起伏/摇摆复合运动的轴向附加力为：

(9)

起伏/摇摆复合运动的横向附加力为：

(10)

2.4 边界条件

研究表明，海洋运动对反应堆系统堆芯入口参数的影响主要体现在堆芯的入口流量[8]，摇摆运动使堆芯入口流量随运动呈瞬态变化，而起伏运动对堆芯入口流量影响较小。摇摆运动下，不同时刻对应堆芯入口流量是不同的，针对入口流量的瞬态变化，建立了通道入口瞬变边界，使堆芯入口流量为时间的函数，满足：

(11)

式中：m摇摆为摇摆运动下堆芯入口流量；m平均值为摇摆运动下堆芯入口流量平均值；m幅值为摇摆运动下的流量波动振幅。

3 海洋运动改造程序简介

基于COBRA-ⅢC改造的海洋条件的子通道热工水力分析程序，通过在程序中建立通道入口参数的连续瞬变边界和附加力瞬变模型，使程序具备模拟海洋条件下堆芯入口边界参数连续瞬态变化与海洋运动导入的附加力瞬态变化功能。因此改造的程序初步具备对海洋条件下反应堆堆芯瞬态热工水力参数的分析能力。

4 海洋条件下反应堆热工水力子通道计算结果分析

4.1 加热实验通道的验证性计算

为验证本文所改造的子通道程序，对文献[7]中实验测量的起伏运动情况下，一加热通道的热工水力特性，采用本文程序进行了相应计算，并将计算结果与实验结果进行对比。实验装置采用电加热，系统压力为0.3 MPa，入口冷却剂的过冷度为1 K，电加热器棒外径为10 mm，长度为1 300 mm，加热功率为130 kW/m2。初始状态，实验装置在静止条件下运行，运行稳定后启动实验装置起伏运动机构，设置起伏运动的重力加速度变化幅度为0.35g0，周期为5.2 s，由此得到起伏加速度随时间的变化如图2所示，待起伏运动进入稳定阶段后，再分析运动对通道热工参数的影响。

图2 起伏加速度随时间的变化Fig.2 Acceleration amplitude vs. time

所计算的该加热通道内出现CHF的时刻如图3所示。图中，本文计算出加热通道出现CHF的时段与文献中实验观测的基本一致，但计算出现CHF的时刻较实验观测的时刻早。对此，本文分析是实验中流体升温与壁面的加热存在一定的延迟，而计算中没有这种影响，造成计算出现CHF的时刻较实验观测的早。

起伏运动对CHFR和出口空泡份额的影响示于图4。图4a中，CHFR=临界热流密度/实际热流密度，图中文献结果是采用程序COBRA Ⅳ改造后计算的。从图4a可看出，本文所计算的CHFR与文献结果的一致性很好。由于起伏加速度的周期性变化，通道的CHFR也呈现相同的周期性变化。图4b中，计算结果有两组：1) 本文用COBRA ⅢC程序改造后计算的结果；2) 文献[7]用COBRA Ⅳ程序改造后计算的结果。可看出：实验测量的空泡份额随时间呈一定的周期性变化，但文献计算和本文计算的出口空泡份额均基本不变。原因是实验装置是一回路系统，在起伏惯性力的周期性作用下，装置中的冷却剂流量发生周期性变化，导致出口空泡份额的周期性变化。而在子通道计算中，通道入口条件参数和加热功率一定，起伏运动时单一通道不会产生横流，因此，通道中冷却剂流量和空泡份额均不会随时间变化。图中的两组计算结果的一致性较好，文献的计算结果较本文计算结果约小4%。

图3 通道中CHF出现的时刻Fig.3 Moment of CHF occurring in channel

4.2 “奥陆”堆起伏运动的子通道验证计算

针对日本核动力舰船采用的“奥陆”堆，文献[7]还研究了在COBRA Ⅳ程序中添加起伏惯性力并修改CHF关系式，以计算“奥陆”堆起伏运动下，反应堆热通道出口的冷却剂流量和CHFR。

图4 起伏运动对CHFR和出口空泡份额的影响Fig.4 Effect of heaving motion on CHFR and outlet void fraction

文献[7]推荐“奥陆”堆起伏运动下的CHF关系式为：

(12)

式中：qCHF0为地面条件下的临界热流密度；Δg为重位方向变化的起伏加速度。

本文采用文献[7]的CHF关系式和基于COBRA ⅢC改造的子通道程序，在相同条件参数下进行了计算，并将本文计算的热通道出口参数与文献结果进行了比较。“奥陆”堆的堆芯子通道如图5所示，所采用的计算参数列于表1。由于文献[7]未给出计算采用的堆芯功率分布，本文根据文献[9]选取“奥陆”堆的功率因子(表2)。

图5 “奥陆”堆堆芯子通道Fig.5 Core sub-channel of MUTSU reactor

表1 “奥陆”堆的子通道计算参数Table 1 Sub-channel calculation parameter for MUTSU reactor

表2 “奥陆”堆的功率因子Table 2 Power factor of MUTSU reactor

起伏运动对反应堆热通道出口冷却剂流量和CHFR的影响示于图6。从图6a可看出，流量随反应堆的起伏运动呈周期性变化，这一点与上述单加热通道在起伏运动下的流量特征不同。分析原因是“奥陆”堆堆芯热功率存在径向差异，子通道间有横流，起伏运动使得横流量呈周期性变化，导致热通道出口的冷却剂流量呈周期性变化。本文计算的热通道出口冷却剂流量随时间的变化与文献结果一致性较好。从图6b可看出，由于受到起伏惯性力和流量周期性变化的影响，导致热通道的CHFR在起伏运动过程中也呈周期性变化。本文计算热通道CHFR的变化周期和幅值与文献的一致，但本文计算的CHFR平均值较文献结果约高25%，原因可能是本文采用的堆芯功率因子与文献[7]有差异，文献[7]中采用的功率因子不详。

图6 起伏运动对反应堆热通道出口冷却剂流量和CHFR的影响Fig.6 Effect of heaving motion on exit coolant flux and CHFR in hot channel

4.3 “奥陆”堆在地面及3种海洋运动条件下的子通道计算分析

为全面研究海洋运动对反应堆堆芯热工水力参数的影响，本文还针对“奥陆”堆，采用表2所列的功率因子，计算了地面条件、起伏运动、摇摆运动及复合运动4种工况下的子通道热工参数；对比不同运动状况下堆芯热通道的CHFR；获得典型子通道：边通道3、角通道10及中心通道14的热工水力参数。计算摇摆运动时，起伏运动的附加力参数与表1相同。起伏/摇摆复合运动参数与单一起伏、摇摆运动参数一致。图7示出4种工况下相对轴向附加力。

图7 4种工况下相对轴向附加力Fig.7 Relative axial additional force in different motion conditions

1) 地面及3种海洋运动条件下的子通道CHFR计算分析

4种工况下热通道最小CHFR如图8所示。地面条件下热通道热点处的CHFR不随时间变化，但在其他3种运动工况下，热点CHFR均随时间周期性地变化。对比图7，运动情况下CHFR变化的周期与附加力相同。起伏运动的CHFR相对地面条件上下波动；摇摆运动的CHFR均比地面条件的小，这是由于摇摆运动加强了通道间的横流，增强了冷却剂对壁面的冲刷，加大了壁液间的换热，从而降低了CHFR；起伏/复合运动的CHFR体现为两种运动的叠加。

图8 4种工况下热通道最小CHFR Fig.8 MCHFR of hot channel in different motion conditions

为比较“奥陆”堆在地面和不同运动工况下热通道的CHFR，研究附加力分别达到峰值和谷值时沿热通道的CHFR。图9为起伏运动热通道CHFR轴向分布。图中热通道的CHFR沿通道高度呈现先降低后升高的趋势，这是由通道内冷却剂的温度上升及反应堆的轴向功率因子决定的；附加力达到最大时，反应堆的流量也达到最大，此时CHFR相比地面工况的高；反之亦然。

图9 起伏运动热通道CHFR轴向分布Fig.9 CHFR axial distribution of hot channel under heaving motion condition

海洋条件下热通道CHFR的变化示于图10。从图10a可看出，随附加力的变化，热通道的CHFR变化很小，这是因为摇摆运动的轴向附加力变化幅度要较起伏运动的小。

图10 海洋条件下热通道CHFR的变化Fig.10 Variation of CHFR for hot channel under ocean condition

比较图10b与图9可看出，在复合运动中，起伏运动对热通道的CHFR影响较大。综合地面及以上3种运动情况的结果，热通道中最小CHFR的位置均在z/H=0.8处，不因运动方式不同而改变。

2) 地面及3种海洋运动不同位置的子通道计算分析

同样，以“奥陆”堆中心通道为例，在地面及3种运动条件下，通道的详细热工水力参数如图11、12所示。从图11a可看出，起伏运动的压降相对地面压降呈上下对称的周期变化；摇摆运动的压降虽也呈周期性的起伏变化，但均较地面压降小；复合运动的压降表现为起伏和摇摆两种运动的压降叠加。从图11b可看出，起伏运动的出口冷却剂流量与地面情况基本相同，摇摆运动和起伏/摇摆复合运动的流量呈上下对称的周期变化，原因是摇摆运动时子通道间有较大横流，而起伏运动时横流很小。复合运动的流量变化曲线的相位较摇摆运动滞后约π/4，这是因为复合运动的轴向附加力改变了相邻通道的压差，在相邻通道间产生横流，横流产生相对压差有一定滞后。

图11 不同运动对压降和出口冷却剂流量的影响Fig.11 Effect of different motions on pressure drop and exit coolant flux

图12 不同运动对出口冷却剂温度和包壳温度的影响Fig.12 Effect of different motions on exit coolant and cladding temperatures

从图12a可看出，受流量变化的影响，起伏运动的出口冷却剂温度与地面情况基本相同；摇摆运动与复合运动的出口冷却剂温度随时间有小幅的周期变化，复合运动的流量曲线相对摇摆运动滞后π/4，两种运动均使出口冷却剂温度有所降低。从图12b可看出，包壳温度也受流量变化影响，起伏运动的出口冷却剂温度与地面情况基本相同；摇摆运动和复合运动的出口包壳温度均对称于地面温度周期地变化，复合运动的包壳温度变化幅度更大；摇摆运动和复合运动中心通道出口的包壳温度变化曲线存在相位差π/4。

5 结论

本文通过建立海洋条件的附加力模型，在子通道程序COBRA ⅢC基础上，开发适用于海洋条件下的反应堆热工水力子通道分析程序。开发的程序能对海洋条件下反应堆堆芯不同工况进行计算，具备对海洋条件反应堆堆芯热工水力计算分析能力。利用程序计算并验证了起伏运动实验通道的CHFR和空泡份额，以及“奥陆”堆在起伏运动情况下热通道的CHFR、出口空泡份额和冷却剂流量，验证计算的结果与文献结果一致。还详细计算了“奥陆”堆在起伏、摇摆和复合运动情况下，热通道的CHFR和不同位置子通道出口的热工水力参数。

[1] 唐锡文. 用FLICA3分析自然循环实验特性[J]. 核动力工程，1997,18(3):254-257.

TAN Xiwen. Analysis on behavior of natural circulation by FLICA3[J]. Nuclear Power Engineering, 1997, 18(3): 254-257(in Chinese).

[2] 白宁,朱元兵,任志豪,等. 子通道分析程序LINDEN的开发与初步验证[J]. 原子能科学技术，2013,47(2)：299-301.

BAI Ning, ZHU Yuanbing, REN Zhihao, et al. Development and verification of subchannel analysis code LINDEN[J]. Atomic Energy Science and Technology, 2013, 47(2): 299-301(in Chinese).

[3] JACKSON J W. COBRA ⅢC/MIT2: A digital computer program for stead state and transient thermal-hydraulic analysis of rod bundle nuclear fuel elements, BNWL-1695[R]. Cambridge: Massachusetts Institute of Technology, 1973.

[4] 钱立波,田文喜,秋穗正,等. 一维冷却剂通道海洋条件下附加力模型研究[J]. 核动力工程，2012,33(2):104-109.

QIAN Libo, TIAN Wenxi, QIU Suizheng, et al. Research on model of additional forces of ocean conditions in one-dimensional coolant channel[J]. Nuclear Power Engineering, 2012, 33(2): 104-109(in Chinese).

[5] 谭长禄,张虹,赵华. 基于RELAP5的海洋条件下反应堆热工水力系统分析程序开发[J].核动力工程，2009,30(6):53-62.

TAN Changlu, ZHANG Hong, ZHAO Hua. Development of ocean-condition code based on RELAP5[J]. Nuclear Power Engineering, 2009, 30(6): 53-62(in Chinese).

[6] KIM J K, PARK G C. Development of RETRAN-03/MOV code for thermal-hydraulic analysis of nuclear reactor under moving conditions[J]. Journal of the Korean Nuclear Society, 1996, 28(6): 542-550.

[7] ISHIDA I, KUSUNOKI T, MURATA H, et al. Thermal-hydraulic behavior of a marine reactor during oscillations[J]. Nuclear Engineering and Design, 1990, 120: 213-225.

[8] 谭思超，庞凤阁，高璞珍. 摇摆运动对自然循环传热特性影响的实验研究[J]. 核动力工程, 2006,27(5):33-37.

TAN Sichao, PANG Fengge, GAO Puzhen. Experimental research of effect of rolling upon heat transfer characteristic of natural circulation[J]. Nuclear Power Engineering, 2006, 27(5): 33-37(in Chinese).

[9] TAKANORI S, MASAFUMI I, MASAO O, et al. Static core performance simulator SCOPERS-2 for light water reactors and its application, NEACRP-L205[R]. Tokai: Japan Atomic Energy Research Institute, 1978.

Sub-channel Analysis on Thermal-hydraulic Characteristic of PWR under Ocean Condition

LI Zhi-wei1, ZHANG Xiao-ying1, CHEN Huan-dong1, BAI Ning2, LI Jing-gang2

(1.SchoolofElectricPower,SouthChinaUniversityofTechnology,Guangzhou510640,China;2.ChinaNuclearPowerTechnologyResearchInstitute,Shenzhen518026,China)

According to the thermal-hydraulic analysis of sub-channel, additional force models of heaving, rolling and complex movement under ocean conditions were established and applied to the sub-channel momentum equation in COBRA ⅢC sub-channel analysis code. The inlet boundary of ocean conditions was also established in this paper. The critical heat flux ratio, exit void fraction and coolant flow of the hot channel were calculated for MUTSU reactor and experiment heating tube, and the results were compared with the experiment data from literature. Calculated thermal-hydraulic characteristics of MUTSU reactor were studied under ocean conditions. The research results show that thermal-hydraulic parameters are varied periodically under fluctuating of ocean conditions. The heaving motion has a large impact on pressure drop of the sub-channel. The sub-channel coolant flow and temperature are greatly influenced by swing motion.

ocean condition; pressurized water reactor; motion additional force; sub-channel code

2014-05-23;

2014-07-11

国家自然科学基金资助项目(51376065,51176052)

李志威(1990—)，男，湖南岳阳人，硕士研究生，工程热物理专业

TL33

A

1000-6931(2015)10-1758-08

10.7538/yzk.2015.49.10.1758