Transmutation results different from FISPACT-II

Validation and Verification
1

Hello all!

I’ve been doing transmutation calculation with Rhodium with OpenMC and comparing them with known results from FISPACT-II. I cannot post my code because I’m a new user but I will add as text at the end. Nevertheless, I’m getting a much lower result for the same neutron energy spectrum and neutron flux.

Isotope OpenMC FISTPACT-II
Rh103 - 1.18 % -9%

In my simulation I use an isotropic neutron source with a particular energy spectrum, this source is at the center of an 1.5 cm rhodium sphere so all the neutrons are generated within it and interact with the material.

Also, the FISTPACT-II says that almost 90% becomes Pd (Pd102+Pd104) and my OpenMC simulation has this value around 54% with most of the rest, 44.6%, becoming Ru102.

Is this an error of how I setup the simulation, should I tally it in another way? Or is this an difference between these codes?

Thanks!

My code

import openmc as op
import openmc.deplete
import numpy as np
import math
import matplotlib.pyplot as plt

chain_filename = "/home/tcsousa99/openmc_projects/Rh-transmutation/chain_endfb80_pwr.xml"
openmc.config['chain_file'] = chain_filename
chain = openmc.deplete.Chain.from_xml('/home/tcsousa99/openmc_projects/Rh-transmutation/chain_endfb80_pwr.xml')

# Define rhodium material
radius = 1.5

rhodium = op.Material(name='Rhodium')
rhodium.add_element('Rh', 1.0, percent_type='ao')
rhodium.set_density('g/cm3', 12.41)
volume_sphere = 4/3 * math.pi * math.pow(radius, 3)
rhodium.volume = volume_sphere
rhodium.depletable = True

materials = op.Materials([rhodium])
materials.export_to_xml()


sph1 = openmc.Sphere(r=radius, boundary_type='vacuum')


# Define cell containing Rhodium
source_cell = op.Cell(fill=rhodium, region=-sph1)
#universe = op.Universe(cells=[source_cell])
geometry = op.Geometry([source_cell])

energy_values = [
    1.2e-01, 1.3e-01, 1.4e-01, 1.5e-01, 1.7e-01, 1.9e-01, 2.1e-01, 2.4e-01, 2.6e-01, 2.9e-01,
    3.2e-01, 3.6e-01, 4.0e-01, 4.6e-01, 5.6e-01, 7.1e-01, 8.8e-01, 1.12e+00, 1.46e+00, 1.88e+00,
    2.4e+00, 3.06e+00, 3.89e+00, 4.88e+00, 6.09e+00, 7.67e+00, 9.71e+00, 1.214e+01, 1.532e+01, 1.913e+01,
    2.357e+01, 2.865e+01, 3.485e+01, 4.397e+01, 5.733e+01, 7.2e+01, 8.911e+01, 1.1044e+02, 1.361e+02, 1.6592e+02,
    2.0674e+02, 2.5577e+02, 3.0997e+02, 3.6947e+02, 4.2646e+02, 5.129e+02, 6.3026e+02, 7.8001e+02, 9.8166e+02, 1.24392e+03,
    1.58295e+03, 1.89701e+03, 2.24678e+03, 2.84736e+03, 3.60268e+03, 4.69139e+03, 6.0916e+03, 7.86558e+03, 9.87871e+03, 1.253548e+04,
    1.654332e+04, 2.058981e+04, 2.35302e+04, 2.533669e+04, 2.919571e+04, 3.514691e+04, 4.284074e+04, 5.182769e+04, 6.01786e+04, 6.789949e+04,
    7.79768e+04, 9.431335e+04, 1.120268e+05, 1.290276e+05, 1.527468e+05, 1.831146e+05, 2.17507e+05, 2.529592e+05, 2.813009e+05, 3.183361e+05,
    3.908826e+05, 4.759756e+05, 5.582925e+05, 6.221372e+05, 6.551776e+05, 7.348871e+05, 8.637569e+05, 1.002912e+06, 1.190408e+06, 1.398369e+06,
    1.649158e+06, 1.895047e+06, 2.119077e+06, 2.400841e+06, 2.670828e+06, 3.137938e+06, 3.843547e+06, 4.635082e+06, 5.365174e+06, 6.106454e+06,
    7.063861e+06, 8.051461e+06, 8.895515e+06, 9.495932e+06, 9.968319e+06, 1.088905e+07, 1.181362e+07, 1.218616e+07, 1.244237e+07, 1.261551e+07,
    1.274771e+07, 1.292576e+07, 1.310631e+07, 1.328937e+07, 1.352221e+07, 1.395315e+07, 1.439316e+07, 1.464451e+07, 1.479733e+07, 1.495234e+07,
    1.510843e+07, 1.521358e+07, 1.537296e+07, 1.564230e+07, 1.597117e+07
]

probability_values = [
    5.75e-12, 6.19e-12, 9.27e-12, 1.39e-11, 2.43e-11, 1.61e-11, 2.86e-11, 4.81e-11, 4.20e-11, 7.36e-11,
    9.97e-11, 1.55e-10, 1.60e-10, 3.50e-10, 6.06e-10, 8.47e-10, 9.64e-10, 1.61e-09, 2.18e-09, 2.81e-09,
    3.50e-09, 4.72e-09, 5.73e-09, 6.71e-09, 9.07e-09, 1.22e-08, 1.63e-08, 1.85e-08, 2.80e-08, 2.85e-08,
    3.85e-08, 3.89e-08, 5.71e-08, 8.61e-08, 1.27e-07, 1.11e-07, 1.72e-07, 1.84e-07, 2.48e-07, 2.56e-07,
    4.38e-07, 4.02e-07, 5.39e-07, 5.08e-07, 5.10e-07, 1.07e-06, 1.12e-06, 1.76e-06, 2.22e-06, 3.11e-06,
    4.03e-06, 2.79e-06, 4.70e-06, 7.84e-06, 8.60e-06, 1.67e-05, 1.63e-05, 2.70e-05, 2.35e-05, 4.45e-05,
    6.32e-05, 5.30e-05, 4.16e-05, 2.26e-05, 1.03e-04, 9.07e-05, 1.75e-04, 1.50e-04, 1.70e-04, 1.44e-04,
    2.94e-04, 4.31e-04, 3.91e-04, 4.88e-04, 8.66e-04, 9.41e-04, 1.22e-03, 1.17e-03, 7.76e-04, 1.78e-03,
    3.97e-03, 3.37e-03, 3.69e-03, 1.31e-03, 1.01e-03, 5.78e-03, 5.37e-03, 6.20e-03, 8.39e-03, 6.95e-03,
    1.06e-02, 5.65e-03, 7.07e-03, 6.90e-03, 4.52e-03, 1.12e-02, 1.11e-02, 1.30e-02, 8.26e-03, 1.07e-02,
    1.25e-02, 1.20e-02, 1.01e-02, 7.69e-03, 9.51e-03, 2.61e-02, 9.63e-03, 6.14e-03, 7.67e-03, 3.46e-03,
    1.20e-02, 1.76e-02, 2.47e-02, 4.06e-02, 1.10e-01, 3.12e-01, 1.35e-01, 6.06e-02, 1.70e-02, 1.21e-02,
    4.44e-03, 2.33e-03, 1.17e-03, 5.40e-04, 2.82e-04
]


# Define the source
source = op.IndependentSource()
source.space = op.stats.Point((0,0,0))
source.angle = op.stats.Isotropic()
source.energy = op.stats.Discrete(energy_values, probability_values)
source.particle = 'neutron'

# SETTINGS
settings = op.Settings()
settings.batches = 4
settings.inactive = 0
settings.particles = 10000
settings.source = source
settings.run_mode = 'fixed source'

model = op.model.Model(geometry, materials, settings)

time_steps = [30]*6
flux_neutrons = 3.9e14 #in cm^-2 s^-1

#flux times surface area of the sphere since this is the considerable area and all neutrons that start isotropically are going outward to the sphere
rate_per_second = flux_neutrons * 4 * math.pi * math.pow(radius,2) 
source_rates = [rate_per_second]*6 #this rate per second is considered constant during the 30 day period

dpl = openmc.deplete.CoupledOperator(model, chain_file= chain_filename, normalization_mode='source-rate', reduce_chain = True, reduce_chain_level = 3)

integrator = openmc.deplete.PredictorIntegrator(dpl, time_steps, source_rates, timestep_units='d')

integrator.integrate()


results = openmc.deplete.Results("depletion_results.h5")

times, number_of_Rh104 = results.get_atoms(rhodium, 'Rh104', time_units="d")
_ ,number_of_Rh103 = results.get_atoms(rhodium, 'Rh103', time_units='d')
_ ,number_of_Rh102 = results.get_atoms(rhodium, 'Rh102', time_units='d')
_ ,number_of_Ru102 = results.get_atoms(rhodium, 'Ru102', time_units='d')
_ ,number_of_Pd102 = results.get_atoms(rhodium, 'Pd102', time_units='d')
_ ,number_of_Pd104 = results.get_atoms(rhodium, 'Pd104', time_units='d')


# Total number of Rh103 atoms lost
initial_Rh103 = number_of_Rh103[0]
final_Rh103 = number_of_Rh103[-1]
lost_Rh103 = initial_Rh103 - final_Rh103

# Number of atoms of each product formed
number_of_Rh104_formed = number_of_Rh104[-1]
number_of_Rh102_formed = number_of_Rh102[-1]
number_of_Ru102_formed = number_of_Ru102[-1]
number_of_Pd102_formed = number_of_Pd102[-1]
number_of_Pd104_formed = number_of_Pd104[-1]

# Calculate the percentage of each product
percentage_Rh103 = (lost_Rh103 / initial_Rh103) * 100
percentage_Rh102 = (number_of_Rh102_formed / lost_Rh103) * 100
percentage_Rh104 = (number_of_Rh104_formed / lost_Rh103) * 100
percentage_Ru102 = (number_of_Ru102_formed / lost_Rh103) * 100
percentage_Pd102 = (number_of_Pd102_formed / lost_Rh103) * 100
percentage_Pd104 = (number_of_Pd104_formed / lost_Rh103) * 100

# Print the results
print(f"Percentage of Rh103: - {percentage_Rh103:.4f}%")
print(f"Percentage of Rh104: - {percentage_Rh104:.4f}%")
print(f"Percentage of Rh102: + {percentage_Rh102:.4f}%")
print(f"Percentage of Ru102: + {percentage_Ru102:.4f}%")
print(f"Percentage of Pd102: + {percentage_Pd102:.4f}%")
print(f"Percentage of Pd104: + {percentage_Pd104:.4f}%")
2

@tcsousa99 What data are you using for FISPACT-II? I suspect this difference is due to different data that is used between OpenMC and FISPACT-II. In your case, it looks like you are using ENDF/B-VIII.0 for OpenMC. I took your model and did a quick reaction rate tally in OpenMC to see what the dominant reactions are under the spectrum imposed and got the following:

   (n,2n)                               0.103019 +/- 0.000349698
   (n,p)                                0.00163331 +/- 4.68323e-06
   (n,a)                                0.000753999 +/- 2.20613e-06
   (n,gamma)                            0.00456995 +/- 9.04538e-05

So you can see that (n,2n) is the dominant reaction mode. Looking at the chain file, we have the following data:

  <nuclide name="Rh103" reactions="6">
    <reaction type="(n,2n)" Q="-9318000.0" target="Rh102"/>
    <reaction type="(n,3n)" Q="-16756000.0" target="Rh101"/>
    <reaction type="(n,p)" Q="0.0" target="Ru103"/>
    <reaction type="(n,a)" Q="0.0" target="Tc100"/>
    <reaction type="(n,gamma)" Q="6999000.0" target="Rh104" branching_ratio="0.924"/>
    <reaction type="(n,gamma)" Q="6999000.0" target="Rh104_m1" branching_ratio="0.076"/>
  </nuclide>

Rh103(n,2n) produces Rh102, which has a 207 d half-life. Rh103(n,gamma) has a reaction rate that is roughly 5% of the (n,2n) rate but its product (Rh104) decays with a 42 s half-life (almost all of which goes to Pd104).

Looking at EAF-2010, which is commonly used for FISPACT-II, it appears that the Rh103 data is missing the (n,2n) channel altogether. Looking at the ENDF file, I can see that they have some information related to radionuclide production for this channel (MF=8 and MF=10 in ENDF lingo) but importantly it is missing cross sections in MF=3. Perhaps @jsublet can comment whether this is an error in the EAF data.

You can also see this notebook I put together quickly looking at some of the cross sections. There are other reactions in EAF that produce 2 neutrons but they have higher thresholds and are likely not relevant here.

1 Like
3

Thank you very much for your quick and insightful response, @paulromano.

Could you please share the code you used for tallying these reactions? Although I didn’t perform the FISPACT-II calculations myself, I can confirm that the EAF-2010 library was used, as you suspected.

Transmuted material

I’m intrigued by the differences between these two libraries, particularly the missing reactions. While this likely explains the variation in dominant reactions, could it also account for the 10% difference in the transmuted Rh-103 material?

If you examine the neutron energy spectrum for D-T reactions, you’ll notice a peak at 14.1 MeV, corresponding to the neutrons generated from the fusion reaction.
image

Looking at the cross-section for (n,2n) for Rh103 in https://www.oecd-nea.org/janisweb/book/neutrons, one can see that it is quite significant for the tens of MeV energy (also that there’s no EAF-2010 as you mentioned).

Looking at the cross-section for (n, gamma) for Rh103, (that both ENDF and EAF-2010 have) there’s resonance curves and high cross-sections at lower energies (also in https://www.oecd-nea.org/janisweb/book/neutrons)

Missing data

Thank you for setting up the notebook; it greatly helped clarify the differences between these reactions. Regarding the missing data in ENDF, could you elaborate on how you identified the missing cross-sections in MF=3? From the notebook, it seemed that most of the important cross-sections were indeed present in ENDF.

Is it possible to manually include the missing data from EAF-2010 to correct this issue? Additionally, could you guide me on how to import the EAF-2010 cross-section database into OpenMC so I can run a comparison and analyze the differences in the results?

4

@paulromano I am amazed how efficient such communication method is… having said that and regarding EAF-2010, I would say the following:

To both: if EAF-2010 (n,2n) was missing and the dominant reactions where the one shown by Paul FISPACT-II could not burn more than OpenMC!

  • Paul; you are using or showing the gxs-709, ENDF-6 translated version of EAF-2010, the Rh103(n,2n) exists but in in MF-10 (not MF-3) as it leads to partial Rh102g and Rh102m, summing those channels will give you a MF-3 mt16 total. In the EAF format all channels are in pseudo MF-3
  • Thomas, you would need to compare what is comparable, if I get it right OpenMc burn Rh103 less than FISPACT-II in a very pronounced fusion spectra !!. To go to the bottom of it you’ll need to compare Paul’s OpenMC tallys with FISPACT-II printlib one group average cross section produced with the above spectra using EAF-2010. FISPACT-II simulation will also gives you the dominant paths.

Hope this helps

1 Like
5

Dear @jsublet,

Thank you for your answer!

Your answer to the first point is very insightful because this was exactly what we tought it was happening, that because EAF-2010 was missing (n, 2n) then FISPACT-II was going to overstimate the other reactions, but that’s not how Monte-Carlo works so your point makes a lot of sense, we should expect low burn in comparison with OpenMC.

Actually I was in communication with some people using FISPACT-II in ITER and they did a simulation closer to mine and the (%) lost of Rh103 was closer. The simulation before was 15 years old and we were a bit unsure of all the conditions.

Isotope OpenMC FISTPACT-II
Rh103 - 1.18 % -3.1%
  • In relation to the point for @paulromano, does this mean then that all cross-sections in MF-3 are for production of non-radioactive nuclei?
  • In relation to the point directed to me, I don’t have access to FISPACT-II but I’m working with colleagues that do, so I will ask of them tomorrow to try this comparison and to consider your advice!

I’ll keep you updated with news such that we can all benefit from understading the differences between the softwares and achieve a comparable result.