Abstract
The effect of BWR fuel assembly 3D model on void reactivity coefficient (VRC) estimation is investigated. VRC values were calculated for different BWR assembly models applying deterministic T-NEWT and Monte Carlo KENO-VI functional modules of SCALE 6.1 code package. The difference between deterministic T-NEWT and Monte Carlo KENO-VI simulations is negligible (0.18 pcm/%). The influence of the assumed more detailed coolant density profile was estimated as well. VRC increases with the application of a larger number of coolant density values across fuel assembly height. It was shown that the coolant density profile described by 6 values per height could be considered sufficient from prospect of VRC estimation, as a more detailed density profile has impact below 1% on total assembly void effects. VRC values were decomposed to values for individual nodes and isotopes, since decomposition provides useful insights to describe the overall behaviour of VRC in detail.
1. Introduction
Many nuclear reactors are characterized by a rather strong coupling between coolant density and neutronic behaviour. For example, in Light Water Reactors (LWRs) cores, where the water coolant also serves as a moderator, a local decrease in water density causes a decrease in moderation and hence a decrease in local power density. As the coolant passes up through the core, it absorbs heat from the fuel elements and eventually initiates subcooled/bulk boiling, and, thus, the coolant void fraction increases. The steam void formation in the core plays a very important role in Boiling Water Reactors (BWRs) [1]. A similar coupling between coolant properties and core neutronic behaviour can arise in sodium fast reactors. Liquid metal coolant causes an appreciable moderation or softening of neutron energy spectrum in fast reactor core (number of neutrons produced per neutron absorbed in fuel increase with average neutron energy), and this leads to a local decrease in reactivity, whereas coolant density decrease leads to an increase in local reactivity (opposite direction as for LWRs). Such thermal-neutronic coupling that arises in LWRs and sodium-cooled fast reactors (SFRs) is much weaker in gas-cooled reactors, because coolant phase change does not occur, and the coolant does not provide appreciable moderation in the core due to its very low density. Thereby, BWR core physics is more complex than physics of other reactors, because the coolant enters the core in a single phase and rapidly develops different two-phase flow regimes with very strong axial dependence. It is worth mentioning that reactor core physics has become even more complicated in recent years, because modern/innovative fuel assemblies have new sophisticated features, such as newly designed inner bypass regions and part-length rods, introducing larger uncertainties in relation to axial void distributions and their impact on integral reactor physical parameters [2–4].
A change of the system reactivity caused by variations in the fuel cross sections due to the moderator density change is defined by the void reactivity coefficient. Since the continuous operational control of BWR core is typically accomplished by the adjustment of the coolant flow rate through the core, which affects the neutron moderation efficiency, understanding the nature of coefficient values is of prime importance for BWR reactor safety assessment. However, the question is whether this effect can be estimated accurately by employing 2D model and using the average coolant density.
The analysis presented in this paper is a continuation of the work on study of void reactivity coefficient in innovative BWR fuel assemblies [5]. The investigation of the void reactivity variation using deterministic T-NEWT functional module of SCALE package in the entire 0–100% void fraction range was performed in the previous study. This voiding range represents all possible moderation conditions appearing during normal operation and accident transients in BWR core. The previous study results indicated that positive values appear for the void fraction above 75% in case of fuel assembly loaded with fresh mixed oxide (MOX) fuel. Positive values may lead to decreased level of core controllability during transients such as a loss of coolant accident. VRC variation analysis was performed only for two-dimensional (2D) geometry in the mentioned study case.
There is a demand for a detailed investigation of voiding effects in fuel assemblies’ axial direction. Thereby, the investigation of the VRC was initiated by employing functional Monte Carlo criticality transport module KENO-VI of SCALE package. In addition, values of void reactivity coefficient were decomposed to separate values in this study in order to interpret VRC differences between various fuel assembly simulation models. The decomposed void reactivity coefficient values provide useful quantitative insights for describing the overall behaviour of VRC in detail. Such analysis allows explaining the overall trend of the void reactivity coefficient by fractional absorption variations.
2. Computational Methodology
2.1. SCALE Code and Nuclear Data
The neutron transport calculations were performed using SCALE 6.1 code package [6]. SCALE provides a framework with 89 computational modules including three deterministic and three Monte Carlo radiation transport solvers that are selected based on the desired solution strategy. In current study, 2D neutron transport calculations in fuel assembly were performed employing T-NEWT module and KENO-VI module (Monte Carlo criticality transport module), which was used for simulations in 3D geometry. Used neutron cross section libraries were generated by utilizing the ENDF/B-VII library based on 238 neutron energy groups.
Validation of the applied SCALE code combined with the technique used for preparation of BWR assembly numerical model was performed in the beginning. For this purpose, T-NEWT module in 2D simulations was involved. NEA BWR-MOX benchmark data [7, 8] were employed for comparison of the results. 30 benchmark solutions were submitted to benchmark by 20 participants using different neutron transport code systems. Previous version of Monte Carlo criticality module KENO-V of SCALE 4.4 code package was used in the benchmark calculations as well.
Two characteristics were analysed: infinitive neutron multiplication and power peaking factor. The first is an important neutronic characteristic that describes fuel assembly criticality. The second indicates how the local power relates to the average power. The comparison of the T-NEWT module results with the ones of the benchmark participants is depicted in Figure 1. It shows that T-NEWT results reflect the average values for all voiding cases, which were derived from presented data of benchmark participants.
(a)
(b)
In all voiding cases, power peaking factors calculated by T-NEWT and KENO-V are almost the same, while variation of eigenvalue evaluated by KENO-V is higher than T-NEWT values. Based on data comparison, it can be stated that the employed fuel assembly model results reflect the data of benchmark participants.
2.2. Fuel Assembly
A modern 10 × 10 BWR design, corresponding to ATRIUM-10 type, with large internal water channel and fuel rods loaded with MOX fuel was chosen for this study. It is the same BWR assembly design as that used in the NEA BWR-MOX benchmark [7, 8].
Figure 2 shows the geometrical layout of the problem and the fuel composition. Fuel assembly consists of 91 fuel rods with 6 different (U-Pu)UO2 and one UO2-Gd2O3 (U-Gd) mixtures. UO2 matrix enrichment by U-235 in the MOX rods was assumed to be 0.20%. The plutonium isotopic vector is given in Table 1. The mean fissile Pu concentration was 3.93 w/o, averaged over both MOX and U-Gd rods.
The radial Pu distribution within the assembly reflects moderation efficiency within the assembly. Higher moderation efficiency exists in peripheral fuel rods since thermal neutron transport from the water between fuel assemblies. Even higher moderation efficiency exists in corner rods since thermal neutron transport from the water comes through two sides of the lattice cell. Thus, there are lower concentrations of plutonium in these rods in order to reduce radial power peaking. Similarly, there are smaller plutonium concentrations around the water channel, though larger than those at the periphery of the assembly. The uranium-gadolinium rods are placed around the water channel to take advantage of the fact that extramoderation in the water channel increases the effectiveness of Gd.
The U-Gd rods have a density of 9.867 g/cm3; meanwhile, MOX fuel rods have a density of 9.921 g/cm3. The fuel temperature was set at 627°C. The density of the zircaloy cladding and water channel was 6.55 g/cm3. The temperature of structure material was set to 327°C.
2.3. Void Reactivity Coefficient Calculation Methods
As it has been mentioned, BWR assemblies are very heterogeneous, since the two-phase flow is presented in the reactor core during the operation. There are the bypass regions in BWR assemblies, which are separated from the two-phase flow: central water channel, gaps between surrounding assemblies. Single-phase flow (pure liquid water) for bypass regions, having density equal to 0.740 g/cm3, was assumed. Meanwhile, the void fraction equals 40% in average for the coolant surrounding the fuel rods in case of a typical BWR core. However, the fuel channel outlet voidage could vary, in general, from 0% to 100% including severe accident conditions. Based on this fact, the VRC values were estimated based on two reactor core boundary conditions: in case the fuel channels outlet voidage equals 40% (water density 0.452 g/cm3) and 100% (water density 0.036 g/cm3). The difference of the reactivity between these two boundary conditions reflects the void reactivity effect in case of the transient from normal operation conditions to a loss of coolant during accident.
During analysis of VRC in 2D geometry, two outlet voidage states 40% and 100% were analysed. The fuel channel outlet regions are marked in light blue colour (Figure 2). Meanwhile, the void fraction for coolant outside channel box and inside the water channel (dark blue colour in Figure 2) was set to 0% in both voiding cases. However, it should be noted that this approach is valid only for 2D simulations. The fuel channel outlet voidage could vary, in general, from 0% to 70%, or even more, across fuel assembly height. The void fraction increases with the coolant flow stream from the bottom to the top of the fuel assembly.
The uneven distribution of two-phase flow was described by dividing the fuel assembly in 3, 6, and 12 nodes per height and assigning appropriate values of coolant density for each separate node. The moderation amount (average water density per height) was preserved the same in all simulation models. The distributions of void fraction and coolant density values in case the average void fraction equals 40% are depicted in Figure 3. Meanwhile, in case average void fraction equals 100%, even distributions across all axial nodes of the fuel assembly were assumed.
The void reactivity coefficient was defined as a ratio of the reactivity change to void fraction change:The decomposition analysis was involved in the study to explain the void reactivity feedback. For this purpose, the total VRC value was decomposed into separate components, which represent contribution of the individual isotope or axial node to the total void reactivity feedback. The VRC through fractional absorption for individual isotopes or nodes () can be expressed as follows:
3. Results
3.1. Void Reactivity Coefficient for 2D and 3D Model
This section is intended to determine whether VRC estimation depends on the chosen neutron transport solver. The comparison of neutronic characteristics calculated using deterministic T-NEWT (2D model) and Monte Carlo KENO-VI (3D model) functional modules was performed (Table 2). The moderation conditions across the fuel assembly height are the same for both simulations, since the same average coolant density values were adopted.
In addition, further specific assumptions were made in case of 3D model having the purpose to preserve similarity between 2D and 3D models. First, mirror boundary conditions on all surfaces were assumed. Second, partial length fuel rods (grey colour in Figure 2) were simulated as full length rods.
The ratio of fast and thermal neutron fluxes is determined by the efficiency of neutron moderation in the fuel assembly. Accordingly, values of neutron multiplication factors are directly proportional to these ratio values. However, multiplication factor calculated using T-NEWT model in case void fraction equals 40% is greater than value calculated using KENO-VI model, even when value of neutron fluxes ratio calculated using KENO-VI model is greater. This indicates simulation differences resulting from the application of different neutron transport models. However, differences of all neutronic characteristics in both void cases are relatively small (below 1%).
3.2. The Dependence of Void Reactivity Coefficient on Fine Water Density Distribution
The single-phase flow in the form of pure liquid exists in the bottom of BWR fuel assembly. Nevertheless, heat produced by nuclear fission reactions causes the coolant to boil, producing the steam. The void fraction increases with the coolant flow stream from the bottom to the top of fuel assembly. This uneven distribution of water density affects the fundamental neutronic characteristics of BWR fuel assembly, including VRC. This section is intended to determine the effect of fine water density distribution on VRC estimation. The comparison of neutronic characteristics calculated by using Monte Carlo simulation in 3D in case of 4 different nodalisation schemes (1, 3, 6, and 12 nodes across fuel assembly height) was carried out (Table 3). Number of nodes means the number of coolant density values used to describe two-phase flow across the fuel assembly height in case void fraction equals 40% (see Figure 4). Meanwhile, only the steam flow was assumed across the fuel assembly height in case void fraction equals 100%.
(a)
(b)
It is seen that large differences occur for fast/thermal neutron ratio values when void fraction is 40%. The ratio decreases (the efficiency of neutron moderation increases) with the increase in number of coolant density values. The difference of ratio estimated between cases where 1 and 6 coolant density values across fuel assembly were applied is 16.2%. In comparison, the difference between cases where 1 and 12 coolant density values were applied is 16.8%. Meanwhile, the differences of ratio estimated in case void fraction equals 100% are not so significant—below 0.01%. Thus, it is evident that the estimation of neutron moderation efficiency in case of normal conditions significantly depends on more realistic coolant density profile.
The system criticality () depends on neutron moderation conditions. Thus, the differences in estimation of fast/thermal neutron ratio, respectively, affect the estimation of system criticality in case void fraction equals 40%. The neutron multiplication factor increases with the application of a more detailed coolant density profile since neutron moderation improves as well. However, the differences increase only from 1.56%, between cases with 1 and 6 applied coolant density values, to 1.62%, between cases with 1 and 12 applied coolant density values. In comparison, the differences in case void fraction equals 100% are below 0.01%.
Since the system criticality of fuel assembly increases in case void fraction equals 40% and does not significantly vary in case void fraction equals 100%, VRC value increases with applying larger number of coolant density values to describe two-phase flow across the fuel assembly height. The VRC increased by 34% between cases where 1 and 6 coolant density values were applied. However, in comparison, the VRC increases further only to 35%, applying even more detailed coolant density profile with 12 values (increasing the difference by 1% only).
The comparison of fractional neutron absorption in individual nodes between analysed nodalisation schemes is depicted in Figure 4. Figure 4(a) displays the comparison between nodalisation scheme with 3 nodes and nodalisation schemes with 6 and 12 nodes, which are condensed to 3 nodes. Figure 4(b) displays the comparison between nodalisation scheme with 6 nodes and nodalisation scheme with 12 nodes which are condensed to 6 nodes. It is seen that there is a significant difference of neutron absorption level across the fuel assembly between 3-node model and nodalisation schemes with 6 and 12 nodes in case void fraction equals 40%. Obviously, the differences depend on more realistic coolant density profile description and nodalisation. Meanwhile, the differences in case void fraction equals 100% depend only on nodalisation, as all nodes have the same coolant density. The differences between 6 and 12 node models are similar in both voidage cases. Thus the indicated differences between 6 and 12 nodes models, which represent more realistic coolant density profile conditions, are comparable to the differences related to nodalisation option.
Neutron moderation conditions are better in the bottom part of fuel assembly in case of 40% fraction. Thus, approximately 75% of total neutron absorption events occur in the first 1/3 part of fuel assembly across flow. Due to equal moderation conditions across the fuel assembly height in case void fraction is 100%, fractional neutron absorption is distributed almost equally. It is seen in Figure 4 that, by the implication of more realistic coolant density profile description, the fractional neutron absorption increases in 1/3 bottom part of fuel assembly and decreases in the 2/3 part.
The decomposition analysis was performed in order to investigate the origin that makes VRC differences for used nodalised schemes. For this purpose, the total VRC value was decomposed into separate components, which represent a contribution of the individual axial node to the total void reactivity feedback. Figures 5 and 6 show the comparison of decomposed VRC to values for individual nodes. Positive contribution appears in the 1/3 bottom part of fuel assembly and negative contribution is estimated in the 2/3 top part. The positive feedback is estimated, since fractional neutron absorption decreases (positive contribution to the system reactivity) in the bottom part in case of 100% fraction. Accordingly, negative contribution to the system reactivity is estimated, since fractional neutron absorption increases in the top of fuel assembly.
The differences of decomposed VRC values between 3-node model and nodalisation schemes with 6 and 12 nodes condensed to 3-node models are depicted in Figure 5. It is seen that VRC increase in the 1/3 bottom part of fuel assembly and decrease in the 2/3 top part. It is related to corresponding phenomenon (related to more realistic coolant density profile description) of fractional absorption in case void fraction equals 40%.
The differences of decomposed VRC values between 6- and 12-node models condensed to 6-node model are depicted in Figure 6. These differences are insignificant in comparison with previous figures. It was previously mentioned that both more realistic coolant density profile description and nodalisation scheme option contribute to these differences. Thus, the coolant density profile description with 6 coolant density values can be considered sufficient for right estimation of voiding effects across the whole fuel assembly.
3.3. Isotopic Decomposition Analysis
This section is intended to determine the effect of fine water density distribution on VRC estimation by comparing system’s fractional absorption data. The same models as in previous section were used (Table 3). Fractional neutron absorption events for separate isotopes in case void fraction equals 40% (for different nodalisation schemes) and 100% (applying 1 node model) are depicted in Figure 7. Only isotopes that absorb more than 1% of total neutrons are included for the analysis. Main isotopes, which contribute to neutron absorption, are Pu-239 (~40%, in case of normal conditions), U-238 (~20%), Pu-240 (~14%), and Pu-241 (~12%). It is seen that fractional absorption of H-1, Gd-157, Pu-239, and Pu-241 isotopes is lower in case of 100% void fraction in comparison to 40%. Meanwhile, fractional absorption of U-238, Pu-240, and Pu-242 isotopes increases. The opposite phenomena are seen in case of application of a more realistic coolant density profile: fractional absorption of H-1, Gd-157, Pu-239, and Pu-241 isotopes is increasing as fractional absorption of U-238, Pu-240, and Pu-242 isotopes is decreasing, since the neutron moderation efficiency is increasing.
Neutron absorption, by H-1, burnable poison Gd-157 and fissile material, Pu-239, and Pu-241, decreases/increases due to moderation efficiency decrease/increase, since these isotopes have larger absorption cross sections in thermal neutron energy range. Meanwhile, absorption by Zr and fertile materials (U-238, Pu-240, and Pu-242) increases/decreases due to moderation efficiency decrease/increase, since these isotopes have larger absorption cross sections in fast neutron energy range.
The results of VRC isotopic decomposition analysis (Figure 8) show that negative VRC value is mainly affected by U-238 isotope. Other isotopes, which influence negative VRC, are Zr, U-235, Pu-240, and Pu-242. Accordingly, H-1, Gd-157, Pu-239, and Pu-241 isotopes have a positive influence on the total void reactivity effect. These positive/negative VRC values are related to the same relationship between absorption cross sections of isotopes and moderation efficiency decrease/increase in case void fraction equals 100% as mentioned above.
The results of isotopic decomposition analysis confirm the previous mentioned remarks. The coolant density profile described by 6 values can be considered sufficient for the right estimation of voiding effects initiated by individual isotopes, since the profile described by 12 values gives smaller VRC differences than 1 pcm/%.
4. Conclusions and Discussions
In this study, the effect of fuel assembly 3D model on void reactivity coefficient estimation was investigated. VRC values were calculated using 2D and 3D models with different nodalisation schemes, which represented different coolant density profiles. Coefficient values were decomposed for individual nodes and isotopes to get insights of the overall behaviour of VRC in detail across the BWR fuel assembly.
It was estimated that total VRC value difference between the deterministic simulation in 2D and the Monte Carlo simulation in 3D is 0.28%. However, the neutron moderation efficiency in case void fraction equals 40% and, accordingly, the neutron multiplication factor are increasing by implementing a more realistic profile description in 3D model simulations. Since the system criticality increases with the application of a larger number of coolant density values in case of 40% void fraction and does not change significantly in case of 100% void fraction, VRC increases accordingly. The VRC difference is 34% between cases where 1 and 6 coolant density values were applied to describe flow across the fuel assembly height. However, the VRC difference is only 35% between cases where 1 and 12 coolant density values were applied (the difference increases only by 1%).
The decomposition of VRC values to separate values for each individual node revealed that the bottom part of the fuel assembly contributes to the positive void effect when negative void effect appears on the top. The isotopic decomposition of VRC confirmed that a more realistic coolant density profile described with 6 values per height could be considered sufficient for accurate VRC estimation. The impact of a more realistic coolant density profile description in case void fraction equals 40% is on about the same level as the impact of nodalisation in case void fraction equals 100%.
Competing Interests
The authors declare that they have no competing interests.