Radially and azimuthally dependent resonance self-shielding treatment for general multi-region geometry based on a unified theory

ABSTRACT

A unified resonance self-shielding method, which can treat general sub-divided fuel regions, is developed for lattice physics calculations in reactor physics field. In a past study, a hybrid resonance treatment has been developed by theoretically integrating equivalence theory and ultra-fine-group slowing-down calculation. It can be applied to a wide range of neutron spectrum conditions including low moderator density ranges in severe accident states, as long as each fuel region is not sub-divided. In order to extend the method for radially and azimuthally sub-divided multi-region geometry, a new resonance treatment is established by incorporating the essence of sub-group method. The present method is composed of two-step flux calculation, i.e. ‘coarse geometry + fine energy’ (first step) and ‘fine geometry + coarse energy’ (second step) calculations. The first step corresponds to a hybrid model of the equivalence theory and the ultra-fine-group calculation, and the second step corresponds to the sub-group method. From the verification results, effective cross-sections by the new method show good agreement with the continuous energy Monte-Carlo results for various multi-region geometries including non-uniform fuel compositions and temperature distributions. The present method can accurately generate effective cross-sections with short computation time in general lattice physics calculations.

KEYWORDS:

• Reactor physics
• lattice physics
• cross-section, neutron spectrum
• radially and azimuthally sub-divided multi-region
• unified resonance self-shielding treatment
• equivalence theory
• ultra-fine-group slowing-down calculation
• sub-group method
• two-step flux calculation

1. Introduction

Resonance self-shielding treatment [1] is an important part for the lattice physics calculations in reactor physics field. Accurate treatment of ²³⁸U wide resonance cross-sections and their influences on flux depressions is a key issue to guarantee a sufficient prediction accuracy of criticality and reactivity coefficients for commercial LWR (light water reactor) core analysis.

In order to establish a practical resonance treatment method, the effective cross-section is required to be generated with short computation time, while keeping its sufficient accuracy. In the past studies conducted by the authors, many of the technical issues have been solved to establish a sophisticated resonance treatment. A brief summary and the development history of the past and present studies are shown in Figure 1 [2–11]. In the present study, a unified resonance treatment method is developed to obtain the sufficient accuracy for effective cross-sections without a direct heterogeneous ultra-fine-group calculation. As shown in Figure 1, treatment of the complicated spatial self-shielding effect for radially and azimuthally sub-divided multi-region geometry within each fuel rod is addressed in the present paper.

Figure 1. Brief summary and development history of the past and present studies.

In the first technical achievement by the authors [4], approximation for geometry and black neutron absorber assumption in the equivalence theory [1] are solved. Several spatial effects in a lattice physics theory (i.e. heterogeneous effect between fuel and non-fuel regions, shadowing effect due to the neighboring fuel rods, and irregular lattice effect by the local water enriched regions) can be incorporated into the pin-by-pin rational coefficients. The flux response as a function of macroscopic total cross-sections obtained from method of characteristics (MOC) [12] fixed-source calculations is directly taken into account for pin-by-pin resolution.

In the second technical achievement by the authors [7], scattering source and resonance interference treatments are improved. The essence of the ultra-fine-group slowing-down calculation [5,6] is incorporated into the equivalence theory, and the hybrid resonance treatment method is established. At this stage, the resonance treatment for both pin-cell and fuel assembly geometries is established to generate accurate effective cross-sections with short computation time. It can be applied to the wide range of neutron spectrum conditions including low moderator density ranges in severe accident states, as long as each fuel region is not sub-divided.

In the current development state, extension of the above hybrid resonance treatment is desired for a complicated spatial self-shielding treatment of each fuel lump. Especially for the cylindrical fuel rod geometry, which is widely adopted for commercial LWR fuel assemblies, generation of radially and azimuthally dependent effective cross-sections is not easy with keeping sufficient accuracy.

From the above background, an objective of the present study is to develop a new resonance treatment method that can accurately take into account the complicated spatial self-shielding effect with short computation time. In order to solve this issue, the essence of sub-group method [10,11] is incorporated into the hybrid resonance treatment method in the second achievement [7] mentioned above. As a result, the unified resonance treatment is established as a multi-hybrid model of the equivalence theory, the ultra-fine-group calculation and the sub-group method. Although another concept of unified approach for the three methods has already been suggested by the authors [13], the present method is based on a more robust scheme from the theoretical point of view.

In the current methodologies related to the resonance treatment, the simpler approaches have already been proposed, e.g. a pin-cell based direct ultra-fine-group slowing-down calculation with pre-tabulated collision probabilities [6]. The present method in this study adopts rather complicated calculation procedures in actual implementation than the direct slowing-down calculation mentioned above. However, the present method has an advantage on accurate and efficient spatial treatment for radially/azimuthally sub-divided fuel regions defined on the detailed core analysis applications. Especially for the azimuthally sub-divided fuel rods adjacent to a large water region, the pin-cell based slowing-down calculation cannot directly be applied. Therefore, the present method has a potential to become an alternative approach against the direct slowing-down calculation, from the view point of its capability to treat complicated calculation conditions.

In this paper, the fundamental theory and the verification results of the new resonance treatment are described in detail. This paper is organized as follows. Section 2 proposes the unified resonance treatment based on a two-step cross-section collapsing scheme. The equivalence theory, the ultra-fine-group calculation, and the sub-group method are also integrated. Section 3 presents a two-step reaction-rate preservation scheme for the reduction of energy discretization error, which is consistent with the unified resonance treatment. Section 4 presents the verification of the new method. Finally, Section 5 concludes the paper.

2. A unified resonance self-shielding treatment

2.1. Concept for two-step resonance calculation

The new method is established by enhancing and integrating the conventional methods, and its concept is based on two-step neutron flux calculations. The first-step calculation is performed to obtain ultra-fine-group fluxes on a simple geometry based on the hybrid treatment of equivalence theory and ultra-fine-group slowing-down calculation. Then, by using the obtained ultra-fine-group fluxes, sub-group cross-sections and sources are generated for discrete energy ranges explicitly defined by the resonance cross-section level. Next, the second-step calculation is performed to obtain sub-group fluxes on an exact geometry based on the sub-group method. Finally, multi-group effective cross-sections are generated by the sub-group cross-sections and the sub-group fluxes.

A concept of the two-step resonance treatment is shown in Table 1. By using a concept of multi-stage cross-section collapsing strategies [14] widely used in reactor physics field, two-step calculation, i.e. ‘coarse geometry + fine energy’ (first step) and ‘fine geometry + coarse energy’ (second step) calculations, is performed. Here, the main calculations in the present method are summarized in Table 2. The corresponding section numbers for each main calculation are also shown in the table.

Table 1. Concept of two-step resonance treatment

	Step
Item	First step	Second step
Theory	Equivalence theory + Ultra-fine-group calculation	Sub-group method
Spatial resolution	or	or or
Energy resolution	Ultra-fine-group (∼100,000)	Sub-group (3–10 for each multi-group)
Collapsed cross-section	${\sigma }_{sg}^i=\frac{{\int }_{sg}dE\sigma \left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}$	${\sigma }_g^i=\frac{\sum _{sg\in g}{\sigma }_{sg}^i{\varphi }_{sg}^i}{\sum _{sg\in g}{\varphi }_{sg}^i}$

Table 2. Summary of the main calculation procedures for the present method

Main calculation	Corresponding sections or	Corresponding steps
procedures	reference documents	in Section 2.6
Generation of pin-by-pin rational coefficients	References 4, 7	(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) )
Ultra-fine-group slowing-down calculation for simplified geometry	Reference 7, Section 2.2	(1(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) )–(5(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) )
Extension of ultra-fine-group flux for each fuel rod to radial multi-region problem	Section 2.3	(6(6) $$\begin{array}{c}\hfill {\Sigma }_t^i\left(E\right)V_i{P}_{i\to j}\left(E\right)={\Sigma }_t^j\left(E\right)V_j{P}_{j\to i}\left(E\right).\end{array}$$(6) )–(9(9) $$\begin{array}{ccc}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)& \cong & {\Sigma }_t^f\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^f\left(E\right)}{S}_f\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i\sum _{j\in f}{P}_{i\to j}\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i·\{1-P_{i\to nf}\left(E\right)\},\hfill \end{array}$$(9) )
Generation of spatially dependent sub-group parameter	Sections 2.4-2.5	(1(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) ), (10(10) $$\begin{array}{c}\hfill \sum _{j\in f}{P}_{i\to j}\left(E\right)+P_{i\to nf}\left(E\right)=1.\end{array}$$(10) )–(13(13) $$\begin{array}{ccc}& & {\varphi }_i\left(E\right)\hfill \\ & & =\frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right].\hfill \end{array}$$(13) )
Sub-group flux calculation	Section 2.5	(14(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) )
Generation of effective cross-section	Section 2.5	(15(15) $$\begin{array}{ccc}& & ({\gamma }_{i,1},{\gamma }_{i,2},{\gamma }_{i,3},{\gamma }_{i,4})\hfill \\ & & \equiv \left(\frac{{\rho }_i{l}_{i,1}}{l_i},-\frac{{\rho }_i{l}_{i,2}}{l_i},-\frac{{\rho }_{i-1}{l}_{i,3}}{l_i},\frac{{\rho }_{i-1}{l}_{i,4}}{l_i}\right).\hfill \end{array}$$(15) )

The main features of the new method are summarized as follows:

(1)	In the first-step calculation, the hybrid resonance treatment based on the equivalence theory and the ultra-fine-group calculation [7] is applied, and the ultra-fine-group flux for each fuel region is accurately obtained.
(2)	In the case of radial multi-region problem with treating the concentric cylinder within a fuel rod, the ultra-fine-group flux for each fuel lump is expanded to each ring region coupled with the basic idea of Stoker–Weiss method [8] and SDDM (Spatially Dependent Dancoff Method) [9]. The direct heterogeneous ultra-fine-group transport calculation is not required at all.
(3)	Because the effective cross-sections generated by direct collapsing with the above ultra-fine-group fluxes for each ring region are not sufficiently accurate for the radial multi-region problem, the ultra-fine-group cross-section is collapsed not for the whole multi-group energy range, but for the partial sub-group energy range defined by the resonance cross-section level. Note that several sub-groups are generally included in each multi-group.
(4)	In the second-step calculation, the sub-group method is applied to the exact geometry by using the sub-group cross-sections obtained from the first-step calculation. In the case of azimuthally dependent self-shielding treatment, MOC is used as a flux calculation method.
(5)	The effective cross-section, i.e. the final product of resonance calculation, is generated by collapsing the sub-group cross-section with the sub-group flux as a weighting function.

2.2. Integration of equivalence theory and ultra-fine-group slowing-down calculation

In the past study, the hybrid resonance self-shielding method has been established by the authors through integrating the equivalence theory and the ultra-fine-group slowing-down calculation [7]. According to Reference [7], a final form of the neutron flux for each fuel region in a fuel assembly is written as (1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) where (2) $$\begin{array}{c}\hfill {\Sigma }_{sd}^f\left(E\right)=E\sum _k\frac{N_k^f}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_f\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(2) (3) $$\begin{array}{c}\hfill \mu \left(E\right)=\theta +(1-\theta )\frac{E}{{\Sigma }_p^{nf}}\sum _k\frac{N_k^{nf}{\sigma }_p^k}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\varphi }_{nf}\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(3) (4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4)

2.3. Extension of neutron flux to the radial multi-region geometry

Although the hybrid resonance self-shielding method described in Section 2.2 can treat any type of fuel lump geometry in an assembly that can be handled by MOC flux calculation, the spatially dependent resonance self-shielding effect within a specific fuel lump cannot be directly treated. One of the important issues in this research field is consideration of the radially dependent resonance self-shielding effect within a cylindrical fuel rod. In this section, the neutron flux in a local ring region within a fuel rod is derived as a natural extension of Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ). Here, the radially sub-divided fuel rod discussed in this section is shown in Figure 2.

Figure 2. Radially sub-divided fuel rod.

The integral form of neutron transport equation for multi-region heterogeneous system consisting of multiple fuel ring regions and a non-fuel region is written as (5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5)

In the following derivation, the non-fuel region is defined as the multiple regions except for fuel region. For a typical LWR unit cell, non-fuel region is composed of cladding and moderator regions.

The reciprocity theorem between regions i and j is written as (6) $$\begin{array}{c}\hfill {\Sigma }_t^i\left(E\right)V_i{P}_{i\to j}\left(E\right)={\Sigma }_t^j\left(E\right)V_j{P}_{j\to i}\left(E\right).\end{array}$$(6)

By applying Equation (6)(6) $$\begin{array}{c}\hfill {\Sigma }_t^i\left(E\right)V_i{P}_{i\to j}\left(E\right)={\Sigma }_t^j\left(E\right)V_j{P}_{j\to i}\left(E\right).\end{array}$$(6) to the neutron source term, the first term in the right-hand side of Equation (5)(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) is rewritten as (7) $$\begin{array}{c}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)={\Sigma }_t^i\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^j\left(E\right)}{S}_j\left(E\right).\end{array}$$(7)

The second term in the right-hand side of Equation (5)(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) is also rewritten as (8) $$\begin{array}{c}\hfill S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right)={\Sigma }_t^i\left(E\right)V_i\frac{P_{i\to nf}\left(E\right)}{{\Sigma }_p^{nf}}{S}_{nf}\left(E\right),\end{array}$$(8) where the total cross-section for the non-fuel region is approximated as the potential scattering cross-section, i.e. Σ^nf_t(E)≅Σ_s^nf(E)≅Σ^nf_p.

By assuming that the macroscopic total cross-section and scattering source are spatially flat within a fuel pellet, Equation (7)(7) $$\begin{array}{c}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)={\Sigma }_t^i\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^j\left(E\right)}{S}_j\left(E\right).\end{array}$$(7) is approximated as (9) $$\begin{array}{ccc}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)& \cong & {\Sigma }_t^f\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^f\left(E\right)}{S}_f\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i\sum _{j\in f}{P}_{i\to j}\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i·\{1-P_{i\to nf}\left(E\right)\},\hfill \end{array}$$(9) where S_f(E) denotes the scattering source for the fuel region, and the following relation for collision probability is used: (10) $$\begin{array}{c}\hfill \sum _{j\in f}{P}_{i\to j}\left(E\right)+P_{i\to nf}\left(E\right)=1.\end{array}$$(10)

When the scattering source for the non-fuel region is written as $S_{nf}\left(E\right)=\mu \left(E\right)\frac{{\Sigma }_p^{nf}}{E}$, which is consistent with the hybrid resonance treatment briefly described in Section 2.2, Equation (8)(8) $$\begin{array}{c}\hfill S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right)={\Sigma }_t^i\left(E\right)V_i\frac{P_{i\to nf}\left(E\right)}{{\Sigma }_p^{nf}}{S}_{nf}\left(E\right),\end{array}$$(8) is approximated as (11) $$\begin{array}{ccc}\hfill S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right)& \cong & \left\{\mu \left(E\right)\frac{{\Sigma }_p^{nf}}{E}\right\}\hfill \\ & & ·{\Sigma }_t^f\left(E\right)V_i\frac{P_{i\to nf}\left(E\right)}{{\Sigma }_p^{nf}}\hfill \\ & =& \frac{1}{E}{\Sigma }_t^f\left(E\right)V_i\mu \left(E\right)P_{i\to nf}\left(E\right),\hfill \end{array}$$(11) where μ(E) is already defined by Equation (3)(3) $$\begin{array}{c}\hfill \mu \left(E\right)=\theta +(1-\theta )\frac{E}{{\Sigma }_p^{nf}}\sum _k\frac{N_k^{nf}{\sigma }_p^k}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\varphi }_{nf}\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(3) in Section 2.2.

Substitution of Equations (9)(9) $$\begin{array}{ccc}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)& \cong & {\Sigma }_t^f\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^f\left(E\right)}{S}_f\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i\sum _{j\in f}{P}_{i\to j}\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i·\{1-P_{i\to nf}\left(E\right)\},\hfill \end{array}$$(9) and (11(11) $$\begin{array}{ccc}\hfill S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right)& \cong & \left\{\mu \left(E\right)\frac{{\Sigma }_p^{nf}}{E}\right\}\hfill \\ & & ·{\Sigma }_t^f\left(E\right)V_i\frac{P_{i\to nf}\left(E\right)}{{\Sigma }_p^{nf}}\hfill \\ & =& \frac{1}{E}{\Sigma }_t^f\left(E\right)V_i\mu \left(E\right)P_{i\to nf}\left(E\right),\hfill \end{array}$$(11) ) into the right-hand side of Equation (5)(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) , and some transformations can yield (12) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\{1-P_{i\to nf}\left(E\right)\}·\frac{S_f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\frac{1}{E}\mu \left(E\right)P_{i\to nf}\left(E\right).\end{array}$$(12)

Since the scattering source for the fuel region is expressed as $S_f\left(E\right)=\frac{{\Sigma }_{sd}^f\left(E\right)}{E}$, which is consistent with the hybrid resonance treatment in Section 2.2, Equation (12)(12) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\{1-P_{i\to nf}\left(E\right)\}·\frac{S_f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\frac{1}{E}\mu \left(E\right)P_{i\to nf}\left(E\right).\end{array}$$(12) is further modified as: (13) $$\begin{array}{ccc}& & {\varphi }_i\left(E\right)\hfill \\ & & =\frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right].\hfill \end{array}$$(13)

In the multi-term rational approximation of the fuel escape probability, P_{i → nf}(E) is formulated as [4]: (14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14)

Equation (14)(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) is identical to the SDGM (Spatially Dependent Gray resonance self-shielding Method) formulation [4], and its derivation idea is based on the Stoker–Weiss method [8] and SDDM [9]. A geometrical treatment of the escape probability, which is composed of four lump components (m = 1, 2, 3, 4), is shown in Figure 3.

Figure 3. Geometrical treatment of spatially dependent fuel escape probability.

Here, γ_{i, m} in Equation (14)(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) denotes the coefficient that captures the geometrical information of the ring i and is defined as (15) $$\begin{array}{ccc}& & ({\gamma }_{i,1},{\gamma }_{i,2},{\gamma }_{i,3},{\gamma }_{i,4})\hfill \\ & & \equiv \left(\frac{{\rho }_i{l}_{i,1}}{l_i},-\frac{{\rho }_i{l}_{i,2}}{l_i},-\frac{{\rho }_{i-1}{l}_{i,3}}{l_i},\frac{{\rho }_{i-1}{l}_{i,4}}{l_i}\right).\hfill \end{array}$$(15)

l_{i, m} denotes the mean chord length of m-th lump derived by the following analytical form: (16) $$\begin{array}{ccc}& & l_{i,m}=\frac{2R}{\pi }·\left(\sqrt{1-{\rho }^2}+\frac{1}{\rho }{sin}^{-1}\rho \pm \frac{\pi }{2}\rho \right),\hfill \\ & & \left\{\begin{array}{c}m=1:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=2:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=3:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \\ m=4:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \end{array}\right..\hfill \end{array}$$(16)

By substituting Equation (16)(16) $$\begin{array}{ccc}& & l_{i,m}=\frac{2R}{\pi }·\left(\sqrt{1-{\rho }^2}+\frac{1}{\rho }{sin}^{-1}\rho \pm \frac{\pi }{2}\rho \right),\hfill \\ & & \left\{\begin{array}{c}m=1:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=2:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=3:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \\ m=4:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \end{array}\right..\hfill \end{array}$$(16) into Equation (15)(15) $$\begin{array}{ccc}& & ({\gamma }_{i,1},{\gamma }_{i,2},{\gamma }_{i,3},{\gamma }_{i,4})\hfill \\ & & \equiv \left(\frac{{\rho }_i{l}_{i,1}}{l_i},-\frac{{\rho }_i{l}_{i,2}}{l_i},-\frac{{\rho }_{i-1}{l}_{i,3}}{l_i},\frac{{\rho }_{i-1}{l}_{i,4}}{l_i}\right).\hfill \end{array}$$(15) , ∑⁴_{m = 1}γ_{i, m} = 1. Considering that ∑^N_{n = 1}β_n = 1, the following relation is derived: (17) $$\begin{array}{c}\hfill \sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n=1.\end{array}$$(17)

By substituting Equations (14)(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) and (17(17) $$\begin{array}{c}\hfill \sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n=1.\end{array}$$(17) ) into Equation (13)(13) $$\begin{array}{ccc}& & {\varphi }_i\left(E\right)\hfill \\ & & =\frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right].\hfill \end{array}$$(13) , the energy-dependent neutron flux for the ring region i within a fuel pellet is derived as (18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18)

The scattering source S_i(E) is also calculated as Equation (19)(19) $$\begin{array}{c}\hfill S_i\left(E\right)=\left\{\begin{array}{c}{\displaystyle \sum _k\frac{N_k^i}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_i\left(E^{\text{'}}\right)}{E^{\text{'}}}(i\in \mathrm{fuel})}\hfill \\ {\displaystyle \mu \left(E\right)\frac{{\Sigma }_p^i}{E}(i\in \mathrm{non}\mathrm{fuel})}\hfill \end{array}.\right.\end{array}$$(19) , which is consistent with the hybrid resonance treatment: (19) $$\begin{array}{c}\hfill S_i\left(E\right)=\left\{\begin{array}{c}{\displaystyle \sum _k\frac{N_k^i}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_i\left(E^{\text{'}}\right)}{E^{\text{'}}}(i\in \mathrm{fuel})}\hfill \\ {\displaystyle \mu \left(E\right)\frac{{\Sigma }_p^i}{E}(i\in \mathrm{non}\mathrm{fuel})}\hfill \end{array}.\right.\end{array}$$(19)

Here, the isotropic and elastic down-scattering approximations in the center-of-mass system are applied to each ring region within a fuel pellet. The scattering source for the non-fuel region is consistent with the hybrid resonance treatment.

By considering the spatially dependent scattering source within a fuel pellet, the neutron flux is modified from Equation (18)(18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18) as (20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) where (21) $$\begin{array}{c}\hfill {\Sigma }_{sd}^i\left(E\right)\equiv E{S}_i\left(E\right),\end{array}$$(21) and η(E) denotes the energy-dependent flux normalization factor in which the fuel-pellet averaged flux obtained from Equation (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) can be preserved. Σⁱ_sd(E) and η(E) are the only differences between Equations (20)(20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) and (18(18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18) ).

2.4. First-step resonance calculation based on equivalence theory and ultra-fine-group slowing-down calculation

As described in Section 2.1, the two-step resonance calculations are performed in the present method. The first-step calculation is based on the equivalence theory and the ultra-fine-group slowing-down calculation, and the space-dependent ultra-fine-group neutron fluxes for each fuel region are appropriately obtained.

As already shown in Section 2.2, the ultra-fine-group fluxes for each fuel region within an assembly are calculated from a set of slowing-down equations shown in Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ). Besides, as derived in Section 2.3, spatially dependent ultra-fine-group fluxes and scattering sources for each ring region within a fuel pellet can be obtained. It should be noted that the fluxes for each ring can be obtained analytically with small computational burdens. Σ^f_sd(E) and μ(E) in Equation (18)(18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18) are the byproducts of the solution of Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ), and thus fluxes for each ring region are easily regenerated.

In principle, the radially dependent effective cross-sections for each ring region can be generated at this single step by using the ultra-fine-group flux solution as a collapsing weight. However, as shown in a derivation scheme of the equation, the following assumptions are applied:

(1)

The ultra-fine-group macroscopic total cross-section and scattering source are spatially flat within a fuel pellet (Equation (9)(9) $$\begin{array}{ccc}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)& \cong & {\Sigma }_t^f\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^f\left(E\right)}{S}_f\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i\sum _{j\in f}{P}_{i\to j}\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i·\{1-P_{i\to nf}\left(E\right)\},\hfill \end{array}$$(9) ).

(2)

The fuel escape probabilities for each fuel ring are expressed by the multi-term rational approximation based on the SDGM formulation (Equation (14)(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) ).

As shown in the verification results of Section 4, these assumptions cause the prediction error of effective cross-sections for each ring region, as a result of the prediction error of the ultra-fine-group flux.

In order to mitigate the influence of prediction error for the flux, the effective cross-sections are not generated at this step. Instead, the sub-group cross-sections are generated by using the ultra-fine-group fluxes obtained from the first-step calculation. A more accurate spatial transport of neutrons within a fuel pellet is taken into account in the second-step calculation by the additional transport calculations based on the sub-group method. The lack of accurate collision probability information for all the region-to-region combinations in a target system, which is due to the assumptions in Equations (9)(9) $$\begin{array}{ccc}\hfill \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)& \cong & {\Sigma }_t^f\left(E\right)V_i\sum _{j\in f}\frac{P_{i\to j}\left(E\right)}{{\Sigma }_t^f\left(E\right)}{S}_f\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i\sum _{j\in f}{P}_{i\to j}\left(E\right)\hfill \\ & =& S_f\left(E\right)V_i·\{1-P_{i\to nf}\left(E\right)\},\hfill \end{array}$$(9) and (14(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) ), is compensated by the second-step calculation. The details are described in the next section.

2.5. Second-step resonance calculation based on sub-group method

The second-step calculation is based on the sub-group method, and the space-dependent sub-group neutron fluxes for each region are obtained.

In the typical resonance energy ranges, the scattering source is isotropic [6] and the fission source can be ignored. Therefore, the integro-differential form of neutron transport equation is written as (22) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_i(E,\Omega )+{\Sigma }_t^i\left(E\right){\psi }_i(E,\Omega )=\frac{S_i\left(E\right)}{4\pi },\end{array}$$(22) where Ω denotes unit vector for direction, and ψ_i(E, Ω) denotes angular flux for region i.

By integrating Equation (22)(22) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_i(E,\Omega )+{\Sigma }_t^i\left(E\right){\psi }_i(E,\Omega )=\frac{S_i\left(E\right)}{4\pi },\end{array}$$(22) against an explicit sub-group energy range sg ( ∈ g), in which the discrete energy range is permitted, the sub-group transport equation is derived as (23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) where ${\Sigma }_{t,sg}^i\equiv \frac{{\int }_{sg}dE{\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}$ : macroscopic sub-group total cross-section,ψⁱ_sg(Ω) ≡ ∫_sgdEψ_i(E, Ω): sub-group angular flux,Sⁱ_sg ≡ ∫_sgdES_i(E): sub-group scattering source.

Here, the angular dependence of energy-dependent flux within each sub-group is ignored.

The typical scheme to solve Equation (23)(23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) is MOC or S_N (discrete ordinate) method. If MOC is adopted as a flux calculation scheme for Equation (23)(23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) , the spatially dependent flux for each sub-group can be obtained for any type of complicated geometry, which can be represented by the MOC. Although the integro-differential form of neutron transport equation is shown as a typical example, the integral form of equation based on the collision probability method can also be utilized.

The microscopic sub-group cross-section for resonance nuclide r and reaction type x is also defined in the same manner as (24) $$\begin{array}{c}\hfill {\sigma }_{x,sg}^{r,i}\equiv \frac{{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}.\end{array}$$(24)

The sub-group cross-sections are generated by the direct energy collapsing of σ^r_x(E) with φ_i(E) weight obtained from the first-step calculation.

As described in Section 2.4, the sub-group transport calculation is performed in order to mitigate the influence of prediction error of φ_i(E) for σ^{r, i}_{x, sg}. From this point of view, sub-group energy structure should be divided not by the neutron energy but by the resonance cross-section level. Based on a reference ultra-fine-group cross-section Σ_ref(E), the sub-group energy structure is determined as (25) $$\begin{array}{c}\hfill sg=\left\{E\left|{\Sigma }_{sg,min}^{\mathrm{bound}}\le {\Sigma }_{\mathrm{ref}}\left(E\right)<{\Sigma }_{sg,max}^{\mathrm{bound}}\right.\right\},\end{array}$$(25) where Σ^bound_{sg, max} and Σ^bound_{sg, min} denote the maximum and minimum resonance cross-section boundaries for sub-group sg ( ∈ g), respectively. Typically, Σ_ref(E) is defined by a volume-averaged ultra-fine-group cross-section for all the fuel regions in a target system such as a fuel assembly. According to the preliminary sensitivity studies, the macroscopic absorption cross-section is a better index than the total cross-section. Though the boundaries Σ^bound_{sg, max} and Σ^bound_{sg, min} are determined for each sub-group sg with arbitrary interval division of resonance cross-section within a group g, the equal interval with logarithmic scale is better than that with linear scale. A concept for determination of sub-group structure is shown in Figure 4.

Figure 4. Concept for determination of sub-group structure.

Here, it should be noted that the defined sub-group energy structure is commonly used for all the spatial regions in sub-group cross-section generation. Even if the fuel temperature distribution within a fuel pellet is taken into account, the only one sub-group structure, e.g. a structure defined by Σ_ref(E) for effective fuel temperature [15] condition, must be applied. If the different sub-group energy structures are used for each region in sub-group cross-section generation, i.e. the different Σ_ref(E)s are applied to each region, the neutron balance among each region and each sub-group is not preserved due to inconsistency of energy group structure for each sub-group, which is an essential issue for the conventional sub-group method.

The sub-group scattering source Sⁱ_sg = ∫_sgdES_i(E) is also generated from S_i(E) based on Equation (19)(19) $$\begin{array}{c}\hfill S_i\left(E\right)=\left\{\begin{array}{c}{\displaystyle \sum _k\frac{N_k^i}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_i\left(E^{\text{'}}\right)}{E^{\text{'}}}(i\in \mathrm{fuel})}\hfill \\ {\displaystyle \mu \left(E\right)\frac{{\Sigma }_p^i}{E}(i\in \mathrm{non}\mathrm{fuel})}\hfill \end{array}.\right.\end{array}$$(19) . From the viewpoint of computational efficiency, the sub-group scattering source calculation in the present method has an advantage compared with the conventional sub-group method. The scattering source can be calculated with fine energy resolution by Equation (19)(19) $$\begin{array}{c}\hfill S_i\left(E\right)=\left\{\begin{array}{c}{\displaystyle \sum _k\frac{N_k^i}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_i\left(E^{\text{'}}\right)}{E^{\text{'}}}(i\in \mathrm{fuel})}\hfill \\ {\displaystyle \mu \left(E\right)\frac{{\Sigma }_p^i}{E}(i\in \mathrm{non}\mathrm{fuel})}\hfill \end{array}.\right.\end{array}$$(19) , and thus the accurate sub-group scattering source is easily obtained in advance by energy integration of S_i(E). This fact enables to treat the scattering source term as a constant in Equation (23)(23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) , as well as the sub-group cross-section. Therefore, Equation (23)(23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) can be solved independently for each sub-group, and any considerations of sub-group to sub-group transmission probability are not required at all.

Finally, the sub-group flux φⁱ_sg is obtained by integrating the angular flux solutions of Equation (23)(23) $$\begin{array}{c}\hfill \Omega ·\nabla {\psi }_{sg}^i\left(\Omega \right)+{\Sigma }_{t,sg}^i{\psi }_{sg}^i\left(\Omega \right)=\frac{S_{sg}^i}{4\pi },\end{array}$$(23) . In addition to the integro-differential form of equation, φⁱ_sg can also be obtained from the integral form of transport equation based on the collision probability method.

By using the sub-group cross-sections and fluxes, the microscopic effective cross-section for resonance nuclide r, reaction type x and energy group g is derived as (26) $$\begin{array}{ccc}\hfill {\sigma }_{x,g}^{r,i}& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}\hfill \\ & =& \frac{\sum _{sg\in g}\frac{{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}{\int }_{sg}dE{\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\sigma }_{x,sg}^{r,i}{\varphi }_{sg}^i}{\sum _{sg\in g}{\varphi }_{sg}^i}.\hfill \end{array}$$(26)

As described before, σ^{r, i}_{x, sg} is generated by the ultra-fine-group flux solution φ_i(E) in the first-step resonance calculation (see Sections 2.2–2.4). σ^{r, i}_{x, sg} has a fuel rod position dependence within a fuel assembly because Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ) are solved with position-dependent rational coefficients. In addition, σ^{r, i}_{x, sg} has a radial dependence within each fuel pellet because the solution of Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ) can be extended to the radial multi-region geometry for cylindrical fuel rod based on a derivation in Section 2.3.

On the other hand, φⁱ_sg is calculated by using σ^{r, i}_{x, sg} and Sⁱ_sg. φⁱ_sg is obtained with a fine spatial resolution which can be treated by the adopted flux solver. Therefore, at least, fuel rod position and intra-pellet radial dependences of effective cross-sections can be taken into account in the present method. Besides, intra-pellet azimuthal dependence of effective cross-sections can also be considered through φⁱ_sg if MOC is adopted for a sub-group flux calculation.

2.6. Calculation flow

A calculation flow of the two-step resonance treatment is shown in Figures 5 and 6. The calculation procedure is described below. Steps (1(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) )-(9) and (10(10) $$\begin{array}{c}\hfill \sum _{j\in f}{P}_{i\to j}\left(E\right)+P_{i\to nf}\left(E\right)=1.\end{array}$$(10) )-(15) correspond to the first- and second-step resonance calculations, respectively. Note that the details for Step (4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ) are described in Reference [7].

• (1) Read the microscopic ultra-fine-group cross-section data σ_fg (fg: ultra-fine-group number) from the library. (Note: In the following description, continuous energy-dependent parameters are converted to ultra-fine-group form, e.g. σ(E) → σ_fg, for numerical treatment.)

• (2) Read the relative atomic weight to the neutron, A, and the microscopic potential scattering cross-section σ_p.

• (3) Read the heterogeneous region volume V and the number density N.

• (4) Perform the hybrid resonance self-shielding calculation [7] based on integration of the equivalence theory and the ultra-fine-group calculation.

• (5) Generate the ultra-fine-group source terms (Σ^f_{sd, fg} and μ_fg) for 1 ⩽ fg ⩽ FG (FG: the number of ultra-fine-group).

• (6) Generate l_{i, m} and γ_{i, m} based on Equations (16)(16) $$\begin{array}{ccc}& & l_{i,m}=\frac{2R}{\pi }·\left(\sqrt{1-{\rho }^2}+\frac{1}{\rho }{sin}^{-1}\rho \pm \frac{\pi }{2}\rho \right),\hfill \\ & & \left\{\begin{array}{c}m=1:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=2:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_i\hfill \\ m=3:+\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \\ m=4:-\pi \rho /\phantom{\pi \rho 2}2,\rho ={\rho }_{i-1}\hfill \end{array}\right..\hfill \end{array}$$(16) and (15(15) $$\begin{array}{ccc}& & ({\gamma }_{i,1},{\gamma }_{i,2},{\gamma }_{i,3},{\gamma }_{i,4})\hfill \\ & & \equiv \left(\frac{{\rho }_i{l}_{i,1}}{l_i},-\frac{{\rho }_i{l}_{i,2}}{l_i},-\frac{{\rho }_{i-1}{l}_{i,3}}{l_i},\frac{{\rho }_{i-1}{l}_{i,4}}{l_i}\right).\hfill \end{array}$$(15) ), respectively.

• (7) Calculate the ring-region-dependent ultra-fine-group flux φ_{i, fg} for 1 ⩽ fg ⩽ FG based on Equation (18)(18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18) , by using Σ^f_{sd, fg} and μ_fg.

• (8) Calculate the ring-region-dependent ultra-fine-group scattering source S_{i, fg} for 1 ⩽ fg ⩽ FG based on Equation (19)(19) $$\begin{array}{c}\hfill S_i\left(E\right)=\left\{\begin{array}{c}{\displaystyle \sum _k\frac{N_k^i}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_i\left(E^{\text{'}}\right)}{E^{\text{'}}}(i\in \mathrm{fuel})}\hfill \\ {\displaystyle \mu \left(E\right)\frac{{\Sigma }_p^i}{E}(i\in \mathrm{non}\mathrm{fuel})}\hfill \end{array}.\right.\end{array}$$(19) , by using φ_{i, fg} and μ_fg.

• (9) Update φ_{i, fg} based on Equations (20)(20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) and (21(21) $$\begin{array}{c}\hfill {\Sigma }_{sd}^i\left(E\right)\equiv E{S}_i\left(E\right),\end{array}$$(21) ).

• (10) Set the reference cross-section Σ_ref and the cross-section boundaries (maximum: Σ^bound_{sg, max}, minimum:Σ^bound_{sg, min}) for each sub-group used in Equation (25)(25) $$\begin{array}{c}\hfill sg=\left\{E\left|{\Sigma }_{sg,min}^{\mathrm{bound}}\le {\Sigma }_{\mathrm{ref}}\left(E\right)<{\Sigma }_{sg,max}^{\mathrm{bound}}\right.\right\},\end{array}$$(25) .

• (11) Define the discrete energy structure for each sub-group based on Equation (25)(25) $$\begin{array}{c}\hfill sg=\left\{E\left|{\Sigma }_{sg,min}^{\mathrm{bound}}\le {\Sigma }_{\mathrm{ref}}\left(E\right)<{\Sigma }_{sg,max}^{\mathrm{bound}}\right.\right\},\end{array}$$(25) .

• (12) Generate the sub-group cross-section σ^{r, i}_{x, sg} for 1 ⩽ sg ⩽ SG (SG: the number of sub-group) based on Equation (24)(24) $$\begin{array}{c}\hfill {\sigma }_{x,sg}^{r,i}\equiv \frac{{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}.\end{array}$$(24) , by using φ_{i, fg}.

• (13) Generate the sub-group scattering source Sⁱ_sg for 1 ⩽ sg ⩽ SG based on Sⁱ_sg = ∫_sgdES_i(E), by using S_{i, fg}.

• (14) Perform the one-group fixed-source calculations for each sub-group, and obtain the fluxes φⁱ_sg for each region.

• (15) Generate the effective cross-section σ^{r, i}_{x, g} based on Equation (26)(26) $$\begin{array}{ccc}\hfill {\sigma }_{x,g}^{r,i}& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}\hfill \\ & =& \frac{\sum _{sg\in g}\frac{{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}{\int }_{sg}dE{\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\sigma }_{x,sg}^{r,i}{\varphi }_{sg}^i}{\sum _{sg\in g}{\varphi }_{sg}^i}.\hfill \end{array}$$(26) , by using σ^{r, i}_{x, sg} and φⁱ_sg.

Figure 5. Calculation flow of unified resonance treatment (first step).

Figure 6. Calculation flow of unified resonance treatment (second step).

2.7. Relation to the conventional methods

As described in Section 1, the present resonance treatment is a unified approach of the conventional three methods, i.e. equivalence theory, ultra-fine-group slowing-down calculation and sub-group method.

The present method is equivalent to the equivalence theory, if (i) NR approximation of the scattering source, i.e. Σ^f_sd(E) ≈ Σ_p^f and μ(E) ≈ 1 (see Equations (2)(2) $$\begin{array}{c}\hfill {\Sigma }_{sd}^f\left(E\right)=E\sum _k\frac{N_k^f}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_f\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(2) and (3(3) $$\begin{array}{c}\hfill \mu \left(E\right)=\theta +(1-\theta )\frac{E}{{\Sigma }_p^{nf}}\sum _k\frac{N_k^{nf}{\sigma }_p^k}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\varphi }_{nf}\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(3) )), is applied to the first-step calculation, (ii) the resonance absorptions of background isotopes are ignored, and (iii) the number of sub-groups is set as one for the second-step calculation. The present method is also equivalent to the ultra-fine-group slowing-down calculation, if the sub-group is the same as the ultra-fine energy group structure for the second-step calculation. If the problem-independent sub-group parameters are prepared in advance, the present method is equivalent to the sub-group method. Therefore, the present resonance treatment is a generalized approach of the conventional three resonance self-shielding methods.

3. Reaction-rate preservation for reduction of energy discretization error

3.1. Concept for two-step reaction-rate preservation scheme

The effective cross-sections obtained from the direct energy collapsing with flux weight cannot be directly applied to the subsequent multi-group flux calculation due to the energy collapsing error. In order to reduce the error, a reaction-rate preservation scheme is necessary for heterogeneous geometry. In the conventional scheme [6,7] based on the SPH method [16], the group-dependent correction factor fⁱ_g is derived as (27) $$\begin{array}{ccc}\hfill f_g^i& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i}\hfill \\ & =& \left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]/\phantom{\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]}\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]\hfill \\ & =& \frac{{\int }_{g}dE{\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i},\hfill \end{array}$$(27) where ${\tilde{\varphi }}_g^i$ denotes the flux obtained with fⁱ_gσ_{x, g}^{r, i}.

The present resonance self-shielding treatment, especially for the first-step calculation, is based on the ultra-fine-group calculation, and thus the conventional scheme can be applied in principle.

Besides, the present method includes the essence of the equivalence theory. The energy-dependent neutron flux can be regenerated by the analytical form of the multi-term rational equation. Considering the fundamental concept of the reaction-rate preservation scheme in the equivalence theory [17], the following equation to solve fⁱ_g is formulated as a natural extension of the hybrid resonance treatment in Reference [7]: (28) $$\begin{array}{c}\hfill f_g^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,g}^i{l}_{i,m}+{\mu }_g{\alpha }_n}{f_g^i{\Sigma }_{t,g}^i{l}_{i,m}+{\alpha }_n}={\int }_{g}dE{\varphi }_i\left(E\right),\end{array}$$(28) where (29) $$\begin{array}{c}\hfill {\Sigma }_{sd,g}^i\equiv {\int }_{g}dE\eta \left(E\right)\frac{1}{E}{\Sigma }_{sd}^i\left(E\right)={\int }_{g}dE\eta \left(E\right)S_i\left(E\right),\end{array}$$(29) (30) $$\begin{array}{c}\hfill {\mu }_g\equiv {\int }_{g}dE\frac{1}{E}\eta \left(E\right)\mu \left(E\right).\end{array}$$(30)

Equations (20)(20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) and (21(21) $$\begin{array}{c}\hfill {\Sigma }_{sd}^i\left(E\right)\equiv E{S}_i\left(E\right),\end{array}$$(21) ) are used in the above derivation. Here, it should be noted that Equations (27)(27) $$\begin{array}{ccc}\hfill f_g^i& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i}\hfill \\ & =& \left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]/\phantom{\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]}\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]\hfill \\ & =& \frac{{\int }_{g}dE{\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i},\hfill \end{array}$$(27) or (28(28) $$\begin{array}{c}\hfill f_g^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,g}^i{l}_{i,m}+{\mu }_g{\alpha }_n}{f_g^i{\Sigma }_{t,g}^i{l}_{i,m}+{\alpha }_n}={\int }_{g}dE{\varphi }_i\left(E\right),\end{array}$$(28) ) is valid only if φ_i(E) is obtained with high accuracy. As described in Section 2.3, some approximations are applied to the derivation process of φ_i(E), and thus φ_i(E) is not necessarily so accurate for direct application to Equations (27)(27) $$\begin{array}{ccc}\hfill f_g^i& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i}\hfill \\ & =& \left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]/\phantom{\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]}\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]\hfill \\ & =& \frac{{\int }_{g}dE{\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i},\hfill \end{array}$$(27) or (28(28) $$\begin{array}{c}\hfill f_g^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,g}^i{l}_{i,m}+{\mu }_g{\alpha }_n}{f_g^i{\Sigma }_{t,g}^i{l}_{i,m}+{\alpha }_n}={\int }_{g}dE{\varphi }_i\left(E\right),\end{array}$$(28) ). Therefore, a new reaction-rate preservation scheme is developed in the present study so that the influence of prediction error of φ_i(E) is mitigated.

The fundamental form of reaction-rate preservation equation is written as: (31) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i.\end{array}$$(31)

From the viewpoint of consistency with the present two-step resonance treatment described in Section 2, Equation (31)(31) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i.\end{array}$$(31) is divided into the following two equations: (32) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=\sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i,\end{array}$$(32) (33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) where ${\tilde{\varphi }}_{sg}^i$ denotes the flux obtained with fⁱ_sgσ_{x, sg}^{r, i}. Equation (32)(32) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=\sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i,\end{array}$$(32) is the reaction-rate preservation equation between ultra-fine-group and sub-group treatments, and Equation (33)(33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) is the reaction-rate preservation equation between sub-group and multi-group treatments. Equation (32)(32) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=\sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i,\end{array}$$(32) is decomposed for each sub-group as: (34) $$\begin{array}{c}\hfill {\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i.\end{array}$$(34)

Note that Equation (32)(32) $$\begin{array}{c}\hfill {\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=\sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i,\end{array}$$(32) is satisfied if Equation (34)(34) $$\begin{array}{c}\hfill {\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i.\end{array}$$(34) is valid.

In the present resonance treatment, the sub-group-wise correction factor fⁱ_sg and the corrected flux ${\tilde{\varphi }}_{sg}^i$ are generated as the first-step calculation by Equation (34)(34) $$\begin{array}{c}\hfill {\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i.\end{array}$$(34) . Then, the multi-group-wise correction factor fⁱ_g is obtained as the second-step calculation by Equation (33)(33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) . If Equations (28)(28) $$\begin{array}{c}\hfill f_g^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,g}^i{l}_{i,m}+{\mu }_g{\alpha }_n}{f_g^i{\Sigma }_{t,g}^i{l}_{i,m}+{\alpha }_n}={\int }_{g}dE{\varphi }_i\left(E\right),\end{array}$$(28) –(30(30) $$\begin{array}{c}\hfill {\mu }_g\equiv {\int }_{g}dE\frac{1}{E}\eta \left(E\right)\mu \left(E\right).\end{array}$$(30) ) are directly used, a significant error for fⁱ_g is induced, as shown in Section 4.3.4. By utilizing the sub-group based information, direct propagation of φ_i(E) prediction error to fⁱ_g is mitigated. This treatment is called the two-step reaction-rate preservation scheme in the present study, and its detail is described in the following sections.

3.2. First-step correction factor calculation for sub-group level

The correction factor calculation for the first step is carried out based on Equation (34)(34) $$\begin{array}{c}\hfill {\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i.\end{array}$$(34) . By the same derivation scheme for the reaction-rate preservation in the hybrid resonance treatment, Equation (34)(34) $$\begin{array}{c}\hfill {\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)=f_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i.\end{array}$$(34) is rewritten as (35) $$\begin{array}{c}\hfill f_{sg}^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i{\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}={\int }_{sg}dE{\varphi }_i\left(E\right),\end{array}$$(35) where (36) $$\begin{array}{c}\hfill {\Sigma }_{sd,sg}^i\equiv {\int }_{sg}dE\eta \left(E\right)S_i\left(E\right),\end{array}$$(36) (37) $$\begin{array}{c}\hfill {\mu }_{sg}\equiv {\int }_{sg}dE\frac{1}{E}\eta \left(E\right)\mu \left(E\right).\end{array}$$(37)

Equations (20)(20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) and (21(21) $$\begin{array}{c}\hfill {\Sigma }_{sd}^i\left(E\right)\equiv E{S}_i\left(E\right),\end{array}$$(21) ) are used in the above derivation.

Equation (35)(35) $$\begin{array}{c}\hfill f_{sg}^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i{\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}={\int }_{sg}dE{\varphi }_i\left(E\right),\end{array}$$(35) is used to calculate a correction factor for reaction-rate preservation, which is consistent with the resonance treatment in Sections 2.2 and 2.3. All the terms except for fⁱ_sg are byproducts of the solution of the methods described in Sections 2.2 and 2.3. fⁱ_sg is easily obtained from iterative calculation until fⁱ_sg is converged. From Equation (35)(35) $$\begin{array}{c}\hfill f_{sg}^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i{\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}={\int }_{sg}dE{\varphi }_i\left(E\right),\end{array}$$(35) , fⁱ_sg for (s + 1)th iteration is generated by using sth result as follows: (38) $$\begin{array}{c}\hfill f_{sg}^i(\mathrm{s}+1)=\frac{{\int }_{sg}dE{\varphi }_i\left(E\right)}{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i\left(\mathrm{s}\right){\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}}.\end{array}$$(38)

The method includes the iteration scheme, but any additional iterative flux calculations are no longer required. The calculation time of fⁱ_sg is negligible compared with the conventional scheme.

3.3. Second-step correction factor calculation for multi-group level

The correction factor calculation for the second step is carried out based on Equation (33)(33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) . By substituting Equation (26)(26) $$\begin{array}{ccc}\hfill {\sigma }_{x,g}^{r,i}& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}\hfill \\ & =& \frac{\sum _{sg\in g}\frac{{\int }_{sg}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{sg}dE{\varphi }_i\left(E\right)}{\int }_{sg}dE{\varphi }_i\left(E\right)}{\sum _{sg\in g}{\int }_{sg}dE{\varphi }_i\left(E\right)}=\frac{\sum _{sg\in g}{\sigma }_{x,sg}^{r,i}{\varphi }_{sg}^i}{\sum _{sg\in g}{\varphi }_{sg}^i}.\hfill \end{array}$$(26) into Equation (33)(33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) , the correction factor is derived as (39) $$\begin{array}{ccc}\hfill f_g^i& =& \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{\frac{{\sum }_{sg\in g}{\sigma }_{x,sg}^{r,i}{\varphi }_{sg}^i}{{\sum }_{sg\in g}{\varphi }_{sg}^i}{\tilde{\varphi }}_g^i}\hfill \\ & =& \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\sigma }_{x,sg}^{r,i}{\varphi }_{sg}^i}\frac{{\sum }_{sg\in g}{\varphi }_{sg}^i}{{\tilde{\varphi }}_g^i}.\hfill \end{array}$$(39)

In the above equation, the microscopic cross-section is included in fⁱ_g, and thus it depends on resonance nuclide and reaction type. This definition of correction factor leads to the large computational burdens.

In order to avoid the nuclide and reaction type dependent correction factor calculations, the modified correction factor ${\tilde{f}}_g^i$ and the effective cross-section ${\tilde{\sigma }}_{x,g}^{r,i}$ are defined so as to satisfy the following relation: (40) $$\begin{array}{c}\hfill {\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}=f_g^i{\sigma }_{x,g}^{r,i}.\end{array}$$(40)

By substituting Equation (40)(40) $$\begin{array}{c}\hfill {\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}=f_g^i{\sigma }_{x,g}^{r,i}.\end{array}$$(40) into the right-hand side of Equation (33)(33) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i=f_g^i{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i,\end{array}$$(33) , the following equation is obtained: (41) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i={\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}{\tilde{\varphi }}_g^i.\end{array}$$(41)

Here, ${\tilde{\sigma }}_{x,g}^{r,i}$ is defined as an effective cross-section in which the energy discretization error associated with the collapsing only from the ultra-fine-group to the sub-group is removed, i.e., (42) $$\begin{array}{c}\hfill {\tilde{\sigma }}_{x,g}^{r,i}\equiv \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}.\end{array}$$(42)

By substituting Equation (42)(42) $$\begin{array}{c}\hfill {\tilde{\sigma }}_{x,g}^{r,i}\equiv \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}.\end{array}$$(42) into Equation (41)(41) $$\begin{array}{c}\hfill \sum _{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i={\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}{\tilde{\varphi }}_g^i.\end{array}$$(41) , the modified multi-group correction factor ${\tilde{f}}_g^i$ is derived as (43) $$\begin{array}{c}\hfill {\tilde{f}}_g^i=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\tilde{\sigma }}_{x,g}^{r,i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{{\tilde{\varphi }}_g^i}.\end{array}$$(43)

As shown in Equation (43)(43) $$\begin{array}{c}\hfill {\tilde{f}}_g^i=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\tilde{\sigma }}_{x,g}^{r,i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{{\tilde{\varphi }}_g^i}.\end{array}$$(43) , ${\tilde{f}}_g^i$ does not depend on resonance nuclide and reaction type like an SPH factor. In other words, ${\tilde{\sigma }}_{x,g}^{r,i}$ is defined in Equation (42)(42) $$\begin{array}{c}\hfill {\tilde{\sigma }}_{x,g}^{r,i}\equiv \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}.\end{array}$$(42) so that the resonance nuclide and reaction type dependence of the correction factor can be removed.

Finally, the corrected effective cross-section, which is provided for the subsequent multi-group flux and depletion calculations in a lattice physics code, is ${\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}$. It should be noted that the effective cross-section, which is comparable to the continuous energy Monte-Carlo results tallied by a simple flux weight scheme, is σ^{r, i}_{x, g}.

3.4. Calculation flow

Calculation flow of the two-step reaction-rate preservation scheme established in Sections 3.1–3.3 is shown in Figure 7. The calculation procedure is described as follows:

• (1) Generate the sub-group-wise correction factor fⁱ_sg based on Equation (35)(35) $$\begin{array}{c}\hfill f_{sg}^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i{\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}={\int }_{sg}dE{\varphi }_i\left(E\right),\end{array}$$(35) with iteration scheme of Equation (38)(38) $$\begin{array}{c}\hfill f_{sg}^i(\mathrm{s}+1)=\frac{{\int }_{sg}dE{\varphi }_i\left(E\right)}{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,sg}^i{l}_{i,m}+{\mu }_{sg}{\alpha }_n}{f_{sg}^i\left(\mathrm{s}\right){\Sigma }_{t,sg}^i{l}_{i,m}+{\alpha }_n}}.\end{array}$$(38) .

• (2) Perform one-group fixed-source calculations for each sub-group by using the corrected sub-group cross-section fⁱ_sgσ_{x, sg}^{r, i} as an input, and obtain the sub-group flux ${\tilde{\varphi }}_{sg}^i$. This calculation can be carried out independently.

• (3) Generate the partially corrected effective cross-section ${\tilde{\sigma }}_{x,g}^{r,i}$, in which the energy discretization error from the ultra-fine-group to the sub-group collapsing is removed, based on Equation (42)(42) $$\begin{array}{c}\hfill {\tilde{\sigma }}_{x,g}^{r,i}\equiv \frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}.\end{array}$$(42) .

• (4) Set the initial value of ${\tilde{f}}_g^i$ as ${\tilde{f}}_g^i=1$.

• (5) Perform one-group fixed-source calculation for each multi-group by using the corrected effective cross-section ${\tilde{f}}_g^i{\tilde{\sigma }}_{x,g}^{r,i}$ as an input, and obtain the multi-group flux ${\tilde{\varphi }}_g^i$. This calculation can be carried out independently.

• (6) Generate the multi-group-wise modified correction factor ${\tilde{f}}_g^i$ based on Equation (43)(43) $$\begin{array}{c}\hfill {\tilde{f}}_g^i=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\tilde{\sigma }}_{x,g}^{r,i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{\frac{{\sum }_{sg\in g}{f}_{sg}^i{\sigma }_{x,sg}^{r,i}{\tilde{\varphi }}_{sg}^i}{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{\tilde{\varphi }}_g^i}=\frac{{\sum }_{sg\in g}{\tilde{\varphi }}_{sg}^i}{{\tilde{\varphi }}_g^i}.\end{array}$$(43) .

• (7) Iterate (5(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) ) and (6) until ${\tilde{f}}_g^i$ is converged.

Figure 7. Calculation flow of two-step reaction-rate preservation scheme.

Steps (1(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) )–(3) and (4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) )–(7) correspond to the first- and second-step correction factor calculations, respectively.

4. Verification

4.1. Concept

In this section, the present resonance treatment is verified by the numerical calculations. The present method is composed of multiple components based on a unified theory, and thus each component is verified first (Section 4.3). Then, the effective cross-sections, i.e. final products of the resonance calculation, are analyzed for various cases as a total verification (Section 4.4). The reference solutions are obtained from the continuous energy Monte-Carlo calculations based on the exact geometry modeling and fine cross-section representation with continuous energy resolution. A part of the references is obtained from the heterogeneous ultra-fine-group slowing-down calculation based on a deterministic method.

The analysis condition is shown in Section 4.2. The verification results for fundamental parameters of the present method are shown in Section 4.3. The application results for various pin-cell and multi-cell problems are shown in Section 4.4. The verification list is shown in Table 3.

Table 3. Verification list

Section	Objective	Geometry	Verified model	Verified parameter	Sub-section
4.3	Confirm prediction accuracy of main products in the present method	Pin-cell	Slowing-down calculation with multi-term rational equation	Ultra-fine-group neutron flux	4.3.1
			Sub-group definition according to resonance cross-section level	Sub-group cross-section	4.3.2
			Final product of resonance calculation	Effective cross-section	4.3.3
			Reaction-rate preservation scheme	Cross-section correction factor	4.3.4
			Final product of lattice physics calculation	Reaction-rate	4.3.5
4.4	Confirm applicability of the present method for various pin-cell and multi-cell problems	Pin-cell	Spatial self-shielding treatment with spatial variation of ultra-fine-group macroscopic cross-section	Radially dependent effective cross-sections with isotope composition distribution	4.4.1
				Radially dependent effective cross-sections with temperature distribution	4.4.2
			Spatial self-shielding treatment	Radially and azimuthally dependent effective cross-sections	4.4.3
		3×3 multi-cell	Spatial self-shielding treatment for irregular lattice	Azimuthally dependent effective cross-sections with large water region	4.4.4

4.2. Analysis condition

Analysis condition of the new resonance treatment is described in this section.

In order to perform each process in the calculation scheme of the present method (see Sections 2.6 and 3.4), a standalone program is developed. The MHI (Mitsubishi Heavy Industries, Ltd.) lattice physics code GALAXY (Geometrically Arbitrary, Lattice physics and Assembly calculation code in X-Y coordinate system) [4,7,18,21] is used for several steps in the process. In concrete, ‘(1)–(5) and (14(14) $$\begin{array}{c}\hfill P_{i\to nf}\left(E\right)=\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\end{array}$$(14) )’ in Figures 5 and 6 and ‘(2) and (5(5) $$\begin{array}{ccc}\hfill {\Sigma }_t^i\left(E\right){\varphi }_i\left(E\right)V_i& =& \sum _{j\in f}{S}_j\left(E\right)V_j{P}_{j\to i}\left(E\right)\hfill \\ & & +S_{nf}\left(E\right)V_{nf}{P}_{nf\to i}\left(E\right).\hfill \end{array}$$(5) )’ in Figure 7 are performed by GALAXY. The GALAXY ultra-fine-group cross-section library (120,000 groups) [7] based on ENDF/B-VII.0 [19] is also utilized.

The number of rational equation terms used for one-region fuel is set to N = 2. Note that the total number of terms used for multi-region fuel is 4N = 8 (see Equation (18)(18) $$\begin{array}{ccc}\hfill {\varphi }_i\left(E\right)& =& \frac{1}{E}·\left[\{1-P_{i\to nf}\left(E\right)\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}+\mu \left(E\right)P_{i\to nf}\left(E\right)\right]\hfill \\ & =& \frac{1}{E}\left[\begin{array}{c}\left\{\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\right.\hfill \\ -\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right\}·\frac{{\Sigma }_{sd}^f\left(E\right)}{{\Sigma }_t^f\left(E\right)}\hfill \\ +\mu \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\hfill \end{array}\right]\hfill \\ & =& \frac{1}{E}\left[\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right.\hfill \\ & & +\left.\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}\right]\hfill \\ & =& \frac{1}{E}\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n}.\hfill \end{array}$$(18) ). The microscopic effective capture cross-sections of ²³⁸U are generated for 88th energy group (6.16–7.52 eV) of the XMAS 172 energy group structure [20], in which the wide resonance cross-section of ²³⁸U is included. The number of sub-groups is 5, and the sub-group energy structure is determined by Equation (25)(25) $$\begin{array}{c}\hfill sg=\left\{E\left|{\Sigma }_{sg,min}^{\mathrm{bound}}\le {\Sigma }_{\mathrm{ref}}\left(E\right)<{\Sigma }_{sg,max}^{\mathrm{bound}}\right.\right\},\end{array}$$(25) based on the ultra-fine-group macroscopic absorption cross-section of the fuel region as Σ_ref(E). The sub-group boundary for the absorption cross-section is determined by a division with equal interval in logarithmic scale. The sub-group transport calculation is performed by the collision probability method based on the equivalent Dancoff method [21]. These conditions are applied to almost all the analyses in this paper. The sensitivities for the above calculation conditions against effective cross-sections are investigated in advance, and it is confirmed that the calculation conditions are sufficient from the view point of prediction accuracy for effective cross-sections.

UO₂ (4.8wt% ²³⁵U) pin-cells with HFP (hot full power) operating conditions are set as the typical LWR neutron spectrum condition. The specifications of the pin-cell problem are shown in Table 4, and the geometrical configuration is shown in Figure 8. The pellet region is sub-divided into 10 equal volume rings for radial direction. The fuel composition and temperature are assumed to be flat within a pellet. Calculation for the one-region pellet model is also performed for comparisons.

Table 4. Specifications of the pin-cell model

Item		Specification
	Fuel	UO2 (4.8wt% ^235U)
Material	Cladding	Zr
	Moderator	Borated water
	Fuel	976 K
Temperature	Cladding	600 K
	Moderator	580 K
Boron concentration		1000 ppm
	Cell pitch	1.26 cm
	Pellet radius	0.4095 cm
Geometry	Cladding outer radius	0.475 cm
	Cladding thickness (gap is omitted)	0.0655 cm

Figure 8. Geometry of pin-cell model.

In Sections 4.3.1, 4.3.3, and 4.3.5, the calculation results by the present method are compared with those by the continuous energy Monte-Carlo code MVP (Monte Carlo code for vector processors) [22]. The ENDF/B-VII.0 nuclear data library is used in all the MVP calculations to be consistent with the calculations by the present method. The total number of neutron sampling for MVP calculations is set to 100 million histories, in which the 1σ statistical error of the effective cross-sections is about 0.1%.

For the scattering kernel treatment, all the analyses including the reference calculations by the continuous energy Monte-Carlo are performed with the conventional asymptotic scattering model, as shown in Equation (2)(2) $$\begin{array}{c}\hfill {\Sigma }_{sd}^f\left(E\right)=E\sum _k\frac{N_k^f}{1-{\alpha }_k}{\int }_E^{E/{\alpha }_k}\frac{d{E}^{\text{'}}{\sigma }_s^k\left(E^{\text{'}}\right){\varphi }_f\left(E^{\text{'}}\right)}{E^{\text{'}}},\end{array}$$(2) . In this paper, discussion for the selection of scattering model, i.e. the asymptotic model or the exact resonance scattering model [23,24], is out of scope. The present resonance treatment is verified with reference calculations based on the same scattering model.

In Sections 4.3.2 and 4.3.4, the calculation results by the present method are compared with those by the direct ultra-fine-group slowing-down calculation with heterogeneous geometry. The heterogeneous slowing-down calculation is performed by the equivalent Dancoff method [21] implemented in GALAXY, which is a different resonance treatment from the hybrid resonance treatment [7] and the present method.

The equivalent Dancoff method is applied to a part of the reference calculations, in which the continuous energy Monte-Carlo calculations cannot be used to obtain reference results since specific conditions dedicated for the deterministic method are used and cannot be reproduced by the Monte-Carlo method. The reliability of the equivalent Dancoff method itself is also shown in Section 4.3.3 through comparisons of effective cross-sections against the MVP results.

4.3. Verification for fundamental parameters of the present method

4.3.1. Ultra-fine-group neutron flux

The ultra-fine-group neutron flux by the first-step calculation of the present method is compared with the continuous energy Monte-Carlo result in this section.

As described in Section 2.2, a set of the slowing-down equations (Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) )) are solved for one-region pellet system, based on the hybrid resonance treatment [7] of the equivalence theory and the ultra-fine-group slowing-down calculation. Σ^f_sd(E) and μ(E) are numerically obtained as byproducts of the solution for Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ). Then by Equations (20)(20) $$\begin{array}{c}\hfill {\varphi }_i\left(E\right)=\frac{1}{E}\eta \left(E\right)\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^i\left(E\right)l_{i,m}+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_{i,m}+{\alpha }_n},\end{array}$$(20) and (21(21) $$\begin{array}{c}\hfill {\Sigma }_{sd}^i\left(E\right)\equiv E{S}_i\left(E\right),\end{array}$$(21) ), which are derived in Section 2.3, the ultra-fine-group flux for each sub-region of a pellet is analytically obtained.

The ultra-fine-group fluxes for each sub-region of a pellet by the present method and those by the continuous energy Monte-Carlo calculation (MVP) are shown in Figure 9 including 88th energy group (6.16–7.52 eV) of XMAS 172 energy group structure. The results for Regions 1, 8, 9, and 10 (sequentially numbered from the center to the surface of a pellet) are plotted for convenience. As shown in Figure 9, the present method (first-step calculation) can roughly predict the ultra-fine-group neutron flux comparable to the continuous energy Monte-Carlo results. However, a slight difference is observed. The flux differences are induced by the approximation of spatially dependent fuel escape probability and scattering source within a pellet, which are applied to the flux derivation scheme in Section 2.3. As shown in Section 4.3.3, the differences of the ultra-fine-group flux directly cause the differences of effective cross-section between the present method ‘without second-step calculation’ and the MVP. The differences of effective cross-section owing to ultra-fine-group flux estimation are mitigated by the second-step calculation based on the sub-group flux, as shown in Section 4.3.3.

Figure 9. Comparison of ultra-fine-group fluxes between the present method (first-step calculation) and the continuous energy Monte-Carlo calculation (MVP).

A breakdown of the calculation time for ultra-fine-group flux calculation in all the ring regions within a pellet is shown in Table 5. As shown in Table 5, the calculation time is very short, and thus the time is negligible against the total computation time of the lattice physics calculations. On the basis of its short computation time, the present method can be easily applied to the large and heterogeneous geometry such as a fuel assembly. In the present method, at first, the fuel-rod-wise ultra-fine-group flux is independently obtained for each fuel cell by Equations (1)(1) $$\begin{array}{c}\hfill {\varphi }_f\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\alpha }_n},\end{array}$$(1) –(4(4) $$\begin{array}{c}\hfill {\varphi }_{nf}\left(E\right)=\frac{1}{E}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd}^f\left(E\right)l_f+\mu \left(E\right){\epsilon }_n{\alpha }_n}{{\Sigma }_t^f\left(E\right)l_f+{\epsilon }_n{\alpha }_n}.\end{array}$$(4) ) and the cell-dependent coefficients α_n, β_n, ϵ_n and θ. Then, the obtained ultra-fine-group data in each fuel rod is expanded for each ring region by using the ring-dependent coefficients γ_m. Utilization of a set of coefficients for the rational equation improves the computational efficiency of the ultra-fine-group calculation, compared with the direct heterogeneous ultra-fine-group transport calculation.

Table 5. Calculation time for the ultra-fine-group flux

Process	CPU time (s)
Ultra-fine-group cross-section library read (only one time)^a	0.31
Ultra-fine-group flux calculation for two region problem^a	0.22
Spatially dependent flux calculation	0.07
Spatially dependent scattering source calculation	0.04
Flux update	0.08

^aReference [7].

4.3.2. Sub-group cross-section

The sub-group cross-sections in fuel regions are generated by the present method. The results are compared with the reference solutions obtained from the direct heterogeneous ultra-fine-group calculation by the equivalent Dancoff method.

The sub-group cross-sections and the corresponding differences from the direct heterogeneous ultra-fine-group calculation results are shown in Figure 10. The differences of the effective cross-section between the present method ‘without second-step calculation’ and the reference solution are also shown in the figure. As shown in Figure 10, the present method can accurately predict the sub-group cross-sections, while the effective cross-section generated by the present method ‘without second-step calculation’ includes some differences.

Figure 10. Sub-group cross-sections and their differences from the direct heterogeneous ultra-fine-group calculation results.

The influence of ultra-fine-group flux error (see Section 4.3.1) on sub-group cross-section is mitigated by defining the sub-group energy structure, not for energy but for macroscopic absorption cross-section level.

4.3.3. Effective cross-section

The effective cross-sections in fuel regions are generated by the present method. The results are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The results based on the equivalence theory (SDGM [4] with the NR approximation), the direct heterogeneous ultra-fine-group calculation (the equivalent Dancoff method [21]) and the present method ‘without second-step calculation’ are also obtained for comparison.

The effective cross-sections and the corresponding differences from the MVP results are shown in Figure 11. As shown in Figure 11, the present method (first + second-step calculation) can accurately predict the effective cross-sections.

Figure 11. Effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP).

In contrast, the results obtained by the equivalence theory with NR approximation and by the present method ‘without second-step calculation’ show the large difference from the reference solutions. For the latter method, the difference is due to the prediction error of ultra-fine-group fluxes in the case of treating radially sub-divided fuel region (see Section 4.3.1).

From the result, the effectiveness of incorporating the sub-group method as the second-step calculation is confirmed. The influence of ultra-fine-group flux error on the final effective cross-section is efficiently reduced by the two-step flux calculation scheme of the present method, i.e. the flux calculation for ‘coarse geometry + fine energy’ (first step) and that for ‘fine geometry + coarse energy’ (second step). The prediction accuracy of effective cross-section is appropriately improved by reflecting the more detailed information of spatial neutron transport within a pellet based on the second-step sub-group calculation.

The differences of the effective cross-sections from those of the MVP results are summarized in Table 6. The results obtained by one-region pellet model are also shown in the table. From Table 6, the present two-step method provides an excellent result for multi-region case, which is almost comparable to the direct heterogeneous ultra-fine-group calculation.

Table 6. Differences of effective cross-section from the continuous energy Monte-Carlo calculation (MVP)

	Relative difference of effective cross-section from MVP (%)
	Pellet one-region	Pellet multi-region problem
Method	problem	Average	Maximum
Continuous energy Monte-Carlo calculation (MVP)^a	0.06	0.11	0.15
Equivalence theory (NR approximation)	1.53	6.23	15.92
Direct heterogeneous ultra-fine-group calculation	0.41	0.20	1.04
Present method (without second-step calculation)	1.08	−4.08	−9.52
Present method (first + second-step calculation)	1.17	0.37	1.01

^aStatistical error.

The calculation efficiency is also an important viewpoint for the total performance of a target method. The required number of one-group fixed-source transport calculations is estimated in Table 7. In the conventional equivalence theory based on the Dancoff method, the one-group calculation is required only one time in the black neutron absorber condition. If the equivalence theory based on the gray resonance treatment [4] is adopted, the one-group calculation is required several times for a wide range of macroscopic total cross-section conditions. After all, the conventional equivalence theory requires only 1-20 times of one-group transport calculations. In contrast to the equivalence theory, the direct heterogeneous ultra-fine-group method generally requires 10,000–200,000 times of one-group transport calculations, which is large computational burden.

Table 7. Estimation for the number of one-group fixed-source transport calculations

	The number of one-group
Method	flux calculations (–)
Equivalence theory	1–20
Direct heterogeneous ultra-fine-group calculation	10,000–200,000
Present method	500–1000

The present method has an intermediate feature between the above two methods. If the number of multi-groups, in which the resonance self-shielding treatment is required, is 100, and the number of sub-groups is set to 5 for each group, the required number of one-group transport calculations is 100×5 = 500. As a result, the required number of one-group transport calculations in the present method is less than 1/10 of that in the direct heterogeneous ultra-fine-group method, as shown in Table 7.

The calculation time for a fuel assembly geometry is an important point for practical lattice physics calculations. The calculation time for a fuel assembly is estimated in Table 8. In the estimation, the steady-state calculations for a typical PWR type 17×17 4.8wt% UO₂ assembly are carried out by the lattice physics code GALAXY.

Table 8. Brief estimation for the calculation time on fuel assembly geometry

		CPU time ratio (–)
Pellet division	Method	Comparison with direct het. ufg. calc.^a	Contribution of resonance calc.^b
	Equivalence theory	0.17	0.15
No	Direct heterogeneous ultra-fine-group calculation	1.00	0.49
	Present method	0.30	0.23
	Equivalence theory	0.08	0.39
Yes	Direct heterogeneous ultra-fine-group calculation	1.00	0.88
	Present method^c	0.12	0.46

^aCPU time ratio between ‘specific method’ and ‘direct ultra-fine-group calculation’ for resonance calculation part.

^bCPU time ratio between ‘resonance calculation’ and ‘resonance + multi-group flux calculation’ for each method.

^cEstimated values from other results.

As for the treatment of fuel pellet, the two cases (‘no-division’ and ‘sub-division for radial direction with 10 equal volume rings’) are analyzed. For each case, the resonance treatment methods based on the equivalence theory [4] and the direct ultra-fine-group calculation [21] are applied for comparison with the present method, and their calculation times are directly measured. Since the present method is not implemented into the GALAXY code at the current development state, the result for the present method is inferred by the results for the other methods and the quantitative information in Tables 5 and 7. The multi-group flux calculation is performed by ‘172 group CCCP (current-coupling collision probability) method + 22 group MOC’ [21] for all cases.

As shown in Table 8, the estimated calculation times by the present method are less than one-third of those by the direct ultra-fine-group calculation for resonance calculation part. The ratios of resonance calculation part against the ‘resonance + flux’ calculations are less than 1/2 in the present method. From the results, the present method can efficiently generate effective cross-sections for large and complicated geometry such as a fuel assembly, comparing with the direct ultra-fine-group calculation.

Finally, the total performance is benchmarked for major resonance treatments, taking into account of the quantitative information in Tables 6–8. The result of qualitative comparison is summarized in Table 9. As shown in Table 9, the present method is better than all the three conventional methods, from the viewpoint of both calculation accuracy and time.

Table 9. Qualitative comparison of overall performance for resonance self-shielding treatments

Method		Calculation accuracy	Calculation time
	Equivalence theory	Acceptable (normal design condition)	Excellent
Conventional	Ultra-fine-group calculation	Excellent	Acceptable (pin-cell)
	Sub-group method	good	good
Present		Excellent	good

4.3.4. Cross-section correction factor

The correction factors for effective cross-section in fuel regions to reduce energy discretization error are generated by the present method. The results are compared with the reference solutions obtained from the direct heterogeneous ultra-fine-group calculation by the equivalent Dancoff method (see Equation (27)(27) $$\begin{array}{ccc}\hfill f_g^i& =& \frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\sigma }_{x,g}^{r,i}{\tilde{\varphi }}_g^i}\hfill \\ & =& \left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]/\phantom{\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i}\right]\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]}\left[\frac{{\int }_{g}dE{\sigma }_x^r\left(E\right){\varphi }_i\left(E\right)}{{\int }_{g}dE{\varphi }_i\left(E\right)}\right]\hfill \\ & =& \frac{{\int }_{g}dE{\varphi }_i\left(E\right)}{{\tilde{\varphi }}_g^i},\hfill \end{array}$$(27) in Section 3.1). The results based on the present method ‘without second-step sub-group calculation’ (see Equation (28)(28) $$\begin{array}{c}\hfill f_g^i\sum _{m=1}^4{\gamma }_{i,m}\sum _{n=1}^N{\beta }_n\frac{{\Sigma }_{sd,g}^i{l}_{i,m}+{\mu }_g{\alpha }_n}{f_g^i{\Sigma }_{t,g}^i{l}_{i,m}+{\alpha }_n}={\int }_{g}dE{\varphi }_i\left(E\right),\end{array}$$(28) in Section 3.1) are also obtained for comparison.

In the two-step reaction-rate preservation scheme, the relative convergence criterion of the correction factors is set as 0.1%. Only 10 iterations are necessary for the second-step calculation in this verification with this convergence criterion.

The correction factors and the corresponding differences from the direct heterogeneous ultra-fine-group calculation results are shown in Figure 12. As shown in Figure 12, the present method (first + second-step calculation) can accurately predict the correction factor.

Figure 12. Correction factors and their differences from the direct heterogeneous ultra-fine-group calculation results.

In contrast, the results obtained by the present method ‘without second-step calculation’ show the large difference from the reference solutions. The prediction error of ultra-fine-group fluxes in the case of treating radially sub-divided fuel region is the cause of the difference (see Section 4.3.1).

From the result, the effectiveness of incorporating the sub-group method as the second-step calculation is confirmed. The influence of ultra-fine-group flux error on the final correction factor is efficiently reduced by the two-step reaction-rate preservation scheme, i.e. the ultra-fine-group to the sub-group collapsing, and the sub-group to the multi-group collapsing.

4.3.5. Reaction-rate

One of the final products of lattice physics calculation is reaction-rate. In order to confirm the influence of effective cross-section differences on reaction-rate, the energy-integrated macroscopic absorption rates in fuel regions are generated from the multi-group transport calculation by the lattice physics code GALAXY. The radially dependent microscopic effective cross-section set, which is generated from a standalone program based on the present method, is directly supplied to the GALAXY.

The results are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The results based on the direct heterogeneous ultra-fine-group calculation by GALAXY are also obtained for comparison. In the reaction-rate calculation for each method, the neutron fluxes are normalized so that the volume and energy integrated neutron generation rate is unity.

The absorption rates and the corresponding differences from the MVP results are shown in Figure 13. In this analysis, the 1σ statistical error of absorption rates by MVP is about 0.02%. The maximum relative differences from the MVP results are −0.4% for the present method, and −0.2% for the direct ultra-fine-group calculation, respectively. The difference is less than ±0.1% for most ring regions in both methods, thus the present method can accurately predict the reaction-rates, which is comparable to the direct ultra-fine-group calculation.

Figure 13. Reaction-rates and their differences from the continuous energy Monte-Carlo calculation (MVP).

4.4. Application for various pin-cell and multi-cell problems

4.4.1. Radially dependent effective cross-sections with non-uniform isotope composition for unit pin-cell

In this section, the pin-cell with non-uniform isotope composition is treated, in which the fuel composition is radially distributed. As a typical application for the LWR lattice calculations, the depleted fuel and the annular fuel are taken into account. For both fuels, the calculation conditions except for the number densities in fuel regions are the same as those given in Section 4.2.

For the depleted fuel case, the GALAXY depletion calculation is performed in advance, and the number densities for each nuclide and each ring region within a pellet are obtained. In this verification, the number densities only for main actinide nuclides, i.e. ²³⁵U, ²³⁸U, ²³⁹Pu, ²⁴⁰Pu, and ²⁴¹Pu, are extracted for simplicity on the burnup condition of 60GWd/t, as shown in Figure 14. The typical Pu build-up effect on a peripheral region of the pellet can be observed from the figure.

Figure 14. Distribution of fuel isotope composition within a pellet.

In this verification, the microscopic effective capture cross-sections of ²³⁸U are generated with the given region-dependent number densities. The effective cross-sections generated by the present method are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The results are shown in Figure 15. As shown in the figure, the present method can accurately predict the effective cross-sections with radially distributed isotope composition. The spatially dependent biases against fresh fuel condition also agree well with the Monte-Carlo results, and therefore the consistency between uniform and non-uniform fuel composition treatments is confirmed for the present method.

Figure 15. Effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) with non-uniform isotope composition.

For the annular fuel case, material for the inner-most region of the pellet, which is sub-divided into 10 equal volume rings for radial direction, is assumed as air. Similar with the above depleted fuel case, the comparison of the effective cross-sections is shown in Figure 16. As shown in the figure, the present method can accurately predict the effective cross-sections of annular fuel. The consistency between solid and annular fuel treatments is confirmed for the present method.

Figure 16. Effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) for annular fuel.

4.4.2. Radially dependent effective cross-sections with non-uniform temperature for unit pin-cell

In this section, the pin-cell with non-uniform temperature is treated, in which the fuel temperature is radially distributed. The calculation conditions except for the fuel temperatures are the same as those given in Section 4.2.

The fuel temperature distribution used for the verification is shown in Figure 17. The corresponding effective fuel temperature which is averaged with a chord weight [15] is the same as that in Table 4.

Figure 17. Distribution of fuel temperature within a pellet.

In this verification, the microscopic effective capture cross-sections of ²³⁸U are generated with the given region-dependent temperatures. The effective cross-sections generated by the present method are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The results are shown in Figure 18. As shown in the figure, the spatially dependent biases against flat temperature condition agree well with the Monte-Carlo results, and therefore the consistency between uniform and non-uniform fuel temperature treatments is confirmed for the present method.

Figure 18. Effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) with non-uniform fuel temperature.

4.4.3. Radially and azimuthally dependent effective cross-sections for unit pin-cell

In this section, the radially and azimuthally divided pin-cell is treated. The calculation conditions except for the sub-division of fuel region in azimuthal direction are the same as those given in Section 4.2. The fuel pellet is sub-divided into 32 sector regions for 2π with equal azimuthal angle interval.

In this verification, the azimuthally dependent microscopic effective capture cross-sections of ²³⁸U are generated. The sub-group cross-sections are generated as the azimuthally independent values in the first step of the present method. The azimuthal dependence against the final effective cross-sections is taken into account by the second-step sub-group flux calculation with direct two-dimensional geometry modeling based on MOC.

The effective cross-sections generated by the present method are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The total number of neutron sampling for MVP calculations is set to 500 million histories, in which the 1σ statistical error of the effective cross-sections is about 0.1-0.2%.

The results are shown in Figure 19. As shown in the figure, the present method can accurately predict the azimuthally dependent effective cross-sections, in which the cross-sections tend to be larger for the diagonal directions than that for the horizontal and vertical directions due to the difference of adjacent moderator region volume. The spatially dependent biases against one-region fuel condition also agree well with the Monte-Carlo results, and therefore the consistency between one-region and azimuthal multi-region treatments is confirmed for the present method.

Figure 19. Azimuthally dependent effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) for unit pin-cell.

Both the radially and azimuthally dependent effective cross-sections are also compared with the reference solutions from MVP. The results for azimuthally -dependent effective cross-section ratios for each ring region are shown in Figure 20. The results for radially divided regions 1, 7, 8, 9, and 10 (sequentially numbered from the center to the surface of a pellet) are plotted for convenience. As shown in the figure, the present method can accurately predict both the radial and azimuthal dependences of effective cross-sections.

Figure 20. Azimuthally dependent effective cross-section ratios in each ring region for unit pin-cell.

As can be observed from Figure 20, the azimuthal dependence of effective cross-sections tends to be larger for the pellet surface (region 10) than that for the pellet center (region 1) due to the difference of distance from moderator region. The heterogeneous effect, which is a main contributor for yielding azimuthal dependence, is small for the pellet center. Thus the effective cross-sections are azimuthally flat for pellet center region.

Though the difference of effective cross-sections between the present method and MVP is not shown in Figure 20, the difference is confirmed to be less than ±2% for all the 320 regions, which is a consistent result with Figure 11 in Section 4.3. Therefore, the present method can accurately predict spatially dependent effective cross-sections for both radially and azimuthally sub-divided fuel condition.

4.4.4. Azimuthally dependent effective cross-sections for 3×3 multi-cell including large water region

Finally, the 3×3 multi-cell divided in azimuthal direction is treated. A large water cell is set in the center region, simply simulating a guide thimble or an instrumentation tube in a typical PWR (Pressurized Water Reactor) fuel assembly. The calculation conditions of fuel regions except for the sub-division of fuel region in azimuthal direction are the same as those given in Section 4.2. The fuel pellet is sub-divided into 32 sector regions for 2π with equal azimuthal angle interval.

The effective cross-sections generated by the present method are compared with the reference solutions obtained from the continuous energy Monte-Carlo code MVP. The total number of neutron sampling for MVP calculations is set to 1 billion histories, in which the 1σ statistical error of the effective cross-sections is about 0.3%.

The results are shown in Figure 21 and Figure 22. As shown in the figures, the present method can accurately predict the azimuthally dependent effective cross-sections. The spatially -dependent biases against one-region fuel also agree well with the Monte-Carlo results, and therefore the consistency between one-region fuel and azimuthal multi-region treatments is confirmed for the present method.

Figure 21. Azimuthally dependent effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) for 3×3 multi-cell with large water region (corner fuel).

Figure 22. Azimuthally dependent effective cross-sections and their differences from the continuous energy Monte-Carlo calculation (MVP) for 3×3 multi-cell with large water region (vertical fuel).

In the conventional methods, only the direct heterogeneous ultra-fine-group slowing-down calculation and the sub-group method with two-dimensional exact geometry can generate the azimuthally dependent effective cross-sections. Therefore, the present method is an efficient alternative to treat azimuthally dependent resonance self-shielding effect, while keeping sufficient accuracy.

5. Conclusion

The unified resonance self-shielding method is developed for general application of lattice physics calculations. The present method can accurately treat the radially and azimuthally dependent resonance self-shielding effect without the direct heterogeneous ultra-fine-group calculation.

A calculation scheme of the present method is based on the concept of multi-stage cross-section collapsing strategies, and is composed of two-step calculation, i.e. ‘coarse geometry + fine energy’ (first step) and ‘fine geometry + coarse energy’ (second step) calculations. The first-step calculation corresponds to the hybrid model of the equivalence theory and the ultra-fine-group calculation, and the second step corresponds to the sub-group method. The two-step reaction-rate preservation scheme is also established to reduce energy discretization error.

From the various verification results, radially and azimuthally dependent effective cross-sections generated by the new method show good agreement with the continuous energy Monte-Carlo results for pin-cell and multi-cell geometries including non-uniform fuel compositions and temperature distributions within a pellet.

The present method can accurately generate effective cross-sections with short computation time in lattice physics calculations, and has a potential for the general application toward the next generation core analysis codes that require high fidelity and sophisticated modeling of reactors.

Nomenclatures
Ak	=	relative atomic weight of nuclide k to the neutron
αk	=	maximum energy loss ratio for nuclide k defined as ${\alpha }_k\equiv {\left(\frac{A_k-1}{A_k+1}\right)}^2$
αn	=	nth rational coefficient considering the heterogeneous and lattice effects
βn	=	weight of nth rational equation normalized as ${\sum }_{n=1}^N{\beta }_n=1$
E	=	neutron energy
ϵn	=	nth rational coefficient for mitigation of self-shielding effect in non-fuel region
f	=	fuel region
g	=	index for energy group (multi-group resolution)
k	=	index for nuclide
lf	=	mean chord length of fuel lump
li	=	mean chord length of ring region i defined as $l_i\equiv 4V_i/\phantom{4V_i{S}_f}{S}_f$ (Sf: pellet surface)
μ	=	self-shielding correction factor for non-fuel region
N	=	the number of rational equation terms
n	=	index for rational equation term
nf	=	non-fuel region
Nfk	=	number density of nuclide k in fuel region
Nnfk	=	number density of nuclide k in non-fuel region
φf	=	neutron flux in fuel region
φi	=	neutron flux in ring region i
ϕnf	=	non-asymptotic component of non-fuel flux
Pj → i	=	collision probability from regions j to i
r	=	index for resonance nuclide
ρi	=	relative outer radius of ring i (normalized by pellet radius R)
ρi − 1	=	relative inner radius of ring i (normalized by pellet radius R)
sg	=	index for energy group (sub-group resolution)
Si	=	neutron source for region i
Σnfp	=	macroscopic potential scattering cross-section in non-fuel region
Σfsd	=	slowing-down cross-section in fuel region
Σft	=	macroscopic total cross-section in fuel region
Σit	=	macroscopic total cross-section in ring region i ∈ f
σks	=	microscopic elastic scattering cross-section for nuclide k
σrx	=	microscopic cross-section for resonance nuclide r and reaction type x
θ	=	ratio of asymptotic spectrum in non-fuel region (0 ⩽ θ ⩽ 1)
ufg	=	index for energy group (ultra-fine-group resolution)
Vi	=	volume for region i
x	=	index for reaction type

Disclosure statement

No potential conflict of interest was reported by the authors.