Abstract
The Chebyshev Rational Approximation Method (CRAM) has become one of the dominant methods for solving the Bateman equations for nuclear fuel depletion analysis. Since its introduction over a decade ago, several improvements have been made to CRAM improving its accuracy and reducing its run time. We analyzed its run time using two previously published methods for solving the CRAM system of equations, direct matrix inversion (DMI) and sparse Gaussian elimination (SGE), for depletion systems of varying numbers of nuclides to see how the two methods perform relative to one another. In addition to these two methods, we introduced the Gauss-Seidel (GS) method for solving the CRAM system of equations and compared its performance relative to DMI and SGE for depletion systems with varying numbers of nuclides. We demonstrated that for practical purposes, GS is faster than SGE and DMI and achieves a practical level of accuracy. All testing was performed using the CRAM implementation in the Griffin reactor physics analysis application.
I. INTRODUCTION
The Bateman equations form the mathematical model that describes the concentrations or number densities of various nuclides in a given system where radioactive decay is occurring.[Citation1] While originally developed for radioactive decay, the Bateman equations can be modified to account for changes in nuclide number densities (NNDs) as a result of neutron transmutation. This is shown in EquationEq. (1)(1) (1) , adapted from Ref. [Citation2]:
where
= | = | time-dependent rate of change in the number density of nuclide |
= | = | time-dependent number density of nuclide |
= | = | probability that nuclide will be produced as a fission product as a result of the neutron-induced nuclear fission of nuclide |
= | = | monoenergetic microscopic fission cross section of nuclide |
= | = | time-dependent monoenergetic neutron flux within the system |
= | = | monoenergetic microscopic nonscattering, nonfission reaction cross section of reaction type of nuclide , which produces nuclide |
= | = | probability of nuclide producing nuclide as a result of radioactive decay |
= | = | decay constant of nuclide |
= | = | monoenergetic microscopic absorption (nonscattering) cross section of nuclide |
= | = | decay constant of nuclide . |
Since real-world nuclear systems may have dozens or hundreds of nuclides, matrix-vector notation is convenient to describe the change in NNDs:
where
= | = | time-dependent vector of the NNDs |
= | = | time-dependent vector of the rate of change of the NNDs |
= | = | time-dependent matrix of the terms from EquationEq. (1)(1) (1) for each nuclide in the system, also known as the depletion matrix. |
When expressed in matrix-vector notation, the NNDs at any point in time can be determined by solving the matrix exponential, where for some length of time is constant (i.e., microscopic cross sections and neutron flux are constant) and the initial NNDs in the system are known:
The definition of a matrix exponential is a power series:
Depending on the nuclides tracked in the system, the depletion matrix may be extraordinarily stiff because of the decay constants of short-lived nuclides.[Citation3] As an example, 4H has a half-life of s, which is an extraordinarily small half-life that will result in a very large decay constant term in the depletion matrix relative to the decay and transmutation terms for other nuclides in the system. While an extreme example, it is not uncommon to consider depletion systems where there may be a difference of 30 or 40 orders of magnitude between the largest and the smallest terms in the depletion matrix, which presents numerical challenges for many traditional matrix exponential solvers, such as solving the power series directly.[Citation4]
Because of the challenges in solving the Bateman equations, numerous methods have been developed to solve the Bateman depletion equations, with one of the most prolific methods being the Chebyshev Rational Approximation Method (CRAM).[Citation3,Citation5]
II. CHEBYSHEV RATIONAL APPROXIMATION METHOD
This section provides a brief overview of CRAM. CRAM was developed by Pusa and Leppänen to solve the Bateman depletion equations and was deployed in Serpent, a Monte Carlo particle transport code developed at the VTT Technical Research Centre of Finland.[Citation6] CRAM has been deployed in other applications to solve the Bateman depletion equations, including ORIGEN, developed at Oak Ridge National Laboratory,[Citation7] and Griffin, developed at Idaho National Laboratory and Argonne National Laboratory.[Citation8,Citation9]
CRAM relies on a min/max approximation theory applied to the scalar exponential, which leads to a partial fraction decomposition (PFD) that is a sum of rational partial fractions of order unity, where is also the polynomial degree. The approximation requires matrix inversions. The large number of required matrix inversions was one of the primary shortcomings of CRAM. Further development work on CRAM yielded two major improvements:
Sparse Gaussian elimination (SGE) applied to the CRAM equations resulted in the ability to avoid performing a direct matrix inversion (DMI), significantly reducing the run time for a depletion system with 1600 nuclides.[Citation10]
Incomplete partial fraction (IPF) decomposition, which involved forming a second-order rational partial fraction by combining two partial fractions, resulted in a reduction in round-off error and capability to perform solves at higher approximation orders (AOs).[Citation11]
These two improvements resulted in four different versions of CRAM that can be deployed, any of which can be referred to as CRAM despite their differences:
PFD with DMI, which was the first form of CRAM developed.
PFD with SGE.
IPF decomposition with DMI.
IPF decomposition with SGE, which is the most modern version of CRAM and the version currently deployed in Serpent.
Note that each of the above versions of CRAM also requires a specified AO. For the remainder of this work, we will concern ourselves with only CRAM-IPF as the comparable performance to PFD at similar AOs with the added ability to use higher AOs with IPF than PFD makes it the superior version. We will consider applying both DMI and SGE to solving the CRAM-IPF equations as well as using Gauss-Seidel (GS) to solve the equations and compare the performance of the various methods for varying numbers of nuclides.
For completeness, the equations used to solve CRAM with PFD and CRAM with IPF are included here:
Partial fraction decomposition:
Incomplete partial fraction:
where
= | = | length the time step depletion occurs over |
= | = | AO of the CRAM solve being used |
= | = | number densities of all nuclides at the beginning of the time step |
= | = | time-dependent vector of the number densities of all nuclides in the system at time |
= | = | depletion matrix, which is constant over time |
= | = | identity matrix with the same dimensions as matrix |
, , = | = | constants provided in the CRAM documentation for the given AO . |
In theory, the convergence of the CRAM solution is completely independent of the method used for solving the matrix inversions present in EquationEqs. (5)(5) (5) and Equation(6)(6) (6) ; however, practically speaking, the use of double-precision floating-point operations will inherently limit the potential accuracy of the CRAM solution. The convergence of PFD and its truncation error is discussed in detail in Ref. [Citation12], and the same general behavior is expected to be applicable to IPF as well. Consequently, we do not repeat the same error analysis here and consider it sufficient to say that in general, CRAM with a large AO, such as 48 for IPF, is expected to theoretically provide results with truncation error, which is less than what can be represented in double-precision floating-point operations.[Citation11]
III. GAUSS-SEIDEL METHOD
Gauss-Seidal is a simple conventional stationary iterative method for solving diagonally dominant sparse linear systems. It is defined by the following iteration:
where , , and = strictly lower triangular, diagonal, and strictly upper triangular components of a matrix, respectively. The error reduction rate of GS is proportional to the spectral radius of the iteration matrix , and the spectral radius becomes lower as the original matrix becomes more strongly diagonally dominant.
Applying GS to CRAM is done by simply solving the following linear system with GS for each order iteratively to invert the matrix:
Gauss-Seidel has a low per iteration cost but is difficult to parallelize and converges slowly compared to Krylov subspace methods, so it is not desirable for most modern large-scale numerical analyses. However, the CRAM depletion calculation is inherently parallelizable across independent depletion zones, meaning the linear system solution for each depletion zone does not need to be parallelized. Thus, the inability to parallelize GS is irrelevant for practical depletion analyses using CRAM. Furthermore, the matrices that appear in CRAM are strongly diagonally dominant due to the subtraction of , whose values are orders of magnitude larger than the typical magnitude of physical quantities constituting the depletion matrices, from the diagonal elements. Owing to the strong diagonal dominance, GS exhibits extremely fast convergence for the CRAM matrices, as illustrated in , where the error is defined as the relative -norm of the solution change between two iterations. As a result, the number of GS iterations required to reach a sufficient convergence level becomes small enough to make GS potentially cheaper than SGE. To summarize, the unique characteristics of CRAM applied to the depletion analysis make GS a highly compelling option.
The convergence of GS can be further improved by introducing the successive over-relaxation method with minimal algorithmic changes and virtually the same computational cost, often by a significant degree for high-spectral-radius problems. However, the extremely fast convergence of GS for the CRAM matrices implies that the spectral radii of all CRAM matrices considered in this work are very close to zero. Therefore, the improvement of the convergence rate by introducing successive over-relaxation will be marginal for CRAM even though optimal relaxation factors are used, which are cumbersome to determine. On the contrary, it can even worsen the spectral radius and deteriorate the convergence of GS in this case if an improper over-relaxation factor is used.
IV. RESULTS
We considered several depletion systems, each with different numbers of nuclides, 1693, 299, and 37. We then determined the run times and made comparisons using the simple standards of maximum nuclide absolute relative difference (MNARD) and average nuclide absolute relative difference (ANARD) as needed. Additionally, CRAM occasionally calculates negative nuclide number densities (NNNDs). When these NNNDs occur, the number of these negative densities is included as this is another metric by which solving CRAM using a given method at a given AO can be judged since negative densities are clearly nonphysical.
Additionally, we considered two different criteria when comparing methods. In one instance, we considered the absolute relative differences (ARDs) of the NNDs of all nuclides in the system, regardless of how small the concentration of the nuclide is. This is intended to investigate the computational robustness of a given method. In another instance, we considered only those nuclides with NNDs greater than atoms/b∙cm and with ARDs greater than . This is intended to investigate the practical applicability of a given method to a real-world analysis where real-world uncertainties and tolerances restrict analyses to two or three significant figures.
The depletion conditions used were a one-group neutron flux of n/cm2∙s for a period of s (100 days). These depletion conditions were chosen to both be somewhat represented of the flux levels experienced in thermal reactors as well as challenge CRAM, which is documented as having difficulties modeling fresh fuel systems compared to depleted fuel systems.[Citation5]
IV.A. 1693 Nuclides
The 1693-nuclide library is generated using decay data from Evaluated Nuclear Data File (ENDF) VIII.0 decay data[Citation13] and pretabulated one-group cross-section data taken from ORIGEN[Citation7] data distributed in SCALE 6.2.X versions for beginning-of-life pressurized water reactor fuel at 5% enrichment for a Westinghouse assembly. This library is meant to represent the maximum anticipated number of nuclides tracked in depletion calculations. While ENDF provides decay data for 3820 nuclides,[Citation13] in reality, the lack of cross-section data for many of these nuclides restricts the number of nuclides that can be produced via transmutation in nuclear fission reactors. When SGE was first applied to CRAM, a system using 1606 nuclides was used, which was considered quite comprehensive.[Citation10]
IV.A.1. Direct Matrix Inversion
For DMI, Griffin uses the LAPACK routine ZGESV to solve linear equations where the matrix inverse is required.[Citation14] We first compared the results of the various AOs of CRAM compared to the highest AO for which there are published coefficients, AO-48, and solved using DMI. The results of this comparison are shown in and and .
Based on these results, we can determine that AO-12 can model nuclides with NNDs greater than atoms/b∙cm within an ARD compared to AO-48, which has been demonstrated previously to achieve near floating-point double precision.[Citation11] We also note that as shown in , there is a numerical instability for several nuclides with NNDs below atoms/b∙cm with the results calculated by higher AOs diverging from the AO-48 result. We will investigate this divergence further in Sec. IV.A.2.
IV.A.2. Sparse Gaussian Elimination
We repeated the same process in Sec. IV.A.1 except now applying SGE to CRAM with the results for various AOs shown in and and .
Based on the results from and , we once again determined that AO-12 is capable of achieving the desired precision for practical modeling purposes.
Unlike DMI, SGE does not exhibit the same numerical instability. However, we cannot conclude that SGE is more accurate than DMI based on this alone. In order to determine which is more accurate, we used the adding and doubling method (ADM) implemented in Griffin[Citation15] to generate a set of results independent of the results calculated by CRAM and compare both the DMI and SGE results to the results calculated by ADM. ADM with approximation power 25 was used, which, based on previous testing, is expected to achieve a 9- or 10-digit agreement with CRAM-AO-48. The results of this comparison are shown in and .
Overall, SGE achieves much better agreement with ADM than DMI does, leading to the conclusion that SGE is more accurate than DMI. It should be noted that the system considered in this work corresponds to a fresh fuel case, which CRAM with DMI has been previously documented as experiencing large relative errors for some nuclides with previous research calculating a 246Cm relative error of .[Citation5] Only one nuclide has an ARD greater than between ADM and SGE, which corresponds to 8Be, an extraordinarily short-lived nuclide that ADM has had difficulty modeling.[Citation15] For the remainder of this work, we will still consider the performance of DMI applied to CRAM for the sake of completeness, but based on the two currently published methods and the results from and , we can already confirm the previously published results that CRAM-SGE is superior to CRAM-DMI in terms of both accuracy and run time.
For completeness, we provide information regarding the nuclides for which CRAM-DMI-48 provides NND ARDs greater than when compared to the ADM-25 NND solution as shown in . With the exception of , which has a remarkably large decay constant as discussed previously, all nuclides that experience large ARDs between CRAM-DMI-48 and ADM-25 can be produced as fission products only from the spontaneous fission of . Curium-242 has a decay constant of and a spontaneous fission branching ratio of . Most of the yield fraction terms are on the order of , which results in the depletion matrix having off-diagonal production terms for the metastable radionuclides in on the order of while the metastable nuclides themselves have diagonal destruction terms between to . This difference of up to 25 orders of magnitude can cause numerical problems when working in double-precision, which occurs when performing the LU decomposition required by DMI. A similar issue is present for which is produced via the reaction. This production term will be many orders of magnitude smaller than the destruction term of , once again presenting numerical challenges when using double-precision.
IV.A.3. Gauss-Seidel
We then considered the application of GS to CRAM. Since GS requires a number of iterations to perform, we first considered the performance of GS with a varying number of iterations for AO-48 and compared the results to that of CRAM-SGE-AO-48. The results of these comparisons are shown in and and .
Based on the results, we can determine that even GS-10 does not achieve the same precision as SGE for nuclides with very small NNDs. When considering that GS-10 has a run time of 56 ms compared to the SGE run time of 58 ms from , SGE clearly seems to be the superior method when it comes to computational precision for these low nuclides with very small NNDs. However, as shown in and , GS-5 is capable of achieving an MNARD of less than for all nuclides with NNDs greater than atoms/b∙cm. Additionally, the run time of CRAM-GS-5-AO-48 is only 29 ms, which is half of the CRAM-SGE-AO-48 run time. This means that for practical modeling purposes, CRAM-GS-5-AO-48 could be used instead of CRAM-SGE-AO-48 and reduce the overall run time of the depletion calculation by 50%. Of course, if we consider practical purposes, also demonstrates that CRAM-SGE-AO-12 is capable of achieving an MNARD of less than for nuclides with NND greater than atoms/b∙cm with a run time of only 20 ms, which is less than that of CRAM-GS-5-AO-48. Consequently, we then tested the application of GS with a varying number of iterations to CRAM-AO-12 and compared the numerical results with those of CRAM-SGE-AO-48, as shown in and and .
The CRAM-GS-AO-12 results do not converge to the CRAM-SGE-AO-48 result with increasing number of iterations to the same extent that the CRAM-GS-AO-48 results did, but this is to be expected as the results would be expected to converge to the CRAM-SGE-AO-12 results, which achieved limited agreement with CRAM-SGE-AO-48. As shown in , CRAM-GS-5-AO-12 achieves an MNARD of less than compared to CRAM-SGE-AO-48 for all nuclides with NNDs greater than atoms/b∙cm. The total run time of CRAM-GS-5-AO-12 is 8 ms, which is less than half the run time of CRAM-SGE-AO-12, once again indicating that for practical modeling purposes, CRAM-GS-5-AO-12 can reduce the depletion calculation run time by more than 50% for the 1693-nuclide case.
IV.B. 299 Nuclides
While we have already verified the performance of GS applied to CRAM for a 1693-nuclide library, which is anticipated to be near the upper limit of the number of nuclides tracked in any depletion calculation, we continued our investigation considering smaller nuclide libraries. The motivation for this is that many depletion tools restrict their analysis to a few hundred nuclides or fewer. There are several motivations for this, namely, the desire to avoid tracking extremely short-lived nuclides, which can introduce stiffness to the depletion matrix while having a minimal impact on the depletion calculation. Consequently, the performance of various methods of solving CRAM with varying depletion matrix sizes has practical interest and is investigated here.
The decay and cross-section data used for the 299-nuclide case comes from the DRAGON5[Citation16] application.
IV.B.1. Direct Matrix Inversion
We repeated the process from Sec. IV.A.1 for the 299-nuclide case, the results of which are shown in and and . Unlike , the 299-nuclide case does not exhibit the same numerical instability for a few nuclides when solving CRAM using DMI when compared to the 1693-nuclide case. Because of this lack of instability, we did not perform the same comparison between CRAM-DMI and CRAM-SGE with ADM since there is not inconsistent behavior between the two versions of CRAM.
Based on the results from and , we can determine that AO-12 is suitable for practical modeling purposes, the same as for the 1693-nuclide case.
IV.B.2. Sparse Gaussian Elimination
We repeated the process from Sec. IV.A.2 for the 299-nuclide case, the results of which are shown in and and .
IV.B.3. Gauss-Seidel
Once again, we considered the performance of GS with varying number of iterations applied to CRAM-AO-48, with the results shown in and and .
Interestingly, unlike in Sec. IV.A.3, CRAM-GS-3-AO-48 is able to achieve an MNARD of less than for nuclides with NNDs greater than atoms/b∙cm. This indicates that the 299-nuclide case is easier to converge, and thus, fewer GS iterations are required to achieve the desired agreement based on the given criteria. However, we show in and that even GS-10 cannot achieve convergence with CRAM-SGE-AO-48 for nuclides with very low NNDs, and based on the trends exhibited in , it would likely take GS-20 or GS-30 to achieve the same level of agreement between CRAM-GS-AO-48 and CRAM-SGE-AO-48 compared to CRAM-SGE-AO-44.
Turning back to the practical case, CRAM-GS-3-AO-48 has a run time of only 3 ms, which is almost one-third of the 8.7-ms run time of CRAM-SGE-AO-48. However, once again, we consider that CRAM-SGE-AO-12 achieves a agreement for nuclides with NNDs greater than atoms/b∙cm with a run time of only 2.6 ms. Thus, we consider CRAM-GS-AO-12 with varying numbers of iterations as shown in and and . The results demonstrated that CRAM-GS-3-AO-12 is able to achieve the desired agreement for nuclides with NNDs greater than atoms/b∙cm with a run time of only 0.96 ms, which is roughly one-third the run time of CRAM-SGE-AO-12.
IV.C. 37 Nuclides
We then considered a 37-nuclide case, representing a very small number of nuclides meant to approximate the fewest nuclides expected to be tracked when performing most depletion calculations.
IV.C.1. Direct Matrix Inversion
Once again, we compared the run times of CRAM-DMI at varying AOs as well as the differences in the results of the various AOs compared to CRAM-DMI-AO-48. The results are shown in and and . No numerical instability is seen for the 37-nuclide case for DMI.
IV.C.2. Sparse Gaussian Elimination
As before, we repeated the process for CRAM-SGE for the 37-nuclide case. The results are shown in and and .
IV.C.3. Gauss-Seidel
We then considered varying the number of iterations for CRAM-GS-AO-48. The results are shown in and and . Once again, we observed that GS-3 is able to satisfy the criterion of an MNARD of less than for NNDs greater than atoms/b∙cm with a run time of 0.204 ms, which is slightly more than half the run time of CRAM-SGE-AO-48, which is 0.367 ms. However, we once again note that CRAM-SGE-AO-12 also satisfies this criterion; thus, we once again applied GS to CRAM-AO-12. The results are shown in and and . As shown, CRAM-GS-3-AO-12 is able to achieve the desired agreement with CRAM-SGE-AO-48 with a run time of 0.063 ms compared to the run time of 0.113 ms for CRAM-SGE-AO-12. This once again demonstrates that applying GS to CRAM can yield faster results than SGE while still achieving reasonable accuracy for practical modeling purposes.
V. CONCLUSION
We have demonstrated that GS is a viable alternative to solving the linear equations in CRAM than the previously published DMI and SGE methods. When compared to the results calculated by CRAM-SGE-AO-48, GS was able to achieve an agreement of for all NNDs greater than atoms/b∙cm with a run time roughly half that of the equivalent SGE solution at the same AO. When considering that full-core depletion calculations may involve solving the CRAM equations potentially billions of times, a run-time reduction of 50% could result in hours or potentially days of saved run time. While we have considered only one specific depletion condition (which CRAM was expected to have some challenge solving) for three different nuclide libraries, we anticipate that the reduction in run time and overall accuracy of the solution would be similar for other nuclide libraries with different depletion conditions; however, analysts seeking to use GS should consider performing a cursory comparison between GS and SGE in order to verify that the observed differences are acceptable for their specific use case for nuclides of interest.
Additionally, throughout the course of this work, we have demonstrated that the performance gap between DMI and SGE decreases as the number of nuclides decreases. This expounds on the original CRAM-SGE publication, which demonstrated that SGE was a factor of 10 faster for a specific nuclide library. This work reproduces similar results for large depletion matrices while also demonstrating that the performance gap between DMI and SGE decreases as the number of nuclides in the system decreases. Consequently, developers of depletion-capable applications should not expect a factor of 10 performance improvement when replacing DMI with SGE in their CRAM implementation for all possible depletion problems. Rather, the factor of 10 improvement is the maximum possible run-time improvement when replacing DMI with SGE for CRAM. It is also possible that for incredibly small matrices, such as , DMI may be faster than SGE, but such depletion systems are unlikely to be encountered when performing practical depletion analysis, and even when encountered, the run time needed to solve the CRAM equations for such a small system would be quite insignificant regardless of the method employed.
Last, we note that the convergence of the 1693-nuclide case with GS is much slower than for the 299- and 37-nuclide cases. This can be clearly seen by comparing to and . We attribute this to the existence of very short-lived radionuclides in the 1693-nuclide case, which are not present in the 299- and 37-nuclide cases with the shortest-lived radionuclides in each library shown in . These short-lived radionuclides result in large decay constant terms on the diagonal of the depletion matrix for the radionuclide itself as well as large negative decay constant terms on the off-diagonals of the decay products of these short-lived radionuclides. This can reduce the diagonal dominance of the depletion matrix, thereby degrading the convergence of GS. However, more research is needed to investigate and quantify the impact of these short-lived radionuclides
The run times of the CRAM-AO-48 and CRAM-AO-12 results for various nuclide libraries are shown for DMI, SGE, and GS-3 in and demonstrate the reduced run time of GS-3 compared to SGE and DMI for all cases considered.
Acronyms
ADM: | = | adding and doubling method |
ANARD: | = | average nuclide absolute relative difference |
AO: | = | approximation order |
ARD: | = | absolute relative difference |
CRAM: | = | Chebyshev Rational Approximation Method |
DMI: | = | direct matrix inversion |
ENDF: | = | Evaluated Nuclear Data File |
GS: | = | Gauss-Seidel |
IPF: | = | incomplete partial fraction |
MNARD: | = | maximum nuclide absolute relative difference |
NND: | = | nuclide number density |
NNND: | = | negative nuclide number density. |
PFD: | = | partial fraction decomposition |
SGE: | = | sparse Gaussian elimination |
Acknowledgments
The authors would like to thank Yaqi Wang for their support in implementing this functionality in Griffin.
Disclosure Statement
No potential conflict of interest was reported by the author(s).
Additional information
Funding
References
- H. BATEMAN, “Solution of a System of Differential Equations Occurring in the Theory of Radioactive Transformations,” Proc. Cambridge Philos. Soc. Math. Phys. Sci., 15, 423 (1910).
- G. I. BELL and S. GLASSTONE, Nuclear Reactor Theory, Litton Educational Publishing, Inc. (1970).
- M. PUSA and J. LEPPÄNEN, “Computing the Matrix Exponential in Burnup Calculations,” Nucl. Sci. Eng., 164, 2, 140 (2010); https://doi.org/10.13182/NSE09-14.
- C. MOLER and C. V. LOAN, “Nineteen Dubious Ways to Compute the Exponential of a Matrix, Twenty-Five Years Later,” Soc. Ind. Appl. Math. Rev., 45, 1, 3 (2003).
- M. PUSA, “Rational Approximations to the Matrix Exponential in Burnup Calculations,” Nucl. Sci. Eng., 169, 2, 155 (2011); https://doi.org/10.13182/NSE10-81.
- J. LEPPÄNEN et al., “The Serpent Monte Carlo Code: Status, Development and Applications in 2013,” Ann. Nucl. Energy, 82, 142 (2015); https://doi.org/10.1016/j.anucene.2014.08.024.
- O. W. HERMANN and R. M. WESTFALL, “ORIGEN-S: SCALE System Module to Calculate Fuel Depletion, Actinide Transmutation, Fission Product Buildup and Decay, and Associated Radiation Source Terms,” in “SCALE: A Modular Code System for Performing Standardized Computer Analyses for Licensing Evaluation,” NUREG/CR-0200, Vol II, Sec. F7, Oak Ridge National Laboratory (1998).
- C. H. LEE et al., “Griffin Software Development Plan,” INL/EXT- 21-63185, ANL/NSE-21/23, Idaho National Laboratory/Argonne National Laboratory (2021).
- O. CALVIN et al., “Implementation of Depletion Architecture in the MAMMOTH Reactor Physics Application,” Proc. 2019 Int. Conf. Mathematics and Computational Methods Applied to Nuclear Science and Engineering, Portland, Oregon, August 25–29, 2019, p. 2806, American Nuclear Society (2019).
- M. PUSA and J. LEPPÄNEN, “Solving Linear Systems with Sparse Gaussian Elimination in the Chebyshev Rational Approximation Method,” Nucl. Sci. Eng., 175, 3, 250 (2013); https://doi.org/10.13182/NSE12-52.
- M. PUSA, “Higher-Order Chebyshev Rational Approximation Method and Application to Burnup Equations,” Nucl. Sci. Eng., 182, 3, 297 (2016); https://doi.org/10.13182/NSE15-26.
- O. CALVIN, S. SCHUNERT, and B. GANAPOL, “Global Error Analysis of the Chebyshev Rational Approximation Method,” Ann. Nucl. Energy, 150, 107828 (2021); https://doi.org/10.1016/j.anucene.2020.107828.
- D. A. BROWN et al., “ENDF/B-VIII.0: The 8th Major Release of the Nuclear Reaction Data Library with CIELO-Project Cross Sections, New Standards and Thermal Scattering Data,” Nucl. Data Sheets, 148, 1 (2018); https://doi.org/10.1016/j.nds.2018.02.001.
- E. ANDERSON et al., LAPACK Users’ Guide, Society for Industrial and Applied Mathematics, Philadelphia, Pennsylvania (1999).
- O. W. CALVIN, B. D. GANAPOL, and R. A. BORRELLI, “Introduction of the Adding and Doubling Method for Solving Bateman Equations for Nuclear Fuel Depletion,” Nucl. Sci. Eng., 197, 558 (2023); https://doi.org/10.1080/00295639.2022.2129950.
- G. MARLEAU et al., “A User Guide for DRAGON Version 5,” User Manual IGE-355, Ecole Polytechnique de Montreal (2018).