Abstract
A high-leakage core has been known to be a challenging problem not only for a two-step homogenization approach but also for a direct heterogeneous approach. In this paper the DIMPLE S06 core, which is a small high-leakage core, has been analyzed by a direct heterogeneous modeling approach and by a two-step homogenization modeling approach, using contemporary code systems developed for reactor core analysis. The focus of this work is a comprehensive comparative analysis of the conventional approaches and codes with a small core design, DIMPLE S06 critical experiment. The calculation procedure for the two approaches is explicitly presented in this paper. Comprehensive comparative analysis is performed by neutronics parameters: multiplication factor and assembly power distribution. Comparison of two-group homogenized cross sections from each lattice physics codes shows that the generated transport cross section has significant difference according to the transport approximation to treat anisotropic scattering effect. The necessity of the ADF to correct the discontinuity at the assembly interfaces is clearly presented by the flux distributions and the result of two-step approach. Finally, the two approaches show consistent results for all codes, while the comparison with the reference generated by MCNP shows significant error except for another Monte Carlo code, SERPENT2.
1. Introduction
Over the years, many different computer code systems have been developed for reactor physics analysis [1], the main goal of which is to obtain the information on the behavior of neutrons in nuclear reactor cores. The analysis of the neutron behaviors is affected by many different factors such as deterministic approach versus probabilistic approach, transport theory versus diffusion theory, and direct one-step approach versus two-step approach [2]. Contemporary computer codes can be categorized as either Monte Carlo (MC) transport codes or deterministic transport codes for the direct one-step approach and nodal diffusion codes for the two-step approach when used together with the transport solver for equivalent cross section data generation. Among these methods, the MC method is usually considered to be most accurate as long as the probability density functions reflect the physical reality, but this method takes a very long computation time or computational resources especially when large amount of tallies is required such as in the case of three-dimensional (3D) whole core PWR depletion analysis. The MC method is usually used to analyze research reactors or benchmark problems which are not in typical shapes or compositions rather than the daily reactor core analysis of commercial power reactors which requires computational efficiency. Instead, a two-step calculation approach has been used for production work of real world commercial power reactor analyses. A two-step approach consists of deterministic two-dimensional (2D) lattice transport calculations and deterministic 3D nodal diffusion calculations. In this paper, comprehensive comparative analyses of two approaches are conducted with well-known conventional reactor analysis codes. The first approach is a direct one-step modeling with Monte Carlo codes such as MCNP6 [3] and SERPENT2 [4] and with deterministic codes such as CASMO-4E [5] and NEWT in SCALE code system [6]. The second approach is a two-step modeling by the nodal calculation with PARCS3.0 [7] together with 2D lattice calculations for the generation of homogenized two-group cross section data and the assembly discontinuity factors (ADFs). This paper will also explicitly introduce the procedure for the generation of those equivalent homogenized diffusion parameters (two-group constants and ADFs) by using color-set models. DIMPLE S06 core was selected because it is a small high-leakage (~10%) core. A high-leakage core has been known as a challenging problem not only for a two-step homogenization approach but also for a direct heterogeneous approach. Since there haven’t been many verification results for small core designs, this paper presents calculation results with comprehensive comparisons between code-to-experiment, code-to-code, and the two approaches.
2. DIMPLE S06 Critical Experiment
The experiments (LEU-COMP-THERM-055) were performed at AEA Technology’s Winfrith site during the late 1980s and early 1990s [8]. The experimental program is comprised of critical experiments with low-enriched uranium dioxide fuel rods containing 3.0 w/o 235U with light water moderation and reflection. These experiments were designed to simulate peripheral regions of a PWR and modeled the equivalent of twelve 16 × 16 PWR fuel bundles, arranged in a cruciform array. There are two configurations of the cruciform array. The first configuration (S06A) incorporates a cruciform array of 3072 3 w/o-enriched rods with a water reflector as shown in Figure 1; the second configuration S06B adds a stainless steel baffle of 2.67 cm thick around the cruciform array with a surrounding water reflector as shown in Figure 2. Ordinary water is used for the reflector and moderator. The pin pitch for both configurations is 1.2507 cm. The fuel pin is composed of fuel, aluminum wrapper, and cladding whose diameters are 1.0130, 1.0398, and 1.0937 cm, respectively. DIMPLE S06 critical experiment analyzed in this paper consists of two cores, with and without a stainless steel baffle. Core S06A has uniform enrichment and has no baffle, absorbers, and water holes; therefore it can be presented in terms of the effect of radial leakage. Core S06B has identical configuration with addition of a baffle; therefore it can provide the information on the effect of baffle. These two configurations provide consistent and comprehensive comparative analyses (e.g., flux distributions, two-group homogenized cross sections, and assembly power distribution).
3. Reactor Core Analysis Codes
3.1. MCNP6
MCNP6 code has a general 3D geometry modeling capability with particles such as neutron, photon, electron, or coupled neutron/photon/electron transports, including the ability to calculate the eigenvalues for critical system [3]. For the comparative neutronic study of DIMPLE S06 critical assembly, the MCNP6 code is used for the reference solutions.
3.2. SERPENT2
SERPENT2 code is a 3D continuous energy MC calculation code developed since 2004 and optimized for reactor physic applications [4]. The development of a new code version, SERPENT2, started in 2010 and is available to licensed users for beta-testing purpose. SERPENT2 code can be used to generate the homogenized few-group constants, and it can also be used for a direct modeling by the flexible description of 2D or 3D geometry.
3.3. CASMO-4E
CASMO-4E code is a 2D multigroup lattice physics code developed by Studsvik Scandpower [5]. The code gives 2D transport solution with MOC solver dealing with full heterogeneous geometries with multigroup cross sections and a quadrature of polar angles, azimuthal angles, and ray tracing. It can handle cylindrical fuel rods surrounded with gap or clad of various compositions in a square pitch array. The CASMO-4E code is used to generate the homogenized two-group cross sections for nodal calculations. The direct modeling of DIMPLE S06 with CASMO-4E is also possible due to its multiassembly capability using MXN card.
3.4. NEWT
NEWT (New ESC-Based Weighting Transport code) is a multigroup discrete-ordinates radiation transport computer code with flexible meshing capabilities that allow two-dimensional neutron transport calculations using complex geometric models [6]. The advantage of the Extended Step Characteristic (ESC) method is that it can easily take higher-order scattering because it has the integrodifferential form of the transport equation. One of the main disadvantages is the computational burden for the direct angular treatment depending on the quadrature order () [9].
3.5. PARCS3.0
PARCS is a 3D reactor core simulator which solves the steady-state and time-dependent, multigroup neutron diffusion and SP3 transport equations in square and hexagonal geometries [7]. PARCS models for DIMPLE S06A and S06B were constructed to perform 2D diffusion calculation with homogenized square nodes generated by SERPENT2, CASMO-4E, and NEWT codes. The nodal expansion method for multigroup is used for the two-group diffusion calculation.
4. Calculational Model
4.1. Direct Modeling in 3D
For the two analysis approaches of DIMPLE S06 benchmark calculation, 2D model is used without using the axial buckling which is calculated through the buckling search with 2D calculation result. This study focuses on the analysis of the two approaches for DIMPLE S06 benchmark calculation. Therefore, the use of the 2D model is enough for comparative study of two approaches. One challenge is that there is no reference data for these infinite 2D calculations. Therefore, a full 3D calculation is preceded to present the validity of the code input data for 2D calculation. The detailed description of the S06A and S06B is same as described in Section 2, and there are some additional geometric descriptions as shown in Figures 3 and 4. First, two cores are placed at the center of a cylinder with the radius of 129.55 cm filled with water. Second, the cruciform array of fuel rods is in between two sets of upper and lower aluminum lattice plates. The lower lattice plates were secured to the aluminum fuel support plates on the fuel support beam base. Lastly, in S06A and S06B configurations, critical moderator heights, measured from the bottom of the fuel stack, are 47.50 cm and 52.49 cm, respectively [8]. The calculation is performed with MCNP6 and SERPENT2 using ENDF/B-VII.0 continuous energy cross section library. The total histories used in the calculation are 1,200 million, in which the number of active cycles is 1,000 and the number of skip cycles is 200 in both codes. The results show good agreements with the experimental values as shown in Table 1. Therefore, it demonstrates that the code input models of the benchmark problems are accurate for MCNP6 code. Then, in this study, the yellow dashed lined 2D cross section as represented in Figures 3 and 4 is taken from the input model of the full 3D calculation. The MCNP6 2D results calculated through this process are used as the reference data.
4.2. Direct Modeling in 2D
The analysis of 2D version of the DIMPLE S06 benchmark is performed with computer codes such as MCNP6, SERPENT2, CASMO-4E, and NEWT. The ENDF/B-VII.0 continuous energy cross section library is used in all calculations using MCNP6 and SERPENT2. The ENDF/B-VI 70 group library is used in all CASMO-4E simulations, and ENDF/B-VII.0 238 group library is used for NEWT calculations. The number of histories used in 3D calculation is also used for the MC calculations of 2D models. The configurations for the 2D direct modeling are shown in Figures 1 and 2. The discretization scheme for 2D DIMPLE S06 NEWT model uses 240 × 240 mesh grid which discretizes the pin cell into 25 sections and 6 quadrature order considering both the accuracy and the execution time as shown in Figure 5. The CASMO-4E calculation with direct modeling was performed using default uniform quadrature (32 azimuthal angles, 3 polar angles, and a ray spacing of 0.1 cm). It should be noted that CASMO-4E cannot use the exact thickness of a baffle for S06B model because it needs to be defined in the multiples of pin pitch. The thickness of a baffle was modified to 2.5014 cm (two pin pitches), compared to the actual baffle thickness of 2.67 cm. A stainless steel baffle has a significant transport effect in S06B, but no special treatment was performed to correct transport cross section of Fe.
4.3. Two-Step Modeling in 2D
Two-step calculation is carried out through the 2D lattice calculation and then the 2D nodal calculation with the data from the 2D lattice calculation. The 2D nodal calculation requires the homogenized cross section sets on each node which should be generated by 2D lattice calculations. These calculations are performed with SERPENT2, CASMO-4E, and NEWT using three unit models. As illustrated in Figure 6, these unit models are named as F, R, and RB, respectively. F is a 16 × 16 fuel assembly with the pitch of 20.0112 cm. R is a reflector of same pitch next to a single fuel assembly. RB is a baffle-reflector of the same pitch next to a single fuel assembly. These two reflector models can catch the important neutron behaviors in the reflector region of the core periphery for the nodal calculations. All the calculations are performed with reflective boundary condition. Then, the nodal calculation is carried out using the homogenized two-group cross section sets and the assembly discontinuity factor (ADF) based on the configuration shown in Figure 7. Furthermore, additional studies are conducted about a two-step modeling. The first one is the comparison of a heterogeneous flux distribution and a homogeneous flux distribution in the reflector region. The second one is the necessity of ADF for the interface of two homogeneous nodes.
5. Results
5.1. 2D Direct Modeling Results
The primary reason for the difference of the eigenvalues from each code in Table 2 is the use of different libraries. MCNP6 and SERPENT2 codes use the same continuous energy cross section library, but CASMO-4E and NEWT codes use multigroup cross section libraries. The different calculation method can also be an additional reason for such differences; for example, the deterministic method may produce exact solution to the modified problem with some approximations from the real problem, and the MC method may have solutions with statistical uncertainties to an exactly specified problem with minimal amount of approximations. Therefore, all codes but the two MC codes cannot treat the problem in the same manner and the result shows that two deterministic codes produce eigenvalues underestimated relative to the results of the two MC codes. S06B results of CASMO-4E also have additional uncertainty from the modified baffle thickness.
5.2. Unit Model Result for DIMPLE S06 Two-Step Calculation
The eigenvalue calculation for three unit models is also carried out with each code, the results of which are summarized in Table 3. Additionally, the homogenized two-group cross sections are calculated with SERPENT2, CASMO-4E, and NEWT with the thermal energy group boundary of 0.625 eV. For the three unit models, the homogenized two-group cross sections are calculated with an infinite medium spectrum. The two-group constants are summarized in Table 4. RB model used in CASMO-4E calculation is specified with a baffle of the modified thickness for the consistency between two approaches, the direct one-step approach and the two-step approach. Furthermore, additional calculations are performed to confirm results to those of other codes, in which the thickness of baffle is specified by 2.67 cm to be consistent with the actual thickness. It is available to specify specific thickness of a baffle in this reflector calculation. The eigenvalue of RB model with the actual thickness of 2.67 cm is 1.05968, the difference of which is 47 pcm compared to the model with the modified baffle thickness.
The calculation results in Table 3 show similar tendency for the discrepancy of eigenvalue to the results of S06A and S06B in Section 5.1. As shown in Table 4, the homogenized cross sections generated by each code can vary depending on the method used to solve the problem. For each configuration, the neutron spectrum calculated in each code can be different due to the difference in the method for calculating the neutron spectrum. In particular, the models of R and RB can produce large variations of neutron spectra at the interface of different media.
In Table 4, it is noted that the transport cross sections of CASMO-4E are significantly different relative to the values of the other two codes in all cases. The primary cause is the difference of the transport approximation in calculating transport cross section to treat the anisotropic scattering effect. There are several approximations available for the transport corrected cross section [10]. Both SERPENT2 and NEWT use a simple method [6, 11] for the transport approximation, which results in similar values of transport cross sections. CASMO-4E uses the inflow transport approximation which gives more accurate solution for high-leakage problems. It should be noted that the out-scatter correction predicts overestimated leakage in high-leakage problems such as R and RB models. Table 4 shows that the transport cross sections in the fast energy region in SERPENT2 and NEWT are smaller than that of CASMO-4E. For the other cross sections, it can also be affected by missed nuclides which are not available in CASMO-4E such as P, S, Ti, and Zn.
5.3. 2D Two-Step Modeling Results without ADF
A two-step calculation is performed using PARCS3.0 code with the nodal parameters in Table 4. The results are compared with those from the direct modeling approach. The results of S06A in Table 5 show the maximum error of 591 pcm. On the other hand, all two-step calculation results for S06B, which contains a baffle at periphery of core, show huge error bigger than 1000 pcm. The primary reason is that the heterogeneity of a baffle and water in reflector node is homogenized simply by flux-volume weighting. The baffle mainly consists of steel, which has a significant absorption cross section for thermal neutrons, so that such a simple homogenization technique will not be able to model accurately the thermal flux dip at the baffle region. It is a critical reason for the significant error in these two-step calculation results for S06B.
5.4. Flux Distribution in Reflector Region and ADF
Figure 8 presents the configuration of a fuel assembly and a reflector model used for generating two-group constants of the reflector node. The flux distribution is calculated to show the homogenization effect in the baffle-reflector. The heterogeneous flux distribution is calculated by the mesh tally card of MCNP6 and the homogeneous flux distribution is calculated by solving the 1D-2G diffusion equations:where is the diffusion coefficient, is the absorption cross section, is the down-scattering cross section, and is the nu-fission cross section; and are the energy group indices. The homogenized two-group cross sections in the reflector region are used to solve the equation, and, additionally, the fuel-reflector interface current is used as the boundary condition [12]. Table 6 shows interface currents, and , from each code. For both MCNP6 and SERPENT2, the interface current is generated with the tally and the detector option. Other codes edit the neutron current data in output files automatically with the generation of group constants. Then, the solution of the diffusion problem for the reflector model provides the homogeneous flux distribution in the reflector region as in Figures 9–12.
In Table 6, all interface currents and fluxes agree very well with each other except one outlier; that is, the thermal interface current of CASMO-4E RB model is . The slight difference from other code values is due to the modified thickness of a baffle in CASMO-4E model. If the actual thickness is used, then the currents become , much closer to other results. The thermal interface currents of R reflector are negative which represents the fact that fast neutrons leak from fuel assembly to reflector and get moderated in the water reflector, and then the moderated neutrons enter back to the fuel assembly. On the other hand, the thermal currents for RB model are positive which means that the amount of neutrons coming back to the fuel assembly is smaller due to the absorption in the baffle. Additionally, the reduced thickness of a baffle in CASMO-4E model gives less positive value which means reduction of thermal neutron absorption in the baffle due to smaller amount of steel.
In Figures 9–12, the comparisons of a heterogeneous flux distribution and a homogeneous flux distribution show the discontinuities of the fluxes at the interface of homogenized nodes. In Figure 12, the heterogeneous flux distribution clearly shows the absorption effect for thermal neutrons by 2.67 cm thick baffle. And the homogeneous flux distribution explicitly demonstrates the homogenization effect of a baffle and a reflector. The RB model has the huge discontinuity of a homogeneous flux in the interface between the fuel assembly and the baffle-reflector. It clearly explains the reason why the result of two-step modeling has larger error in the RB model. The errors can be significantly reduced by applying ADF defined as the ratio of a heterogeneous interface flux to a homogeneous interface flux in the following:where ; is the interface heterogeneous flux and is the interface homogeneous flux.
Table 7 summarizes the ADF values obtained from the homogeneous solutions of the homogenized two-group diffusion problem and the heterogeneous solutions from each code. The ADF for single fuel assembly is the ratio of the surface averaged flux to the cell averaged flux [13] and is calculated by each code. The reflector DFs are calculated separately using the data from the codes.
5.5. 2D Two-Step Modeling Results with ADF
Table 8 summarizes the results of a two-step modeling with ADF, which show considerable improvements around 1000 pcm for S06B. The effect of ADF for S06A is not as large as that for S06B since there is less flux discontinuity in the interface of FA-R model as in Figures 9 and 10. In a direct modeling, the corner reflector with baffles on two sides is explicitly modelled as in Figure 5, while a 1D RB model with a baffle on one side is used in two-step calculations. A single set of ADFs generated from the 1D RB model is used for both of the two sides. This 1D approximation in the corner reflector can cause additional errors compared against the direct modeling approach. As shown in CASMO-4E results of S06A, it underestimates the reactivity slightly more relative to the result without ADF.
5.6. Assembly Power Distribution Comparison
The assembly power distributions are compared with each other for a direct modeling approach and a two-step modeling approach. The DIMPLE S06 nodal calculation model is illustrated in Figure 7 by its symmetric structure; only two types of assemblies are used to compare the assembly power distribution. The first one is a fuel assembly (F) in center region and the second one is a fuel assembly (F) in edge as represented by Figure 13. The results are shown in Tables 9 and 10.
Assembly power comparison also shows similar tendency to that of eigenvalue calculation results; that is, for S06A results, the accuracy is already good even without ADF and therefore the improvement by ADF is marginal, while, for S06B results, the accuracy without ADF is low and the use of ADF brings huge improvements in assembly power distribution.
6. Conclusions
The focus of this work is a comprehensive comparative analysis of different modeling approaches and codes for a small high-leakage core, DIMPLE S06 critical experiment. The direct heterogeneous modeling analysis was performed with four different neutron transport solvers: MCNP6, SERPENT2, CASMO-4E, and NEWT. For this analysis, full-core 2D reference data (axially infinite) was generated which has good agreement with measurements. The direct heterogeneous modeling results were compared with these reference data. The two-step homogenization modeling was performed using the nodal diffusion code PARCS 3.0. The homogenized diffusion parameters were generated using the lattice code calculation, and ADFs and homogenized cross sections were implemented in the nodal diffusion code PARCS 3.0. For a consistent comparison, ADFs were calculated from 1D solution of diffusion equation with the generated two-group homogenized constants. In addition, the flux distribution was calculated to show the discontinuity at the fuel-reflector interface. Furthermore, it was clearly shown that the error due to the fuel-reflector discontinuity is significantly reduced by the implementation of ADFs. Comparison for two-group homogenized cross sections from each lattice physics code shows that the generated transport cross section has significant difference according to transport approximation to treat anisotropic scattering effect. Finally, the two modeling approaches for DIMPLE S06 benchmark show a reasonable agreement in both the eigenvalue and the assembly power distributions, while the comparison with the reference generated by MCNP shows a significant error in the eigenvalue prediction except for another Monte Carlo code, SERPENT2. This shows that a small core problem with high leakage is still challenging except for Monte Carlo simulation which has no spatial approximation in the calculation process.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgment
This work was partially supported by a National Research Foundation of Korea (NRF) grant funded by the Korean government (MSIP).