一种核反应堆堆芯棒束通道内硼扩散过程的计算方法

技术领域

本发明属于核反应堆安全分析技术领域，具体涉及一种应用于核反应堆堆芯棒束通道内冷却剂中硼扩散过程及硼浓度分布特性的计算方法。

背景技术

在压水反应堆中，硼由于其中子吸收截面较大，作为溶质溶解在水中可以用于调节堆芯反应性，进而调节堆芯功率，是压水堆的重要组成部分。在核反应堆启动阶段，稀释堆芯硼浓度就是一个非常重要的步骤。在堆芯硼浓度稀释的过程中，必然会导致堆芯硼浓度分布不均匀，然而由于核反应堆内各部分硼浓度难以监测，在硼稀释过程中的不当操作可能导致核反应性瞬态事故或堆芯局部功率过大而导致的堆芯损毁事故。因此，进行核反应堆堆芯硼浓度稀释过程的模拟分析对其设计和运行具有十分重要的意义。

目前，学界对核反应堆安全分析的方法一般有三种：1.系统分析，如RELAP5程序，TRACE程序等一维系统分析程序；2.计算流体动力学(CFD)模拟，如ANSYS Fluent软件等；3.子通道分析程序，如COBRA 程序等。其中，一维系统分析程序无法对核反应堆堆芯棒束通道内横向通道内冷却剂流动换热现象进行分析计算；CFD方法建模复杂，瞬态计算耗时长，两相流模型不完善，用于核反应堆安全分析有较大局限性；子通道分析方法建模较简单，对计算要求低，可以实现核反应堆堆芯棒束通道内(准)三维热工水力计算，是一种重要的核反应堆安全分析方法。

现有的子通道分析程序在进行硼扩散过程计算时，只考虑了各控制体间由于对流导致的硼输运过程。这种简单的输运模型只适用于计算层流条件下溶质的输运扩散。在核反应堆正常运行阶段，一般堆芯棒束通道内冷却剂流速较大，呈现湍流特征，这种情况下溶质湍流扩散现象非常明显。而且湍流扩散主要影响堆芯棒束通道内冷却剂溶质横向输运过程，这一过程对硼浓度在堆芯棒束通道内均匀分布非常重要。因此，现有的子通道分析方法在复杂的湍流流动条件下并不能准确计算硼溶质在各个棒束通道内的实际扩散过程。

2011年Ozdemir对核反应堆LOCA事故下可能发生的冷却剂稀释事故基于商用子通道程序COBRA‐TF建立了子通道硼扩散模型，使用子通道程序对堆芯棒束通道内硼浓度分布变化过程进行计算模拟。基于子通道坐标系，硼扩散方程表示如下：

上式中：

c_b——硼质量分数，ppm；

A——子通道截面积，m²；

v_f——流体轴向流速，m/s；

w_f——流体横向流速，m/s；

L_g——横向流面积，m²。

以上方程按图2所示控制体分布进行方程离散如下：

上式中：

e_ik——用于判定横向流动方向的函数；

k——与i通道相邻的通道轴向编号；

A——流通面积，m²；

ΔX——轴向控制体长度，m；

ρ——冷却剂密度，kg/m²；

c——硼浓度，ppm；

m——冷却剂轴向质量流速，kg/s；

w——冷却剂横向质量流速，kg/s；

r——液相中硼含量(使用均匀流模型时，r＝1)；

上述方程基于子通道流场进行求解，使子通道程序具有分析硼扩散过程的能力。但是，上述方程只考虑了相邻通道间净质量交换导致的硼质量交换过程，而没有考虑湍流的影响。在雷诺数较高时，由于湍流的存在，流体的流线不再清晰可辨，流场中有许多小漩涡，相邻流层间不但有滑动，还有混合，这种混合在宏观上并没有流体的净质量交换，但是对溶质的浓度分布影响很大。

发明内容

为了解决上述问题，本发明提供了一种核反应堆堆芯棒束通道内硼扩散过程的计算方法，该方法基于现有子通道方法，在流场计算的基础上添加硼组分输运方程，能更精确模拟核反应堆内硼扩散过程。

为了达到上述目的，本发明的采用如下技术方案：

一种核反应堆堆芯棒束通道内硼扩散过程的计算方法，该计算方法包括如下步骤：

步骤1：根据模拟对象确定堆芯棒束通道几何尺寸，进行子通道划分，并给定边界条件，包括进口流体温度、流量、硼浓度和出口压力、加热功率，给定初始条件，包括各控制体内流体压力、温度、流速、硼浓度；

步骤2：每一个时间步长使用子通道分析方法求解子通道内流体的质量守恒方程、动量守恒方程和能量守恒方程，以得到堆芯棒束通道内的流场状态参数；

步骤3：建立考虑硼扩散过程的组分输运方程

式中：

i——通道编号；

j——通道层数编号；

A_i——通道i流通截面积，m²；

ΔX_j——第j层控制体长度，m；

ρ_ij——通道i第j层控制体中流体密度，kg·m^‐3；

c_ij——通道i第j层控制体中硼浓度，ppm；

m_ij-1——i通道第j‐1层控制体轴向质量流速，kg/s；

——上一个时间步长通道i第j‐1层控制体中硼浓度，ppm；

w_k——横流流量，kg/s；

——与通道i第j层控制体有质量交换的控制体上一个时间步长的硼浓度，ppm；

w'_k——湍流交混系数，通过湍流交混模型获得；

上标*——表示上一个时间步长的状态参数；

上式中等号左边项为控制体内硼质量；右边第一项为从上一控制体轴向流入的硼质量；右边第二项为从控制体轴向流出的硼质量；右边第三项为横向净流量，当w_k为正即流体从i通道j层控制体横向流出时，k＝i，当w_k为负时即流体从周边向i通道j层控制体横向流入，为与通道i第j层控制体有质量交换的控制体上一个时间步长的硼浓度；第四项是湍流交混项，用于计算由于湍流交混导致的硼质量扩散，其中c_kj取为i通道j层控制体及其周边控制体内硼浓度；

根据步骤2中计算结果获得流场状态参数，即ρ_ij，m_ij-1，m_ij为已知量，结合上一时间步长的状态参数及周边控制体硼浓度求解考虑湍流扩散现象的组分扩散方程，以获得组分浓度分布状态及其变化过程；假设横向流向为i通道向n通道流动，其组分扩散方程为：

求解得t时刻i通道j层控制体内硼浓度为：

式中：

Δt——时间步长；

从上式能够看出，控制体内硼浓度计算式需要本时间步长内的参数，在求解时，经多次迭代求解得到计算结果，在该时间步长内经过大量迭代直到计算收敛，收敛条件为：

[c_ij(m)-c_ij(m-1)]/c_ij(m-1)<δ

上式中：

m——迭代次数；

c_ij(m)——第m次迭代计算的通道i第j层控制体内硼浓度，ppm；

c_ij(m-1)——第m‐1次迭代计算的通道i第j层控制体内硼浓度，>

δ——收敛因子，设定为0.0001；

在核反应堆中，硼在冷却剂中所占的质量分数小，对冷却剂的物性影响小，因此在求解组分输运方程时，为简化计算，忽略了组分浓度对冷却剂物性参数的影响；组分输运方程在考虑了组分的湍流扩散现象后，能更精确的模拟硼在核反应堆堆芯棒束通道内的准三维扩散过程；

步骤4：开展下一时间步长计算，重复步骤(2)和步骤(3)，直至达到要求的计算时间跨度。

本发明与现有技术相比，具有如下优点：

1.克服了现有子通道计算方法在计算核反应堆堆芯棒束通道内硼溶质扩散过程的缺陷，考虑了流体中实际存在的湍流扩散现象对硼扩散过程的影响；

2.子通道划分方式简单，易操作，容易实现核反应堆堆芯棒束通道内硼扩散的过程模拟；

3.能够实现不同运行工况下核反应堆堆芯棒束通道内硼扩散过程的模拟；

4.该方法可以适用于多种溶质在核反应堆堆芯棒束通道内的流动搅浑过程的模拟；

5.模型独立，方法通用性强，可以适用于多种子通道分析程序；

6.模型建立周期短、计算资源消耗少且计算结果精度高，用于各种硼浓度在堆芯棒束通道内的分布特性计算。

本发明已通过实践证明，该方法能够成功模拟压水堆堆芯棒束通道内高雷诺数棒束通道内的硼扩散过程，本发明中提出的棒束通道内的湍流扩散方程在进行棒束通道内的硼扩散过程计算时，不但可以避免传统子通道分析方法无法模拟因湍流扩散导致的硼横向扩散过程，同时可以保证采用该模型获得数值结果的精确性，完全能够克服已有方法对堆芯棒束通道内硼扩散分析计算的缺陷。

附图说明

图1是本发明方法的流程图。

图2是棒束子通道模型控制体划分图。

图3是帮束子通道模型计算示意图

图4是棒束子通道模型横向流动示意图。

具体实施方式

下面结合附图和具体实施方式对本发明方法进行详细的说明。

如图1所示，本发明提供了一种核反应堆堆芯棒束通道内硼扩散过程的计算方法，具体方法如下：

步骤1：根据模拟对象确定堆芯棒束通道几何尺寸，进行子通道划分，如图2、图3所示，并给定边界条件，包括进口流体温度、流量、硼浓度和出口压力、加热功率，给定初始条件，包括各控制体内流体压力、温度、流速、硼浓度；

步骤2：每一个时间步长使用子通道分析方法求解子通道内流体的质量守恒方程、动量守恒方程和能量守恒方程，以得到堆芯棒束通道内的流场状态参数；

子通道是棒束之间流道的自然几何划分，它以燃料棒本身和燃料棒之间的假想连线所包围的流动面积定义为一个子通道的横截面积。如图2所示，流体在这样的流道中流动，一面与周围的加热棒进行能量和动量交换，一面通过假想边界与相邻通道进行质量、能量和动量交换。为了描述这种多分支的流动，不采用常规的二维或三维坐标系统，而是采用局部坐标系，规定子通道间横向流动的流出为正，其正方向为间隙的法向方向，同时，子通道沿轴向也将划分若干控制体，对于通道n与通道i之间的间隙k，径向引入通道间隙宽度s_k、交混长度l_k的概念以方便计算，如图4所示，通道间隙宽度s_k为两相邻通道间隙的长度，间隙的高度为控制体高度，交混长度l_k为两通道中心点距离。子通道内流体的控制方程组如下：

(1)质量守恒方程

式中：

i,n——相邻通道编号；

A_i——通道i流通截面积，m²；

ρ_i——通道i流体密度，kg/m³；

m_i——通道i轴向质量流量，kg/s；

w_in——通道i与n间横向流量，kg/s；

(2)能量守恒方程

式中：

h_i——通道i的流体比焓，J/kg；

h_n——通道n的流体比焓，J/kg；

h——通道间横流交混的比焓，J/kg；

——通道i内燃料线功率，W/m；

w′_in——湍流交混质量流量，kg/s；

s_in——通道i与n间隙宽度，m；

l_in——通道i与n间交混长度，m；

T_i——通道i的温度，℃；

T_n——通道n的温度，℃；

K——流体的导热率，W/(m·℃)。

(3)轴向动量守恒方程

式中：

p_i——通道i的压力，Pa；

g——重力加速度，m/s²；

θ——通道倾斜角度，°；

f_T——横向动量因子；

u_i——通道i的轴向流速，m/s；

u_n——通道n的轴向速度，m/s；

u——通道i与n间横向交混供体流速，m/s；

f_i——通道i的摩擦系数；

D_i——通道i的等效水力直径，m；

K_s——局部阻力系数；

Δx——轴向控制体高度，m。

(4)横向动量守恒方程

式中：

s_in——通道i与n间隙宽度，m；

l_in——通道i与n间交混长度，m；

K_G——通道i与n间横流阻力系数；

ρ^*——通道i与n间横流交混供体密度，kg/m³。

在子通道分析计算区域对上述方程进行离散化处理，将轴向流动流量m与压力p定义在控制体单元的边界上，而将状态参数如比焓h、密度ρ与空泡份额α等定义在控制体单元的中心。采用有限差分法，控制方程组中各个方程的离散方程如下：

质量守恒方程离散：

式中：

上标*——上一时间步长的数值；

A_i——通道i流通截面积，m²；

ρ_ij——通道i第j层控制体中流体密度，kg/m³；

w_in,j——第j层中通道i与n间横向流动；

Δx_j——第j层控制体高度，m。

能量守恒方程离散：

式中：

h_ij——通道i第j层控制体中流体比焓，，J/kg；

——通道间横流交混的比焓，J/kg，

C_in——表征通道间热传导率,W/(m²·℃)，C_in＝K/l_in。

轴向动量守恒方程离散：

式中：

——上一时间步长通道i与n间隙第j层的横向交混供体流速， m/s；

——通道i第j层控制体的两相倍增因子；

v′_ij——通道i第j层控制体的动量比体积，m³/kg，

u_ij——通道i第j层控制体轴向速度，m/s，；

u_nj——通道n第j层控制体轴向速度，m/s；

横向动量守恒方程离散：

式中：

下标j——控制体层数编号；

ρ——两通道间横流交混的供体密度,kg/m³；

K_G——通道间隙横流阻力系数。

所有网格的方程最后可以组成一个压力的三对角矩阵，它可以利用直接消去法求解，也可以利用Gauss‐Seidel迭代法求解。总体计算思路是：①用动量方程求解出速度；②利用上一时层求解压力方程；③用解出的压力回代求解速度和焓值。

步骤3：建立考虑硼扩散过程的组分输运方程

式中：

i——通道编号

j——通道层数编号

A_i——通道i流通截面积，m²；

ΔX_j——通道i第j层控制体长度，m；

ρ_ij——通道i第j层控制体中流体密度，kg·m^‐3；

c_ij——通道i第j层控制体中硼浓度，ppm；

m_ij-1——i通道第j层控制体轴向质量流速，kg/s；

w_k——横流流量，kg/s；

——上一个时间步长通道i第j‐1层控制体中硼浓度，ppm；

——与通道i第j层控制体有质量交换的控制体上一个时间步长硼浓度，ppm；

w'_k——湍流交混系数，通过湍流交混模型获得；

上标*——表示上一个时间步长的状态参数；

上式中第一项为控制体内硼质量；右边第一项为从上一控制体轴向流入的硼质量；右边第二项为控制体轴向流出项；右边第三项为横向净流量，当w_k为正(即流体从i通道j层控制体横向流出)时，k＝i，当w_k为负时(即流体从周边向i通道j层控制体横向流入)，为与通道i第j层控制体有质量交换的控制体上一个时间步长硼浓度；第四项是湍流交混项，用于计算由于湍流交混导致的硼质量扩散，其中c_kj可取为i通道j层控制体及其周边控制体内硼浓度，湍流交混系数w'_kj由湍流交混模型获得。在子通道中，湍流横流用于描述由于湍流导致的质量交混，目前使用较多的湍流交混模型有：

w_k′＝as_kG_k

式中：

w_k′——湍流横流，kg/(m·s)；

s_k——间隙宽度，m；

l_k——交混长度，m；

G_k——间隙的质量流速，kg/(m²·s)；

D_k——通道平均水力直径，m；

Re_k——雷诺数，其中参数为通道l和l′的平均值。

湍流交混模型是经验模型，在使用时需根据实际几何来选择湍流交混模型以获得准确结果。在硼扩散方程中，通过步骤2获得的流场参数可以直接求解方程，获得各控制体内浓度分布计算结果。

根据步骤2中计算结果获得流场状态参数，即ρ_ij，m_ij-1，m_ij为已知量，结合上一时间步长的状态参数及周边控制体硼浓度求解考虑湍流扩散现象的组分扩散方程，以获得组分浓度分布状态及其变化过程。对于如图3所示的双通道几何结构,假设横向流向为i通道向n通道流动，其组分扩散方程为：

求解得t时刻i通道j层控制体内硼浓度为：

式中：

Δt——时间步长；

从上式可以看出，控制体内硼浓度计算式需要本时间步长内的参数，在求解时，经多次迭代求解得到计算结果，在该时间步长内经过大量迭代直到计算收敛，收敛条件为：

[c_ij(m)-c_ij(m-1)]/c_ij(m-1)<δ

上式中：

m——迭代次数；

c_ij(m)——第m次迭代计算的通道i第j层控制体内硼浓度，ppm；

c_ij(m-1)——第m‐1次迭代计算的通道i第j层控制体内硼浓度，>

δ——收敛因子，设定为0.0001；

在核反应堆中，硼在冷却剂中所占的质量分数较小，对冷却剂的物性影响很小，因此在求解组分输运方程时，为简化计算，忽略了组分浓度对冷却剂物性参数的影响；组分输运方程在考虑了组分的湍流扩散现象后，能更精确的模拟硼在核反应堆堆芯棒束通道内的准三维扩散过程。

步骤4：开展下一时间步长计算，重复步骤(2)和步骤(3)，直至达到要求的计算时间跨度。