堆芯氧化腐蚀产物沉积层内两相沸腾传热传质特性

刘 延，刘晓晶，何 辉

(上海交通大学 核科学与工程学院，上海 200240)

反应堆堆芯长期处于高温、高压、高放射性的环境下，燃料性能是评价堆芯安全性和经济性的关键性指标。反应堆运行时，蒸汽发生器等持续受到冷却剂的冲刷腐蚀，形成大量以铁、镍离子及其氧化物为主的氧化腐蚀产物。这些腐蚀产物在过冷沸腾(SNB)的驱动下，沉积在堆芯上部的燃料包壳表面，形成一层薄的疏松多孔的结垢层，被称为氧化腐蚀产物沉积层(chalk river unidentified deposit, CRUD)[1]。一方面，由于腐蚀产物的沉积，燃料棒径向热阻增大，传热被恶化。另一方面，沉积层内呈现疏松多孔的结构，SNB被强化，使得硼元素不断浓缩直至析出并被沉积层吸附，导致硼元素在燃料棒轴向分布不均匀。由于10B具有显著的吸中子能力，因此沉积层内的硼吸附会改变堆芯轴向的中子分布，使得堆芯功率分布向底部畸变，引发堆芯功率漂移现象(CRUD-induced power shift, CIPS)[2]。

为了评估CRUD内的传热传质现象，许多学者开发了一系列软件，对在压水堆运行条件下CRUD内的传热传质进行了模拟研究[3-11]。例如，EPRI的BOA(boron-induced offset anomaly)[12]和CASL的MAMBA-3D[13]。但关于CRUD内烟囱沸腾和形貌结构影响的研究仍相对较少。为此本文建立一套二维的传热传质预测方法，实现CRUD内的传热、压降、流速、溶质扩散和化学反应等多物理现象的耦合，并获得冷却剂在烟囱表面蒸发对传热和硼吸附的影响，以及传热传质随沉积层孔隙率、厚度、烟囱半径、烟囱直径等形貌参数的变化规律。

1 建模方法

图1 CRUD内传热、流动、传质和化学反应等多物理现象(a)和计算区域及边界(b)示意图Fig.1 Schematics of multi-physical phenomena within CRUD depositions (a) and simulation domain and boundary (b)

各种溶质之间保持化学平衡。为区分CRUD内蒸汽和水的不同位置，将烟囱区域称为干区，而环绕它周围的圆柱壳体称为湿区[8,15]。烟囱表面的沸腾增强了溶质的浓缩，不断浓缩的溶质反过来影响着溶质间的化学平衡和饱和温度。选择1个湿区作为计算区域，以沉积层厚度方向z(燃料棒径向方向)和垂直于z的r方向(燃料棒轴向方向)建立二维圆柱坐标，如图1b所示。

1.1 传热模块

CRUD内部的传热可分为多孔介质的导热、燃料包壳的传热和烟囱边界沸腾带走的热量3个部分，由能量守恒得到的控制方程及边界条件如式(1)～(5)所示。

(1)

(2)

(3)

(4)

(5)

式中：kCRUD为CRUD湿区内固体和液体加权后的总导热系数，受孔隙率和温度等参数的影响；qclad为包壳与CRUD交界面的热流密度；d为CRUD厚度；hc为冷却剂与CRUD表面的对流换热系数，在计算中取定值12 000 W/(m·K)[10]；Tf为冷却剂主流区温度；rc和Rc分别为烟囱半径和烟囱中心到湿区对称边界的距离；he为烟囱表面水的蒸发换热系数，由Pan[9]提出，用于计算烟囱表面发生SNB带走的能量，如式(6)所示；Tsat为饱和温度，受溶质浓度等参数的影响。

(6)

1.2 毛细流动模块

CRUD毛细热管内的流速较小，因此内部的流动符合达西定律。假设CRUD内冷却剂是不可压缩流体，控制方程及边界条件如式(7)～(11)所示。

(7)

(8)

p|z=d=pf

(9)

(10)

(11)

式中：ε为CRUD孔隙率；ρw为水的密度；κ为CRUD渗透率；μw为冷却剂动力黏度；pf为冷却剂主流区压力。

根据达西定律，水在多孔介质中的速度与压力梯度呈正比，如式(12)所示。

(12)

1.3 溶质扩散与化学反应模块

冷却剂在毛细作用下以ul的速度进入CRUD，同时将H3BO3和Li+等溶质带入CRUD中。在冷却剂流动带来的对流、化学势驱动下的扩散和电场对带电粒子作用下的迁移三者的共同作用下，溶质不断被浓缩，同时维持着溶质之间的化学平衡[11,16]。各溶质的浓度控制方程和边界条件如式(13)～(17)所示。

(13)

(14)

C|z=d=Cf

(15)

(16)

(17)

式中：控制方程右侧第1项代表化学反应引起的浓度变化，第2项代表输运引起的浓度变化；对于稳态情况下，式(13)左右侧均等于0；J为溶质摩尔流量；Cf为冷却剂主流区溶质浓度；nR为径向方向单位向量。

输运过程由能斯特-普朗克方程控制，如式(18)所示。

(18)

式中：D为扩散系数，需要基于分形理论对纯水中的扩散系数进行修正[7]，如式(19)所示，其中τ为曲折度；z为溶质电荷数；F为法拉第常数；Φ为电势。

(19)

(20)

式中：kf和kr分别为正、逆反应速率常数；对于稳态，式(20)左右侧等于0，此时反应物和生成物满足化学平衡，如式(21)所示。

(21)

式中：m为质量摩尔浓度；γ为活度系数，受到温度和浓度的影响；ζ为溶液体积修正项，是浓度的函数。因此，化学平衡状态下的正逆反应速率常数K可以表示为式(21)。

本文所考虑的化学反应主要包括水的电离、硼酸的形式转化和沉淀的析出和溶解。

水的电离：

(22)

硼酸的形式转化：

(23)

(24)

(25)

沉淀的析出和溶解：

4H3BO3+2Li+

(26)

1.4 数值计算方法

采用有限体积法进行柱坐标下的离散，使用温度、压力和各溶质浓度中的最大相对误差作为收敛条件，允许的最大相对误差小于10-7。上述三模块之间互相影响，在耦合计算过程中互相反馈直至最终计算结果达到收敛。CRUD传热传质计算流程图如图2所示。传热模块获得温度分布，并为毛细流动模块和溶质扩散与化学反应模块更新物性。毛细流动模块得到的压力和流速分别用于更新蒸汽物性和确定溶质传输的流速。溶质扩散与化学反应模块得到的浓度更新饱和温度进而影响温度和压力分布。

图2 CRUD传热传质计算流程图Fig.2 Flow chart of simulation process of CRUD heat and mass transfer

2 对比验证

Haq等[10]在Cohen多孔介质模型和Pan模型的基础上，实现了CRUD内的温度场、流动场和浓度场的耦合计算。使用Haq等相同的数据(表1)，获得的结果与Haq等的结果进行比较，如图3、4所示。对比显示，本文获得的结果在趋势上与Haq等的计算结果一致，但在数值上略有差别，最大相对误差为0.83%。这主要是由于导热系数上的差别，算例基于分形理论确定了CRUD的导热系数，而Haq等在算例中的CRUD导热系数采用了定值0.506 W/(m·K)。

表1 对比验证计算参数Table 1 Value of parameter used in validation analysis

图3 多孔壳体中心轴向方向温度变化Fig.3 Axial temperature in shell center

图4 CRUD与包壳表面径向方向温度变化Fig.4 Radial temperature at CRUD and cladding interface

计算结果显示，CRUD内温度最大位置出现在多孔沉积层中心靠近包壳的区域，略高于烟囱表面区域的温度，这主要是由于毛细流在烟囱表面发生沸腾带走了部分热量。整个CRUD区域内，烟囱表面靠近包壳区域的压力最小，但硼酸和锂浓度都在此区域内达到最大值。

3 结果讨论

为评估CRUD形貌参数对CRUD的传热传质的影响，使用多组CRUD形貌参数和典型压水堆热工参数进行了计算，如表2所列。

表2 输入参数Table 2 Value of input parameter

采用表中标准输入参数得到的温度、压力、硼酸浓度分布如图5所示。CRUD的温度在与包壳交界面处达到最大值(634.60 K)，与相同r位置CRUD与冷却剂交界面的温差达到18.74 K。r方向上温差均较小，最大值仅为1.20 K。最大SNB热流密度发生在烟囱表面靠近包壳区域内，达到0.31 MW/m2。同时这一区域也是最大压降的位置，达到1.8 kPa，使得r方向流速达到0.52 cm/s。将总能量分为通过烟囱表面SNB带走的能量和通过CRUD传导至CRUD与冷却剂表面通过冷却剂对流换热带走的能量两部分，将两者占总能量的比例分别称为SNB占比和热传导占比，在分析单元内，SNB占比达到75.07%。

图5 温度(a)、压力(b)和硼酸浓度(c)分布Fig.5 Distributions of temperature (a), pressure (b), and boric acid concentration (c)

3.1 厚度的影响

图6 硼质量和SNB占比随厚度的变化Fig.6 Boron mass and SNB ratio vs. CRUD thickness

压水堆包壳表面发现的CRUD厚度变化范围为20～90 μm，与压水堆所处的热工和水化学条件、运行时间和冷却剂杂质去除情况有关。随CRUD厚度变化的单位面积内硼质量和SNB占比如图6所示。随着CRUD厚度的增加，硼也不断得到积累，质量不断增加，在厚度达到40 μm时开始析出Li2B4O7，这与EPRI[2]报道的硼析出大致发生在35～42 μm相符合。在发生硼析出后，可溶性硼的增加趋势变缓，更多硼以沉淀的形式析出。同时，厚度达到40 μm后SNB占比有微弱增大。

3.2 孔隙率的影响

针对压水堆内CRUD的表征分析表明CRUD的孔隙率约为0.4～0.8，并且会随着硼和金属离子的沉淀反应、金属氧化物和氢气的氧化还原反应等变化，但在现有模拟中一般不考虑这些反应对孔隙率的影响。随孔隙率变化的最大温差、SNB占比和最大硼酸浓度如图7所示。由于水的导热系数小于构成CRUD多孔骨架的金属氧化物，因此孔隙率越高，导热系数越低，温差越大。SNB占总能量的比例与孔隙率呈近似线型比例关系，这可能是由于孔隙率越大，毛细热管内的冷却剂越多，相应地通过SNB带走的热量也更多。采用分形理论修正后的扩散系数随孔隙率变大而变大，在冷却剂毛细流动流速较小且变化不大的情况下，更大的扩散系数仅需较小的扩散梯度即可与硼酸的对流建立平衡，因此孔隙率越高硼酸浓度越低。

3.3 烟囱内径的影响

沉积层内的烟囱对于CRUD传热具有极其重要的意义，因为绝大部分的蒸发相变都发生在烟囱表面。研究表明，CRUD烟囱的直径范围一般为2～5 μm。随烟囱内径变化的最大温差、SNB占比、最大压降和最大硼酸浓度如图8所示。对于固定的烟囱密度，烟囱内径的增大意味着“湿区”所占面积份额的减小和烟囱表面相变面积的增大，因此相应的SNB占比也相应增大，热传导占比减小，最大温差也减小。SNB的增加，也引起了最大压降的增大，毛细流动流速增大，这最终引起了最大硼酸浓度的增大。

图7 最大温差(a)、SNB占比(b)和最大硼酸浓度(c)随孔隙率的变化Fig.7 Maximum temperature difference (a), SNB ratio (b), and boric acid concentration maximum (c) vs. porosity

图8 最大温差(a)、SNB占比(b)、最大压降(c)和最大硼酸浓度(d)随烟囱内径的变化Fig.8 Maximum temperature difference (a), SNB ratio (b), maximum pressure drop (c), and boric acid concentration maximum (d) vs. chimney diameter

图9 最大温差(a)、SNB占比(b)、最大压降(c)和最大硼酸浓度(d)随烟囱密度的变化Fig.9 Maximum temperature difference (a), SNB ratio (b), pressure drop maximum (c), and boric acid concentration maximum (d) vs. chimney density

3.4 烟囱密度的影响

CRUD内烟囱密度的变化范围较大，在不同反应堆、堆芯不同位置甚至同一燃料棒不同方位上都有所不同。一般来说，烟囱的密度范围为330～4 500 mm-2。随烟囱密度变化的最大温差、SNB占比、最大压降和最大硼酸浓度如图9所示。随着单位面积内烟囱数量的增大，CRUD的最大温差减小，SNB占比、最大压降、最大硼酸浓度都相应增大，这一趋势也与烟囱内径增大所引起的趋势相同，均与单位面积内SNB增强有关。

4 结论

本研究建立了一种用于CRUD内传热传质现象的预测方法，实现了CRUD内传热、流动、溶质输运和化学反应等多物理现象的耦合计算。此方法合理地预测了CRUD内的温度、压力、流速和浓度分布。

在所分析的多物理现象中，在CRUD独特结构烟囱表面发生的沸腾相变对传热传质均具有重要影响，SNB越强，CRUD的最大温差越小，最大压降越大，硼酸浓度也越大。通过参数敏感性分析，明晰了随厚度的增大，硼的存在形式从可溶性硼向Li2B4O7转化的趋势；解明了影响CRUD内传热传质的孔隙率、烟囱内径和烟囱密度对SNB均具有正相关性，而对硼浓度分别具有负、正、正相关性。该研究有利于进一步理解CRUD内硼吸附现象，可为堆芯功率漂移现象CIPS的预测提供技术支持。