4.8. Factors influencing simulation time

Though DCHAIN normally runs quite quickly (at least relative to PHITS), one may notice that the run time of DCHAIN is variable depending on a variety of factors including the neutron reaction cross section library used and various other settings.

To more discretely understand what determines the length of a DCHAIN simulation, one must understand what actions occur and at what frequency during a simulation. At the most basic level, the DCHAIN simulation is subdivided into \(N_{\text{reg}}\) sub-calculations where \(N_{\text{reg}}\) is the number of regions in the calculation. Within each of these \(N_{\text{reg}}\) independent sub-simulations, the cross section data library is first loaded (taking a constant time \(T_{\text{rd }}\)\(_\sigma\)\(_{\text{ lib}}\)) and the neutron flux spectrum and yield data are read in, taking \(T_{\text{rd flux}}\) and \(T_{\text{rd yield}}\), respectively. The decay data is only read in once regardless of region count, taking \(T_{\text{rd decay lib}}\). The flux spectrum file is also only read partially in each region, stopping once the flux of the current region has been reached and read.

Then, the inventory of each nuclide is propagated through all calculation time steps, one step at a time, taking a variable amount of time \(T_{\text{calc}}\) in each. This calculation time is dependent on, first, the number of calculation time steps \(N_{\text{calc steps}}\), explained below.

Essentially, the user provides a number of variables which control the number of calculation time steps. At its most basic level, the simulation is divided into periods of irradiation and periods of cooling; ISTEP is the sum of the number of these separate irradiation and cooling steps with distinct beam power levels. Each cooling step is calculated as a single step. Each irradiation step is actually subdivided into IDIVS (D=50) equally spaced sub-time steps. Furthermore, the user may request that data be outputted at any time during irradiation or cooling; ITOUT is the number of total times at which the user is requesting this output. One additional calculation time step is performed for each time requested during (in the middle of) a cooling period. For requests for output during the middle of an irradiation time step, the existing irradiation time step is divided into smaller irradiation time steps for each output time requested during that irradiation period. An example of this relationship is shown in Figure 4.8.1.

Fig. 4.8.1 Determining the number of calculation steps from beam power schedule, requested output times ITOUT, and number of subdivisions per irradiation step IDIVS

This ultimately means that requesting output in the middle of an irradiation time step has a significantly higher cost on the run-time of the simulation than a request during a cooling period. Note that requests for output at the exact end of a cooling or irradiation time step, which is already defined by the irradiation/cooling step section, do not result in an additional calculated time step; the code simply chooses to output the values for the end of that time step rather than assuming the default behavior of suppressing output.

Additionally, each calculation time step’s length is dependent on the complexity of the problem, as in the number of nuclides present and how they can decay. During irradiation time steps, the calculation is considerably more complex since all of the neutron activation channels are opened as being possible and the code explores all possible routes for decay and production of the present nuclides in each calculation time step. These relations using DCHAIN’s default settings are summarized in Equations (4.8.1), (4.8.2), and (4.8.3).

(4.8.1)\[T_{\text{total}} \propto T_{\text{rd decay lib}} + \sum_{i=1}^{N_{\text{reg}}} \left( T_{\text{rd $\sigma$ lib}} + \frac{i}{N_{\text{reg}}}T_{\text{rd flux}} +T_{\text{rd yield}} +N_{\text{calc steps}} \cdot T_{\text{calc}} \right)\]

(4.8.2)\[N_{\text{calc steps}} \propto \textcolor{BlueViolet}{\texttt{IDIVS}}\cdot(N_{\text{irradiation}}+N_{\text{output,beam-on}}) + N_{\text{cooling}}+N_{\text{output,beam-off}}\]

(4.8.3)\[T_{\text{calc}} \propto N_{\text{nuclides-present}} \cdot (N_{\text{decay-pathways}} + \delta_{\text{beam on/off}} \cdot N_{\text{activation pathways}})\]

With these old default DCHAIN settings, the decay library is read in only once, the neutron reaction cross section library and yield data files are read completely in each region, and the neutron flux file is read partially (up to the current region) for every region in a simulation. Recent performance improvements have sought to dramatically decrease these times.

First, IXSRALL= 1 is an option which allows for the cross section library to be read and stored into memory (RAM) once at the beginning of a calculation rather than being reread in every region. Due to the approximately 100 MB size of the current default library, this clearly increases RAM usage significantly but also dramatically reduces the amount of repetitive reading done for problems with numerous regions.

Second, IREDUFMT= 1 utilizes a new syntax of the the T-Track and T-Yield tally outputs from PHITS. One may notice that in the old default format that there are many zero neutron flux bins in the T-Track tally and many zeros in the T-Yield tally for cases where a certain nuclide is only produced in one of numerous regions. Especially for xyz mesh and tetrahedral mesh problems, this could cause the file sizes of those tally outputs to be quite large. This new reduced format causes T-Dchain’s T-Track and T-Yield tallies to only print the absolute bare minimum of information necessary for DCHAIN, sacrificing some human readability for minimization of size and read time. Especially for problems of many regions, this dramatically reduces file sizes and improves performance in DCHAIN.

Third, IRDONCE= 1 adjusts the logic for when the T-Track and T-Yield files are closed and reopened. Previously, the T-Track output file (which is divided into sections ordered by region) would be reopened for every region and then parsed until the spectrum for the current region of interest was reached. The flux values would then be extracted, and the file would be closed and reopened from the start upon reaching the next region in the calculation. The T-Yield tally output (which is divided into sections ordered by nuclide atomic number) would be opened and parsed in its entirety for each region in the problem. The IRDONCE variable allows for these files to instead be opened only once and stepped through, where possible, eliminating redundant rereading of lines which can be significant for problems with many regions. With either the old or new reduced format of the T-Track tally, it is possible to just step through the output in this manner. However, for the T-Yield output, IRDONCE= 1 only has an effect if using IREDUFMT= 1 which orders the yields by region rather than by atomic number.

All of these performance-enhancing settings are enabled by default in the T-Dchain tally in the current distribution of PHITS. Equation (4.8.4) shows an updated version of Equation (4.8.1) when all of these new parameters are enabled. These settings remain disabled by default within DCHAIN to maintain backwards compatibility with older simulations; however, IXSRALL= 1 and IRDONCE= 1 should still be usable and provide performance improvements even for older DCHAIN calculations.

(4.8.4)\[T_{\text{total}} \propto T_{\text{rd decay lib}} + T_{\text{rd $\sigma$ lib}} + T_{\text{rd flux}} +T_{\text{rd yield}} + \sum_{i=1}^{N_{\text{reg}}} \left( N_{\text{calc steps}} \cdot T_{\text{calc}} \right) %T_{total} \propto T_{load-dcy} + T_{load-xs} + T_{rd-flx} +T_{rd-yld} + \sum_{i=1}^{N_{reg}} \left( N_{calc-steps} \cdot T_{calc} \right)\]

Understanding these relationships is important for scenarios where many DCHAIN simulations need to be conducted and optimization of the problem’s design could bear significant time savings. Otherwise, in the general use-case the simulation time is on the order of seconds (and, at most, minutes) and is typically insignificant relative to the time required for the precursory PHITS simulation.

4.8.1. Other considerations for many-regioned problems

For smaller mesh problems, these following settings do not really matter regarding file size or performance. Concerns do begin to arise with larger geometries (many regions), which can be quite common for tetrahedral mesh and xyz mesh problems. One should be very careful in designing their power cycling schedule and requested output times. Referring to the discussion in Section 4.8 on how the run time of DCHAIN scales, requesting more output times than is necessary, especially during the middle of an irradiation period, can significantly slow a calculation.

As for other performance-related or file size-related settings, here are a few in addition to the already discussed IXSRALL, IREDUFMT, and IRDONCE:

• One may wish to alter IDIVS depending on how long their irradiation time is. Maintaining the default of 50 subdivisions per irradiation time step is suggested, but if a cruder approximation is acceptable, then it can be reduced to significantly shorten run time. As noted in Equation (4.8.2) from earlier, run time almost scales linearly with IDIVS, when large enough, outside of the independent loading of libraries and parsing of PHITS output.

• JMODE is strongly recommended to be left at its default value of 2 (considering both nuclides produced directly in the T-Yield tally as well as those produced through neutron reactions with the T-Track flux spectra). However, if you are very confident that the activation due to secondary neutron reactions is negligible, JMODE can be set to 0 to then only consider nuclides produced directly through beam reactions (from T-Yield). This can significantly reduce the run time of DCHAIN.

• Setting IYILD= 0 prevents the writing of a summarized nuclide yield file which may or may not be useful to you.

• IGGRP can be changed to a gamma spectrum structure with fewer bins if you are trying to absolutely minimize the size of the .act file.

• ACMIN should almost certainly be changed from it’s PHITS default of 10\(^{-20}\) or DCHAIN default of 0 (threshold = 10\(^{-10}\) times the total current activity), but knowing what is a reasonable value to set this at can be difficult, especially if activity spans many orders of magnitude across the whole geometry.

• Setting IGSDEF= 0 prevents the writing of a MCNP-formatted gamma spectrum file which may or may not be useful to you.

• IWRTCHN should almost certainly be set to = 0 since that file is usually quite large. Only leave this enabled if you really think you need the information from this file for every single region and major time step. Use of IWRCHNUC to specify the exact nuclides of concern (rather than just listing everything) is encouraged.

• CHRLVTH should be very carefully set if you choose to leave IWRTCHN= 1.

• If you want the additional nicely formatted regionwise decay heat/photon dose output for tetrahedral problems, be sure to set FOAMOUT and FOAMVALS accordingly. These are not enabled by default.

• Likewise, if you want DCHAIN to automatically generate ANGEL plots of xyz mesh problems, be sure to set IPLOTMODE and IPLOTAXIS (ipltmode and ipltaxis in PHITS, also recognized in DCHAIN) accordingly.