CN107066751B - 针对非均匀几何变分节块方法的平源加速方法 - Google Patents

Abstract

一种针对非均匀几何变分节块方法的平源加速方法，包括如下步骤：1、将有限元划分为平源区，在各个平源区内建立平均中子通量密度展开矩和中子通量密度分布的转换关系；2、求得平源区内部的平均中子通量密度展开矩和有限元节点上中子通量密度展开矩之间的关系；3、将转换关系代入原有的响应矩阵方程，将求解对象从各个有限元节点上的未知量转换为各个平源区的未知量；4、利用通用的迭代方法，实现非均匀几何变分节块方法的求解；本发明可以大幅度地减少计算时间和计算内存。

Description

针对非均匀几何变分节块方法的平源加速方法

技术领域

本发明针对核反应堆堆芯中子学计算领域，提出了一种针对非均匀几何变分节块方法的平源加速方法。

背景技术

核反应堆中子学计算研究以核反应堆堆芯为应用对象，其堆芯由许多不同种类的组件构成。根据堆型的不同，组件内部的几何结构和材料布置复杂多变。因此，实际的反应堆中子学问题是一个三维非均匀几何的中子学问题。对核反应堆进行快速、精确的中子学计算，是反应堆设计和校核的基本保障。目前在堆芯物理设计过程中，主要针对均匀几何的进行扩散方程的求解。

变分节块法是核反应堆物理设计中常用的方法之一。它以二阶偶宇称形式的中子扩散方程为出发点，方程呈现椭圆方程的形式，有利于Garlerkin方法的应用，更适合有限元方法的空间离散。变分节块法的计算思想是：首先通过变分方法在非均匀求解区域建立包含二阶中子输运方程和自然边界条件的泛函；然后采用标准正交多项式进行Ritz离散，同时利用球谐函数实现角度展开，并构造响应矩阵；最后分别在三维堆芯的各个节块内分别求解响应矩阵方程；节块之间以流和其高阶矩耦合，最终得到问题区域的中子通量密度分布。变分节块法能够达到较高的精度，但仅能处理均匀节块的问题。在基于均匀节块的中子学计算过程中，无法避免均匀化过程，以及随之而来的误差。随着科学技术的发展和计算机水平的提高，人们已经越来越开始重视减少近似和假设，追求更高精度的中子学计算方法，消除均匀化过程是中子学发展的必然趋势。因此，研究具备非均匀几何处理能力的计算方法对于中子学计算具有十分重要的意义。

基于变分节块法，后续研究在变分节块法的框架下提出了非均匀几何变分节块方法，通过引入有限元方法以精确描述栅元内部的非均匀结构。在此工作中，各栅元被单独处理为节块，节块内部采用有限元网格进行细化，精确描述栅元和冷却剂的材料、几何分布。为了能够实际描述栅元中的曲边几何结构，且保证足够的计算精度，往往需要大量的有限元网格剖分，因此计算量很大。以图1中所示的一个压水堆栅元的有限元剖分为例，栅元包含燃料和冷却剂两区，分别以红色和蓝色标示，黑点代表有限元网格的节点。采用二阶等参有限元对其进行非结构网格剖分，共包含32个有限元，对应97个有限元节点。在实际计算过程中，需要存储和求解反应堆堆芯中每一个栅元中97个有限元节点上的中子通量密度展开矩，会花费很大的计算代价。

发明内容

为了有效地减小非均匀几何变分节块方法的计算代价，本发明提出了一种针对非均匀几何变分节块方法的平源加速方法，可以大幅度地减少计算时间和计算内存。

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

一种针对非均匀几何变分节块方法的平源加速方法，通过将有限元划分为平源区，以提高计算效率和减小计算内存，步骤如下：

步骤1：在各平源区中，将各个有限元节点的源项近似为平源区内的平均源项，通过此近似减小非均匀求解区域中的自由度数目；首先给出公式(1)和公式(2)响应矩阵方程的表达式；根据公式(4)中的中子通量密度分布φ(x,y)的表达式，求得公式(5)中平源区的平均中子通量密度展开矩

和中子通量密度分布表达式φ(x,y)的关系：

在非均匀几何变分节块方法中，节块内部中子通量密度展开矩的求解方程为

其中：

-1—矩阵的求逆；

φ—节块内部中子通量密度展开矩向量，其中展开矩代表展开系数的值；

j—节块表面净中子流密度展开矩向量；q—中子源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关；

响应矩阵方程为：

式中：

j ⁺—出射中子流密度展开矩向量；

j ^-—入射中子流密度展开矩向量；

u—中子流源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关；

针对二维情况，已知中子源的空间分布表示为：

式中：

q(x,y)—中子源的空间分布；

Σ_sgg′—第g群到第g’群的中子散射截面；

g,g′—两个不同的能群标识；

φ_g′(x,y)—第g^/群的中子通量密度分布；

k_eff—有效增殖因子；

νΣ_fg′—第g’群的中子产生截面；

χ_g—第g群的中子裂变谱；

二维情况下，节块内部中子角通量密度的离散表达式与轴向无关：

φ(x,y)≈g ^T(x,y)φ 公式(4)其中：

φ(x,y)—节块内部中子通量密度的分布；

g(x,y)—x-y方向有限元形状函数向量；

T—向量或矩阵的转置；

φ—节块内部中子通量密度的展开矩向量；

在平源加速假设下，平源区的平均中子通量密度展开矩

能够表示为中子通量密度分布在本平源区内部的积分平均值：

式中：

—各平源区平均中子通量的矩向量；

h_t (x,y)—对应于每个有限元网格的分片常量：h_t(x,y)＝δ_e，e为每个有限元网格；

—对角线元素为每个有限元网格面积的对角矩阵；

步骤2：将公式(4)中节块内部中子角通量密度的离散表达式代入公式(5)中，就能够得到公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ的关系；

将公式(4)代入公式(5)得：

其中：

—平源区的平均中子通量密度展开矩；

同时

步骤3：利用公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ之间的转换关系，代入传统变分节块方法中响应矩阵方程的表达式公式(1)和公式(2)，将求解对象从各个有限元节点上的未知量转换为各个平源区的未知量，即公式(13)和公式(14)；

将公式(1)代入到公式(6)得：

此外，群内源项能够表示为：

假设每个平源区内有限元节点上的源项与平源区的源项相等：

式中：

Σ_x(x,y)—平源区内部不同位置处的x反应的截面，x可以为吸收，总截面，散射；

—各平源区x反应截面所组成的对角矩阵；

因此，将公式(10)代入公式(9)，得：

其中：

将公式(11)代入公式(8)得：

同理，公式(2)变为：

步骤4：利用通用的迭代方法，对公式(13)和公式(14)进行求解，从而得到整个非均匀求解区域的中子通量密度分布和中子流密度分布，实现非均匀几何变分节块方法的求解。

本发明通过将计算对象从各个有限元节点上中子通量密度展开矩向量φ，转变为求解各个平源区内平均中子通量密度展开矩向量展开矩向量，能够降低存储内存和计算时间，实现针对非均匀几何变分节块方法的平源加速方法。

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

1.本发明通过采用平源加速方法划分每个有限元为一个平源区，能够显著减少栅元中自由度的数目。以图1为例，若采用本发明中的方法，栅元中的自由度可以从97减少到32，从而显著降低中子通量密度展开矩的存储计算内存。

2.本发明通过利用公式(14)中的表达式，能够将计算对象从各个有限元节点上中子通量密度展开矩向量φ，转变为求解各个平源区内平均中子通量密度展开矩向量

减小解向量的长度，降低响应矩阵的规模，从而节省求解过程中的浮点数操作数目，显著提高计算效率。

附图说明

图1等参有限元描述下的压水堆栅元非均匀几何。

具体实施方式

下面结合具体实施方式对本发明作进一步详细说明：

本发明一种针对非均匀几何变分节块方法的平源加速方法，通过将有限元划分为平源区，以提高计算效率和减小计算内存，步骤如下：

步骤1：在各平源区中，将各个有限元节点的源项近似为平源区内的平均源项，通过此近似减小如图1所示的非均匀求解区域中的自由度数目；首先给出公式(1)和公式(2)中的响应矩阵方程的表达式；根据公式(4)中的中子通量密度分布φ(x,y)的表达式，求得公式(5)中平源区的平均中子通量密度展开矩

和中子通量密度分布表达式φ(x,y)的关系：

在非均匀几何变分节块方法中，节块内部中子通量密度展开矩的求解方程为

其中：

-1—矩阵的求逆；

φ—节块内部中子通量密度展开矩向量，其中展开矩代表展开系数的值；

j—节块表面净中子流密度展开矩向量；

q—中子源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关；

响应矩阵方程为：

式中：

j ⁺—出射中子流密度展开矩向量；

j ^-—入射中子流密度展开矩向量；

u—中子流源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关。

以二维情况为例，已知中子源的空间分布可以表示为：

式中：

q(x,y)—中子源的空间分布；

Σ_sgg′—第g群到第g’群的中子散射截面；

g,g′—两个不同的能群标识；

φ_g′(x,y)—第g’群的中子通量密度分布；

k_eff—有效增殖因子；

νΣ_fg′—第g’群的中子产生截面；

χ_g—第g群的中子裂变谱。

二维情况下，节块内部中子角通量密度的离散表达式与轴向无关：

φ(x,y)≈g ^T(x,y)φ 公式(4)其中：

φ(x,y)—节块内部中子通量密度的分布；

g(x,y)—x-y方向有限元形状函数向量；

T—向量或矩阵的转置；

φ—节块内部中子通量密度的展开矩向量；

在平源加速假设下，平源区的平均中子通量密度展开矩

能够表示为中子通量密度分布在本平源区内部的积分平均值：

式中：

—各平源区平均中子通量的矩向量；

h_t (x,y)—对应于每个有限元网格的分片常量：h_t(x,y)＝δ_e，e为每个有限元网格；

—对角线元素为每个有限元网格面积的对角矩阵。

步骤2：将公式(4)中节块内部中子角通量密度的离散表达式代入公式(5)中，就能够得到公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ的关系；

将公式(4)代入公式(5)得：

其中：

—平源区的平均中子通量密度展开矩；

同时

步骤3：利用公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ之间的转换关系，代入传统变分节块方法中响应矩阵方程的表达式公式(1)和公式(2)，将求解对象从各个有限元节点上的未知量转换为各个平源区的未知量，即公式(13)和公式(14)；

将公式(1)代入到公式(6)得：

此外，群内源项能够表示为：

假设每个平源区内有限元节点上的源项与平源区的源项相等：

式中：

Σ_x(x,y)—平源区内部不同位置处的x反应的截面，x可以为吸收，总截面，散射等；

—各平源区x反应截面所组成的对角矩阵；

因此，将公式(10)代入公式(9)，得：

其中：

将公式(11)代入公式(8)得：

同理，公式(2)变为：

步骤4：对于非均匀几何求解区域内不同种类的节块，分别计算各节块内部的空间相关响应矩阵；

步骤5：针对特定能群，利用红-黑扫描的方式对响应矩阵方程公式(14)进行迭代求解，得到出、入射偏中子流密度展开矩j ⁺、j ^-；

步骤6：在求得出、入射偏中子流密度展开矩j ⁺、j ^-后，代入公式(13)中，即可解得节块内部对应各个平源区的平均中子通量密度展开矩向量

从而得知节块内部的中子通量密度分布；

步骤7：进入下一能群的计算，最终求解整个非均匀几何求解区域的中子通量密度分布和中子流密度分布。通过将计算对象从各个有限元节点上中子通量密度展开矩向量φ，转变为求解各个平源区内平均中子通量密度展开矩向量

能够降低存储内存和计算时间，实现针对非均匀几何变分节块方法的平源加速方法。

验证结果显示，在保证精度的前提下，平源加速方法可将计算速度提高4倍以上，内存减少70％。

Claims (1)

1.一种针对非均匀几何变分节块方法的平源加速方法，其特征在于：通过将有限元划分为平源区，以提高计算效率和减小计算内存，步骤如下：

步骤1：在各平源区中，将各个有限元节点的源项近似为平源区内的平均源项，通过此近似减小非均匀求解区域中的自由度数目；首先给出公式(1)和公式(2)响应矩阵方程的表达式；根据公式(4)中的中子通量密度分布φ(x,y)的表达式，求得公式(5)中平源区的平均中子通量密度展开矩

和中子通量密度分布表达式φ(x,y)的关系：

在非均匀几何变分节块方法中，节块内部中子通量密度展开矩的求解方程为

其中：

-1—矩阵的求逆；

φ—节块内部中子通量密度展开矩向量，其中展开矩代表展开系数的值；

j—节块表面净中子流密度展开矩向量；q—中子源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关；

响应矩阵方程为：

式中：

j ⁺—出射中子流密度展开矩向量；

j ^-—入射中子流密度展开矩向量；

u—中子流源项展开矩向量；

—响应矩阵，二者表达形式不同，仅与节块内部的材料布置、几何形状有关；

针对二维情况，已知中子源的空间分布表示为：

式中：

q(x,y)—中子源的空间分布；

Σ_sgg′—第g群到第g’群的中子散射截面；

g,g′—两个不同的能群标识；

φ_g′(x,y)—第g^/群的中子通量密度分布；

k_eff—有效增殖因子；

νΣ_fg′—第g’群的中子产生截面；

χ_g—第g群的中子裂变谱；

二维情况下，节块内部中子角通量密度的离散表达式与轴向无关：

φ(x,y)≈g ^T(x,y)φ 公式(4)

其中：

φ(x,y)—节块内部中子通量密度的分布；

g(x,y)—x-y方向有限元形状函数向量；

T—向量或矩阵的转置；

φ—节块内部中子通量密度的展开矩向量；

在平源加速假设下，平源区的平均中子通量密度展开矩

能够表示为中子通量密度分布在本平源区内部的积分平均值：

式中：

—各平源区平均中子通量的矩向量；

h_t (x,y)—对应于每个有限元网格的分片常量：h_t(x,y)＝δ_e，e为每个有限元网格；

—对角线元素为每个有限元网格面积的对角矩阵；

步骤2：将公式(4)中节块内部中子角通量密度的离散表达式代入公式(5)中，就能够得到公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ的关系；

将公式(4)代入公式(5)得：

其中：

—平源区的平均中子通量密度展开矩；

同时

步骤3：利用公式(6)中平源区的平均中子通量密度展开矩

和有限元节点上中子通量密度展开矩φ之间的转换关系，代入传统变分节块方法中响应矩阵方程的表达式公式(1)和公式(2)，将求解对象从各个有限元节点上的未知量转换为各个平源区的未知量，即公式(13)和公式(14)；

将公式(1)代入到公式(6)得：

此外，群内源项能够表示为：

假设每个平源区内有限元节点上的源项与平源区的源项相等：

式中：

Σ_x(x,y)—平源区内部不同位置处的x反应的截面，x可以为吸收，总截面，散射；

—各平源区x反应截面所组成的对角矩阵；

因此，将公式(10)代入公式(9)，得：

其中：

将公式(11)代入公式(8)得：

同理，公式(2)变为：

步骤4：利用通用的迭代方法，对公式(13)和公式(14)进行求解，从而得到整个非均匀求解区域的中子通量密度分布和中子流密度分布，实现非均匀几何变分节块方法的求解。

Images