I was recently trying to reproduce the fairly standard 7% immediate decay heat after shutdown number using the new get_decay_heat function. The results seem a bit lower than I was expecting; it seems we’re missing some large component of the heating.
To check this, I started with the standard pincell_depletion problem that ships with OpenMC. To get a realistic nuclide vector coming out of it, I used the CASL chain. I also tried depleting with the full ENDF chain, and the results are within 5% of each other. I cranked up the particle count by a factor of 50, too.
I changed up the depletion only slightly, as follows:
bu_steps = [0.1, 0.5, 1.0, 2.0]
bu_steps.extend([2.0] * 30)
power = 174 # W/cm, for 2D simulations only (use W for 3D)
integrator = openmc.deplete.PredictorIntegrator(op, bu_steps, power, timestep_units='MWd/kg')
If I plot the decay heat ratio compared to full power, my results are a little disappointing:
This should be moving towards 7%. I think the decay heat function is perhaps incorrect at the moment. What could we possibly be missing? Or have I missed something in computing the power?
I have a full endf chain file with lots of reactions tracked (instead of the default 6) and n,gamma branching ratios added here if you are looking for other decay chains to try, also I would be super keen to see how it compares
I’d be curious to see some number density or atoms plots of the isotopes you expect most to contribute to decay heat to make sure there’s no funny business going on with any of those.
Just double checking as well, according to the documentation
power (float or iterable of float, optional) – Power of the reactor in [W]. A single value indicates that the power is constant over all timesteps. An iterable indicates potentially different power levels for each timestep. For a 2D problem, the power can be given in [W/cm] as long as the “volume” assigned to a depletion material is actually an area in [cm^2]. Either power, power_density, or source_rates must be specified.
You’re sure you’ve assigned the materials the correct “volume” with the correct units? If this script is on GH, might be good to link it so we can make sure there’s no silly typos (I believe in your ability to code this correctly, but sanity checks are a good start).
I am also wondering if some high decay-heat contributing reaction is missing for some reason from your chain.
The heating is not calculated correctly. It’s obvious in hindsight that the CASL chain won’t have all the nuclides to get correct decay heat values. I didn’t notice that the decay chain file still matters when calculating decay heats; it seemed like we should be looping over nuclides in the material file that we pulled from depletion rather than from the chain file. Guess that’s not the case, somewhere in there!
It sure would be nice if we could warn about scenarios where the openmc.config decay chain appears to differ from what was used in the original depletion calculation.
If I ensure that my chain is the full ENDF one for both the depletion calc and the decay heat calc, I get a lovely number hanging right around the classic 7% estimate:
