Particle and Heavy Ion Transport code System, PHITS, version 2.52

Abstract

An upgraded version of the Particle and Heavy Ion Transport code System, PHITS2.52, was developed and released to the public. The new version has been greatly improved from the previously released version, PHITS2.24, in terms of not only the code itself but also the contents of its package, such as the attached data libraries. In the new version, a higher accuracy of simulation was achieved by implementing several latest nuclear reaction models. The reliability of the simulation was improved by modifying both the algorithms for the electron-, positron-, and photon-transport simulations and the procedure for calculating the statistical uncertainties of the tally results. Estimation of the time evolution of radioactivity became feasible by incorporating the activation calculation program DCHAIN-SP into the new package. The efficiency of the simulation was also improved as a result of the implementation of shared-memory parallelization and the optimization of several time-consuming algorithms. Furthermore, a number of new user-support tools and functions that help users to intuitively and effectively perform PHITS simulations were developed and incorporated. Due to these improvements, PHITS is now a more powerful tool for particle transport simulation applicable to various research and development fields, such as nuclear technology, accelerator design, medical physics, and cosmic-ray research.

Keywords:

• PHITS
• particle transport simulation
• Monte Carlo
• JENDL-4.0
• INCL4.6
• INC-ELF
• statistical multi-fragmentation model
• DCHAIN-SP
• Kurotama model
• OpenMP

1. Introduction

Currently, Monte Carlo codes for particle transport simulation in three-dimensional (3D) matter are indispensable in various research and development fields, such as nuclear technology, accelerator design, medical physics, and cosmic-ray research. Therefore, the general purpose Monte Carlo Particle and Heavy Ion Transport code System (PHITS) [1,2] is being developed through a collaboration of several institutes in Japan and Europe. It is written in the Fortran language, and derived from the NMTC/JAM code [3] in combination with HETC-CYRIC [4].

PHITS can deal with the transport of nearly all particles, including neutrons, protons, heavy ions, photons, and electrons, over wide energy ranges using various nuclear reaction models and data libraries. The geometrical configuration of the simulation must be set with either general geometry (GG) or combinatorial geometry (CG). The interactive solid modeler SimpleGeo [5] can be used for generating PHITS readable geometries in the GG format. Various quantities, such as heat deposition, track length, and production yields, can be deduced from the PHITS simulation using implemented “tally” estimator functions. The code also has a function to draw 2D and 3D figures of the calculated results and the setup geometries using an original graphic tool named ANGEL. The platforms on which PHITS can be executed are Windows, Mac, Linux, and Unix.

PHITS has several important features, such as an event-generator mode for low-energy neutron interaction [6], beam transport functions [7–9], a function for calculating the displacement per atom (DPA) [10], and a microdosimetric tally function [11]. Due to these features, it has been widely used for various applications. For example, PHITS was extensively used in the design of the shielding, target, and neutron beam lines for the J-PARC project [12–16]. Both the beam transport functions for simulating particle trajectories in magnet fields, gravity, T0-choppers, and neutron super mirrors and the DPA calculation function are very useful in the design phase. Calculation of the dose and dose equivalents in human bodies irradiated by various particles was carried out using PHITS in order to determine radiological protection needs [17–19] and medical physics issues [20]. The event-generator mode was indispensable for these calculations. The mode is also useful in the estimate of the soft error rates of semi-conductor devices [21,22]. The microdosimetric function of PHITS was used in the development of a new computational model for calculating the relative biological effectiveness (RBE)-weighted dose for charged particle therapy [23].

Recently, we upgraded the code and released the new version as PHITS2.52. It is improved from the previously released version, PHITS2.24, in terms of not only the code itself but also the contents of its package, such as the attached data libraries. The history of the initial development of PHITS has been reported previously [1,2,24–26]. Thus, this paper focuses on describing the important features of PHITS2.52 and the details of its upgraded profiles by showing the results of several benchmark calculations.

2. Physics models in PHITS2.52

Figure 1 summarizes the physics models recommended for use in PHITS2.52 for simulating nuclear and atomic collisions. The intra-nuclear cascade models JAM [27] and INCL4.6 [28] are employed for simulating the dynamic stage of hadron-induced nuclear reactions in the high and intermediate energy regions, respectively. Other intra-nuclear cascade models, such as INC-ELF [29] and a modified version of BERTINI [30], can be also selected in this region. On the other hand, the quantum molecular dynamics model JQMD [31] is generally employed for nucleus-induced reactions, because consideration of the nuclear forces between every combination of nucleons is very important in the simulation. Simulations for deuteron-, triton-, ³He-, and α-induced reactions at intermediate energies are the exceptional cases where INCL4.6 is recommended to be employed, because INCL4.6 can deal with light-ion induced reactions and is much faster than JQMD. The evaporation and fission model GEM [32] is adopted for simulating the static stage for both hadron- and nucleus-induced reactions. The energy losses of charged particles, except for electrons, are calculated using the SPAR [33] or ATIMA [34] codes with the continuous slowing down approximation.

Figure 1 Physics models recommended for use in PHITS2.52 for simulating nuclear and atomic collisions

Nuclear and atomic data libraries are generally used for simulating low-energy neutron-induced nuclear reactions and photo- and electro-atomic interactions, respectively. These data libraries are written in the ACE format, which is the same format as that adopted in MCNP [35], but different from that of original evaluated-data files, such as JENDL-4.0 [36] and ENDF/B-VII.0 [37]. Thus, it is necessary to develop PHITS-readable data libraries by compiling the original evaluated-data files using nuclear data processing systems such as NJOY99 [38]. The neutron and photon data libraries included in the PHITS2.52 package were compiled from JENDL-4.0 [36] and its sub-library based on EPDL97 [39], respectively. The electron library was developed based on EEDL [40] and ITS3.0 [41], in addition to EPDL97. Photo-nuclear reactions can be treated only for the lower energy regions in which the giant resonance is the dominant reaction mechanism.

The highest energies of particles recommended to be simulated using PHITS2.52 are 200 GeV for hadrons, 100 GeV/n for nucleus, 100 MeV for electrons and positrons, and 100 GeV for photons. Note that PHITS2.52 can deal with the transport of particles above these energies, but the accuracy of simulations for such high-energy particles has not been validated. On the other hand, the lowest particle energy that can be simulated using PHITS2.52 is 1 keV, except for neutrons, because the continuous slowing down approximation cannot be applied to the transport simulation of charged particles below this energy. Event-by-event analysis of each ionization is necessary for such low-energy particle transport simulation, but this process is not incorporated in PHITS2.52 because of its long computational time. Thus, the minimum spatial resolution reliably achieved by PHITS simulations is approximately 10⁻⁴ g/cm², e.g. 1 μm in liquid water. For neutrons, PHITS2.52 can handle transport down to 10⁻⁵ eV using the neutron data library.

3. Upgraded features of PHITS2.52

Below are listed the important upgraded features of PHITS2.52 as compared to PHITS2.24:

1.	Two intra-nuclear cascade models, Liège intra-nuclear cascade version 4.6 (INCL4.6) [28] and Intra-Nuclear Cascade with Emission of Light Fragment (INC-ELF) [29], were implemented for simulating the dynamic stage of nuclear reactions.
2.	The Statistical Multi-fragmentation Model, SMM [42], was implemented for simulating the static stage of nuclear reactions before the evaporation process [43].
3.	The Kurotama model [44] was implemented for calculating the total reaction cross sections of nucleon–nucleus and nucleus–nucleus interactions.
4.	Shared-memory parallelization using the Open-MP directives became available [45].
5.	The procedures for calculating the statistical uncertainties of tally results were modified.
6.	A function for calculating the time-evolution of radioactivity during and after irradiation was implemented by incorporating DCHAIN-SP [46] into the PHITS package.
7.	The Kerma factors included in the neutron nuclear data library were revised. Photo- and electro-atomic data libraries were newly developed. These libraries are included in the PHITS package.
8.	The electron-, positron-, and photon-transport algorithms were modified to improve the accuracy of the electromagnetic transport simulation.
9.	The photon-nuclear reaction cross sections were revised based on the JENDL photo-nuclear data file 2004, JENDL/PD-2004 [47].
10.	A function for restarting a PHITS simulation using the results obtained from a previous PHITS simulation was implemented in order to reduce the statistical uncertainties when they are not small enough.
11.	Several time-consuming algorithms, such as that for JQMD, were optimized to reduce computational time.
12.	A tally that can be arbitrarily defined by users was supplied in order to deduce user specific quantities from PHITS simulations.
13.	A function for generating knocked-out electrons, so-called δ-rays, around the trajectory of a charged particle was implemented.
14.	The radiation damage model for calculating the DPA was improved using screened Coulomb scattering [10].
15.	Various user-support tools were developed and included in the PHITS package, such as installers for Windows and Mac PCs, shell scripts for executing PHITS that can be used with drag-and-drop of icons, programs for creating an animation from tally results, and sample input files with recommended parameter settings for various calculation situations.

The details of the first eight items are described in the next sub-sections.

3.1. Implementation of the latest intra-nuclear cascade models

In PHITS2.24, a modified version of BERTINI, JAM, and JQMD can be employed for simulating hadron-induced nuclear reactions at intermediate energies. However, none of these models can reproduce double-differential cross sections for the production of high-energy light clusters, because they do not consider the formation of clusters during the cascade stage.

We therefore implemented two of the latest intra-nuclear cascade models, INCL4.6 [28] and INC-ELF [29], into PHITS2.52. Both models can consider the formation of clusters by introducing coalescence models [48,49]. They also incorporate various new features that are not taken into account in the former models, e.g., INCL4.6 incorporates the Coulomb deflection for entering and outgoing charged particles, while INC-ELF incorporates the knock-out process. In addition, INCL4.6 can be used for simulating light-ion induced reactions, even at lower energies, because it takes care of both the effects of the binding energy of the nucleons inside the incident cluster and the possibility of a fusion reaction. Considering these features, we decided to change PHITS to employ INCL4.6 as the default model for simulating nuclear reactions induced by not only hadrons but also ions lighter than α particles with energies less than 3.0 GeV/n. Note that this change occasionally induces significant increase of the computational time, because the default model previously adopted for simulating hadron-induced nuclear reaction in PHITS, the modified version of BERTINI, is extremely faster than the other intra-nuclear cascade models due to its simplified and optimized algorithm.

Figure 2 shows the double differential cross sections of Pb(p,n) and Fe(p,d) reactions calculated using PHITS2.52 employing JAM, INCL4.6, or INC-ELF. The corresponding experimental data [50,51] are also drawn in the graphs. It can be seen in panel (A) that the three models can reproduce the measured cross sections for neutron production very well. On the other hand, it is evident from panel (B) that JAM fails to reproduce the measured cross sections for high-energy deuteron productions. This tendency is due to the lack of a cluster formation mechanism in the model. The higher estimate of the deuteron-production cross sections above 350 MeV obtained using INC-ELF (compared to that obtained using INCL4.6) is probably attributed to the consideration of the knock-out process in INC-ELF.

Figure 2 Double differential cross sections of Pb(p,n) and Fe(p,d) reactions calculated using PHITS2.52 employing JAM, INCL4.6, or INC-ELF. The incident proton energies were 258 and 558 MeV for panels (A) and (B), respectively. The experimental data were taken from [50,51]

3.2. Implementation of the statistical multi-fragmentation model

After the dynamic stage of nuclear reaction, PHITS2.24 simulates the static stage using the GEM model [32], which assumes that residual nuclei release their entire excitation energies via evaporation and fission processes. However, highly excited residual nuclei occasionally release a large portion of their excitation energies by breaking into more than two pre-fragments via a multi-fragmentation process. Thus, PHITS2.24 tends to underestimate the production cross sections of intermediate-mass fragments created by high-energy nucleus–nucleus interactions, which usually results in high excitation energies for the residual nuclei. Note that the intermediate-mass fragments are those having the mass number approximately more than 10 and less than a half of the original nuclei in this paper.

We therefore implemented the Statistical Multi-fragmentation Model, SMM [42], in PHITS2.52. When the excitation energy of a residual nucleus after the dynamic stage is more than 2 MeV per nucleon, PHITS2.52 activates the SMM and determines the charges, masses, and excitation energies of each pre-fragment by referring to the entropies of 10,000 possible final states using the Monte Carlo technique. The de-excitation of the created pre-fragments down to their ground states is then simulated using the GEM model. Note that PHITS2.52 with its default setting does not activate SMM in the simulation of nuclear reactions, because the use of SMM occasionally increases its computational time by a factor of 2 or 3.

Figure 3 shows the incident-energy dependence of the production cross sections of ²⁴Na and ⁷⁵Se from Pb bombarded by carbon ions calculated using PHITS2.52 with and without the SMM. The corresponding experimental data [52] are also depicted in the graphs. It is evident from the figure that the use of the SMM dramatically improves the accuracy of PHITS simulations when calculating such intermediate-mass fragment yields. This improvement is because of the fact that such fragments are barely produced by evaporation and/or fission processes, and the multi-fragmentation process is the dominant mechanism for creating them. On the other hand, the introduction of the SMM does not change the accuracy of PHITS simulations for calculation of the production cross sections of neutrons and light fragments very much. The details of this implementation are described in [43].

Figure 3 Incident-energy dependence of the production cross sections of ²⁴Na and ⁷⁵Se from Pb bombarded by carbon ions calculated using PHITS2.52 with and without SMM. The experimental data were taken from [52]

3.3. Implementation of the Kurotama model

The total reaction cross sections of nucleon–nucleus and nucleus–nucleus interactions are used in PHITS to determine the flight path of each particle before colliding with a nucleus. In PHITS2.24, a model proposed by Pearlstein [53,54] and modified by Niita et al. [3] was employed for calculating nucleon–nucleus reaction cross sections, while models proposed by Shen et al. [55] and Tripathi et al. [56,57] were alternatively adopted for nucleus–nucleus reaction cross sections. In addition to these models, the Kurotama model was implemented into PHITS2.52 for calculating both nucleon–nucleus and nucleus–nucleus interaction cross sections. In this implementation, we extended the model to be capable of calculating the cross sections for low-energy particles on the basis of the Tripathi model.

As an example, the total reaction cross sections for ¹²C-¹²C interactions calculated using the Shen, Tripathi, and Kurotama models implemented in PHITS2.52 are drawn in Figure 4 , together with the corresponding measured data for ¹²C-C interactions [58,59]. It can be seen in the figure that the Kurotama model in PHITS2.52 can reproduce the experimental data very well, particularly for energies near a few hundred MeV/n – the energy region of importance from the viewpoint of treatment planning of charged particle therapy. It should be noted that the Kurotama model is an optional selection, and the modified Pearlstein and Tripathi models are still employed in the default settings of PHITS2.52, as was the case in PHITS2.24, for calculating the nucleon–nucleus and nucleus–nucleus reaction cross sections, respectively.

Figure 4 Total reaction cross sections for ¹²C-¹²C interactions calculated using the Shen, Tripathi, and Kurotama models implemented in PHITS2.52. Experimental data for the ¹²C-C interactions were taken from [58,59]

3.4. Implementation of memory-shared parallelization using OpenMP

There are two types of parallel computing methods; one is distributed-memory parallelization using the message passing interface (MPI) protocol, and the other is shared-memory parallelization using OpenMP. In PHITS2.24, only distributed-memory parallelization is available, but it is not suited for simulations requiring a large amount of memory, such as dose distribution calculations inside a fine-mesh voxel phantom, because the total RAM memory in the computer system must be greater than that used for single computing multiplied by the number of parallel processing elements (PE).

We therefore implemented shared-memory parallelization into PHITS2.52 by introducing OpenMP directives. The strong scaling efficiency1 for the implemented shared-memory parallelization achieved 92.4% when a typical Linux system with eight cores was used [45]. This value is sufficiently high for practical use of supercomputers, such as the “K computer” at RIKEN, because they usually consist of clusters of nodes that contain eight cores with a shared memory among them. Hybrid parallelization using both MPI and OpenMP is also available in PHITS2.52.

Another advantage of the parallelization using OpenMP is that it does not require the installation of any special software before parallel computing, while the use of MPI requires the installation of its protocol. Thus, with OpenMP, parallel computing is enabled, even on a conventional Windows or Mac PC with factory settings. Figure 5 shows the scalability of OpenMP parallelization measured using a Windows PC equipped with an Intel Core2 Quad CPU and 4 GB RAM memory. A sample input file for calculating nuclear reaction cross sections included in the PHITS2.52 package was employed in this examination, although its history number was increased to 10⁷. It is evident from the figure that the scalability reached a maximum when the number of threads was equal to the number of CPU cores in the PC, i.e., four. For higher thread numbers, the scalability decreases with an increase in the number of threads, rather than remaining constant, due to competition with accession of the shared memory. The scaling efficiencies obtained in this study are generally smaller than those obtained in a previous study using a Linux system [45], because Windows PCs cannot concentrate only on executing PHITS due to the fact that several background tasks continue to run during the execution. Note that the computational time for OpenMP parallelization using only 1 thread is approximately 1.8 times longer than that for single computing. Thus, the introduction of OpenMP parallelization is effective only when the number of cores installed in a PC is greater than 2.

Figure 5 Example of the scalability of OpenMP parallelization measured using a Windows PC equipped with an Intel Core2 Quad CPU and 4 GB RAM memory

3.5. Modification of calculation procedures for statistical uncertainties of tally results

In PHITS2.24, the statistical uncertainties of tally results were deduced from the number of events contributing to the tally. However, this method is suitable for estimating the statistical uncertainties for measurements, not for Monte Carlo simulation results. We therefore modified the code to estimate the statistical uncertainties of tally results by calculating their standard deviations. For this purpose, two alternative calculation modes were implemented in PHITS2.52: a “batch variance mode” and a “history variance mode”, which consider the batch and history as a sample of the Monte Carlo simulation, respectively. Note that PHITS bundles a number of histories as a batch to summarize the tally results, where the history indicates simulation for a source particle.

Using these modes, the standard deviation, σ, of tally results is calculated by

where N is the number of samples, i.e., the total batch and history numbers for “batch variance” and “history variance” modes, respectively, x_i and w_i are the tally results and the source weights of each sample, respectively, and and are the mean values of the tally results and source weights for all of the samples, respectively. The selection of the history variance mode generally gives better results due to its larger N, but occasionally causes an increase in the computational time. The ratio of σ to is outputted as the relative error of the PHITS simulation. This uncertainty estimation is also carried out for geometrical mesh tallies, such as the xyz-mesh tally, and the obtained relative errors in each mesh can be drawn in a 2D plot in the same manner as the tally results themselves. This function enables a visual check on the regions where the statistical uncertainties are not small enough.

3.6. Incorporation of DCHAIN-SP

The time evolution of residual radioactivity must be evaluated in the target and shielding designs of accelerator facilities, as well as in the estimates of the secondary exposure of medical workers. PHITS2.24 can be used for calculating the production yields of radioactive nuclides, but not for their time evolution during or after irradiation. Thus, the simulation code DCHAIN-SP [46] was adopted for evaluating the time evolution of the residual radioactivity from the production yields of radioactive nuclides, as well as the fluxes of neutrons below 20 MeV, calculated by PHITS. DCHAIN-SP, in combination with PHITS, was employed in the target design of the J-PARC project [12,16]. However, some special knowledge is required to prepare the input files of DCHAIN-SP from the output files of PHITS, which is not easy for PHITS users to learn.

We therefore implemented a new tally named [T-DCHAIN] into PHITS2.52 in order to automatically create the input files for DCHAIN-SP from PHITS simulation results. The irradiation and cooling times used in the DCHAIN-SP calculation can also be specified in the tally. In addition, the source, executable, and library files for DCHAIN-SP, together with the shell scripts for successively executing PHITS and DCHAIN-SP, are included in the PHITS2.52 package. Using these items, PHITS users can easily calculate the time evolution of residual radioactivity. As an example, Figure 6 graphically presents the time evolution of the residual radioactivity produced in a water target irradiated by 150 MeV proton beam for 6 minutes as calculated using PHITS2.52 in combination with DCHAIN-SP and the [T-DCHAIN] tally. It can be seen from the graph that the total radioactivity decreased by an order of magnitude 10 minutes after irradiation. This type of information is very useful for avoiding unnecessary secondary exposure of medical workers in proton therapy facilities.

Figure 6 Calculated time evolution of the residual radioactivities produced in a water target irradiated by a 150 MeV proton beam for 6 minutes. This calculation was performed using PHITS2.52 in combination with DCHAIN-SP and the [T-DCHAIN] tally. The lines are guide for the eye

3.7. Improvements of nuclear and atomic data libraries

Neutron-nuclear and photo- and electro-atomic data libraries are indispensable for transport simulations of low-energy neutrons, as well as photons, electrons, and positrons, in PHITS. However, no atomic data library was included in the PHITS2.24 package. Thus, PHITS users had to obtain those libraries by themselves from other code packages, such as MCNP [35]. We therefore developed new atomic data libraries and included them in the PHITS2.52 package. The photon data library was compiled from a sub-library of JENDL-4.0 [36] that is based on EPDL97 [39] for 1 keV up to 100 GeV. The electron data library was developed based on EEDL [40] and ITS3.0 [41], in addition to EPDL97, for 1 keV up to 10 GeV. However, the accuracy for calculating the bremsstrahlung spectra generated by electrons over 100 MeV may not be precise enough, because EEDL does not contain data between 10 MeV and 10 GeV, and we simply interpolated those data on a logarithmic scale in the high-energy region. Thus, the highest energy of electrons and positrons recommended to be simulated using this electron data library is 100 MeV.

The Kerma factors for most nuclei included in the neutron data library attached to PHITS2.24 were determined by processing JENDL-4.0 using the energy balance method [60] in NJOY99 [38]. In the energy balance method, the Kerma factors are calculated by subtracting the average energy of secondary neutrons and photons from the sum of the incident neutron energy and the mass-difference Q value of the reaction. Thus, the rigorous energy balance must be established in the nuclear data library applied to this method, i.e., the evaluation of the complete dataset of neutron and photon production cross sections are needed. However, JENDL-4.0 does not always contain the all data required for this method; e.g. the photon production cross sections for ³⁵Cl are not equipped in JENDL-4.0. Hence, the Kerma factors for some nuclei included in our neutron data library were doubtful and needed to be revised. We therefore modified such doubtful Kerma factors using the kinematic limits method [60] in NJOY99 instead of the energy balance method. The kinematic limits method can give a reasonable approximation of the Kerma factor from its maximum limit, which can be calculated from the available cross-section data, considering the conservation of momentum and energy.

As an example, the Kerma factors for ³⁵Cl included in the neutron data libraries attached to PHITS2.52 and PHITS2.24 are depicted in Figure 7 , together with the corresponding data obtained from ENDF/B-VII.0 [37]. The figure clearly indicates the overestimation of the Kerma factors in our previous library, because the library does not contain all the data required for the energy balance method as described before. Hence, absorbed doses inside a human body exposed to low-energy neutrons calculated by PHITS2.24 significantly depend on the density of ³⁵Cl, although it comprises only approximately 0.1% of the human body by weight. However, ³⁵Cl is an exceptional nucleus, and the Kerma factors for most nuclei were not dramatically changed by this revision.

Figure 7 Kerma factors for ³⁵Cl included in the neutron data libraries attached to PHITS2.52 and PHITS2.24, together with the corresponding data obtained from ENDF/B-VII.0 [37]

We also revised the cross-section data for 19 nuclei based on the updated files for JENDL-4.0, which was released in 2012. Furthermore, thermal neutron scattering law data libraries were developed for 15 frequently used materials, such as hydrogen in water, using the corresponding sub-library of JENDL-4.0, whose numerical values are nearly the same as those of ENDF/B-VI [61]. These libraries are also included in the PHITS2.52 package.

3.8. Modification of electron, positron, and photon transport algorithms

In PHITS2.24, the range of an electron largely depends on its cut-off energy, because the energy losses due to production of secondary electrons with energies above its cut off are double-counted in the calculation of the effective stopping power used in the electron transport simulation. We therefore modified the code to adjust the effective stopping power by excluding the average energy losses due to the production of the higher-energy electrons, in the same manner as that used in EGS5 [62]. Due to this modification, the range of an electron calculated using PHITS is now independent of its cut-off energy. We also modified the code to conserve the energy in the event of photon-atomic interactions, such as the photo-electric effect. This modification is very important for the calculation of the response of a scintillator irradiated by low-energy photons.

Figure 8 shows calculated depth-dose distributions at the central axis of a water cylinder uniformly irradiated by a 10 MeV electron beam using PHITS2.52 and PHITS2.24. In the calculations, the cut-off energy of the electrons and positrons was set to 10 keV or 1 MeV. The corresponding experimental data [63] are also plotted in the graph. It is evident from the graph that the results obtained using PHITS2.52 are in excellent agreement with the experimental data, regardless of the electron cut-off energy. On the other hand, PHITS2.24 reproduced the experimental data only when the electron cut-off energy was set to 1 MeV. This difference arises because the energy losses due to production of secondary electrons with energies above 1 MeV are negligibly small in comparison to the total stopping power, and the double-count problem in PHITS2.24 does not cause a serious influence on the electron transport simulation for that case. Thus, the cut-off energy of the electrons and positrons must be set to 1 MeV or higher in PHITS2.24, while it can be decreased down to 1 keV in PHITS2.52. Note that the recommended value for the cut-off energy is 100 keV in PHITS2.52, except for simulations requiring very high spatial resolution, because the computation time becomes significantly longer with a decrease in the cut-off energy.

Figure 8 Calculated depth-dose distributions at the central axis of a water cylinder uniformly irradiated by a 10 MeV electron beam using PHITS2.52 and PHITS2.24. The cut-off energy of electrons and positrons, E _C, was set to 10 keV or 1 MeV. The experimental data are taken from [63]

4. Conclusions

We upgraded many aspects of the PHITS code and released the new version as PHITS2.52. In the new version, higher accuracy of the simulation was achieved by implementing the latest nuclear reaction models (INCL4.6, INC-ELF, SMM, and the Kurotama model) and revising the attached neutron data library. The reliability of the simulation was also improved by modifying both the algorithms for electron-, positron-, and photon-transport simulations and the procedure for calculating the statistical uncertainties of the tally results. In addition, incorporation of the original atomic data libraries and DCHAIN-SP into the new package enables PHITS users to simulate the electromagnetic transport and time evolution of radioactivity without the need to obtain other code packages. Furthermore, the efficiency of the simulation was improved due to the implementation of shared-memory parallelization together with optimization of several time-consuming algorithms. PHITS2.52 has been distributed to many countries via the Research organization for Information Science and Technology (RIST), the Data Bank of the Organization for Economic Co-operation and Development's Nuclear Energy Agency (OECD/NEA DB), and the Radiation Safety Information Computational Center (RSICC).

Several tasks necessary for further PHITS development still remain. The incorporation of the EGS5 code and higher-energy photo-nuclear reaction mechanisms, such as quasi-deuteron photo-disintegration, are currently in progress for improving the accuracy of electron-, positron-, and photon-transport simulations, particularly for higher energies. A preliminary version of PHITS in combination with EGS5 is included in the PHITS2.52 package, although many new features implemented in PHITS2.52 are not available in that version. Mechanisms for electron- and muon-induced nuclear reactions are also planned to be implemented. Improvements in the nuclear reaction models, particularly for lighter systems such as p-Li reactions, are also necessary, because the models currently implemented in PHITS cannot consider the shell effect, which plays an important role in simulating nuclear reactions involving lighter nuclei. For this purpose, development of a new nuclear reaction model has been initiated by combining an intra-nuclear cascade model and a distorted-wave Born approximation calculation [64]. This model will be incorporated in a future version of PHITS.

Acknowledgements

Apart from our own efforts, the development of PHITS2.52 was extensively supported by many researchers and technicians from all over the world. We wish to thank Dr. J. Cugnon of University of Liège, and Dr. D. Mancusi, Dr. A. Boudard, Dr. J.C. David, and Dr. S. Leray of CEA/Saclay for their help in implementing INCL4.6. We also wish to thank Dr. Y. Uozumi, Mr. Y. Sawada, and Mr. S. Nogamine of Kyushu University for their support in implementing INC-ELF. We would like to express our gratitude to Dr. A. Kohama of RIKEN, Dr. K. Iida of Kochi University, and Dr. K. Oyamatsu of Aichi Shukutoku University for their efforts on implementing the Kurotama model. The implementation of the shared-memory parallelization was supported by the Next-Generation Integrated Simulation of Living Matter, Strategic Programs for R&D of RIKEN, and RIKEN Special Postdoctoral Researchers (SPDR) Program. For this improvement, we used the K computer and RIKEN Integrated Cluster of Clusters (RICC). Programming for the modification of the procedures for calculating the statistical uncertainties and the implementation of the restart calculation function were mainly performed by Mr. D. Obinata of Fujitsu Systems East Limited and Mr. K. Sakamoto of JAEA with the support of Center for Computational Science & e-Systems of JAEA. We wish to thank Dr. H. Hagura of JAEA for his efforts on developing the electro-atomic data library. We also wish to thank Dr. Y. Namito and Dr. H. Hirayama of KEK for their advice on modifying the electron and photon transport algorithms and developing the preliminary version of PHITS in combination with EGS5. We would like to express our gratitude to Dr. K. Batkov of the European Spallation Source for his efforts on the development of programs for creating animations from PHITS tally results, and to Dr. C. Theis of CERN for his efforts in implementing a function for generating PHITS-readable geometry into SimpleGEO.

Notes

1. The ratio of total computational time for each core to execute the same job by single and parallel computing.