- Email: [email protected]
Tort solutions to the three-dimensional MOX benchmark, 3-D Extension C5G7MOX
Progress in Nuclear Energy 48 (2006) 445e455 www.elsevier.com/locate/pnucene
Tort solutions to the three-dimensional MOX benchmark, 3-D Extension C5G7MOX Jesse J. Klingensmith a,*, Yousry Y. Azmy b, Jess C. Gehin c, Roberto Orsi d a
The Pennsylvania State University, 126 Hammond Building, University Park, PA 16802, USA b The Pennsylvania State University, 229 Reber Building, University Park, PA 16802, USA c Oak Ridge National Laboratory, P. O. Box 2008, MS 6363, Oak Ridge, TN 37831, USA d ENEA Centro Dati Nucleari, Via Martiri di Monte Sole, 4, 40129 Bologna, Italy
Abstract Solutions to the OECD/NEA three-dimensional fuel assembly benchmark problem, C5G7MOX, obtained using the discreteordinates neutron transport code TORT are presented. The accuracy and convergence of the effective multiplication factor and maximum pin power are examined while varying the fidelity of the geometric model and the order of the angular quadrature, specifically the Square LegendreeChebychev quadrature. The calculated value of the effective multiplication factor is shown to converge asymptotically as the order of the quadrature is increased. The behavior of that asymptotic limit of keff with decreasing computational cell size and increasing staircase resolution of the curved rod-moderator is less converged. Moreover, rapidly increasing computational cost with geometric configuration refinement prevented us from achieving asymptotic convergence with the spatial approximation. The maximum pin power results are generally shown to asymptotically converge to within the statistical error of the reference solution for the problem slices where control rods are not present. For the zones that contain control rods, the solution does not appear to be tightly converged for the meshes and quadrature orders employed. Ó 2006 Elsevier Ltd. All rights reserved. Keywords: TORT; Discrete-ordinates; C5G7MOX benchmark; Multiplication factor; Pin power distribution
1. Introduction The OECD/NEA Experts Group on Three-Dimensional Radiation Transport has previously conducted a benchmark exercise in which a small core configuration consisting of UO2 and MOX PWR fuel assemblies was modeled in three dimensions without spatial homogenizations (Smith et al., 2002). This benchmark proved to be a stringent test for radiation transport codes modeling the local details of fuel rods in two dimensions, but was relatively uniform in the axial direction. A three-dimensional extension of the C5G7 MOX benchmark was proposed as a more realistic test of radiation transport methods’ ability to accurately model three-dimensional nuclear reactor configurations (Smith et al., 2003). This was done by introducing stronger spatial heterogeneity in the third dimension, i.e. along * Corresponding author. E-mail addresses: [email protected] (J.J. Klingensmith), [email protected] (Y.Y. Azmy), [email protected] (J.C. Gehin), [email protected] (R. Orsi). 0149-1970/$ - see front matter Ó 2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.pnucene.2006.01.011
446
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
the z-axis, in the form of control rods. In addition to the problem geometry, detailed elsewhere in this special issue, the material composition, abbreviated (seven groups) cross-section library, and Monte Carlo reference solution were provided for comparison with deterministic methods’ results. The results detailed in this paper were obtained using the DOORS (Discrete Ordinates Oak Ridge System) codes, in particular TORT (Three-Dimensional Discrete Ordinates Neutron/Photon Transport Code) (Rhoades and Simpson, 1997). The DOORS system was developed at the Oak Ridge National Laboratory (ORNL) for use in radiation transport for shielding applications and for modeling nuclear systems. While an MPI-enabled version of TORT exists, all results presented here were obtained using the serial version of the code. All reported computations were performed on the CPILE Beowulf cluster operated by the Nuclear Science and Technology Division of ORNL. The CPILE is comprised of 48 AMD Athlon 2000 MPþ processors (1.667 GHz and 512 MB RAM) running the Linux operating system. The geometry input preparation and output visualization were done using codes in the BOT3P package (Orsi, 2005). 2. Specification of material cross-sections, problem geometry, and solution methodology The complex geometric configurations for this benchmark exercise challenge the users of most deterministic transport codes that use rectangular cells. Not only is the stair-casing of the curved rod-moderator interfaces tedious, the large number of such rods with different material content, e.g. fuel types and control rods among others, increases the potential for error. The DTM3 code, a member of the BOT3P package, simplifies the process by separating the description of the physical configuration from the mesh-based approximate description. The former is constructed via a set of combinatorial instructions, while the latter is automated within some user-defined constraints then written to a file in FIDO format as required by TORT (Rhoades and Simpson, 1997). As a demonstration of the powerful visualization capabilities of the BOT3P package, various 3-D views of the core arrangement are depicted in Fig. 1. Note that these illustrations are generated directly from the TORT input, enabling the user to debug the geometry setup exactly as ‘‘seen’’ by TORT. Fig. 1 is plotted using the Kdiv0 stair-casing approximation detailed later that, due to the model simplicity, gives the clearest illustration. Fig. 2 shows another powerful feature of BOT3P more explicitly by plotting a low-angle view of the rodded B configuration with all of the fuel pins removed. The fission chamber has been retained to give the viewer a sense of the scale of the core. This feature allows the user to check for the correct positioning of any material of interest and to get an intuitive, three-dimensional ‘‘feel’’ for the complicated core geometry. The results presented here were obtained using the xeyez linear nodal method in TORT. The convergence criteria were set to 105 for the flux outer iteration convergence criteria and 106 for the eigenvalue. Generally 50e60 outer UO2 Fuel
7.0% MOX
4.3% MOX
8.7% MOX
Control Rod
OECD / NEA c5g7 MOX benchmark geometry Rodded B Configuration
Fig. 1. Core model with moderator and reflector removed.
Fission Ch.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455 Fission Ch.
CR - Outer
CR - Inner
447
CR - MOX
Fig. 2. Rodded B configuration showing only control rods and fission chambers.
iterations were required for convergence, with four inner iterations per outer. The inner iterations flux convergence criterion was set to 1010 to ensure that all four inners would be completed per outer iteration, as recommended by the TORT manual. 3. Mesh refinement studies Two separate types of mesh refinement were attempted to determine the effect on the accuracy of the effective multiplication factor and pin powers. For the purposes of this report, the variation of the meshing in the x- and y-coordinate directions will be referred to as the radial mesh refinement, and the variation along the z-coordinate direction will be referred to as the axial mesh. For all variations of mesh refinement, the volume of all bodies is exactly conserved by the mesh generation capability of the DTM3 code. The degree of refinement of the staircase approximation is controlled in the BOT3P program GGTM by the single parameter Kdiv. Fig. 3a shows the effect of increasing radial mesh refinement on a single fuel pin cell. For the simplest case, the pin cell is divided into a 3 3 mesh, in which the fuel pin, in this mesh geometry, takes on a square shape. This is the ‘‘Kdiv’’ ¼ 0 case. The refinement in Fig. 2 is shown to be increasing from left to right, or Kdiv of 0, 1, 2, and 4. Kdiv ¼ 1 is the second approximation, which corresponds to a 5 5 mesh and a cruciform-shaped rod. The Kdiv ¼ 2 case corresponds to a 7 7 mesh, and the Kdiv ¼ 4 to an 11 11 mesh. Increasing the value of Kdiv reduces the computational cell size and improves the staircase approximation of the curved surface of the fuel rod, but does so at the expense of an increased number of computational cells, and thus larger memory and hard disk storage requirements, and a longer computing time. Due to these restrictions, results obtained with the Square LegendreeChebychev (SLC) quadrature were restricted to Kdiv ¼ 0, Kdiv ¼ 1 and Kdiv ¼ 2. However, preliminary results using fully symmetric quadratures and Kdiv ¼ 1, Kdiv ¼ 2, and Kdiv ¼ 4 were conducted (as shown in Fig. 4 below) but this level of radial mesh refinement has not been repeated using the SLC quadrature. Kdiv ¼ 4 may still be executed in the future using the SLC quadrature, but is now only included here to help visualize how increasing the mesh refinement appears to increase the fidelity of the model. The axial refinement was also varied to determine its influence on the accuracy of the effective multiplication factor and pin power distribution. The initial configuration consisted of 14 computational cells along the z-axis, divided into three cells per core slice and five cells in the upper moderator. A slice is defined in the benchmark specification to be one-third of the core height, i.e. 14.28 cm. The alternative primary configuration with axial mesh refinement consisted of 36 cells. The mesh consisting of 14 computational cells is referred to as the ‘‘coarse Z’’ mesh, while the configuration using 36 cells is referred to as the ‘‘fine Z’’ mesh. Whether either mesh is strictly coarse or fine is irrelevant, the terms are simply used to distinguish between the two approximations. The fine Z mesh divides the two slices with little or no control rod insertion, slice 1 and slice 2, into four equal cells each. Slice 3 is divided into 12 equal sections and the upper moderator is divided into 16 equal sections. Fig. 3b shows this concept graphically. The smaller cell sizes were
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
448
a
b Upper Moderator
Slice # 3
Slice # 2
Slice # 1
Coarse Mesh: Z1
Fine Mesh: Z2
Fig. 3. (a) Effect of increasing radial mesh refinement corresponding to Kdiv ¼ 0, 1, 2, 4. (b) Comparison of coarse (Z1) and fine (Z2) axial meshing schemes.
used in the hope of increasing the accuracy of the transport solution in the presence of the strongly absorbing control rods. There were several additional configurations tested that included especially finer meshing at the tips of control elements, however, the results of these cases did not show notable improvement over finer or the above mesh options. In order to provide a sense for the relative effect of changing the axial and radial meshing, Table 1 shows the number of bodies and the number of computational cells that comprise each of the six mesh configurations. A cursory 1.1436
S6 S12
1.1434
Eigenvalue
S16
1.1432
Referen
1.143 1.1428 1.1426 1.1424 150
200
250
300
350
400
Cells per Radial Dimension Fig. 4. Effective multiplication factor, unrodded case, with fully symmetric quadratures.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
449
Table 1 Overview of problem configurations and the effect on problem size Number of bodies Kdiv ¼ 0 Kdiv ¼ 1 Kdiv ¼ 2
Coarse Z Fine Z Coarse Z Fine Z Coarse Z Fine Z
4924 14,668 24,412
Number of computational cells 340,704 876,096 505,400 1,299,600 931,896 2,396,304
glance shows that the increased axial refinement comes at the cost of greatly increasing the number of computational cells in the problem, and thus the storage and time required for solution. The nearly 2.4 million cells that comprise the most refined geometry (with coarse Z meshing and a w400 angle quadrature set) consumed approximately two weeks of wall-clock time with the serial version of TORT, hence they were considered to be the most that could be realistically used due to constraints on computational times. While the numerical results show that there is a benefit associated with increasing the fidelity of the model geometry, the number of angles required and computational cost must be taken into account when choosing a model for actual production runs. Table 2 shows the computational time, in days, for the rodded A case with the Kdiv1 level of radial meshing and both the coarse and fine axial meshing. In many cases the number of outer iterations required for convergence varied between the three rodding configurations, but the difference was generally in a narrow range, i.e. one or two iterations, and no statistically significant trends were observed. 4. Angular quadrature refinement studies The two types of angular quadratures employed in this study are the fully symmetric, Sn, and the Square Legendree Chebychev (SLC) quadratures, denoted here by Qn. For both cases, n is the order of the quadrature set. The Sn quadrature provides a perfectly symmetric, hence rotationally invariant, scattering source in the radial plane given the configuration’s symmetry. However, the order of the fully symmetric quadratures is limited by the inherent problem of producing negative angular weights for quadratures of order greater than 20. Since negative angular weights are non-physical in nature, the order of the quadrature set is limited to 16 in this study. The S6, S12, and S16 sets are used in the form that is distributed with the DOORS code system. The Qn, or SLC, quadrature sets do not strictly enforce full symmetry to ensure that the multi-dimensional moments of the angular flux are better conserved. The angular approximations are considered to be separable in the two angular components, using Legendre moments in the m-direction and Chebychev moments in the f-direction. SLC quadratures of orders 6, 8, 10, 12, 14, 16, 18, and 20 were used. Due to the properties of the different numerical expansions used for integration, the SLC quadrature does not produce negative angular weights even for the high orders employed in the present study (Lathrop and Carlson, 1965; Smith, 2003). Table 2 Computational time, in days, required for convergence for the rodded A cases employing Kdiv1 radial meshing SLC quadrature order
Axial refinement Z1
Z2
6 8 10 12 14 16 18 20
0.6 1.0 1.5 2.7 3.9 5.1 6.8 8.3
1.8 3.1 4.5 7.2 10.8 12.4 19.3 23.4
450
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
Table 3 gives a comparison between the two quadrature sets used in terms of the order and the number of non-zero weight angles. Zero-weight angles are included in quadrature sets supplied to TORT as a separator for h-levels. This is necessary to initialize the angular flux in curvilinear geometry, but involves practically no computation in Cartesian geometry, and clearly does not contribute to the scalar flux or other quantities of interest here, namely the eigenvalue or pin powers. Table 3 illustrates one of the limitations of traditional discrete-ordinates calculations, especially in three dimensions. The fully symmetric, three-dimensional quadratures distributed with the DOORS package include only orders 2e16. Fig. 4 shows some of the preliminary results for the unrodded case obtained using the fully symmetric quadrature set with the coarse Z mesh and increasing refinement of the stair-casing approximation. Increasing number of cells per dimension corresponds to Kdiv1, Kdiv2, and Kdiv4. Note that Kdiv4 was used only with the fully symmetric quadrature set, while the final results were obtained using the non-symmetric SLC quadrature. The plot shows that there is no clear convergence with either increasing radial or angular quadrature refinement. This, in keeping with the two-dimensional results, is motivated using a different quadrature set.
5. Results and discussion Figs. 5 through 7 show the convergence of the effective multiplication factor for the three control rod configurations provided in the benchmark on the six meshes described above with increased order of the quadrature. The zero-weight angles have been discounted because they do not contribute to the scalar neutron flux. Again, Kdiv0, Kdiv1, and Kdiv2 are used to denote cases of increasing radial mesh refinement, increasing from zero to two. Coarse and fine Z denote the level of axial mesh refinement. Additionally, the provided reference value for each case is plotted as a horizontal line, with the statistical error in the reference solution represented by the lesser weight lines on either side of the reference. The first thing to note is that there are fewer plotted points for Kdiv ¼ 2 when compared with the Kdiv ¼ 1 and Kdiv ¼ 0 results. The missing cases are not included due to the long computational time they require to converge the outer iterations to a solution. As a general rule, two weeks of wall-clock time per run was the longest time investment allowed; all runs requiring more time were not executed. This study is focused on the three quantities that can have an effect on the accuracy of the solution: the number of angles (angular quadrature order), the radial mesh refinement (increased stair-casing of the fuel rod and generally small cell size), and axial mesh refinement. While these are not completely independent variables, the results will be discussed with respect to these parameters individually, then collectively. The results of the two-dimensional benchmark suggested that very high order angular quadratures (containing more than 500 non-zero weight angles) are necessary for convergence to the asymptotic limit (Azmy, 2004). Figs. 5 through 7 indicate that this trend continues in the three-dimensional extension of the benchmark. For the highest order quadrature (containing w800 angles) the results appear to have nearly converged. This result is consistent with the expectation that using more discrete ordinates yields a more accurate result, but this is done at the cost of increasing the computational cost of the simulation.
Table 3 Order and number of non-zero weight angles for Sn and Qn quadratures Sn
Qn
Order
# of angles
Order
# of angles
6 12 16
48 168 288
6 8 10 12 14 16 18 20
72 128 200 288 392 512 648 800
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
451
effective multiplication factor
1.1435
1.143
1.1425
1.142
reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, fine Z
1.1415
1.141
1.1405 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 5. Convergence of the effective multiplication factor using the SLC quadrature, unrodded case.
The other notable trend is that for all three radial meshes (Kdiv0, Kdiv1, and Kdiv2) the effective multiplication factor accuracy is improved with increased axial mesh refinement. On the plots above, this is denoted as the change from ‘‘coarse Z’’ to ‘‘fine Z’’. The final point to note is that the results do not monotonically approach the reference when the radial meshing is refined. When moving from Kdiv0 to Kdiv1, the solution accuracy actually decreases, but then improves when moving to Kdiv2. This is indicative of the solution failing to reach an asymptotic convergence regime for the Kdiv0 case. This failure could be indicative of the Kdiv0 geometry inaccurately computing leakage from the pin cells. Ideally, Kdiv4 results are necessary to test the convergence as the radial meshing is refined beyond Kdiv2, but this is not practical given the large problem size and long execution time required. Fig. 8 provides a measure of the relative error in the effective multiplication factor with increasing angular quadrature. The results are plotted for the Kdiv1, fine Z configuration, which is chosen for being generally representative of the convergence trends. Clearly Kdiv0 and Kdiv2 relative errors are actually smaller than the corresponding Kdiv1
effective multiplication factor
1.1285
1.128
1.1275
1.127 reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, Fine Z
1.1265
1.126
1.1255 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 6. Convergence of the effective multiplication factor using the SLC quadrature, rodded A configuration.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
452
1.078
effective multiplication factor
1.0775
1.077
1.0765
1.076
reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, fine Z
1.0755
1.075
1.0745 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 7. Convergence of the effective multiplication factor using the SLC quadrature, rodded B configuration.
values shown in Fig. 8. However, this is believed to be a coincidence indicating that the Kdiv0 mesh has not entered the asymptotic regime yet, and therefore its accuracy is not fully reliable. The figure shows that for the lowest quadrature orders, the solution does not appear to have reached the asymptotic regime. Fortunately, the solution does begin to converge smoothly for quadrature sets with greater than w250 angles. This lack of convergence for few-angle approximations was also observed in the radial mesh refinement studies, though computational times limited the radial refinement that were executed. Also notable is that the error for the unrodded and rodded A configurations appears to overlap, while the error for the rodded B configuration is substantially worse. This is directly attributed to the level of control rod insertion in the problems. When moving from the unrodded case to the rodded A case, control rods are inserted in only one-twelfth of the core. Moving to the rodded B case involves control rod insertion in one-third of the core. The error does increase slightly between the unrodded and rodded A case, but only marginally. Conversely, the error increases approximately by one-third between the rodded A and rodded B cases. This is indicative of shortcomings of the method in the presence of control elements that cause large thermal flux gradients. These effects can be observed more clearly by examining the maximum pin power, as shown below.
Effective Multiplication Factor Percent Error
0.3 Unrodded Rodded A Rodded B
0.25 0.2 0.15 0.1 0.05 0 0
200
400
600
800
1000
Number of Non-Zero Angles (Q-type Quad) Fig. 8. Relative errors for all Kdiv1 configurations and fine Z with increasing quadrature order.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
453
The following plots show the maximum calculated pin power as a function of the number of non-zero weight angles. The heavier weight solid line represents the provided reference solution, with the one standard deviation uncertainty bands above and below the reference value shown as dashed lines. The points denoted with an ‘‘x’’ are from the fine Z meshing, while the ‘‘þ’’s denote coarse Z results. Again, the Kdiv1 results were chosen as representative of the three radial meshing refinements attempted. The shape of the curves here indicates that the solution has reached the asymptotic regime for most of the plots. Also, the fine Z meshing is seen to increase the accuracy of the results, often to within the uncertainty band of the reference calculation. More important details can be seen when looking at the differences between the three cases for each of the three axial slices and the overall case (Figs. 9 and 10). When looking at the ‘‘Overall’’ results (overall here implies the entire core height), there are many similarities to both the effective multiplication factor results and to the two-dimensional results. In these cases, the calculated values actually cross the reference and then converge asymptotically with increased number of angles. Again the refinement of the axial mesh leads to increased pin power accuracy. Comparing the results from the three rodding configurations, one can see that the under-prediction of the maximum pin power becomes more pronounced with increased control rod insertion. When comparing the ‘‘Slice #1’’ results, the desired asymptotic behavior is still present as well as the increasing under-prediction of the maximum pin power. A notable difference is the effect of the axial meshing. For the unrodded and rodded A cases, there is little improvement when using the fine Z meshing, but the improvement is much more pronounced in the rodded B case. This is likely due to the axial meshing being more important in the presence of control rods and the associated large flux gradients (Fig. 11). The ‘‘Slice #2’’ results are particularly informative. The results for the unrodded case show moderate improvement with increased axial meshing, though both curves fall within the uncertainty bands of the reference. The rodded A results improve dramatically, with the fine Z meshing giving values within the reference uncertainty. Finally, the rodded B results show a similar trend to the previous curves, but the change in meshing does not improve the accuracy of the solution. This behavior, in some ways, conflicts with the previous conjecture that the axial meshing is more
Unrodded - Slice #1 1.117
2.494
1.116 1.115
2.492
Max. Pin Power
Max. Pin Power
Unrodded - Overall 2.496
2.490 2.488 2.486 2.484 2.482
1.114 1.113 1.112 1.111 1.110 1.109 1.108
2.480
1.107 1.106
2.478 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Unrodded - Slice #2
Unrodded - Slice #3
0.888
0.493
0.887
0.492
Max. Pin Power
Max. Pin Power
300
0.886 0.885 0.884 0.883
0.491 0.490 0.489 0.488
0.882 0.881
0.487 0
100
200
300
400
500
600
Number of angles
700
800
900
0
100
200
300
400
500
600
Number of angles
Fig. 9. Maximum pin power convergence for the unrodded configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
454
Rodded A - Overall
Rodded A - Slice #1
2.262
1.205 1.204
2.260
Max. Pin Power
Max. Pin Power
1.203 2.258 2.256 2.254 2.252
1.202 1.201 1.200 1.199 1.198 1.197
2.250
1.196 1.195
2.248 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Rodded A - Slice #3
Rodded A - Slice #2 0.835
0.306
0.834
0.306
Max. Pin Power
Max. Pin Power
300
0.833 0.832 0.831 0.830
0.305 0.305 0.304 0.304
0.829 0.828
0.303 0
100
200
300
400
500
600
700
800
0
900
100
200
300
400
500
600
Number of angles
Number of angles
Fig. 10. Maximum pin power convergence for the rodded A configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
Rodded B - Slice #1
Rodded B - Overall 1.206
1.839
1.204
1.837
Max. Pin Power
Max. Pin Power
1.838
1.836 1.835 1.834 1.833 1.832
1.202 1.200 1.198 1.196 1.194
1.831
1.192
1.830 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Rodded B - Slice #3
Rodded B - Slice #2 0.555
0.218
0.555
0.218
Max. Pin Power
Max. Pin Power
300
0.554 0.554 0.553 0.553 0.552
0.217 0.217 0.216 0.216 0.215
0.552
0.215 0
100
200
300
400
500
600
Number of angles
700
800
900
0
100
200
300
400
500
600
Number of angles
Fig. 11. Maximum pin power convergence for the rodded B configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
455
important in the presence of control elements. Slice #2 in the rodded B configuration actually contains control elements in one of the four assemblies and while the plotted points appear to overlap, they are actually numerically slightly different. It is possible, however, that the strong neutron absorption from the control rods requires an even finer meshing to ensure that the solution is asymptotically converged. While the maximum pin power would not occur in one of the assemblies containing control rods, the results indicate that proximity to strongly absorbing control elements requires a higher fidelity model in order to achieve convergent results. The cases that showed notable improvement previously were slices that had control rods positioned above them, but not actually within the slice. ‘‘Slice #3’’ continues the trends discussed above. The unrodded case results show great improvement with increased axial meshing, while the rodded A and rodded B configuration results both indicate that the solution may not be fully converged. It should be noted that the lack of convergence in the presence of control elements does not appear in the ‘‘Overall’’ results because the power in those pins is much less than for the rest of the assembly, making the effect on the ‘‘Overall’’ less appreciable. 6. Conclusion In summary, the results show that TORT, when using a sufficiently refined model, can be used to accurately model complex three-dimensional geometries that include strongly absorbing control elements, such as LWR fuel assemblies without homogenization. Convergence of the effective multiplication factor has been observed with increased spatial and angular refinement, though there is evidence of a threshold for point-wise convergence. When increasing the number of non-zero angles in the angular quadrature of the SLC type, the solution displays asymptotic convergence. The improvement of the solution with increased stair-casing of the curved rod-moderator interface has been demonstrated, but conclusive proof has not been obtained due to computational limitations. Both of these results are consistent with those from the earlier two-dimensional exercise using the DORT code. The additional area of investigation is the effect of the axial refinement of the problem geometry. In general, increasing the refinement of the axial mesh was also shown to increase the accuracy of the TORT calculation with respect to the MCNP reference solution. The pin power results exhibited very similar trends as observed in the effective multiplication factor results, especially in areas that are not strongly affected by the presence of control elements. For both the ‘‘Slice #1’’ and ‘‘Overall’’ regions, the maximum pin power shows asymptotic convergence with increased angular refinement and increased accuracy with increased axial refinement. For the regions that contained control elements, the results were not as predictable or well behaved, indicating that the solution is not likely converged in those regions. While this observation is not encouraging, these regions are comparatively less important due to their low fission densities. References Azmy, Y.Y., et al., 2004. DORT solutions to the two-dimensional C5G7MOX benchmark problem. Progress in Nuclear Energy 45 (2e4), 215e 231. Lathrop, K.D., Carlson, B.G., 1965. Discrete Ordinates Angular Quadrate of the Neutron Transport Equation, Los Alamos Scientific Laboratory, LA-3186. Orsi, R., 2005. Bologna Radiation Transport Analysis Pre-Post-Processors Version 4.0. Nuclear Science and Engineering 150, 368e373. Rhoades, W.A., Simpson, D.B., 1997. The TORT Three-Dimensional Discrete Ordinates Neutron/Photon Transport Code, ORNL/TM-13221, Oak Ridge National Laboratory. Smith, M.A., et al., 2002. Benchmark on Deterministic 2-D/3-D MOX Fuel Assembly Transport Calculations without Spatial Homogenization: C5G7MOX Benchmark, NEA/NSC/DOC. Smith, M.A., et al., 2003. Benchmark Specification for Deterministic MOX Fuel Assembly Transport Calculations without Spatial Homogenization (3-D Extension C5G7MOX), NEA/NSC/DOC. Smith, M.A., 2003. Personal communications.
Tort solutions to the three-dimensional MOX benchmark, 3-D Extension C5G7MOX Jesse J. Klingensmith a,*, Yousry Y. Azmy b, Jess C. Gehin c, Roberto Orsi d a
The Pennsylvania State University, 126 Hammond Building, University Park, PA 16802, USA b The Pennsylvania State University, 229 Reber Building, University Park, PA 16802, USA c Oak Ridge National Laboratory, P. O. Box 2008, MS 6363, Oak Ridge, TN 37831, USA d ENEA Centro Dati Nucleari, Via Martiri di Monte Sole, 4, 40129 Bologna, Italy
Abstract Solutions to the OECD/NEA three-dimensional fuel assembly benchmark problem, C5G7MOX, obtained using the discreteordinates neutron transport code TORT are presented. The accuracy and convergence of the effective multiplication factor and maximum pin power are examined while varying the fidelity of the geometric model and the order of the angular quadrature, specifically the Square LegendreeChebychev quadrature. The calculated value of the effective multiplication factor is shown to converge asymptotically as the order of the quadrature is increased. The behavior of that asymptotic limit of keff with decreasing computational cell size and increasing staircase resolution of the curved rod-moderator is less converged. Moreover, rapidly increasing computational cost with geometric configuration refinement prevented us from achieving asymptotic convergence with the spatial approximation. The maximum pin power results are generally shown to asymptotically converge to within the statistical error of the reference solution for the problem slices where control rods are not present. For the zones that contain control rods, the solution does not appear to be tightly converged for the meshes and quadrature orders employed. Ó 2006 Elsevier Ltd. All rights reserved. Keywords: TORT; Discrete-ordinates; C5G7MOX benchmark; Multiplication factor; Pin power distribution
1. Introduction The OECD/NEA Experts Group on Three-Dimensional Radiation Transport has previously conducted a benchmark exercise in which a small core configuration consisting of UO2 and MOX PWR fuel assemblies was modeled in three dimensions without spatial homogenizations (Smith et al., 2002). This benchmark proved to be a stringent test for radiation transport codes modeling the local details of fuel rods in two dimensions, but was relatively uniform in the axial direction. A three-dimensional extension of the C5G7 MOX benchmark was proposed as a more realistic test of radiation transport methods’ ability to accurately model three-dimensional nuclear reactor configurations (Smith et al., 2003). This was done by introducing stronger spatial heterogeneity in the third dimension, i.e. along * Corresponding author. E-mail addresses: [email protected] (J.J. Klingensmith), [email protected] (Y.Y. Azmy), [email protected] (J.C. Gehin), [email protected] (R. Orsi). 0149-1970/$ - see front matter Ó 2006 Elsevier Ltd. All rights reserved. doi:10.1016/j.pnucene.2006.01.011
446
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
the z-axis, in the form of control rods. In addition to the problem geometry, detailed elsewhere in this special issue, the material composition, abbreviated (seven groups) cross-section library, and Monte Carlo reference solution were provided for comparison with deterministic methods’ results. The results detailed in this paper were obtained using the DOORS (Discrete Ordinates Oak Ridge System) codes, in particular TORT (Three-Dimensional Discrete Ordinates Neutron/Photon Transport Code) (Rhoades and Simpson, 1997). The DOORS system was developed at the Oak Ridge National Laboratory (ORNL) for use in radiation transport for shielding applications and for modeling nuclear systems. While an MPI-enabled version of TORT exists, all results presented here were obtained using the serial version of the code. All reported computations were performed on the CPILE Beowulf cluster operated by the Nuclear Science and Technology Division of ORNL. The CPILE is comprised of 48 AMD Athlon 2000 MPþ processors (1.667 GHz and 512 MB RAM) running the Linux operating system. The geometry input preparation and output visualization were done using codes in the BOT3P package (Orsi, 2005). 2. Specification of material cross-sections, problem geometry, and solution methodology The complex geometric configurations for this benchmark exercise challenge the users of most deterministic transport codes that use rectangular cells. Not only is the stair-casing of the curved rod-moderator interfaces tedious, the large number of such rods with different material content, e.g. fuel types and control rods among others, increases the potential for error. The DTM3 code, a member of the BOT3P package, simplifies the process by separating the description of the physical configuration from the mesh-based approximate description. The former is constructed via a set of combinatorial instructions, while the latter is automated within some user-defined constraints then written to a file in FIDO format as required by TORT (Rhoades and Simpson, 1997). As a demonstration of the powerful visualization capabilities of the BOT3P package, various 3-D views of the core arrangement are depicted in Fig. 1. Note that these illustrations are generated directly from the TORT input, enabling the user to debug the geometry setup exactly as ‘‘seen’’ by TORT. Fig. 1 is plotted using the Kdiv0 stair-casing approximation detailed later that, due to the model simplicity, gives the clearest illustration. Fig. 2 shows another powerful feature of BOT3P more explicitly by plotting a low-angle view of the rodded B configuration with all of the fuel pins removed. The fission chamber has been retained to give the viewer a sense of the scale of the core. This feature allows the user to check for the correct positioning of any material of interest and to get an intuitive, three-dimensional ‘‘feel’’ for the complicated core geometry. The results presented here were obtained using the xeyez linear nodal method in TORT. The convergence criteria were set to 105 for the flux outer iteration convergence criteria and 106 for the eigenvalue. Generally 50e60 outer UO2 Fuel
7.0% MOX
4.3% MOX
8.7% MOX
Control Rod
OECD / NEA c5g7 MOX benchmark geometry Rodded B Configuration
Fig. 1. Core model with moderator and reflector removed.
Fission Ch.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455 Fission Ch.
CR - Outer
CR - Inner
447
CR - MOX
Fig. 2. Rodded B configuration showing only control rods and fission chambers.
iterations were required for convergence, with four inner iterations per outer. The inner iterations flux convergence criterion was set to 1010 to ensure that all four inners would be completed per outer iteration, as recommended by the TORT manual. 3. Mesh refinement studies Two separate types of mesh refinement were attempted to determine the effect on the accuracy of the effective multiplication factor and pin powers. For the purposes of this report, the variation of the meshing in the x- and y-coordinate directions will be referred to as the radial mesh refinement, and the variation along the z-coordinate direction will be referred to as the axial mesh. For all variations of mesh refinement, the volume of all bodies is exactly conserved by the mesh generation capability of the DTM3 code. The degree of refinement of the staircase approximation is controlled in the BOT3P program GGTM by the single parameter Kdiv. Fig. 3a shows the effect of increasing radial mesh refinement on a single fuel pin cell. For the simplest case, the pin cell is divided into a 3 3 mesh, in which the fuel pin, in this mesh geometry, takes on a square shape. This is the ‘‘Kdiv’’ ¼ 0 case. The refinement in Fig. 2 is shown to be increasing from left to right, or Kdiv of 0, 1, 2, and 4. Kdiv ¼ 1 is the second approximation, which corresponds to a 5 5 mesh and a cruciform-shaped rod. The Kdiv ¼ 2 case corresponds to a 7 7 mesh, and the Kdiv ¼ 4 to an 11 11 mesh. Increasing the value of Kdiv reduces the computational cell size and improves the staircase approximation of the curved surface of the fuel rod, but does so at the expense of an increased number of computational cells, and thus larger memory and hard disk storage requirements, and a longer computing time. Due to these restrictions, results obtained with the Square LegendreeChebychev (SLC) quadrature were restricted to Kdiv ¼ 0, Kdiv ¼ 1 and Kdiv ¼ 2. However, preliminary results using fully symmetric quadratures and Kdiv ¼ 1, Kdiv ¼ 2, and Kdiv ¼ 4 were conducted (as shown in Fig. 4 below) but this level of radial mesh refinement has not been repeated using the SLC quadrature. Kdiv ¼ 4 may still be executed in the future using the SLC quadrature, but is now only included here to help visualize how increasing the mesh refinement appears to increase the fidelity of the model. The axial refinement was also varied to determine its influence on the accuracy of the effective multiplication factor and pin power distribution. The initial configuration consisted of 14 computational cells along the z-axis, divided into three cells per core slice and five cells in the upper moderator. A slice is defined in the benchmark specification to be one-third of the core height, i.e. 14.28 cm. The alternative primary configuration with axial mesh refinement consisted of 36 cells. The mesh consisting of 14 computational cells is referred to as the ‘‘coarse Z’’ mesh, while the configuration using 36 cells is referred to as the ‘‘fine Z’’ mesh. Whether either mesh is strictly coarse or fine is irrelevant, the terms are simply used to distinguish between the two approximations. The fine Z mesh divides the two slices with little or no control rod insertion, slice 1 and slice 2, into four equal cells each. Slice 3 is divided into 12 equal sections and the upper moderator is divided into 16 equal sections. Fig. 3b shows this concept graphically. The smaller cell sizes were
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
448
a
b Upper Moderator
Slice # 3
Slice # 2
Slice # 1
Coarse Mesh: Z1
Fine Mesh: Z2
Fig. 3. (a) Effect of increasing radial mesh refinement corresponding to Kdiv ¼ 0, 1, 2, 4. (b) Comparison of coarse (Z1) and fine (Z2) axial meshing schemes.
used in the hope of increasing the accuracy of the transport solution in the presence of the strongly absorbing control rods. There were several additional configurations tested that included especially finer meshing at the tips of control elements, however, the results of these cases did not show notable improvement over finer or the above mesh options. In order to provide a sense for the relative effect of changing the axial and radial meshing, Table 1 shows the number of bodies and the number of computational cells that comprise each of the six mesh configurations. A cursory 1.1436
S6 S12
1.1434
Eigenvalue
S16
1.1432
Referen
1.143 1.1428 1.1426 1.1424 150
200
250
300
350
400
Cells per Radial Dimension Fig. 4. Effective multiplication factor, unrodded case, with fully symmetric quadratures.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
449
Table 1 Overview of problem configurations and the effect on problem size Number of bodies Kdiv ¼ 0 Kdiv ¼ 1 Kdiv ¼ 2
Coarse Z Fine Z Coarse Z Fine Z Coarse Z Fine Z
4924 14,668 24,412
Number of computational cells 340,704 876,096 505,400 1,299,600 931,896 2,396,304
glance shows that the increased axial refinement comes at the cost of greatly increasing the number of computational cells in the problem, and thus the storage and time required for solution. The nearly 2.4 million cells that comprise the most refined geometry (with coarse Z meshing and a w400 angle quadrature set) consumed approximately two weeks of wall-clock time with the serial version of TORT, hence they were considered to be the most that could be realistically used due to constraints on computational times. While the numerical results show that there is a benefit associated with increasing the fidelity of the model geometry, the number of angles required and computational cost must be taken into account when choosing a model for actual production runs. Table 2 shows the computational time, in days, for the rodded A case with the Kdiv1 level of radial meshing and both the coarse and fine axial meshing. In many cases the number of outer iterations required for convergence varied between the three rodding configurations, but the difference was generally in a narrow range, i.e. one or two iterations, and no statistically significant trends were observed. 4. Angular quadrature refinement studies The two types of angular quadratures employed in this study are the fully symmetric, Sn, and the Square Legendree Chebychev (SLC) quadratures, denoted here by Qn. For both cases, n is the order of the quadrature set. The Sn quadrature provides a perfectly symmetric, hence rotationally invariant, scattering source in the radial plane given the configuration’s symmetry. However, the order of the fully symmetric quadratures is limited by the inherent problem of producing negative angular weights for quadratures of order greater than 20. Since negative angular weights are non-physical in nature, the order of the quadrature set is limited to 16 in this study. The S6, S12, and S16 sets are used in the form that is distributed with the DOORS code system. The Qn, or SLC, quadrature sets do not strictly enforce full symmetry to ensure that the multi-dimensional moments of the angular flux are better conserved. The angular approximations are considered to be separable in the two angular components, using Legendre moments in the m-direction and Chebychev moments in the f-direction. SLC quadratures of orders 6, 8, 10, 12, 14, 16, 18, and 20 were used. Due to the properties of the different numerical expansions used for integration, the SLC quadrature does not produce negative angular weights even for the high orders employed in the present study (Lathrop and Carlson, 1965; Smith, 2003). Table 2 Computational time, in days, required for convergence for the rodded A cases employing Kdiv1 radial meshing SLC quadrature order
Axial refinement Z1
Z2
6 8 10 12 14 16 18 20
0.6 1.0 1.5 2.7 3.9 5.1 6.8 8.3
1.8 3.1 4.5 7.2 10.8 12.4 19.3 23.4
450
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
Table 3 gives a comparison between the two quadrature sets used in terms of the order and the number of non-zero weight angles. Zero-weight angles are included in quadrature sets supplied to TORT as a separator for h-levels. This is necessary to initialize the angular flux in curvilinear geometry, but involves practically no computation in Cartesian geometry, and clearly does not contribute to the scalar flux or other quantities of interest here, namely the eigenvalue or pin powers. Table 3 illustrates one of the limitations of traditional discrete-ordinates calculations, especially in three dimensions. The fully symmetric, three-dimensional quadratures distributed with the DOORS package include only orders 2e16. Fig. 4 shows some of the preliminary results for the unrodded case obtained using the fully symmetric quadrature set with the coarse Z mesh and increasing refinement of the stair-casing approximation. Increasing number of cells per dimension corresponds to Kdiv1, Kdiv2, and Kdiv4. Note that Kdiv4 was used only with the fully symmetric quadrature set, while the final results were obtained using the non-symmetric SLC quadrature. The plot shows that there is no clear convergence with either increasing radial or angular quadrature refinement. This, in keeping with the two-dimensional results, is motivated using a different quadrature set.
5. Results and discussion Figs. 5 through 7 show the convergence of the effective multiplication factor for the three control rod configurations provided in the benchmark on the six meshes described above with increased order of the quadrature. The zero-weight angles have been discounted because they do not contribute to the scalar neutron flux. Again, Kdiv0, Kdiv1, and Kdiv2 are used to denote cases of increasing radial mesh refinement, increasing from zero to two. Coarse and fine Z denote the level of axial mesh refinement. Additionally, the provided reference value for each case is plotted as a horizontal line, with the statistical error in the reference solution represented by the lesser weight lines on either side of the reference. The first thing to note is that there are fewer plotted points for Kdiv ¼ 2 when compared with the Kdiv ¼ 1 and Kdiv ¼ 0 results. The missing cases are not included due to the long computational time they require to converge the outer iterations to a solution. As a general rule, two weeks of wall-clock time per run was the longest time investment allowed; all runs requiring more time were not executed. This study is focused on the three quantities that can have an effect on the accuracy of the solution: the number of angles (angular quadrature order), the radial mesh refinement (increased stair-casing of the fuel rod and generally small cell size), and axial mesh refinement. While these are not completely independent variables, the results will be discussed with respect to these parameters individually, then collectively. The results of the two-dimensional benchmark suggested that very high order angular quadratures (containing more than 500 non-zero weight angles) are necessary for convergence to the asymptotic limit (Azmy, 2004). Figs. 5 through 7 indicate that this trend continues in the three-dimensional extension of the benchmark. For the highest order quadrature (containing w800 angles) the results appear to have nearly converged. This result is consistent with the expectation that using more discrete ordinates yields a more accurate result, but this is done at the cost of increasing the computational cost of the simulation.
Table 3 Order and number of non-zero weight angles for Sn and Qn quadratures Sn
Qn
Order
# of angles
Order
# of angles
6 12 16
48 168 288
6 8 10 12 14 16 18 20
72 128 200 288 392 512 648 800
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
451
effective multiplication factor
1.1435
1.143
1.1425
1.142
reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, fine Z
1.1415
1.141
1.1405 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 5. Convergence of the effective multiplication factor using the SLC quadrature, unrodded case.
The other notable trend is that for all three radial meshes (Kdiv0, Kdiv1, and Kdiv2) the effective multiplication factor accuracy is improved with increased axial mesh refinement. On the plots above, this is denoted as the change from ‘‘coarse Z’’ to ‘‘fine Z’’. The final point to note is that the results do not monotonically approach the reference when the radial meshing is refined. When moving from Kdiv0 to Kdiv1, the solution accuracy actually decreases, but then improves when moving to Kdiv2. This is indicative of the solution failing to reach an asymptotic convergence regime for the Kdiv0 case. This failure could be indicative of the Kdiv0 geometry inaccurately computing leakage from the pin cells. Ideally, Kdiv4 results are necessary to test the convergence as the radial meshing is refined beyond Kdiv2, but this is not practical given the large problem size and long execution time required. Fig. 8 provides a measure of the relative error in the effective multiplication factor with increasing angular quadrature. The results are plotted for the Kdiv1, fine Z configuration, which is chosen for being generally representative of the convergence trends. Clearly Kdiv0 and Kdiv2 relative errors are actually smaller than the corresponding Kdiv1
effective multiplication factor
1.1285
1.128
1.1275
1.127 reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, Fine Z
1.1265
1.126
1.1255 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 6. Convergence of the effective multiplication factor using the SLC quadrature, rodded A configuration.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
452
1.078
effective multiplication factor
1.0775
1.077
1.0765
1.076
reference Kdiv0, coarse Z Kdiv0, fine Z Kdiv1, coarse Z Kdiv1, fine Z Kdiv2, coarse Z Kdiv2, fine Z
1.0755
1.075
1.0745 0
100
200
300
400
500
600
700
800
900
number of non-zero weight angles Fig. 7. Convergence of the effective multiplication factor using the SLC quadrature, rodded B configuration.
values shown in Fig. 8. However, this is believed to be a coincidence indicating that the Kdiv0 mesh has not entered the asymptotic regime yet, and therefore its accuracy is not fully reliable. The figure shows that for the lowest quadrature orders, the solution does not appear to have reached the asymptotic regime. Fortunately, the solution does begin to converge smoothly for quadrature sets with greater than w250 angles. This lack of convergence for few-angle approximations was also observed in the radial mesh refinement studies, though computational times limited the radial refinement that were executed. Also notable is that the error for the unrodded and rodded A configurations appears to overlap, while the error for the rodded B configuration is substantially worse. This is directly attributed to the level of control rod insertion in the problems. When moving from the unrodded case to the rodded A case, control rods are inserted in only one-twelfth of the core. Moving to the rodded B case involves control rod insertion in one-third of the core. The error does increase slightly between the unrodded and rodded A case, but only marginally. Conversely, the error increases approximately by one-third between the rodded A and rodded B cases. This is indicative of shortcomings of the method in the presence of control elements that cause large thermal flux gradients. These effects can be observed more clearly by examining the maximum pin power, as shown below.
Effective Multiplication Factor Percent Error
0.3 Unrodded Rodded A Rodded B
0.25 0.2 0.15 0.1 0.05 0 0
200
400
600
800
1000
Number of Non-Zero Angles (Q-type Quad) Fig. 8. Relative errors for all Kdiv1 configurations and fine Z with increasing quadrature order.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
453
The following plots show the maximum calculated pin power as a function of the number of non-zero weight angles. The heavier weight solid line represents the provided reference solution, with the one standard deviation uncertainty bands above and below the reference value shown as dashed lines. The points denoted with an ‘‘x’’ are from the fine Z meshing, while the ‘‘þ’’s denote coarse Z results. Again, the Kdiv1 results were chosen as representative of the three radial meshing refinements attempted. The shape of the curves here indicates that the solution has reached the asymptotic regime for most of the plots. Also, the fine Z meshing is seen to increase the accuracy of the results, often to within the uncertainty band of the reference calculation. More important details can be seen when looking at the differences between the three cases for each of the three axial slices and the overall case (Figs. 9 and 10). When looking at the ‘‘Overall’’ results (overall here implies the entire core height), there are many similarities to both the effective multiplication factor results and to the two-dimensional results. In these cases, the calculated values actually cross the reference and then converge asymptotically with increased number of angles. Again the refinement of the axial mesh leads to increased pin power accuracy. Comparing the results from the three rodding configurations, one can see that the under-prediction of the maximum pin power becomes more pronounced with increased control rod insertion. When comparing the ‘‘Slice #1’’ results, the desired asymptotic behavior is still present as well as the increasing under-prediction of the maximum pin power. A notable difference is the effect of the axial meshing. For the unrodded and rodded A cases, there is little improvement when using the fine Z meshing, but the improvement is much more pronounced in the rodded B case. This is likely due to the axial meshing being more important in the presence of control rods and the associated large flux gradients (Fig. 11). The ‘‘Slice #2’’ results are particularly informative. The results for the unrodded case show moderate improvement with increased axial meshing, though both curves fall within the uncertainty bands of the reference. The rodded A results improve dramatically, with the fine Z meshing giving values within the reference uncertainty. Finally, the rodded B results show a similar trend to the previous curves, but the change in meshing does not improve the accuracy of the solution. This behavior, in some ways, conflicts with the previous conjecture that the axial meshing is more
Unrodded - Slice #1 1.117
2.494
1.116 1.115
2.492
Max. Pin Power
Max. Pin Power
Unrodded - Overall 2.496
2.490 2.488 2.486 2.484 2.482
1.114 1.113 1.112 1.111 1.110 1.109 1.108
2.480
1.107 1.106
2.478 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Unrodded - Slice #2
Unrodded - Slice #3
0.888
0.493
0.887
0.492
Max. Pin Power
Max. Pin Power
300
0.886 0.885 0.884 0.883
0.491 0.490 0.489 0.488
0.882 0.881
0.487 0
100
200
300
400
500
600
Number of angles
700
800
900
0
100
200
300
400
500
600
Number of angles
Fig. 9. Maximum pin power convergence for the unrodded configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
454
Rodded A - Overall
Rodded A - Slice #1
2.262
1.205 1.204
2.260
Max. Pin Power
Max. Pin Power
1.203 2.258 2.256 2.254 2.252
1.202 1.201 1.200 1.199 1.198 1.197
2.250
1.196 1.195
2.248 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Rodded A - Slice #3
Rodded A - Slice #2 0.835
0.306
0.834
0.306
Max. Pin Power
Max. Pin Power
300
0.833 0.832 0.831 0.830
0.305 0.305 0.304 0.304
0.829 0.828
0.303 0
100
200
300
400
500
600
700
800
0
900
100
200
300
400
500
600
Number of angles
Number of angles
Fig. 10. Maximum pin power convergence for the rodded A configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
Rodded B - Slice #1
Rodded B - Overall 1.206
1.839
1.204
1.837
Max. Pin Power
Max. Pin Power
1.838
1.836 1.835 1.834 1.833 1.832
1.202 1.200 1.198 1.196 1.194
1.831
1.192
1.830 0
100
200
300
400
500
600
700
800
900
0
100
200
Number of angles
400
500
600
700
800
900
700
800
900
Number of angles
Rodded B - Slice #3
Rodded B - Slice #2 0.555
0.218
0.555
0.218
Max. Pin Power
Max. Pin Power
300
0.554 0.554 0.553 0.553 0.552
0.217 0.217 0.216 0.216 0.215
0.552
0.215 0
100
200
300
400
500
600
Number of angles
700
800
900
0
100
200
300
400
500
600
Number of angles
Fig. 11. Maximum pin power convergence for the rodded B configuration with Kdiv1, coarse Z (þ) and fine Z (x) meshing.
J.J. Klingensmith et al. / Progress in Nuclear Energy 48 (2006) 445e455
455
important in the presence of control elements. Slice #2 in the rodded B configuration actually contains control elements in one of the four assemblies and while the plotted points appear to overlap, they are actually numerically slightly different. It is possible, however, that the strong neutron absorption from the control rods requires an even finer meshing to ensure that the solution is asymptotically converged. While the maximum pin power would not occur in one of the assemblies containing control rods, the results indicate that proximity to strongly absorbing control elements requires a higher fidelity model in order to achieve convergent results. The cases that showed notable improvement previously were slices that had control rods positioned above them, but not actually within the slice. ‘‘Slice #3’’ continues the trends discussed above. The unrodded case results show great improvement with increased axial meshing, while the rodded A and rodded B configuration results both indicate that the solution may not be fully converged. It should be noted that the lack of convergence in the presence of control elements does not appear in the ‘‘Overall’’ results because the power in those pins is much less than for the rest of the assembly, making the effect on the ‘‘Overall’’ less appreciable. 6. Conclusion In summary, the results show that TORT, when using a sufficiently refined model, can be used to accurately model complex three-dimensional geometries that include strongly absorbing control elements, such as LWR fuel assemblies without homogenization. Convergence of the effective multiplication factor has been observed with increased spatial and angular refinement, though there is evidence of a threshold for point-wise convergence. When increasing the number of non-zero angles in the angular quadrature of the SLC type, the solution displays asymptotic convergence. The improvement of the solution with increased stair-casing of the curved rod-moderator interface has been demonstrated, but conclusive proof has not been obtained due to computational limitations. Both of these results are consistent with those from the earlier two-dimensional exercise using the DORT code. The additional area of investigation is the effect of the axial refinement of the problem geometry. In general, increasing the refinement of the axial mesh was also shown to increase the accuracy of the TORT calculation with respect to the MCNP reference solution. The pin power results exhibited very similar trends as observed in the effective multiplication factor results, especially in areas that are not strongly affected by the presence of control elements. For both the ‘‘Slice #1’’ and ‘‘Overall’’ regions, the maximum pin power shows asymptotic convergence with increased angular refinement and increased accuracy with increased axial refinement. For the regions that contained control elements, the results were not as predictable or well behaved, indicating that the solution is not likely converged in those regions. While this observation is not encouraging, these regions are comparatively less important due to their low fission densities. References Azmy, Y.Y., et al., 2004. DORT solutions to the two-dimensional C5G7MOX benchmark problem. Progress in Nuclear Energy 45 (2e4), 215e 231. Lathrop, K.D., Carlson, B.G., 1965. Discrete Ordinates Angular Quadrate of the Neutron Transport Equation, Los Alamos Scientific Laboratory, LA-3186. Orsi, R., 2005. Bologna Radiation Transport Analysis Pre-Post-Processors Version 4.0. Nuclear Science and Engineering 150, 368e373. Rhoades, W.A., Simpson, D.B., 1997. The TORT Three-Dimensional Discrete Ordinates Neutron/Photon Transport Code, ORNL/TM-13221, Oak Ridge National Laboratory. Smith, M.A., et al., 2002. Benchmark on Deterministic 2-D/3-D MOX Fuel Assembly Transport Calculations without Spatial Homogenization: C5G7MOX Benchmark, NEA/NSC/DOC. Smith, M.A., et al., 2003. Benchmark Specification for Deterministic MOX Fuel Assembly Transport Calculations without Spatial Homogenization (3-D Extension C5G7MOX), NEA/NSC/DOC. Smith, M.A., 2003. Personal communications.