CN110717275B - 一种针对压水堆堆芯的三维中子通量数值模拟方法 - Google Patents

Abstract

本发明公开了一种针对压水堆堆芯的三维中子通量数值模拟方法。该方法首先将所要模拟的三维压水堆堆芯沿轴向划分若干层，针对每一层基于特征线方法建立二维的中子输运模型；然后将所要模拟的三维压水堆堆芯基于径向的栅元几何划分为若干长条，针对每一长条基于离散纵标方法建立一维的中子输运模型；通过JFNK方法迭代二维中子输运模型和一维中子输运模型构成的残差模型至收敛，获得三维压水堆堆芯的中子通量分布。本发明相比现有技术，将二维中子输运模型和一维中子输运模型转化为残差模型同时求解，具有二阶收敛速度，迭代求解稳定性好，可用于数值反应堆的输运模块计算，提高数值反应堆输运计算效率，增加计算稳定性，节约数值计算产生的核时成本。

Description

一种针对压水堆堆芯的三维中子通量数值模拟方法

技术领域

本发明涉及核反应堆堆芯设计和安全领域，具体涉及一种针对压水堆堆芯进行三维中子通量数值模拟方法。

背景技术

随着核电行业的不断发展，为了应对核电厂延寿及安全分析，对压水堆堆芯的数值模拟的精度要求和效率要求越来越高。为了满足压水堆数值模拟的精度要求，数值反应堆这类高保真方法越来越多的应用在压水堆的数值模拟当中。

高保真方法采用全堆芯直接求解，计算的网格数目较多，其输运模块的计算三维中子输运模型的负担较大。目前常采用横向积分的方式，将三维中子输运模型转化为二维中子输运模型和一维中子输运模型分别计算，并通过泄漏项对两个模型进行耦合。

由于泄漏项的存在，使得二维一维数值模拟过程中的中子源为负，分别计算一维中子输运模型和二维中子输运模型的过程中中子通量会出现不符合物理规律的负值，导致数值模拟计算发散；且分别计算一维中子输运模型和二维中子输运模型，再通过泄漏项进行耦合，使得两个模型的耦合过程不够紧密，导致迭代次数较多，计算效率较低。

综上所述，为了稳定高效的对压水堆堆芯进行数值模拟，有必要将一维中子输运模型和二维中子输运模型紧密耦合的计算。基于无需显式构造雅可比矩阵的牛顿克罗诺夫方法(JFNK方法)，可以构造两个模型的残差模型，通过迭代计算残差模型将两个模型紧密耦合，避免数值模拟过程中的负中子源和负中子通量的问题，增加计算过程的稳定性和计算效率。

发明内容

为了克服上述现有技术中存在的问题，本发明提供一种针对压水堆堆芯进行三维中子通量数值模拟方法。该方法可用于传统的一维中子输运模型和二维中子输运模型分别数值模拟时难以收敛的问题，与传统方法相比，该方法采用JFNK方法计算两个输运模型构成的残差模型，可同时计算两个输运模型，具有二阶收敛速度，计算效率高，且只要初始值选取再接近真解的附近，就一定能够收敛，稳定性好。

为了实现上述目的，本发明采用了以下技术方案予以实施：

一种针对压水堆堆芯的三维中子通量数值模拟方法，包括如下步骤：

步骤1：读取所要模拟的压水堆堆芯的几何信息、材料信息和边界条件；

步骤2：根据步骤1获取的压水堆堆芯的几何信息和材料信息，将压水堆堆芯沿轴向划分为若干层，针对每一层基于特征线MOC方法建立二维平面中子输运模型；将压水堆堆芯基于径向栅元划分为长条，针对每一长条基于离散纵标S_N方法建立一维中子输运模型；如式(1)；

式中：

g——一维/二维中子输运模型的能群编号；

i——二维中子输运模型的平源区编号；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

p——一维中子输运模型的栅元编号；

φ^k _g,i——第g群，第k层，平源区i的二维中子通量；

ψ^k _g,m,i——第g群，第k层，平源区i，角度m方向的二维中子通量；

MOC_2D——特征线方法计算二维中子通量的迭代格式；

φ^p _g,k——栅元p，第g群，第k层的一维中子通量；

ψ^p _g,m,k——栅元p，第g群，第k层，角度m方向的一维中子角通量；

S_N_1D——离散纵标方法计算一维中子通量的迭代格式；

步骤3：将二维平面中子输运模型和一维中子输运模型作为黑箱，构建关于中子通量的残差模型R(Φ)，如式(2)所示；

式中：

R(Φ)——两个中子输运模型的残差模型

Φ——中子通量构成的向量形式；

g——一维/二维中子输运模型的能群编号；

G——一维/二维中子输运模型的能群总数；

i——二维中子输运模型的平源区编号；

I——二维中子输运模型的平源区总数；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

K——一维/二维中子输运模型的层数总数；

p——一维中子输运模型的栅元编号；

P——一维中子输运模型的栅元总数；

步骤4：根据步骤1获取的压水堆堆芯的边界条件，通过无需显式构造雅可比矩阵的牛顿克罗诺夫方法计算残差方程组(3)，获得第l个牛顿迭代步的收敛的中子通量残差ΔΦ^l*；通过式(4)获得第l+1个牛顿迭代步的中子通量Φ^l+1；

R(Φ^*)＝R(Φ^l)+J_RΔΦ^l*＝0 (3)

n→∞,ΔΦ^l,n→ΔΦ^l*

Φ^l+1＝Φ^l+ΔΦ^l* (4)

l→∞,Φ^l→Φ^*

式中:

l——牛顿迭代步，l＝0,1,2,3,……；

n——线性迭代步，n＝0,1,2,3,……；

J_R——残差方程组R(Φ^*)的雅可比矩阵；

Φ^*——中子通量的收敛值；

Φ^l——第l个牛顿迭代步中的中子通量(l＝0是人为假设的中子通量的初始值)；

Φ^l+1——第l+1个牛顿迭代步中的中子通量(如式(4)所示)；

ΔΦ^l,n——第l个牛顿迭代步中的第n个线性迭代步的中子通量残差(n＝0是人为假设的中子通量残差的初始值)；

ΔΦ^l*——线性迭代收敛后第l个牛顿迭代步的中子通量残差；

ε——差分步长；

步骤5：通过式(5)判断中子通量Φ^l+1是否满足收敛条件，若满足收敛条件，则获得中子通量的收敛值Φ^*，若不满足收敛条件，则重复步骤4至步骤5；

式中：

R(Φ^l+1)——第l+1牛顿迭代步的中子通量的残差模型；

R(Φ^l)——第l+1牛顿迭代步的中子通量的残差模型；

ξ——压水堆堆芯内中子通量的相对收敛限值。

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

本发明采用JFNK方法将传统的分别计算一维中子输运模型和二维中子输运模型的过程转化迭代计算两个输运模型构成的残差模型，可以同时计算两个输运模型。该方法具有二阶收敛速度，迭代次数少，计算效率高，且只要初始值选取再接近真解的附近，就一定能够收敛，稳定性好。

附图说明

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

图2是压水堆堆芯二维一维计算对象示意图。

具体实施方式

下面结合附图和具体实施方式对本发明做进一步详细说明。

具体步骤如图1所示。本发明将所要模拟的三维压水堆堆芯沿轴向划分建立二维平面中子输运模型；沿径向栅元划分为长条，建立一维中子输运模型；基于两个中子输运模型建立残差模型，通过JFNK方法，计算残差方程组，同时计算两个中子输运模型，获得压水堆的中子通量。具体步骤如下：

步骤1：读取所要模拟的压水堆堆芯的几何信息、材料信息和边界条件；

步骤2：根据步骤1获取的压水堆堆芯的几何信息和材料信息，将压水堆堆芯沿轴向划分为若干层，针对每一层基于特征线MOC方法建立二维平面中子输运模型；将压水堆堆芯基于径向栅元划分为长条，针对每一长条基于离散纵标S_N方法建立一维中子输运模型；如式(1)；

式中：

g——一维/二维中子输运模型的能群编号；

i——二维中子输运模型的平源区编号；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

p——一维中子输运模型的栅元编号；

φ^k _g,i——第g群，第k层，平源区i的二维中子通量；

ψ^k _g,m,i——第g群，第k层，平源区i，角度m方向的二维中子通量；

MOC_2D——特征线方法计算二维中子通量的迭代格式；

φ^p _g,k——栅元p，第g群，第k层的一维中子通量；

ψ^p _g,m,k——栅元p，第g群，第k层，角度m方向的一维中子角通量；

S_N_1D——离散纵标方法计算一维中子通量的迭代格式；

步骤3：将二维平面中子输运模型和一维中子输运模型作为黑箱，构建关于中子通量的残差模型R(Φ)，如式(2)所示；

式中：

R(Φ)——两个中子输运模型的残差模型

Φ——中子通量构成的向量形式；

g——一维/二维中子输运模型的能群编号；

G——一维/二维中子输运模型的能群总数；

i——二维中子输运模型的平源区编号；

I——二维中子输运模型的平源区总数；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

K——一维/二维中子输运模型的层数总数；

p——一维中子输运模型的栅元编号；

P——一维中子输运模型的栅元总数；

步骤4：根据步骤1获取的压水堆堆芯的边界条件，通过无需显式构造雅可比矩阵的牛顿克罗诺夫方法计算残差方程组(3)，获得第l个牛顿迭代步的收敛的中子通量残差ΔΦ^l*；通过式(4)获得第l+1个牛顿迭代步的中子通量Φ^l+1；

R(Φ^*)＝R(Φ^l)+J_RΔΦ^l*＝0 (3)

n→∞,ΔΦ^l,n→ΔΦ^l*

Φ^l+1＝Φ^l+ΔΦ^l* (4)

l→∞,Φ^l→Φ^*

式中:

l——牛顿迭代步，l＝0,1,2,3,……；

n——线性迭代步，n＝0,1,2,3,……；

J_R——残差方程组R(Φ^*)的雅可比矩阵；

Φ^*——中子通量的收敛值；

Φ^l——第l个牛顿迭代步中的中子通量(l＝0是人为假设的中子通量的初始值)；

Φ^l+1——第l+1个牛顿迭代步中的中子通量(如式(4)所示)；

ΔΦ^l,n——第l个牛顿迭代步中的第n个线性迭代步的中子通量残差(n＝0是人为假设的中子通量残差的初始值)；

ΔΦ^l*——线性迭代收敛后第l个牛顿迭代步的中子通量残差；

ε——差分步长；

步骤5：通过式(5)判断中子通量Φ^l+1是否满足收敛条件，若满足收敛条件，则获得中子通量的收敛值Φ^*，若不满足收敛条件，则重复步骤4至步骤5；

式中：

R(Φ^l+1)——第l+1牛顿迭代步的中子通量的残差模型；

R(Φ^l)——第l+1牛顿迭代步的中子通量的残差模型；

ξ——压水堆堆芯内中子通量的相对收敛限值。

图2展示的是，压水堆堆芯采用二维一维耦合计算时，一维中子输运模型和二维平面中子输运模型的几何示意图及栅元的平源区划分示意图。

Claims (1)

1.一种针对压水堆堆芯的三维中子通量数值模拟方法，其特征在于：包括如下步骤：

步骤1：读取所要模拟的压水堆堆芯的几何信息、材料信息和边界条件；

步骤2：根据步骤1获取的压水堆堆芯的几何信息和材料信息，将压水堆堆芯沿轴向划分为若干层，针对每一层基于特征线MOC方法建立二维平面中子输运模型；将压水堆堆芯基于径向栅元划分为长条，针对每一长条基于离散纵标S_N方法建立一维中子输运模型；如式(1)；

式中：

g——一维/二维中子输运模型的能群编号；

i——二维中子输运模型的平源区编号；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

p——一维中子输运模型的栅元编号；

φ^k _g,i——第g群，第k层，平源区i的二维中子通量；

ψ^k _g,m,i——第g群，第k层，平源区i，角度m方向的二维中子通量；

MOC_2D——特征线方法计算二维中子通量的迭代格式；

φ^p _g,k——栅元p，第g群，第k层的一维中子通量；

ψ^p _g,m,k——栅元p，第g群，第k层，角度m方向的一维中子角通量；

S_N_1D——离散纵标方法计算一维中子通量的迭代格式；

步骤3：将二维平面中子输运模型和一维中子输运模型作为黑箱，构建关于中子通量的残差模型R(Φ)，如式(2)所示；

式中：

R(Φ)——两个中子输运模型的残差模型

Φ——中子通量构成的向量形式；

g——一维/二维中子输运模型的能群编号；

G——一维/二维中子输运模型的能群总数；

i——二维中子输运模型的平源区编号；

I——二维中子输运模型的平源区总数；

m——一维/二维中子输运模型的特征线的方向编号；

k——一维/二维中子输运模型的层数编号；

K——一维/二维中子输运模型的层数总数；

p——一维中子输运模型的栅元编号；

P——一维中子输运模型的栅元总数；

步骤4：根据步骤1获取的压水堆堆芯的边界条件，通过无需显式构造雅可比矩阵的牛顿克罗诺夫方法计算残差方程组(3)，获得第l个牛顿迭代步的收敛的中子通量残差ΔΦ^l*；通过式(4)获得第l+1个牛顿迭代步的中子通量Φ^l+1；

Φ^l+1＝Φ^l+ΔΦ^l*

l→∞,Φ^l→Φ^* (4)

式中:

l——牛顿迭代步，l＝0,1,2,3,……；

n——线性迭代步，n＝0,1,2,3,……；

J_R——残差方程组R(Φ^*)的雅可比矩阵；

Φ^*——中子通量的收敛值；

Φ^l——第l个牛顿迭代步中的中子通量，l＝0是人为假设的中子通量的初始值；

Φ^l+1——第l+1个牛顿迭代步中的中子通量；

ΔΦ^l,n——第l个牛顿迭代步中的第n个线性迭代步的中子通量残差，n＝0是人为假设的中子通量残差的初始值；

ΔΦ^l*——线性迭代收敛后第l个牛顿迭代步的中子通量残差；

ε——差分步长；

步骤5：通过式(5)判断中子通量Φ^l+1是否满足收敛条件，若满足收敛条件，则获得中子通量的收敛值Φ^*，若不满足收敛条件，则重复步骤4至步骤5；

式中：

R(Φ^l+1)——第l+1牛顿迭代步的中子通量的残差模型；

R(Φ^l)——第l+1牛顿迭代步的中子通量的残差模型；

ξ——压水堆堆芯内中子通量的相对收敛限值。

Images