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?