© A. Laureau et al., published by EDP Sciences,

1 Introduction

Low void effect sodium fast reactors provide an improved behavior during accidents thanks to a negative feedback coefficient due to sodium expansion. This effect is provided by a large sodium plenum that increases the neutron leakage if the sodium density decreases. An optimization of the core geometry [1] leads to the CFV (low sodium void worth) design applied to the ASTRID (Advanced Sodium Technological Reactor for Industrial Demonstration [2]) reactor. This concept includes axial heterogeneities that increase the neutron gradient and consequently the neutron leakage due to a sodium density reduction.

The heterogeneities built into the design of this kind of reactor may induce a spatial decoupling between fissile zones during different accident scenarios. For this reason, spatial neutron kinetic models must be developed to verify that the flux redistribution remains limited in such situations. To this end, the spatial kinetic Transient Fission Matrix (called TFM) approach previously developed and presented in [3,4] has been adapted to model the effect of local medium perturbations [5]. The TFM approach is based on a conversion to discretized Green functions of the Monte Carlo response of the system in order to perform kinetic calculations without new reference calculation during the transient and thus with a reduced computation time. This approach is used to model neutron kinetics during transient coupled calculations; in order to take into account the evolution of the system during the transient, TFM is associated to specific interpolation models. The innovative interpolation model produces on the fly estimations of the matrix variations during the transient without requiring new Monte Carlo calculations. Such an interpolation model based on the correlated sampling technique and suitable for heterogeneous cores with a fast neutron spectrum has been developed. The resulting neutronic approach called ‘perturbative TFM’ is applied in this paper on the ASTRID concept.

After a brief introduction on the fission matrices and of this perturbative TFM approach in Section 2, this paper presents a first study and a validation of the developed approach on low void effect sodium fast reactors. Direct Monte Carlo and S _N calculations have been performed and the results used for the comparison and validation. Two application cases are introduced in Section 3: a simple one highlighting the neutron leakage phenomena associated to the low void effect and a second one corresponding to a representative assembly of the internal core of the ASTRID sodium cooled reactor. In Section 4, we describe the generation of the perturbed fission matrices together with the influence of different calculation parameters. The validation of the interpolation model that reconstructs the core-perturbation as a sum of individual local-perturbations is presented in Section 5 with comparisons to direct Monte Carlo and deterministic S _N calculations. Finally, in Section 6, we show how this approach can be applied to the calculation of point kinetics local feedback parameters.

2 Perturbative Transient Fission Matrix approach

The Transient Fission Matrix approach is based on a time dependent version of the fission matrices. The objective is to precalculate the time dependent transport characteristics of the neutrons in order to perform transient coupled calculations with a reduced calculation time. The raw information contained in the fission matrices is the probability that a fission neutron created in a volume j produces a new fission neutron in a volume i where i and j are volumes of the discretized geometry. This information is summarized in fission matrices (line i column j), that are derived according to the emission spectrum χ (prompt or delayed), the neutron multiplicity ν (prompt or delayed), e.g. . The average prompt propagation time from j to i is stored in the time matrix (line i column j).

Recent developments [5] extend the application field of this TFM approach to fast reactors with a heterogeneous core such as the sodium fast reactor. For sodium fast reactors with a low void effect based on a sodium plenum, if the density of the sodium decreases, the probability to produce a new fission decreases due to the increased neutron leakage. Because the core geometry of such reactors is very heterogeneous, not only must the neutron departure and arrival volumes be taken into account, but also the intermediate volumes crossed by the neutron. Finally, the raw information used is the probability that a fission neutron created in j produces a new fission neutron in i considering a perturbation in the crossed volumes k. The theory and the methodology of this perturbative approach using the correlated sampling technique is detailed in [5], the present paper focuses on the application of this approach to the study of the low void effect sodium fast reactor ASTRID [1] and to its validation with direct Monte Carlo and deterministic S _N calculations. The effect of the perturbation in the crossed volume k on the fission matrix is named for the density effect and for the Doppler effect. Considering a perturbation of −1% for the density and +300 K for the Doppler effect during the correlated sampling process, the matrix interpolation is the following:(1)

3 Presentation of the application cases

We describe below two application cases, A and B, that we used to illustrate these developments on one dimensional geometries. Case A is a simplified version (geometry and composition) of case B. With case A, we focus on the effect of a local perturbation in the sodium placed between the fuel and the B₄C. Case B corresponds to a beginning of life representative assembly of a sodium cooled reactor with a negative sodium void effect [1]. The objective here is to test the ability of the approach developed here to accurately represent this type of system with its complex feedback effects.

In order to compare our results with those obtained by S _N calculations, the axial composition of the assembly is considered homogeneous and the Doppler perturbation is performed using a uniform temperature variation of all the material components (both fuel and coolant) although these two assumptions are not required for the TFM perturbed matrix calculation.

3.1 Case geometry

As mentioned above, case A, detailed in Figure 1, is a simplified core configuration with three distinct areas: fissile, sodium (with no structure), and B₄C. The axial boundary condition is a neutron leakage and the radial condition is a boundary reflection.

The second case called B is also detailed in Figure 1. It is a representative average assembly of the ASTRID [6] reactor even if the material compositions have been simplified for modeling purposes. This configuration is very heterogeneous with fertile areas in the core, and a sodium plenum that is optimized to ensure a negative sodium void effect. The axial boundary condition is a neutron leakage. The radial boundary condition is also a neutron leakage, and the dimension of the system is adjusted to ensure criticality. With a hexagonal representation of the core-assembly, the distance between flats is equal to 125 cm corresponding to a multiplication factor of k _eff = 0.99980 ± 0.00002.

3.2 Composition

The material temperatures and isotopic reference compositions of cases A and B are given in Table 1. Note that the sodium plenum (or “Na Plenum”) for case B is called “Na” for case A since its composition is different, it comprises only sodium, but with the same content as in the plenum of case B. These compositions are considered radially homogeneous so that, for example, the B₄C zone contains sodium and steel. In the rest of this paper, unless otherwise specified, a variation of −1% of the sodium density is applied on all the areas for the density perturbation. Likewise, a step modification of +300 K is applied to estimate the contribution of the temperature perturbation (Doppler effect).

3.3 Notation

Different configurations of cases A and B with different sodium densities and temperatures will be used in this paper. These configurations will be referred to as follows:

• A ^init and B ^init refer to the initial configuration, without perturbation.

• and refer to a global temperature variation of +300 K on the geometry as a whole.

• and refer to a global sodium density variation of −1% on the geometry as a whole.

• refers to a −300 K fuel temperature local variation at the top of the fuel area (25–50 cm) for case A.

• refers to a −2% sodium density local variation in the Na area (50–75 cm) for case A.

3.4 Calculation parameters

3.4.1 Monte Carlo – Serpent

The TFM matrices are calculated with a modified version of Serpent 2.1.21 [7]. The modifications concern the calculation of the fission matrices, including the distinction between prompt and delayed neutrons, the fission to fission time matrix, and the correlated sampling technique to generate the locally perturbed matrices.

The nuclear database used is JEFF 3.1 [8]. The number of simulated neutrons for each calculation is one billion. The system boundary conditions for both geometries are leakage on axial boundaries. Concerning the radial boundary a reflection is used for case A and a leakage for case B.

3.4.2 Deterministic code ERANOS

The ERANOS [9] calculations are based on the traditional two level lattice/core scheme and the JEFF 3.1 nuclear database. First, self shielded cross sections are computed by the ECCO code cell, using the fundamental mode assumption for each kind of material of the 1D core description. For fissile material, a buckling search algorithm is used to obtain the critical flux for the cross section collapsing to a 33 energy group mesh. For the subcritical materials such as fertile or structural parts of the 1D sub-assembly, the process is based on source calculations using the spectrum coming from previous fissile calculations.

For core calculations, the subassembly is modeled by a 1D core and a discrete ordinate S _N method with n = 16 and with 33 energy groups corresponding the self shielded cross sections prepared as described above. The boundary conditions are a flux leakage on the axial boundaries for both cases, while a critical buckling models the radial leakage for case B.

3.5 Neutron flux and multiplication factor at steady state

Reference calculations on cases A ^init and B ^init provide a multiplication factor value of k _eff = 1.02951 ± 0.00003 and k _eff = 0.99980 ± 0.00002 respectively with the SERPENT code. Considering this cross section preparation scheme, the values obtained with ERANOS can be considered in quite good agreement with respectively k _eff = 1.02603 and k _eff = 1.00088.

The source neutron and neutron flux distributions in case A ^init are given in Figure 2. Observe that all the neutrons are created on the left of the geometry, in the fuel area. Some of the neutrons are reflected from the sodium area (between 50 and 75 cm) and some of the neutrons are stopped in the B₄C on the right. The results obtained with the ERANOS and SERPENT codes are in very good agreement.

The source neutron and neutron flux distribution in case B ^init are given in Figure 3. The fissile areas are directly visible on the source neutron shape (bottom) with the maxima between 120/145 cm and 165/200 cm. The fertile areas correspond to the non-zero production areas with a low fission and as a consequence a low production rate. A good agreement is obtained, except at the top of the fissile area (between 190 and 200 cm). Serpent estimates a skin effect due to a spectrum shift of the neutrons back scattered at lower energy from the sodium plenum to the fissile area. The flux distribution (bottom) shows that the neutron flux is slightly different in the gas plenum (between 200 and 210 cm); the fission rate difference near the sodium plenum is thus explained not only by the spectrum effect but also by a difference in the way the neutron transport is modeled in this medium.

4 Application and interpretation of the perturbative TFM approach

In this section we present the results obtained with the correlated sampling technique applied to the TFM approach. Different aspects are studied in parallel:

• Different matrices of the TFM approach are presented in Section 4.1 to discuss the physics behind the matrices, the influence of the emission spectra, and the structure of the time matrix.

• In Section 4.2 we present the effect of a global perturbation of the core on the fission matrices.

• In Section 4.3 we focus on the effect of a local perturbation in a volume k.

• In Section 4.4 we compare the estimation of the matrices using the correlated sampling technique and direct independent Monte Carlo calculations.

• Finally the effect of the neutron source weighting by the previous generations in the correlated sampling process during the Monte Carlo calculation is studied in Section 4.5.

4.1 Raw Transient Fission Matrices

The fission matrices , and of case A ^init and case B ^init are computed with a discretization of respectively 60 and 120 bins. They are shown in Figure 4.

The neutron propagation is directly visible on these matrices. Each emission position corresponds to a column, and for this column the position of the neutrons produced by fission corresponds to the different lines. Concerning the matrices (left), we can see on case A that all the fissions come from and occur in the fissile zone with an index between 0 and 30. The probability of generating a new source neutron is reduced near to the small values of index j and i due to the leakage. On the contrary, around bin 30, the sodium is a neutron reflector so that the source neutron production is less impacted by the end of the fuel area. For case B the fissile and the fertile zones are in the middle of the geometry, we can see the impact of the fertile zones: the fission probability is reduced in the lines i (target cell) that belong to the fertile areas. The structure of the matrix depends more on the target area than on the position j of the neutron emission. We can see that the statistics is better for the columns corresponding to a fissile area due to the larger number of neutrons emitted there.

We can see on case A (top) that the delayed production shape (middle) is somewhat sharper than that of the prompt production (left). Indeed the delayed production multiplicity is higher for fast neutrons due to the threshold fission reaction of ²³⁸U that produces a lot of delayed neutrons. Then if the creation occurs close to the production (j ≃ i), the neutron is “younger” and its energy higher. For case B, the fertile-fissile distinction is directly visible due to the increase of ν_d in the fertile area. On the last plot on the right, the propagation time increases if the neutron production by a fission in i is far from the neutron creation in j. The neutron fission to fission time is smaller in the fertile than in the fissile area: since the ²³⁸U fission is a threshold reaction, only a neutron at the beginning of its life (before losing energy by scattering) can induce a fission. Note that this time matrix is not the neutron lifetime but the neutron fission to fission time.

4.2 Global perturbation of the fission matrices

Figure 5 presents the matrices for cases A (top) and B (bottom), together with their variation due to a uniform modification of the density (−1% – middle) and of the temperature (+300 K – right). They correspond to the estimations of cases A ^init, , and B ^init, , .

Concerning the density effect (middle), the physical effect is similar between cases A and B. The neutron production is reduced on the diagonal of the matrix (for a target volume i close to the origin volume j), and the production is relocated far from the neutron emission position. This effect is due to a larger mean free path resulting from the decreased sodium density. Note the reduced effect when the neutron targets a fertile area (horizontal strips) on case B, without impacting the fissile area near to the fertile. Near the boundary between the fuel area and the plenum sodium, the strong negative feedback is explained by more neutron leakage to the B₄C.

Concerning the Doppler effect (right), the impact on the neutron propagation is not a relocalisation such as with the density change, but a negative global feedback due to a modification of the fission-absorption ratio and spectrum. The effect is larger close to the sodium area on case A. A strong local effect may be noticed for case B in the fuel close to the fertile zone. The Doppler effect in the fertile area impacts the neutron spectrum and results in a significant skin effect when the neutrons return to the fuel areas.

4.3 Local perturbation of the fission matrices

The fission matrix variations are calculated for each local perturbation position. To illustrate this, Figure 6 shows the local contributions of a density perturbation for the configurations and , and Figure 7 presents the same contributions due to the Doppler effect for the configurations and .

The main effect for the sodium density feedback (Fig. 6) is a relocalisation of the source neutron production to the other side of the perturbation position in k: reduction of the production if (i and j) >k or (i and j)<k. We can see on case B that the shape is identical as case A, but with a local reduction of the fission rate in the fertile area without impacting the future interactions of neutrons crossing the whole area. The neutron spectrum also becomes harder, increasing the neutron production (positive feedback); this will be quantified in Section 6. For the figures in the right panel with a perturbation k in the sodium area, the leakage from the fuel to the B₄C increases, reducing the source neutron production in the fuel close to the sodium (production in j or target in i close to bin 30 for case A).

The impact of the Doppler effect (Fig. 7) is more global and mainly depends on the i target position. Case B illustrates the strong effect of the fertile area on the fuel (bottom-left): a temperature variation in the fertile does not result in a variation of the fission rate in the fertile but in the fissile zone. This effect is due to the important spectrum variation at the interface between the two areas.

4.4 Evaluation of reference perturbed fission matrices with direct Monte Carlo calculations

Figure 8 illustrates the benefit obtained with the correlated sampling technique for the estimation of the variation of the fission matrices due to a density perturbation on the example of the configuration .

On the top panel are represented the variations of the matrix estimated using only one Monte Carlo calculation with the correlated sampling technique (neutron weight modification), and on the right the estimation of the same matrix using two distinct Monte Carlo calculations (one for the reference case, and one with the perturbed sodium density, their subtraction providing the result represented). Each Monte Carlo calculation uses the same number of neutrons (one billion). The global evolution is the same but the statistical error is much larger with the two independent calculations. Indeed the sodium variation is very small (1%) and, for the estimation with distinct calculations, each element of the matrix corresponds to the difference between two quantities with independent statistical errors and a very slight difference. The same effect can be observed with local variations (bottom) and is exacerbated by the smaller perturbation amplitude. To conclude, Monte Carlo perturbed calculations are a very powerful tool for the precise estimation of the variation of fission matrices in the event of a global or local modification.

4.5 Propagation of the perturbed weight through generations

The correlated sampling process requires to propagate the neutron weight modification to the neutron produced per fission of the next generation. In this way, the effect of the perturbation on the neutron source perturbation is correctly taken into account. As detailed in [5], this perturbed weight propagation is not required for fission matrix generation if the mesh is fine enough since the perturbation of the neutron source is directly taken into account through the eigen vector.

This assumption is checked in this section: the neutron source is perturbed in order to see the impact of the matrices. Figures 9 and 10 represent respectively global and local perturbed matrices with zero memorized generations (left), eight memorized generations (middle) and the difference (right) in percent of the maximum absolute value of the left panel-matrix. The study has been performed for a number of memorized generations from 2 to 16 to verify that the observations are unchanged.

We can see that, for a global variation (Fig. 9), the estimated matrices are almost identical. The difference between the estimations is around 2% for the density and 10% for the Doppler. These differences are due to the statistical errors since no specific pattern is visible. Note that the estimation with zero memorized generations (left) is smoother than with eight generations (middle) because of the increasing dispersion of the neutron weight in the latter; this dispersion decreases the contribution of the low weight neutrons and thus the statistics.

In Figure 10, the same behavior is visible on the local contributions of volume k = 15 (in the fissile area). The error is amplified due to the reduced statistics.

We can conclude here that, as expected, thanks to the fission matrix properties, the weighting of the perturbed source over generations is not required since the matrices are identical with and without perturbed weight propagation. The source redistribution does not affect the local neutron propagation since the origin volume j is assumed to be small enough. If a coarse mesh is used, the neutron weigh propagation through generations would be required.

5 Validation of the TFM interpolation

During the Monte Carlo calculation, different kinds of matrices are required. Considering the sodium density variation for example, the following matrices are estimated (“x” standing for “p” or “d”):

• the reference matrix;

• the globally perturbed ;

• the locally perturbed matrices .

The matrices of interest for the interpolation using equation (2) of Section 5.2 are the locally perturbed matrices using and , and the equivalent one on case B. The estimation of is optional but useful to check the amplitude of the cross effects as described below in Section 5.1. Note that we estimate each set (sodium density and Doppler effect) of fission matrices for each emission spectrum (χ_p and χ_d) and neutron multiplicity (ν_p and ν_d), together with the matrix of the prompt fission to fission propagation times.

The second validation process (Sect. 5.2) concerns variations of density/temperature with a different shape and amplitude from the matrices estimated for the interpolation model. This validation is based on cases and .

5.1 Globally versus sum of locally perturbed matrices

represents the variation of for a global modification of the sodium density. The first step of this validation is to compare this matrix to the sum of the local variations . Indeed, the sum of the local contributions does not take into account the crossed contributions while does. This comparison provides an estimation of the bias inherent to the approach that consists in reproducing the global perturbation with the sum of the local contributions.

5.1.1 Case A

Figure 11 shows the difference obtained between the sum of the local contributions and the global perturbation for both configurations and . The discrepancy is calculated as the difference between the two matrices, normalized by the maximum value of the perturbation (matrix on the right) to avoid a division by a very small value.

We can see that the discrepancy associated to the sum of the local variations of the sodium density is very small (<1%). Concerning the Doppler effect the difference is a bit larger, around −3% corresponding to an over estimation of the effect since the crossed contributions between the local perturbations are not taken into account in the sum of the local individual contributions.

In our application case, an error of less than 5% is small enough for coupled calculations. However, if a better accuracy is required, a multi-level scheme can be defined to interpolate complex distributions in the reactor. For example the perturbation can be reconstructed as a global perturbation plus a sum of local perturbations. In this way, the global perturbation correctly takes into account the cross effects, and the sum of the local perturbations fits the interpolation on the real distribution. This approach could be generalized on different simplified meshes, based on a coarse mesh subscripted on k and a finer one subscripted on K, for example using the geometry discretization between fuel and fertile areas as a coarse mesh. The sum over K contributions is used to obtain the area-averaged perturbation, and then the residual fine perturbation over k is superimposed to fit the exact perturbation shape.

The comparison of the global versus sum of local perturbation results can also be done on the eigen vector and eigen values of these matrices. Figure 12 represents the variation of the normalized neutron source in the core due to the global sodium density (red) and Doppler (blue) variation, the reference calculation (computed with the global perturbation) is drawn with a solid line, and the interpolation (sum of local perturbations) with a dashed line. Additionally, Table 2 presents the reactivity variation due to the density and the Doppler effects. Note that the statistical error on the eigen values (neutron multiplication factor) is of the same order of magnitude than that of the Monte Carlo calculation used to generate the matrices. The matrices corresponding to the global and to the local perturbations are estimated using the same Monte Carlo calculation, and then the difference directly corresponds to the bias that consists in considering that a sum of local contributions is equivalent to a global contribution that takes into account the crossed effects.

The agreement obtained between local and global perturbations on the flux redistribution and reactivity prediction is very good. A discrepancy of a few percent is obtained on the reactivity prediction of the interpolation, a result of the same order of magnitude as the discrepancy-matrix of Figure 11.

Slightly different results on the reactivity variation are obtained with ERANOS, in particular on the Doppler effect with a difference of 14.4%. This difference comes from the calculation scheme since the fundamental mode of the lattice calculation is not really matched in this configuration with significant leakage on the left and reflections on the right. The cross section self-shielding is not representative of the local leakage and of the spectrum evolution in this radially infinite reactor as will be illustrated in Section 6 with the local feedback estimation.

5.1.2 Case B

The same analysis has been performed on case B. Figure 13 represents the discrepancy obtained between the sum of the local contributions and the global perturbation for both configurations and . Figure 14 presents the redistribution of the normalized neutron source in the core, and Table 3 provides the associated reactivity variation together with ERANOS calculations.

Similar results are obtained concerning the matrix variations with a very good agreement on the density effect (Fig. 13 – top) and a small discrepancy on the Doppler effect (bottom). The neutron source redistribution (Fig. 14) is perfectly reproduced with TFM, and a good agreement is obtained with ERANOS. The only slight difference concerns the Doppler effect associated redistribution in the lower fertile area, an area with a reduced importance since less power is released there.

A large difference is obtained concerning the reactivity variation between TFM and ERANOS on the density effect in Table 3. The origin of this difference will be explained in Section 6 where we compare the TFM and ERANOS feedback effect distribution in the reactor.

5.2 Interpolation of a local perturbation

The validation detailed in the previous section is based on the same perturbation amplitude as the fission matrix estimations (global variation of −1% sodium density and +300 K for the temperature). This section deals with a second validation focusing on a different amplitude perturbation (−2% and −300 K) on case A in order to verify the validity of the linear interpolation chosen for the density effect and the logarithmic interpolation for the Doppler effect. Furthermore, local perturbations (on a fraction of the geometry) are also used to validate the capability of this approach to model the full-scale core with a sum of local contributions. The validation is thus based on a different perturbation shape and amplitude than the one used to compute the matrices. It corresponds to the configurations and that require a supplementary calculation compared to and .

The perturbed interpolated matrices are calculated using equation (2) applied to a density reduction in the sodium plenum (bins 30–45) and a reduced temperature to the right of the fuel (bins 15–30):(2)

Three items are used to validate this calculation:

• k _eff estimation: it corresponds to the eigen value and it is compared to a distinct Monte Carlo estimation performed directly with the modified composition – Section 5.2.1.

• νΣ_fψ distribution: it corresponds to the eigen vector and it is compared to a distinct Monte Carlo calculation (the same one as for the k _eff validation) – Section 5.2.2.

• Interpolated matrix: it is compared to a distinct Monte Carlo calculation (the same one as for the k _eff validation) and also to another calculation using the local modification of the medium to generate the perturbed matrices with the correlated sampling technique – Section 5.2.3.

5.2.1 k _eff validation

For each perturbation, the reactivity variation is calculated with the eigen value from the interpolation and with a direct Monte Carlo estimation, using . The results obtained are summarized in Table 4.

The discrepancy is less than a few pcm of reactivity variation (within the statistical error) between the TFM interpolation and classic Monte Carlo calculations. With such a small discrepancy, transient coupled calculations where, e.g. the density and Doppler effects have to be updated at each time step, can be carried out validly.

Note that a larger difference is observed with ERANOS concerning the Doppler effect. The origin of this difference will be explained in Section 6.2, the point kinetics feedback estimation providing an information on the feedback distribution.

5.2.2 νΣ_fψ validation

The second step is a comparison of the νΣ_fψ redistribution in the core due to the perturbation. Figure 15 displays the density and Doppler perturbations. The dashed line represents the results predicted by the TFM interpolation (matrix eigen vector variation), and the continuous line the difference between νΣ_fψ scores from classic Serpent perturbed and reference calculations with the associated uncertainties.

The TFM interpolation is able to predict very efficiently the neutron redistribution in the core due to the local perturbation, for both localized sodium density and fuel temperature variations. The statistical error of the Serpent calculation is larger than that of the TFM interpolation, since this distribution is calculated using two distinct estimations with their own statistical errors. Finally, a very good agreement is also obtained with ERANOS (brown and turquoise lines superimposed with the corresponding results of the other codes).

5.2.3 Validation of the interpolated matrices

Even if the important issue for the neutron kinetics studies is the power redistribution validated in 5.2.2, an interesting aspect is the validation of the matrices themselves since they contain the information on the redistribution of local neutron propagation. We focus here on the sum of the local contributionsthat represent the matrix global variation due to the perturbation considered. Excluding statistical errors, the comparison of these matrices with a reference will verify that the biases of the model are negligible: the limit of the linear and logarithmic dependency of the density and of the temperature, and secondly the cross effects between the local contributions. A reference of these global variation matrices can be estimated in two ways (detailed below in the next two paragraphs):

• Direct calculation: the global variation matrices are estimated as the difference of two “classic” fission matrices estimated with two distinct (initial and perturbed) calculations. The advantage is to have a direct comparison with a result obtained without any “correlated sampling” technique. The drawback is to combine two independent calculations with their independent statistical errors.

• Using the correlated sampling technique: the global variation matrices are estimated with the correlated sampling technique where the neutron perturbed weights are calculated using the final shape and amplitude (without summing the individual contributions and without the lin-log interpolation). Then, even if we are comparing two independent Monte Carlo calculations, each estimation is performed using the correlated sampling technique so that a low statistical error is expected.

Note that all these calculations use the same amount of simulated neutrons (one billion) to have a comparable computation cost.

Comparison using usual Monte Carlo calculation as reference

For this comparison, the reference matrix (middle panel in Fig. 16) is estimated with two independent calculations where the matrices are evaluated without the correlated sampling technique. Since the calculations are independent, the statistical error is directly visible.

The behavior of the reference matrix is globally reproduced by the interpolated matrix. The main difference comes from the statistical noise, which is larger on the Doppler effect (bottom). The good agreement observed on the global trend is very important since the reference calculation does not use the “correlated sampling” technique implemented in this work.

Comparison using a correlated sampling calculation as reference

The second comparison uses the “correlated sampling” technique to obtain a reference calculation (Fig. 17 – middle) where the correct shape and amplitude of the perturbation are provided. No interpolation is required for these reference matrices, thus they cannot be estimated on the fly but they take into account cross effects between the local perturbations. The interpolated matrices (left) are calculated using a sum of local contributions, without cross effects between the different k volume contributions.

The results show a very small discrepancy between the interpolation and the reference calculation. Thanks to the correlated sampling technique, the statistical errors on the direct Monte Carlo calculations are removed. The residual statistical error comes from the difference of statistical errors on the perturbed estimations. In order to avoid this effect, all the matrices should be estimated using only one Monte Carlo calculation, with multiple matrix variation estimations (−1% global, −2% local density, +300 K global and −300 K local temperature perturbations) and not two distinct calculations (one for the global perturbation and one for the local ones). In this way the statistical errors due to the differences in the neutron histories would be suppressed, but our objective is already achieved here, i.e. to validate the capability of the interpolation to predict an arbitrary perturbation shape.

6 Calculation of point kinetics parameters

As presented in [5], the perturbative TFM approach may also be used to calculate the local feedback effects for point kinetics applications. These local feedbacks correspond to the eigen values of the local contributions such as . Figures 18 and 19 show the results obtained for cases A and B together with the results of ERANOS. As previously mentioned, the density feedback coefficient is simply calculated with a linear dependency: , and the Doppler dependence is assumed to be logarithmic: .

6.1 Case A

A good agreement is obtained on the density effect. The blue curve (ERANOS) is slightly below the red one (TFM), and since the feedback is negative, the global reactivity variation is a bit larger. This result is consistent with the reactivity variation of a global perturbation shown in Table 2 (+3.6% larger with ERANOS). Note that the feedback is constant in the sodium area since there is no radial leakage and the sodium is the only component (no steel). Indeed wherever the sodium density variation occurs, the only effect on the reactivity is the leakage and this is related only to the integrated amount of sodium in the whole plenum.

A difference is observed on the Doppler effect. As already mentioned, the difference is a bit larger for ERANOS in Table 2: +14.4%. We can see here that the distribution shape is not the same between TFM and ERANOS, the effect is exacerbated in the right region of the fuel (25–50 cm). During the group condensation process of the fuel cross sections, only one area of fuel is considered at the grid level, without radial leakages. However, at the “core” level, the neutron spectrum close to the interfaces is different from the average spectrum in the fuel. It appears that a possible recommendation for the calculation scheme is to go beyond the fundamental mode assumption and to consider cross exchanges between the grid and the core resolution. In this way a core-representative distribution of the leakage and source current between areas could be shared with the lattice calculation.

6.2 Case B

Figure 19 shows the feedback coefficient distribution of case B. A good global agreement is obtained.

6.2.1 Sodium density feedback analysis

Concerning the density feedback, a very good agreement is obtained, except in the gas plenum (between 200 and 210 cm). This difference is the origin of the large discrepancy of +38% on the global reactivity variation presented in Table 3. Indeed the global variation is the sum of a large positive (in the fuel) and an equivalent large negative (in the sodium plenum) component, resulting in a value close to zero. For this reason, the small difference in the gas plenum leads to a large difference on the global feedback.

This discrepancy has several origins. The scheme employed for the cross section condensation is a possible explanation. The multi-group cross section generation for the heterogeneous 1D core calculation is performed using the fundamental mode for each area, the self-shielding assuming that the neutron spectrum is spatially converged. Since the thickness of the fissile area is only 10 cm, a lattice calculation is not fully representative of the neutron behavior in this region.

Another discrepancy related to the cross section generation concerns the representativeness of the spectrum used for the homogenization. The neutron source term used for non-fissile areas such as the gas plenum area in the ERANOS calculations is the fundamental mode neutron spectrum of the fuel. This modeling implies a bias because the neutrons coming back from the sodium plenum to the fuel have a more thermalized spectrum than the ones outgoing from the fissile area. Moreover, the steel and sodium densities are constant between the fuel and the gas plenum, while the steel density is multiplied by 2.8 and the sodium density by 0.3 between the sodium plenum and the gas plenum. Then the neutron source is not representative for the neutrons coming from the sodium plenum to the gas plenum.

6.2.2 Doppler feedback analysis

The Doppler feedback effect is slightly overestimated in the upper fissile area and underestimated in the middle fertile area with ERANOS. We can see a different behavior between each fissile-fertile area interface with a local amplification of the feedback estimated with TFM. The calculation scheme and the small fissile thickness (between 25 and 35 cm) explain these discrepancies due to a neutron spectrum modification inside the areas. A behavior similar to the sodium density feedback is observed in the gas plenum. The same demonstration as previously discussed on the density feedback applies to the Doppler effect, the neutron spectrum being different between the neutron coming from the fissile area and the one coming from the sodium plenum.

6.2.3 Spectral analysis

In order to confirm the influence of the neutron spectrum variation in small thickness areas, different maps of the core are presented in Figure 20 using a double discretization with 25 000 bins for the energy (abscissa) and 300 for the axial position (ordinate). The ‘A’ map represents the neutron flux in the reactor. This map illustrates the neutron energetic and spatial propagation in the reactor. The ‘B’ and ‘C’ maps represent the neutron spectra corresponding respectively to the positive and negative axial neutron velocities. Compared to the ‘A’ map they show a normalisation of the flux for each axial position to correct the scalar flux axial variation. The ‘D’ map presents the difference between the ‘B’ and ‘C’ maps, illustrating the spectrum variation due to the neutron propagation direction.

We can see on the ‘A’ map that the neutrons are produced at high energy in the fissile area and are slowed down in the sodium plenum before being absorbed in the B₄C or reflected to the fissile matter. Note the effect of the ²³Na elastic scattering resonance at 2.8 keV that depreciates the flux in the entire reactor.

The ‘B’ and ‘C’ maps show that, at high energy, in the fissile areas and away from their boundaries (more than 5 cm), the neutron spectrum is constant, which is in agreement with the assumption of the fundamental mode: a single spectrum is representative of the fissile areas. However, at energies less than 10 keV in the fissile areas and in the fertile areas, the assumption of the fundamental mode is not satisfied due to the spectrum variation. For each zone, a single neutron spectrum is not fully representative. A neutron condensation with a spatial discretization can be a good way to improve the deterministic calculation schemes.

We can see on the ‘C’ map that the spectrum also depends on the neutron direction. The neutron spectrum is harder in the “fissile to sodium plenum” direction (positive difference at high energy and negative at low energy for positions around 200 cm), meaning that the leakage current spectrum is harder. Indeed the neutrons that go back to the fissile area after a few interactions on the sodium have lost some of their energy. The same effect occurs at each transition zone with, for the high energies, a positive component at the top of the fissile areas and a negative one at the bottom. This effect shows that, due to the spectrum anisotropy, adding an angular dependency to the condensed cross sections may also improve the deterministic calculation schemes.

Finally, we can notice a very strongly anisotropic spectrum at the interface between the gas plenum and the sodium plenum (see the zoom on the ‘D’ map). These spots with a large spectrum variation between the up and the down components correspond to the ⁵⁶Fe resonances (see cross sections displayed under the zoom box). This effect can explain the larger difference between TFM and ERANOS on the local feedback estimation in the gas plenum area shown in Figure 19 since ERANOS does not consider the angular dependency of the condensed cross sections.

7 Conclusions

The correlated sampling technique associated to the TFM neutron kinetic approach discussed in this paper proves to be a powerful tool that can provide perturbed fission matrices.

This approach provides a new way to perform spatial kinetic calculation. It requires a unique Monte Carlo calculation prior to transient calculations, once per reactor configuration, and finally the perturbed fission matrices are used to provide an on-the-fly prediction of the reactivity and of the source neutron redistribution.

This approach is quantitatively validated on direct Monte Carlo and reference correlated sampling calculations: the flux redistribution and the reactivity variation associated to an arbitrary perturbation shape are predicted with a discrepancy limited to a few percent. The calculation parameters are discussed. It appears that the neutron source perturbation usually used with the correlated sampling technique, viz. the neutron weight propagation, is not required. The capability of the fission matrices to reconstruct the neutron source shape and their perturbations avoid the statistical convergence reduction. Feedback effects from perturbations on the sodium density and the temperature, assuming respectively a linear and a logarithmic dependency, is quantitatively validated on direct Monte Carlo calculations. Neglecting the cross contributions by summing the local individual contributions prove to be of little consequence.

This tool is also compared to deterministic S _N calculations with the ERANOS code on spatial neutronics modeling and point kinetics parameter calculations. A good global agreement is obtained and a need for further developments is identified to improve deterministic calculation schemes concerning the modeling of heterogeneities in small thickness zones.

Different perspectives can already be identified such as the coupling of this innovative neutronics approach with the thermal hydraulics. For more complex systems, a key to improve the accuracy of the interpolation could consist in associating a multi-scale scheme to the correlated sampling TFM approach: different coarse and fine meshes can be used in parallel to model the crossed volume contribution in order to provide an information on the cross effects between the individual local contributions. Finally another perspective concerns full core scale calculations based on a fine mesh and quasistatic resolution with regular TFM based estimations of the flux shape and of the local feedback coefficients.

Acknowledgments

The authors wish to thank the IN2P3 department of the CNRS (National Center for Scientific Research) for its support during the initial development of the TFM approach. We are also very thankful to our colleague Elisabeth Huffer for her help with the rereading.

References

All Tables

All Figures