Discrepancy between RESNUCL and Fluence Based Activity

Versions

Please provide the used software versions.

FLUKA:4-5.2
Flair:3.5-0a

Description

Please describe the issue. In case of reporting a bug, please provide step-by-step instructions for reproducing the error.

Dear FLUKA-Experts,

I am currently investigating a discrepancy between the results of the induced radioactivity calculated from FLUKA and some offline-fluence based calculations methodes. The main workflow involves the calculation of the particle fluences within the components of the irradiation station and then calculating the induced activity with ActiWiz3-Creator. These results are then compared to FLUKA RESNUCL scoring, to compare the two calculation methodes. Due to the large number of components and nuclides scored, the results were first compared based on the ratios of total induced activity:

Ratio = \frac{A_{tot, ActiWiz}}{A_{tot, FLUKA}}

During the comparison of the results, the overall agreement between the codes was quite good for most nuclides and components within the expected variance of the nuclear data and the specific methodes used for the calculation. Isomeres were a bit harder, but I do know that FLUKA handles isomeres by populating each isomeric state equally and then decaying those states, while ActiWiz uses branching ratios from nuclear data.

However for some components the Ratio between the total induced activity was quite large. Specifically for the natural Molybdenum Plate (density = ) in our irradiation setup, the ratio grew to about 1.93 for the overall induced activity. For the dominant nuclide, Mo-99, the induced activity was:

ActiWiz3: 5.91e13 Bq

FLUKA: 2.68e13 Bq

So a ratio of almost 2.2. Therefore suggesting an over/underestimation of the neutron capture within the codes. This is also reinforced by the overestimation if Tc-101 by ActiWiz compared to FLUKA:

The code for the comparison had to be a bit longer, because the RESNUCL results are printed out at 1s after irradiation and the ActiWiz results were calculated for 1h after irradiation. The decay step was done with the python package radioactivedecay. The image above shows the decay corrected data, so both results are shown after 1h of decay. Furthermore are both results scaled to the entire plate with the dimensions of 13\times24\times3\text{ cm}^3 and a density of 10.22 \frac{g}{cm^3}, since the USRTRACK scorings are based on the results with this density assigned to the Molybdenum.

I know that FLUKA by default uses 10.22 \frac{g}{cm^3} for the density of Molybdenum, while ActiWiz uses 10.28 \frac{g}{cm^3}.

Due to this discrepency I tried to manually calculate the expected activity based on the available data and the particle fluences calculated with FLUKA. The particle fluences (attached below) are postprocessed to yield particles/s/\text{cm}^2.
Within FLUKA I use the JEFF-3.3 pointwise library, and ActiWiz was also build on the JEFF-3.3 data (and EAF and ENDFB/VIII0), as far as I know. The dominant production pathways for Mo-99 should be the neutron capture of the stable Mo-isotopes and the (n,2n)-reaction of Mo-100. The (\gamma/n)-reaction only has neglectable contributions within our setup (<0.00086%). Therefore I wrote small code solving the matrix exponential for the Bateman-equations only considering the Mo-92 to Mo-100. The effective cross sections for the transition matrix were calculated with the rebinned neutron fluence taken from the FLUKA USRTRACK scoring and then folded with the 709 Groups CCFE JEFF-3.3 data. Additionally I calculated the same for the other libraries JEFF-3.2, JEFF-4.0 (1102 Groups CCFE) and EAF. Furthermore I did the same with JEFF-3.3 pointwise data and linear interpolation between the energies, but I would say this is the one I am least confident with. I also got a license for FISPACT, which should be a more complex version of my groupwise matrix exponential calculations, to check, if I missed something:

Manually (JEFF-3.3 709 CCFE): 5.65e13
Manually (JEFF-3.3 Pointwise): 4.9e13
Manually (JEFF-3.2 709 CCFE): 5.40e13
Manually (EAF2010): 5.78e13
Manually (JEFF-4.0 1102 CCEF): 5.21e13
Manually (ENDFB/VIII.0 709 CCFE): 5.24e13
FISPACT (JEFF-3.3 709 CCFE): 5.72e13

All these calculations rather agree with ActiWiz3 than with FLUKA. I also checked the production contribution by energy within ActiWiz and most of the production happens below 20 MeV and predominantly by neutrons.

So FLUKA should be using data from JEFF-3.3, just as ActiWiz does and should not be based on the FLUKA models used above 20 MeV. The manual calculations above suggest, even though i wouldn’t entirely bet on them, that the discrepancy is also not entirely based on the pointwise vs groupwise data.

My question would be, if there is any idea on what else might be the cause for this disparity?

Thank you in advance!

Input files

Please upload all relevant files. (FLUKA input file, Flair project file, user routines, and data files)

Flux_NEUTRON_MoPlate.txt (51.8 KB)
Flux_PROTON_MoPlate.txt (51.8 KB)
Flux_PHOTON_MoPlate.txt (51.8 KB)

99MoBestV12.flair (280.8 KB)
99MoBestv12.inp (107.9 KB)

Dear @David_Ohm

Thank you for your report and for the careful analysis. The comparison you set up is useful, but the fluence-folding approaches (both the ActiWiz estimate and the simplified folding of the USRTRACK spectrum with σ(n,γ)) miss one important effect: resonance self-shielding in the thick Mo target. Your own second plot shows it, since most of the Mo-99 production comes from the resonance region.

The capture probability in a step of length x_i​ is 1−exp(−σN⋅x_i), and it is approximated to σN⋅x_i only when the exponent is small. A fold of the cumulative track length with a cross section, Σx_i⋅σN, implicitly makes this replacement, so it is valid only for thin targets or small cross sections. At the resonance peaks the exponent is not small and the true capture probability saturates leading to errors. On thick targets the accumulation of these tiny errors is amplified, leading to a production overestimation with the fold.
In addition, the binning effect of USRTRACK returns bin-averaged flux which cannot fully resolve the depression of the flux inside inside a resonance 1/Σ_t(E) , strongly diminishing its effect.
The folding has to be done at the step level, as the MC does, with σ(E) at the energy of each step and not on the energy-binned cumulative track length.

I checked this on a Mo-98-only target, comparing the point-wise (PW) and group-wise (GW) treatments:

  • PW: the cross sections used at the collision energies match the evaluated data. The capture sampling is consistent with σ_c/σ_t and plain fold is about a factor 2 too high.
  • GW: the group-averaged cross sections carry no self-shielding, so GW matches the plain fold and overestimates Mo-99 by about a factor 1.9 at full density.
  • Dilution test: I repeated the runs with only Mo-98 as target lowering its density ×100 and by ×10⁴. The GW/PW ratio fell from 1.9 to 1.55 and then to about 1.0. At ×10⁴ dilution, PW, GW and the plain fold agree within the statistical error (about 15%). So the discrepancy follows the self-shielding and is not a difference in the cross-section data.

The FLUKA result with PW is the one to trust here. The fluence-based estimates, and GW, are approximations that overestimate the production in this case. Resonance self-shielding, among other things, was one of the reasons for implementing the point-wise treatment in FLUKA.

For fluence-based activity estimates in thick, strongly absorbing targets, please use RESNUCLE with PW. A plain fold of the flux with σ(n,γ) is only valid for thin or dilute targets

Dear @vasilis,

thank you very much for this clear and detailed explanation!

I was thinking about self-shielding at one point, but came to the wrong conclusion after looking at the plot of the neutron spectrum. Sometimes the first thought is the right one.

I falsely assumed in that moment, that, since the dips at resonance energies were clearly visible, they would be accounted for when using this spectrum in the folding method. But its the other way around.

Maybe a quick followup question towards your suggestion in the end. When assessing the effects of material-composition changes, such as variations in alloying additions or trace impurities, is it advisable to calculate production yields “online” using separate FLUKA simulations with RESNUCLE, specifically for targets containing a strongly absorbing major element? Since the self-shielding effects towards the neutron spectrum, should also affect the minor elements.

Kind regards
David

Hi @David_Ohm

for accurate calculations, yes MC is the best way to go. The other methods can provide you an approximate result which could be always be a good starting point.

However with the MC you might struggle to collect statistics.
One option could be that instead of relying on the RESNUCLEi scoring you do online weighting per each step, (routine fluscw.f on a USRTRACK activated with the USERWEIGht) with the expected reaction probability on the isotope you are interested and the output USRTRACK spectra will be the prompt production yield. You could even use the neutron cross sections directly from FLUKA to do the weighting.

Some small comments on your input:

  1. You have too many scoring cards USRBIN, USRTRACK. Using too many can significantly penalize the CPU tracking time.
  2. You could be more harsh on the EMF thresholds so save some CPU time

Dear @vasilis,

thank you very much once again for the clarification and the suggestions. I will look into the online method with the user routine, I think there was some demonstration in the slides of the topical or advanced course.

Cheers
David