APPLICATION OF THE JSIR2S CODE PACKAGE FOR SHUTDOWN DOSE RATE CALCULATIONS ON JET

In this paper we present a computational exercise for shut-down dose rate calculations for the JET tokamak using the in-house developed JSIR2S code package as part of its validation. The computation is performed in two parts: neutron transport and transport of secondary gamma radiation. In order to calculate neutron activation reaction rates with sufficiently low variance, hybrid variance reduction techniques using the ADVANTG code have been utilized. Probability based sampling of secondary source particles was performed. Calculated gamma dose rates after shut down are compared with dose rate measurements performed on site using ionization chambers. The C/E agreement for 1 octant is between 0.8 to 1 while statistically meaningfull results for the 2 octant are yet to be obtained.


INTRODUCTION
The JSIR2S code package for simulation of radiation due to radioactive decay of neutron activated and fission produced radioactive nuclides has been under development at the "Jožef Stefan" Institute for the last 2 years. It is based on the Rigorous-2-step computational framework (R2S) [1] on a rectangular mesh. It couples together Monte Carlo particle transport code MCNPv6.2 [2], isotopic inventory calculation code FISPACT-II [3] and a code for calculating electron and positron spectra from β decay BetaShape [4] with a series a Python scripts. A typical JSIR2S workflow is displayed in Figure 1 and discussed in the following points: • A standard MCNP geometry model is used as a basis computational model for neutron and secondary radiation particle transport. An additional user supplied input file is used where the following is defined: -Rectangular mesh parameters for division of the MCNP model geometry and list of cells included into calculation (all by default).
-Irradiation history timeline and neutron flux intensities.
-Decay radiation evaluation times.
-Selection of nuclear data libraries (separately for transport and depletion).
-Parameters on the calculation precision and utilization of computer resources.
• Generation of internal data library generated from FISPACT-II library summaries on decay radiation spectra and intensities using exact line energies. Continuous energy spectra are generated for electron and positron energy spectra by beta decay using the BetaShape code.
• A parser extracts the initial isotopic inventory and geometry cell densities from the MCNP model geometry. MCNP model cells are divided by the rectangular mesh. Due to irregular shapes being generated by intersection of a superimposed mesh and geometry cells, their volumes are calculated stochastically up to an absolute precision by setting all model materials to vacuum, cell volumes to 1, setting up a source plane with source particles emitted in parallel and calculating the volume averaged particle flux. These vales are multiplied by the average ray surface area to calculate cell volumes of computational units. The average ray surface area is given by S ray = S vol−surf /nps, where S vol−surf is the plane source area used for volume calculations, and nps number of source particles.
• Each MCNP model cell inside a mesh voxel is treated as a separate computational unit in terms of both activation and as a secondary radiation source [5], commonly known as "Cell under voxel approach". Previously determined volumes are used for calculation of volume averaged neutron flux calculations. Transport of activating particles (usually neutrons) is performed, using several inputs, provided activating cells remain the same throughout. Primary spectra, intensities and consequently reaction rates are calculated using volume-averaged estimators. An example of separate treatment of cells under voxel as compared to the regular mesh is presented in Figure 2.
• Inventory calculations from which decay spectra are calculated from previously generated libraries. Secondary particle source files are generated.
• Transport of secondary particles using uniform sampling with weight adjustment or by probability.   The code package has been thus far validated on gamma field measurements at the JSI TRIGA reactor [6,7] and on a ITER Port-plug computational benchmark [8].
In this paper we present computational analysis and the work flow of shut-down dose rate calculations (SDR) on JET tokamak using the JSIR2S following the the entire operational history, compared to measurements performed in the framework of the Work Package JET3 NEXP [9,10].

JSIR2S WORKFLOW ON JET SDR CALCULATIONS
For this particular case 2 reference computational MCNP models of 1 8 of the JET tokamak were provided by EUROfusion [11] (Figures 3a and 3b). Each octant was divided by a 100×120×160 uniform mesh for volume and delayed radiation source term calculations. Additional empty box surrounding the model was included and reflective surfaces removed in order to perform stochastic volume calculations. Volume calculation precision was set to 1 × 10 −2 cm 3 .  Due to outer torus radius being in excess of 10 m, a hybrid variance reduction technique FW-CADIS using the ADVANTG code [12] was used in order to have sufficiently low neutron flux variance on the above mentioned mesh throughout the model in VITAMIN-J 175 energy group structure [13]. A rectangular weight window mesh with 103×170×309 voxels was used.

Ionization chamber
Since MCNP formatted fixed-source description is required for the ADVANTG calculation, a simplified neutron source description using a collection of point sources was modelled in for weight window generation. These weight window parameters were calculated using the 27n19g nuclear data library based on the ENDF/B-VII.0 [14] and were used for the actual transport using the JET plasma neutron source routine. A FENDL 2.1 [15] library, made specifically for fusion applications was used for particle transport calculations and JEFF 3.2 [16] was used for isotope evaluation missing in FENDL 2.1. Comparison of neutron fluxes with weight windows versus analogue Monte Carlo calculation is displayed in Figure 4 as described by Equation 1, where R is the relative comparison, V analog and V V R are analogue and FW-CADIS based weight-window variance reduction accelerated results. The agreement is within the 1σ statistical uncertainty of the calculation, observing no bias in our approach. The figure of merit (FOM) using the weight windows was about 1.6.
Neutron fluxes were calculated in ≈2 × 10 5 units (cells under voxels) and normalized according to JET neutron emission intensities. EAF-2010 [17] nuclear data libraries were used for calculations of isotopic inventory. Dose rates were calculated for 6 h and 12 d after shut down inside ionization chambers using the ICRP74 [18] flux to H * 10 biological dose equivalent and air kerma conversion factors. Secondary particle source term can be sampled in two ways: • Probability based sampling: Emission position and emission energy are sampled by their probability, using Vose aliasing algorithm. This approach is useful for calculations of dose, total flux and other integral quantities close to the most active materials.
• Uniform based sampling: The computational units and secondary particle emission energies are sampled uniformly and the particle weight is adjusted by their emission probability. This approach is useful for spectrum calculations and for calculating integral values some distance from most active materials.
Choosing the most appropriate secondary source particle sampling technique is at this time decided by the code package user. In the future, various statistical tests will be included as a part of the code package, on basis of which a more informed selection of secondary source sampling technique can be done. At the present time only probability based shut down dose rate calculation results are presented.

JSIR2S RESULTS AND COMPARISON WITH MEASUREMENTS
In this section a comparison of dose rates calculated by the JSIR2S code package as compared to the measurements and calculations by other codes [19] is presented. An attempt on using ADVANTG code for transport of the decay γ radiation using the point source description with first collision estimates was performed but was deemed to computationally expensive for practical use. Hence a limited number of analogue calculations is presented and possibilities on accelerating the calculations discussed in the next section. Results for the 1 st octant model are presented in Table 1 with dose rate C/E ratios ranging from 0.8 to 1.0, comparable to other code packages [19] but with C/E agreement having a different trend. This discrepancy should be investigated further.
Although calculated dose uncertainties are still high even for long simulations, the detector is relatively exposed to the most activated components and reasonable results were obtained using analogue Monte Carlo transport of the decay radiation. Decay gamma air kerma field, 6 h after shut-down is displayed in Figure 5.   Results for the 2 nd octant model are presented in Table 2 with dose rate C/E ratios ranging from 0.4 to 0.2, which is further off compared to other R2S methodology codes [19]. These results are preliminary, and discussion on their improvement by applying variance reduction techniques is addressed in the following section. Modeling discrepancies of the 2 nd octant model as described by [20] contribute to poor agreement in general, even when using variance reduction techniques for the decay gamma field.

FUTURE DEVELOPMENTS
Since only probability based sampling was used, we suspect that the less intense sources in closer proximity to the ionization chamber were under-sampled. The probability sampling also neglects the less intense spectral lines. In present work decay gamma transport was not accelerated by the ADVANTG code due to computational requirements implied by the use of custom MCNP source routines. Currently developments are under way on using novel variance reduction techniques such as responsematrix [21] and MAGIC method [22], with FOM comparable to the FW-CADIS method but are not limited to the use of discretized geometry, multi-group nuclear data libraries and standard source descriptions. Additional calculations for both octant models will be performed using these methods. Several modeling inaccuracies [20] have been identified as main contributors to poor computational agreement for the 2 nd octant model nd are being updated and corrected by the JET community. The entire shutdown dose rate calculation will therefore be performed again on an updated model. In current calculations only EAF2010 nuclear data libraries were used for calculations of activation and fission reaction rates and isotopic inventory calculation. The JSIR2S supports the use of additional nuclear data libraries such as JEFF 3.3 [23] and ENDF/B VII.1 [24] up to 1 GeV. Analysis using the above mentioned libraries for isotopic inventory calculations if also planned for the future.

CONCLUSIONS
The JSIR2S code has been applied for calculations of the JET shut down dose rates. The results for the 1 st octant of the model are in relatively good agreement with the measurements. The agreement of calculated dose rates with measurements is not as good for the 2 nd octant of the JET tokamak model and further investigation is required.

ACKNOWLEDGMENT
The authors would like to thank the JET and EUROfusion communities for providing the necessary data for testing and further developments of the JSIR2S code package for fusion applications.