CN106960090A - 一种反应堆组件几何变形反应性的计算方法 - Google Patents

Abstract

一种反应堆组件几何变形反应性的计算方法，该方法假设组件变形时材料的微观截面不发生改变，并使用材料核子密度的变化来等效组件的几何变化，包括如下步骤：步骤1：将近似1应用到一阶微扰方程中，推导得到反应堆内组件几何变形反应性的计算表达式；步骤2：基于近似2，使用材料核子密度的变化来等效几何变形带来的反应性变化，进而使用步骤1中得到的反应性变化表达式对几何变形问题进行求解；本发明方法具有很强的几何适应性，很好的计算精度和较高的计算效率。

Description

一种反应堆组件几何变形反应性的计算方法

技术领域

本发明属于核反应堆堆芯设计和安全领域，具体涉及一种反应堆组件几何变形反应性的计算方法。

背景技术

反应堆的安全性是核能发展的重要前提，组件在堆内由于温度变化、辐照或者外力作用(例如地震)会造成几何形状的改变，这种几何变化往往是不规则的，而且导致的反应性变化对于堆芯安全是不能忽视的，需要准确快速的计算组件几何变化引入的反应性。

目前的堆芯计算方法可以分为蒙卡方法和确定论方法两类。其中蒙卡方法几何适应性强，可以描述变形的组件，但计算效率低下。确定论方法计算效率相对较高，但大多基于规则几何，不能描述弯曲的组件。目前虽然已经存在一些组件变形的近似处理方法，但仍然存在计算效率低，几何适应性弱以及计算精度低的问题。

发明内容

为了解决上述问题，本发明提出了一种能够处理组件几何变形反应性的计算方法。这种方法基于改进的一阶微扰方程，能够快速准确的得到堆内组件几何变化带来的反应性改变。

本发明和现有方案相比，在实施过程中使用了如下近似：

1.假设反应堆内组件变形时，反应堆内组件材料的微观截面不变；

2.假设组件几何变形带来的反应性变化，可以被组件材料密度变化引起的反应性变化抵消。也就是说，可以通过等效的材料密度变化来计算组件几何变形带来的反应性变化。

该方法的具体实现主要包括如下步骤：

步骤1：将近似1应用到一阶微扰方程中，推导得到反应堆内组件几何变形反应性的计算表达式，该步骤具体包含以下内容：

中子输运方程的扰动方程表示为

其中

<φ*|δSφ>＝∫∫∫φ*dEdΩ∫∫δ(Σ_s(r；E',Ω'→E,Ω))φ(r,E',Ω')dE'dΩ'dV (3)

<φ*|δAφ>＝∫∫∫φ*φδΣ_tdEdΩdV (4)

式中Δρ为变化的反应性，F、S、A分别为裂变项、散射项和吸收项，k_eff为有效增殖因子，φ*为中子通量密度，φ为共轭中子通量密度，χ为裂变谱，υ为每次裂变中子数，Σ_f、Σ_s、Σ_t分别为裂变截面、散射截面和总截面，r为空间位置，E、Ω、V分别为能量、角度和体积，δ为表示参数变化的运算符，例δF表示裂变项的变化；

如果令

φ*|Sφ>＝∫∫∫φ*dEdΩ∫∫Σ_s(r；E',Ω'→E,Ω)φ(r,E',Ω')dE'dΩ'dV (6)

φ*|Aφ>＝∫∫∫φ*φΣ_tdEdΩdV (7)

根据近似1，组件变形时忽略其材料微观截面的变化，那么当材料的核子密度变化率时，有

<φ*|δFφ>＝ε<φ*|Fφ> (8)

<φ*|δSφ>＝ε<φ*|Sφ> (9)

<φ*|δAφ>＝ε<φ*|Aφ> (10)

将公式(8)—(10)代入式(1)，便可得反应性变化的表达式

又因为中子平衡的表达式为

其中Lφ为泄漏项，所以

其中

为描述方便，将式(13)的中间项和右端项的分数部分表示为算子和

那么组件变形后的反应性变化就能够写成

这就是组件变形时反应性变化的计算表达式，但此表达式只适用于各向同性问题，即在所有方向上密度和尺寸的改变量都一致；对于各向异性问题，在不同方向上密度和尺寸的变化会存在差异，因此，可以将公式(17)改写成如下形式

其中i表示不同的方向；

步骤2：基于近似2，使用材料核子密度的变化来等效几何变形带来的反应性变化，进而使用步骤1中得到的反应性变化表达式(18))对几何变形问题进行求解；具体包含以下内容：

如果用N₁表示初始核子密度，N₂表示等效核子密度，那么等效核子密度变化率ε便可以表示为

其中κ_d为密度变化系数；κ_l为线性尺度变化系数；

利用式(19))就能够考虑各向异性问题不同方向上的几何变化，进而计算出总的反应性变化。

本发明和现有技术相比，具有以下优点：

1.不需要重新计算几何变形后的材料截面参数，能够节省大量时间；

2.不直接描述组件的几何形变，具有很强的几何适应性，现有的大部分程序都可以很方便的使用这套理论。

附图说明

图1组件变形问题示意图，其中图1(a)为均匀各向同性问题，其中图1(b)为均匀各向异性问题，其中图1(c)为非均匀各向异性问题。

图2组件偏移问题示意图(处理方法一)，其中图2(a)为真实组件结构在初始位置时和堆芯物理计算程序划分的网格的对应关系，其中图2(b)为组件向内发生小范围偏移时的情况，粗线表示发生了材料密度变化的网格，其中图2(c)为组件向内发生了较大偏移时的情况，粗线同样表示发生了材料密度变化的网格。

图3组件偏移问题示意图(处理方法二)，其中图3(a)中组件偏移情况和图2(b)中相同，其中图3(b)中组件偏移情况和图2(c)中相同，其中图3(c)为组件向外发生了较大范围偏移时的情况。其中图3(a)‐(c)中的粗线表示发生了材料密度变化或者面积变化的网格。

具体实施方式

本发明一种能够计算组件几何变形反应性的方法，该方法假设组件变形时材料的微观截面不变，并使用材料核子密度的变化来等效组件几何变形带来的反应性变化。下面将介绍此方法在一阶微扰方程中的应用。

中子输运方程的扰动方程可以表示为

其中

<φ*|δSφ>＝∫∫∫φ*dEdΩ∫∫δ(Σ_s(r；E',Ω'→E,Ω))φ(r,E',Ω')dE'dΩ'dV (3)

<φ*|δAφ>＝∫∫∫φ*φδΣ_tdEdΩdV (4)

式中Δρ为变化的反应性，F、S、A分别为裂变项、散射项和吸收项，k_eff为有效增殖因子，φ*为中子通量密度，φ为共轭中子通量密度，χ为裂变谱，υ为每次裂变中子数，Σ_f、Σ_s、Σ_t分别为裂变截面、散射截面和总截面，r为空间位置，E、Ω、V分别为能量、角度和体积，δ为表示参数变化的运算符，例δF表示裂变项的变化。

如果令

<φ*|Sφ>＝∫∫∫φ*dEdΩ∫∫Σ_s(r；E',Ω'→E,Ω)φ(r,E',Ω')dE'dΩ'dV (6)

<φ*|Aφ>＝∫∫∫φ*φΣ_tdEdΩdV (7)

组件变形时，忽略微观截面的变化，那么当材料的核子密度变化率时，有

<φ*|δFφ>＝ε<φ*|Fφ> (8)

<φ*|δSφ>＝ε<φ*|Sφ> (9)

<φ*|δAφ>＝ε<φ*|Aφ> (10)

将公式(8)—(10)代入式(1)，可得反应性变化的表达式

因为根据中子平衡有

其中Lφ为泄漏项，所以

其中

为描述方便，将式(13)的中间项和右端项的分数部分表示为算子和

那么扰动后的反应性变化就可以写成

这就是虚拟密度理论的基本表达式。如果用N₁表示初始核子密度，N₂表示等效核子密度，那么等效系数ε可以表示为

其中κ_d为密度变化系数；κ_l为线性尺度变化系数。

对于球体的各向同性膨胀问题，假如半径变化系数为f(半径增加1％，则f＝1.01)，那么密度变化系数和线性尺度变化系数分别为

κ_l＝f

所以等效系数

以上讨论的都是各向同性问题，在所有方向上密度和尺寸的改变量都一致。对于各向异性问题，在不同方向上密度和尺寸的变化会存在差异，因此，可以将公式(17)改写成如下形式

其中i表示不同的方向。

利用式(19)就可以考虑各向异性问题不同方向上的几何变化，进而计算出总的反应性变化。

如图1中图1(a)、图1(b)、图1(c)所示，堆内的组件变形问题可以分为均匀各向同性，均匀各向异性和非均匀各向异性三类。其中均匀各向同性问题可直接使用式(17)计算；组件的均匀轴向膨胀和均匀径向膨胀都属于均匀各向异性问题；而组件弯曲则属于非均匀各向异性问题。下面分别介绍轴向膨胀、径向膨胀和组件弯曲这三种情况下的虚拟密度理论表达式。

1)轴向膨胀

轴向膨胀和径向膨胀都只需要关注轴向和径向的变化，所以可以将式(19)改写为

轴向膨胀时，轴向尺寸增大(κ_d＝f)，密度减小轴向等效系数ε_z＝0；径向尺寸不变(κ_d＝1)，密度减小径向等效系数因此，由公式(20)可得轴向膨胀的反应性变化表达式

2)径向膨胀

径向膨胀时，轴向尺寸不变(κ_d＝1)，密度减小轴向等效系数径向尺寸增大(κ_d＝f)，密度减小径向等效系数同理，可由公式(20)得到径向膨胀的反应性变化表达式

3)组件弯曲

组件弯曲在径向上可以转换为组件偏移问题。而对于偏移问题又可以采用两种方法来进行处理：一种是网格不变，只关心网格内的核子密度变化；另一种是网格随组件偏移而改变，同时关心网格大小和核子密度的变化。下面分别介绍这两种处理方法：

(1)方法一

图2为组件偏移问题第一种处理方法的示意图。其中图2(a)为真实组件结构和堆芯物理计算程序划分的网格的对应关系，一个六边形组件划分为六个等边三角形(图中只画出了部分网格)。对于最右边的组件，六个三角形中燃料及冷却剂的核子密度都是相同的。

图2(b)为组件发生偏移之后的情况。假设组件的偏移没有超过原先六个三角形的范围，那么如图中粗线所示，有六个三角形内的核子密度发生了变化，只需要知道每个三角形节块内的密度变化情况，就可以使用虚拟密度理论进行求解。

图2(c)为组件发生偏移的另一种情况。在这种情况下，组件的偏移量超过了原先六个三角形节块的范围。同样如图中粗线所示，现在有七个三角形节块内的核子密度发生了变化，对于这种情况同样可以根据偏移量确认每个节块内的密度变化情况，然后使用虚拟密度理论得到偏移后的反应性变化。

(2)方法二

图3是组件偏移问题第二种处理方法的示意图。其中图3(a)和图3(b)的组件偏移情况分别对应图2中的图2(b)和图2(c)。

在方法二中，将一个组件中六个三角形节块的交点固定在组件中心不变，组件偏移时网格交点也随之改变。在这种处理方式下，组件偏移后需要同时考虑发生了网格变化或核子密度变化的区域。如图3中粗线所示，图3(a)和图3(c)需要考虑六个三角形区域的变化，图3(b)需要考虑七个区域的变化。确定了各区域的网格变化和核子密度变化之后，就可以等到该区域的等效系数，进而使用虚拟密度理论进行求解。

这两种方法相比，方法一中网格不变的处理方式显然更简单，但是适用范围有限。假设图3(c)中是堆芯外围组件向外发生了较大的偏移，超过了原先的网格范围，这种情况固定网格显然不再适用。而且，对于上文提到的堆芯径向膨胀和轴向膨胀，由于堆芯整体尺寸的增加，固定网格也明显不适用。所以，固定网格的方法只适用于局部的小扰动，对于这类问题，使用固定网格能够使问题简化，方便处理。而考虑网格变化的适用性更广，能够处理堆芯整体尺寸的改变，但是在局部非均匀变形问题的处理上会比固定网格更复杂。

虽然上述两种方法的处理方法和适用范围存在差别。但是，最终都能够得到各个区域内的等效系数。现假设某个三角形节块内只存在燃料和冷却剂两种材料，分开考虑燃料和冷却剂提供的反应性，分别使用ε_f和ε_c表示组件偏移后燃料和冷却剂的等效系数，那么此网格内的反应性变化可以写为

推广到全堆，当堆内有组件弯曲时，考虑发生密度变化的节块，然后对每个节块单独计算，最后求和就得到了组件弯曲的反应性变化，表达式可写为

因为上述表达式只考虑组件偏移对网格大小和核子密度的影响，所以即使组件在偏移的同时还发生了不规则的形状变化，上面的表达式依然适用。

Claims (1)

1.一种反应堆组件几何变形反应性的计算方法，其特征在于：该方法在计算过程中应用了如下近似：

近似1：反应堆内组件变形时，反应堆内组件材料的微观截面不发生改变；

近似2：使用材料的密度变化来等效计算组件几何变形引入的反应性变化；

该方法的具体实现主要包括如下步骤：

步骤1：将近似1应用到一阶微扰方程中，推导得到反应堆内组件几何变形反应性的计算表达式，该步骤具体包含以下内容：

中子输运方程的扰动方程表示为

$\mathrm{\&Delta;}\mathrm{\&rho;}=\frac{(\mathrm{\&phi;}*|(\frac{1}{k_{eff}}\mathrm{\&delta;}F+\mathrm{\&delta;}S-\mathrm{\&delta;}A)\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}---\left(1\right)$

其中

$<\mathrm{\&phi;}*|\mathrm{\&delta;}F\mathrm{\&phi;}>=\&Integral;\&Integral;\&Integral;\frac{\mathrm{\&chi;}\left(E\right)}{4\mathrm{\&pi;}}\mathrm{\&phi;}*dEd\mathrm{\&Omega;}\&Integral;\&Integral;\mathrm{\&delta;}\left({\mathrm{\&upsi;\&Sigma;}}_f\right(r,E^{\&prime;}\left)\right)\mathrm{\&phi;}(r,E^{\&prime;},{\mathrm{\&Omega;}}^{\&prime;}){\mathrm{dE}}^{\&prime;}{\mathrm{d\&Omega;}}^{\&prime;}dV---\left(2\right)$

<φ*|δSφ>＝∫∫∫φ*dEdΩ∫∫δ(Σ_s(r；E',Ω'→E,Ω))φ(r,E',Ω')dE'dΩ'dV (3)

<φ*|δAφ>＝∫∫∫φ*φδΣ_tdEdΩdV (4)

式中Δρ为变化的反应性，F、S、A分别为裂变项、散射项和吸收项，k_eff为有效增殖因子，φ*为中子通量密度，φ为共轭中子通量密度，χ为裂变谱，υ为每次裂变中子数，Σ_f、Σ_s、Σ_t分别为裂变截面、散射截面和总截面，r为空间位置，E、Ω、V分别为能量、角度和体积，δ为表示参数变化的运算符，例δF表示裂变项的变化；

如果令

$<\mathrm{\&phi;}*|F\mathrm{\&phi;}>=\&Integral;\&Integral;\&Integral;\frac{\mathrm{\&chi;}\left(E\right)}{4\mathrm{\&pi;}}\mathrm{\&phi;}*dEd\mathrm{\&Omega;}\&Integral;\&Integral;{\mathrm{\&upsi;\&Sigma;}}_f(r,E^{\&prime;})\mathrm{\&phi;}(r,E^{\&prime;},{\mathrm{\&Omega;}}^{\&prime;}){\mathrm{dE}}^{\&prime;}{\mathrm{d\&Omega;}}^{\&prime;}dV---\left(5\right)$

<φ*|Sφ>＝∫∫∫φ*dEdΩ∫∫Σ_s(r；E',Ω'→E,Ω)φ(r,E',Ω')dE'dΩ'dV (6)

<φ*|Aφ>＝∫∫∫φ*φΣ_tdEdΩdV (7)

根据近似1，组件变形时忽略其材料微观截面的变化，那么当材料的核子密度变化率时，有

<φ*|δFφ>＝ε<φ*|Fφ> (8)

<φ*|δSφ>＝ε<φ*|Sφ> (9)

<φ*|δAφ>＝ε<φ*|Aφ> (10)

将公式(8)—(10)代入式(1)，便可得反应性变化的表达式

$\mathrm{\&Delta;}\mathrm{\&rho;}=\mathrm{\&epsiv;}\frac{<\mathrm{\&phi;}*|(\frac{1}{k_{eff}}F+S-A)\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}---\left(11\right)$

又因为中子平衡的表达式为

$(\frac{1}{k_{eff}}F+S-A)\mathrm{\&phi;}=L\mathrm{\&phi;}---\left(12\right)$

其中Lφ为泄漏项，所以

$\mathrm{\&Delta;}\mathrm{\&rho;}=\mathrm{\&epsiv;}\frac{<\mathrm{\&phi;}*|(\frac{1}{k_{eff}}F+S-A)\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}=\mathrm{\&epsiv;}\frac{<\mathrm{\&phi;}*|L\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}---\left(13\right)$

其中

$<\mathrm{\&phi;}*|L\mathrm{\&phi;}>=\&Integral;\&Integral;\&Integral;{\mathrm{\&phi;}}^{*}\mathrm{\&Omega;}\&CenterDot;\&dtri;\mathrm{\&phi;}dEd\mathrm{\&Omega;}dV---\left(14\right)$

为描述方便，将式(13)的中间项和右端项的分数部分表示为算子和

$\tilde{P}=\frac{<\mathrm{\&phi;}*|(\frac{1}{k_{eff}}F+S-A)\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}---\left(15\right)$

$\tilde{L}=\frac{<\mathrm{\&phi;}*|L\mathrm{\&phi;}>}{<\mathrm{\&phi;}*|F\mathrm{\&phi;}>}---\left(16\right)$

那么组件变形后的反应性变化就能够写成

$\mathrm{\&Delta;}\mathrm{\&rho;}=\mathrm{\&epsiv;}\tilde{P}=\mathrm{\&epsiv;}\tilde{L}---\left(17\right)$

这就是组件变形时反应性变化的计算表达式，但此表达式只适用于各向同性问题，即在所有方向上密度和尺寸的改变量都一致；对于各向异性问题，在不同方向上密度和尺寸的变化会存在差异，因此，将公式(17)改写成如下形式

$\mathrm{\&Delta;}\mathrm{\&rho;}=\underset{i}{\&Sigma;}{\mathrm{\&epsiv;}}_i{\tilde{L}}_i---\left(18\right)$

其中i表示不同的方向；

步骤2：基于近似2，使用材料核子密度的变化来等效几何变形带来的反应性变化，进而使用步骤1中得到的反应性变化表达式(18)对几何变形问题进行求解；具体包含以下内容：

如果用N₁表示初始核子密度，N₂表示等效核子密度，那么等效核子密度变化率ε便表示为

$\mathrm{\&epsiv;}=\frac{N_2-N_1}{N_1}=\frac{{\mathrm{\&kappa;}}_d{\mathrm{\&kappa;}}_l{N}_1-N_1}{N_1}={\mathrm{\&kappa;}}_d{\mathrm{\&kappa;}}_l-1---\left(19\right)$

其中κ_d为密度变化系数；κ_l为线性尺度变化系数；

利用式(19)就能够考虑各向异性问题不同方向上的几何变化，进而计算出总的反应性变化。