ON THE FLY PREDICTION OF TH-DEPENDENT SPATIAL MACROSCOPIC CROSS-SECTIONS USING FFT

Monte Carlo (MC) codes can accurately model neutron transport in nuclear reactors. However, the efficient inclusion of thermal-hydraulic (TH) feedback within the MC calculation sequence is still an open problem, particularly when burnup’s time-evolution must be included in the analysis. For this reason, deterministic codes, leveraging the use of macroscopic cross-sections generated with higher order methods from 2D lattice calculations, are still widely used to perform reduced-order multiphysics analyses. However, traditional cross-sections generation procedures typically decompose the large core problem into multiple assembly-level problems; thus not having the ability to capture inter-nodal effects. Moreover, the pre-generation procedure requires additional pre-computational time to perturb/branch the problem for various operational conditions (e.g. fuel temperature), which, again, is decoupled from the core. In this paper, we propose a new method leveraging the use of Fourier transfer functions to predict the cross-sections distribution due to a variation in TH conditions. The method was tested against a 3D BWR unit-cell problem with realistic density profile and axial fuel heterogeneity. The method was able to compute the mono-energetic cross-sections distribution with maximum error lower than 2%. Insights on the influence of the statistics used to generate the cross-sections on the accuracy of the results is also provided.


INTRODUCTION
Continuous energy Monte Carlo codes allow to accurately calculate the steady-state fission power distribution for a pre-determined temperature distribution. However, in the frame of fuel cycle analysis, the effect of the thermal hydraulic (TH) feedback is usually ignored due to the high computational cost associated to the coupling of MC and TH solvers. Traditionally, to overcome the problem of the computational cost, nodal methods, coupled with low-fidelity thermal hydraulic solvers, are utilized. The latter computational sequence relies on the use of multigroup macroscopic cross-sections generated with Monte Carlo codes or a high order deterministic method from 2D lattice calculations. The main issue of this sequence is that it does not capture three dimensional and inter-nodal effects. Moreover, it requires a pre-generation phase for the computation of the cross-sections at different operating conditions [1].
In this paper, we present a new method to predict the macroscopic cross-sections as a function of thermal hydraulic conditions on-the-fly [2]. The method relies on the use of transfer functions to map the effect of the temperature and density spatial profiles on the macroscopic cross-sections. The transfer functions are computed at a reference TH condition through a Fast Fourier Transform (FFT) algorithm and no additional branching calculations are needed. The method has three main advantages. First, it is extremely simple to implement, and it does not require any modifications at the source code level. Second, it allows to capture three dimensional effects neglected by the traditional 2D lattice calculations. Finally, the use of the FFT, in conjunction with a frequency filter, may be used to efficiently denoise the outputs of the MC codes (examined in future works). The method is here tested against an axially reflected BWR unit-cell to predict the change in nodal cross-sections and power shape due to the TH feedbacks. The analysis was limited to the feedback between neutronics and void distribution (i.e. coolant density). Agreement within 2% was obtained between the predicted one group macroscopic cross-sections and the reference Monte Carlo results. This agreement was also preserved near the boundary, where the spatial flux gradients are maximum due to the proximity to the axial reflectors.

SPATIAL TRANSFER FUNCTIONS FOR CROSS-SECTIONS
This section, the theory is explained. A general definition of linear system response is provided in Section 2.1, while the notion of transfer function for cross-sections is introduced in Section 2.2. Finally, the simple Fourier approximation used to obtain the numerical results is derived in Section 2.3. The derivations are performed without considering the uncertainty associated to the MCgenerated cross-sections. Qualitative observations on the stochastic effects and overall sensitivity will be provided in Section 3.

Definition of Linear System Response
Let the relationship between the cross-section σ(x x x) and the density profile ρ(x x x) be expressed as: where O denotes a non-linear operator. If there exists a linear mapping L that locally approximates O, the equation can be rewritten as: The linear operator is related to the functional derivative of O if the latter is differentiable. After expanding ρ(x x x) in a suitable basis {ψ n } n=0,..,N and exploiting the linearity of the operator L, Eq. 2 becomes: In Eq. 3, ρ (n) is the n th expansion coefficient, corresponding to the n th basis function ψ n (x x x). Therefore, the L [ψ n (x x x)] denotes the response of the system to a density with spatial distribution ψ n (x x x), i.e., the cross-section profile σ(x x x) resulting from a density profile ρ (n) ψ n (x x x). For the sake of simple notation, this response will be denoted as σ n (x x x) in the remainder of the paper and the ≈ symbol will be substituted by the equality sign. The use of this convention leads to: Therefore, the cumulative response is a superposition of the partial responses σ n (x x x) weighted by the expansion coefficients ρ (n) . The latter can be easily computed through projection if the basis set is orthogonal or bi-orthogonal. For a non-orthogonal basis set composed by linearly independent vectors/functions, a system of equation must be inverted.

The Fourier Transfer Functions
The application of the Discrete Fourier Transform (DFT), here denoted by F, to Eq. 4 leads to the following identity in the frequency domain: Then, for F{ρ} = 0, we can define the transfer function: The choice of using the DFT to define the transfer function is based on two considerations. First, the computation of the expansion coefficients is extremely efficient using Fast Fourier Transform (FFT) algorithms [4]. In this paper, the fft function implemented in MATLAB was utilized. Secondly, the introduction of the Fourier transform enables the use of common signal processing techniques for denoising the output of MC codes when used for cross-sections generation. Considerations on using a hard threshold for denoising purposes will be provided in Section 3.

A Simple Approximation
A remarkably simple expression for the transfer function is obtained when ρ (0) >> ρ (n) and the expansion coefficient are normalized such that ρ (0) = 1. In this case, Eq. 6 can be simplified as: In Eq. 7, F{σ 0 } is the system response to a density profile of spatial shape ψ 0 (x x x), projected in the frequency domain. As noticeable, the transfer function is formally identical to the one com-monly defined for Linear Time Independent (LTI) systems [4]. Therefore, the cross-section can be calculated by using the Fourier inverse-transform F −1 : This simplified formulation can be leveraged when either H is a slow varying function of the density through the weights ρ (n) , or the density changes are modest. In this case, H can be calculated only once and utilized to predict the cross-sections' distribution for several density profiles ρ (x x x). It must also be observed that one transfer function must be generated for each cross-section type. A thorough discussion on the validity of Eqs. 7-8 can be found in Ref. [2].

NUMERICAL RESULTS
To test the performance of the proposed approach, we use an axially-reflected 3D BWR unitcell, i.e., pin radially surrounded by water. The fuel lattice pitch is 1.87 cm and radial reflective boundary conditions are imposed. The length of the assembly is 366 cm and 20 cm of water reflector surrounds the active core at the top and bottom. The other geometrical specification are reported in Table 1. The 3% fresh fuel was burned for 3.0 MWd/kgU. The axial burnup distribution is reported in Fig. 1. The unit-cell was discretized in 36 equally-high axial layers for the purpose of tallying the cross-sections' nodal values. The TH feedback was limited to the void fraction and only one-group constants are considered. The problem was chosen to show the main strengths and the limits of the method with the minimum computational effort. The reference results were generated using 300,000 particle per 400 criticality cycles of which 100 inactive, leading to maximum standard deviation below 0.1% for all the quantities of interest. The sequence to compute the transfer function was implemented in a MATLAB 2018b linkage code, while Serpent 2.1.30 was used to compute the cross-sections spatial distribution and the associated uncertainty [3].
The five density profiles reported in Fig. 2, {ρ j } j=1,...,5 , were extracted from a fixed point MC-TH calculation at different iterations. A density variation from ρ 1 to ρ 5 leads to an average power 1 . It is noticeable that good agreement is found in all cases, although the accuracy decreases with increasing distance between initial guess and actual profile. The Mean Absolute Percentage Error (MAPE) for Σ a passes from 0.76% to 1.39% going from ρ 2 (x x x) to ρ 5 (x x x), while it passes from 0.386% to 0.658% for νΣ f . The MAPE of Σ a for all the simulated cases is reported in Table 1.It can be noticed that the MAPE, although always lower than 2%, increases with the distance between ρ 1 (x x x) and ρ j (x x x), as expected based on the theoretical description in Section 2. The macroscopic fission neutrons production follows the same pattern, with MAPE bound by 1%.
The accuracy of the cross-sections' prediction is caused by the low sensitivity of the transfer functions to the density variation. This can be seen from Fig. 5, where the Bode diagram of the transfer functions generated at ρ 1 (x x x) and ρ 5 (x x x) for Σ a are reported. It is noticeable, that all the expansion coefficients are not null. There are two causes for this phenomenon. First, the cross-sections present steep local gradients at the domain's boundaries due to the presence of reflectors. Second, the DFT is applied to a non-periodic signal in a finite space-domain, giving origin to spectral leakage, i.e., higher modes are exited because the periodic extension of the signal outside the observation window is not continuous. Both these phenomena excite high order coefficients. This may constitute a problem if the FFT with hard threshold is utilized to denoise the MC output.  An interesting consideration can be also drawn from Fig. 6.a. The latter represents the transfer function H Σa 1 generated with poor statistics: 1,000 particles/cycle for 100 cycles (50 active). It is evident that the noise undermines the quality of the TFs' estimation. This stochastic fluctuation propagates through the calculations leading to highly inaccurate results as shown in Fig. 6.b, where the prediction performed with high statistics and low statistics are compared.

CONCLUSIONS
In this paper, a new method for the prediction of the macroscopic cross-sections was presented. The methodology relies on transfer functions to map the effect of thermal-hydraulic fields on the spatial cross-sections' distributions. The method was here applied to analyze a BWR 3D unit-cell with realistic axial density profile and non-homogeneous axial fuel composition. It was shown that the method was able to predict the one-group macroscopic cross-sections within a few percent even near the boundaries. According to our analysis, two caveats must be considered when utilizing this technique. First, the accuracy of the solution decreases for increasing distance between the initial guess and the actual distribution of the quantity of interest. Second, the statistics can undermine the quality of the results, by exciting higher components of the Fourier spectrum. In future work, higher dimensional problems will be considered, i.e., full three dimensional, in order to verify the capability of the method to deal with full-core calculations. Additionally, the use of other basis set alternative to Fourier, e.g., wavelets, will be examined for efficiently denoising the MC-generated cross-sections.