Depletion failure

I am running a U233/Th232 molten salt reactor simulation. I am using the following depletion code: Power is 3 Gigawatts

time_steps_day=[246060]*10
model = openmc.Model(geometry=geometry,materials=materials, settings=settings)
operator = openmc.deplete.CoupledOperator(model,chainfile)
integrator = openmc.deplete.PredictorIntegrator(operator, time_steps_day, power)
integrator.integrate()
results = openmc.deplete.ResultsList.from_hdf5(“./depletion_results.h5”)

I get the following warning after 6 depletion iterations:
Creating state point openmc_simulation_n6.h5…

/home/brian/.local/lib/python3.10/site-packages/openmc/deplete/helpers.py:427: RuntimeWarning: overflow encountered in scalar divide
return source_rate / energy
/home/brian/.local/lib/python3.10/site-packages/openmc/deplete/openmc_operator.py:548: RuntimeWarning: invalid value encountered in multiply
rates *= self._normalization_helper.factor(source_rate)
/home/brian/.local/lib/python3.10/site-packages/scipy/sparse/linalg/_dsolve/linsolve.py:276: MatrixRankWarning: Matrix is exactly singular
warn(“Matrix is exactly singular”, MatrixRankWarning)
/home/brian/.local/lib/python3.10/site-packages/scipy/sparse/linalg/_dsolve/linsolve.py:276: MatrixRankWarning: Matrix is exactly singular
warn(“Matrix is exactly singular”, MatrixRankWarning)

And the following failure report:
AllocationError Traceback (most recent call last)
Cell In[16], line 5
3 operator = openmc.deplete.CoupledOperator(model,chainfile)
4 integrator = openmc.deplete.PredictorIntegrator(operator, time_steps_day, power)
----> 5 integrator.integrate()
6 results = openmc.deplete.ResultsList.from_hdf5(“./depletion_results.h5”)

File ~/.local/lib/python3.10/site-packages/openmc/deplete/abc.py:806, in Integrator.integrate(self, final_step, output)
804 # Solve transport equation (or obtain result from restart)
805 if i > 0 or self.operator.prev_res is None:
→ 806 conc, res = self._get_bos_data_from_operator(i, source_rate, conc)
807 else:
808 conc, res = self._get_bos_data_from_restart(i, source_rate, conc)

File ~/.local/lib/python3.10/site-packages/openmc/deplete/abc.py:759, in Integrator._get_bos_data_from_operator(self, step_index, source_rate, bos_conc)
756 “”“Get beginning of step concentrations, reaction rates from Operator
757 “””
758 x = deepcopy(bos_conc)
→ 759 res = self.operator(x, source_rate)
760 self.operator.write_bos_data(step_index + self._i_res)
761 return x, res

File ~/.local/lib/python3.10/site-packages/openmc/deplete/coupled_operator.py:465, in CoupledOperator.call(self, vec, source_rate)
462 openmc.lib.run()
464 # Extract results
→ 465 rates = self._calculate_reaction_rates(source_rate)
467 # Get k and uncertainty
468 keff = ufloat(*openmc.lib.keff())

File ~/.local/lib/python3.10/site-packages/openmc/deplete/openmc_operator.py:533, in OpenMCOperator._calculate_reaction_rates(self, source_rate)
530 for nuc, i_nuc_results in zip(rxn_nuclides, nuc_ind):
531 number[i_nuc_results] = self.number[mat, nuc]
→ 533 tally_rates = self._rate_helper.get_material_rates(
534 mat_index, nuc_ind, react_ind)
536 # Compute fission yields for this material
537 fission_yields.append(self._yield_helper.weighted_yields(i))

File ~/.local/lib/python3.10/site-packages/openmc/deplete/helpers.py:217, in DirectReactionRateHelper.get_material_rates(self, mat_id, nuc_index, react_index)
198 “”“Return an array of reaction rates for a material
199
200 Parameters
(…)
214 reaction rates in this material
215 “””
216 self._results_cache.fill(0.0)
→ 217 full_tally_res = self.rate_tally_means[mat_id]
218 for i_tally, (i_nuc, i_react) in enumerate(
219 product(nuc_index, react_index)):
220 self._results_cache[i_nuc, i_react] = full_tally_res[i_tally]

File ~/.local/lib/python3.10/site-packages/openmc/deplete/helpers.py:186, in DirectReactionRateHelper.rate_tally_means(self)
184 # If the mean cache is empty, fill it once with this transport cycle’s results
185 if self._rate_tally_means_cache is None:
→ 186 self._rate_tally_means_cache = self._rate_tally.mean
187 return self._rate_tally_means_cache

File ~/.local/lib/python3.10/site-packages/openmc/lib/tally.py:276, in Tally.mean(self)
273 @property
274 def mean(self):
275 n = self.num_realizations
→ 276 sum_ = self.results[:, :, 1]
277 if n > 0:
278 return sum_ / n

File ~/.local/lib/python3.10/site-packages/openmc/lib/tally.py:306, in Tally.results(self)
304 data = POINTER(c_double)()
305 shape = (c_size_t*3)()
→ 306 _dll.openmc_tally_results(self._index, data, shape)
307 return as_array(data, tuple(shape))

File ~/.local/lib/python3.10/site-packages/openmc/lib/error.py:21, in _error_handler(err, func, args)
19 # Raise exception type corresponding to error code
20 if err == errcode(‘OPENMC_E_ALLOCATE’):
—> 21 raise exc.AllocationError(msg)
22 elif err == errcode(‘OPENMC_E_OUT_OF_BOUNDS’):
23 raise exc.OutOfBoundsError(msg)

AllocationError: Tally results have not been allocated yet.

I have tried varying source particle settings, varying the duration of the time steps, neither of which made any difference except when I changed the time step duration to a week, it failed one second iteration. Keff is well over 1.0 (1.03) when the failure occurs. Any suggestions?

Have you figured out what caused your “AllocationError: Tally results have not been allocated yet.”?

This error appeared for me when I tried to run a pincell depletion example with an erroneous chain file (jeff33 chain file I generated myself with missing fission reactions and yields)…
Without fission reactions/yields in the chain file I guess it is not surprising that it will fail after the the first transport step. With the corrected chain file the depletion example runs.

Error in pincell depletion example with erroneous jeff33 chain:
Traceback (most recent call last):
File “~/openmc/calculations/AP1/pincell_depletion_07/pincell_depletion.py”, line 360, in
integrator.integrate()
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/abc.py”, line 801, in integrate
n, res = self._get_bos_data_from_operator(i, source_rate, n)
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/abc.py”, line 745, in _get_bos_data_from_operator
res = self.operator(x, source_rate)
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/coupled_operator.py”, line 456, in call
rates = self._calculate_reaction_rates(source_rate)
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/openmc_operator.py”, line 527, in _calculate_reaction_rates
tally_rates = self._rate_helper.get_material_rates(
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/helpers.py”, line 219, in get_material_rates
full_tally_res = self.rate_tally_means[mat_index]
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/deplete/helpers.py”, line 188, in rate_tally_means
self._rate_tally_means_cache = self.rate_tally.mean
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/lib/tally.py”, line 310, in mean
sum
= self.results[:, :, 1]
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/lib/tally.py”, line 340, in results
_dll.openmc_tally_results(self._index, data, shape)
File “~/openmc/venv-openmc/lib/python3.10/site-packages/openmc/lib/error.py”, line 21, in _error_handler
raise exc.AllocationError(msg)
openmc.exceptions.AllocationError: Tally results have not been allocated yet.

The script I used to generate the jeff33 chain file is an adapted one of this one here: data/depletion/generate_tendl_chain.py at master · openmc-dev/data · GitHub

Update: My jeff33 chain file was erroneous (fission product yield and fission reactions were missing since I pointed incorrectly to the fpy_files).

The warnings/messages from openmc.deplete.Chain.from_endf for jeff33 were:
[…]
[Various warnings from decay.py about continuous spectra with log-linear interpolations e.g. for Th230, Th232, … Es 253]
[…]
Decay mode [‘beta-’] of parent In123_m1 has a zero branching ratio.
[…]
The following fissionable nuclides have no fission product yields:
Ra223, replaced with U235
Ra226, replaced with U235
Ac225, replaced with U235
[…]
Fm255, replaced with U235