From 43437e941027470bb5162c61a066ca53490c5461 Mon Sep 17 00:00:00 2001 From: Lorenzo Chierici Date: Fri, 28 Apr 2023 05:11:17 +0200 Subject: [PATCH] TransferRates (old msr continuous capabilities) (#2358) --------- Co-authored-by: Jonathan Shimwell Co-authored-by: Olek <45364492+yardasol@users.noreply.github.com> Co-authored-by: yardasol Co-authored-by: Paul Romano Co-authored-by: Gavin Ridley --- docs/source/methods/depletion.rst | 82 ++++++++ docs/source/pythonapi/deplete.rst | 16 +- docs/source/usersguide/depletion.rst | 70 ++++++- openmc/deplete/__init__.py | 1 + openmc/deplete/abc.py | 37 +++- openmc/deplete/chain.py | 58 ++++++ openmc/deplete/pool.py | 56 +++++- openmc/deplete/transfer_rates.py | 188 ++++++++++++++++++ .../deplete_with_transfer_rates/__init__.py | 0 .../ref_depletion_with_feed.h5 | Bin 0 -> 37328 bytes .../ref_depletion_with_removal.h5 | Bin 0 -> 37328 bytes .../ref_depletion_with_transfer.h5 | Bin 0 -> 37328 bytes .../ref_no_depletion_only_feed.h5 | Bin 0 -> 37328 bytes .../ref_no_depletion_only_removal.h5 | Bin 0 -> 37328 bytes .../ref_no_depletion_with_transfer.h5 | Bin 0 -> 37328 bytes .../deplete_with_transfer_rates/test.py | 80 ++++++++ .../unit_tests/test_deplete_transfer_rates.py | 118 +++++++++++ 17 files changed, 695 insertions(+), 11 deletions(-) create mode 100644 openmc/deplete/transfer_rates.py create mode 100644 tests/regression_tests/deplete_with_transfer_rates/__init__.py create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_feed.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_removal.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_transfer.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_no_depletion_only_feed.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_no_depletion_only_removal.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/ref_no_depletion_with_transfer.h5 create mode 100644 tests/regression_tests/deplete_with_transfer_rates/test.py create mode 100644 tests/unit_tests/test_deplete_transfer_rates.py diff --git a/docs/source/methods/depletion.rst b/docs/source/methods/depletion.rst index 1b131bed5..87a9976dd 100644 --- a/docs/source/methods/depletion.rst +++ b/docs/source/methods/depletion.rst @@ -257,3 +257,85 @@ choose one of two methods for estimating the heating rate, including: The method for normalization can be chosen through the ``normalization_mode`` argument to the :class:`openmc.deplete.CoupledOperator` class. + +-------------- +Transfer Rates +-------------- + +OpenMC allows continuous removal or feed of nuclides by adding an +extra transfer rate term to the depletion matrix. An application of this feature +is the chemical processing of Molten Salt Reactors (MSRs), where one can +model the removal of fission products or feeding fresh fuel into the system. + +A transfer rate as defined here is the rate at which nuclides are +continuously removed/fed from/to a material. + +.. note:: + + A transfer rate can be positive or negative, indicating removal or feed + respectively. + +Mathematically, it can be thought of as an additional term :math:`\mathbf{T}` +in the depletion equation that is proportional to the nuclide density, which can be written as: + +.. math:: + + \begin{aligned}\frac{dN_i(t)}{dt} = &\underbrace{\sum\limits_j f_{j\rightarrow i} + \int_0^\infty dE \; \sigma_j (E,t) \phi(E,t) N_j(t) - \int_0^\infty dE \; \sigma_i(E,t) + \phi(E,t) N_i(t)}_\textbf{R} \\ + &+ \underbrace{\sum_j \left [ \lambda_{j\rightarrow i} N_j(t) - \lambda_{i\rightarrow j} N_i(t) \right ]}_\textbf{D} \\ + &- \underbrace{t_i N_i(t)}_\textbf{T} \end{aligned} + +where the reaction term :math:`\mathbf{R}`, the decay term :math:`\mathbf{D}` +and the new transfer term :math:`\mathbf{T}` have been grouped together so that +:math:`\mathbf{A} = \mathbf{R}+\mathbf{D}-\mathbf{T}`. +The transfer rate coefficient :math:`t_i` defines the continuous transfer of the +nuclide :math:`i`, which behaves similar to radioactive decay. +:math:`t_i` can also be defined as the reciprocal of a cycle time +:math:`T_{cyc}`, intended as the time needed to process the whole inventory. + +Note that this formulation assumes homogeneous distribution of nuclide +:math:`i` throughout the material. + +A more rigorous description of removal rate and its implementation can be found +in the paper by `Hombourger +`_. + +The resulting burnup matrix can be solved with the same integration algorithms +that are used in the absence of the transfer term. + +.. note:: + + If no ``destination_material`` is specified, nuclides that are removed + or fed will not be tracked afterwards. + +Coupling materials +------------------ + +To keep track of removed nuclides or to feed nuclides from one depletable material +to another, the respective depletion equations have to be coupled. This can be +achieved by defining one block matrix, with diagonal blocks corresponding to +depletion matrices :math:`\mathbf{A_{ii}}`, where the index :math:`i` indicates +the depletable material id, and off-diagonal blocks corresponding to inter-material +coupling matrices :math:`\mathbf{T_{ij}}`, positioned so that that the indices :math:`i` and +:math:`j` indicate the nuclides receiving and losing materials, respectively. +The nuclide vectors are assembled together in one single vector and the resulting +system is solved with the same integration algorithms seen before. + +As an example, consider the case of two depletable materials and one +transfer defined from material 1 to material 2. The final system will look like: + +.. math:: + + \begin{aligned}\frac{d}{dt}\begin{pmatrix}\vec{N_1}\\ \vec{N_2}\end{pmatrix} &= + \begin{pmatrix}\mathbf{A_{11}} & \mathbf{0}\\ \mathbf{T_{21}} & \mathbf{A_{22 }} + \end{pmatrix} \begin{pmatrix}\vec{N_1}\\ \vec{N_2}\end{pmatrix} \end{aligned} + +where: + +:math:`\mathbf{A_{11}} = \mathbf{R_{11}}+\mathbf{D_{11}}-\mathbf{T_{21}}`, and + +:math:`\mathbf{A_{22}} = \mathbf{R_{22}}+\mathbf{D_{22}}`. + +Note that mass conservation is guaranteed by transferring the number +of atoms directly. diff --git a/docs/source/pythonapi/deplete.rst b/docs/source/pythonapi/deplete.rst index dd679a40a..dc9370ab5 100644 --- a/docs/source/pythonapi/deplete.rst +++ b/docs/source/pythonapi/deplete.rst @@ -54,7 +54,7 @@ provides the following transport operator classes: IndependentOperator The :class:`CoupledOperator` and :class:`IndependentOperator` classes must also -have some knowledge of how nuclides transmute and decay. This is handled by the +have some knowledge of how nuclides transmute and decay. This is handled by the :class:`Chain` class. Minimal Example @@ -195,14 +195,24 @@ total system energy. helpers.FissionYieldCutoffHelper helpers.FluxCollapseHelper -The :class:`openmc.deplete.IndependentOperator` uses inner classes subclassed +The :class:`openmc.deplete.IndependentOperator` uses inner classes subclassed from those listed above to perform similar calculations. +The following classes are used to define transfer rates to model continuous +removal or feed of nuclides during depletion. + +.. autosummary:: + :toctree: generated + :nosignatures: + :template: myclass.rst + + transfer_rates.TransferRates + Intermediate Classes -------------------- Specific implementations of abstract base classes may utilize some of -the same methods and data structures. These methods and data are stored +the same methods and data structures. These methods and data are stored in intermediate classes. Methods common to tally-based implementation of :class:`FissionYieldHelper` diff --git a/docs/source/usersguide/depletion.rst b/docs/source/usersguide/depletion.rst index a16e98ca7..439e30305 100644 --- a/docs/source/usersguide/depletion.rst +++ b/docs/source/usersguide/depletion.rst @@ -262,7 +262,7 @@ Loading and Generating Microscopic Cross Sections ------------------------------------------------- As mentioned earlier, any transport code could be used to calculate one-group -microscopic cross sections. The :mod:`openmc.deplete` module provides the +microscopic cross sections. The :mod:`openmc.deplete` module provides the :class:`~openmc.deplete.MicroXS` class, which contains methods to read in pre-calculated cross sections from a ``.csv`` file or from data arrays:: @@ -355,9 +355,9 @@ Multiple Materials A transport-independent depletion simulation using ``source-rate`` normalization will calculate reaction rates for each material independently. This can be -useful for running many different cases of a particular scenario. A +useful for running many different cases of a particular scenario. A transport-independent depletion simulation using ``fission-q`` normalization -will sum the fission energy values across all materials into :math:`Q_i` in +will sum the fission energy values across all materials into :math:`Q_i` in Equation :math:numref:`fission-q`, and Equation :math:numref:`fission-q` provides the flux we use to calculate the reaction rates in each material. This can be useful for running a scenario with multiple depletable materials @@ -370,3 +370,67 @@ The values of the one-group microscopic cross sections passed to :class:`openmc.deplete.IndependentOperator` are fixed for the entire depletion simulation. This implicit assumption may produce inaccurate results for certain scenarios. + +Transfer Rates +============== + +Transfer rates define removal or feed of nuclides to or from one or more +depletable materials. This can be useful to model continuous fuel reprocessing, +online fission products separation, etc. + +Transfer rates are defined by calling the +:meth:`~openmc.deplete.abc.Integrator.add_transfer_rate()` method directly from +one of the Integrator classes:: + + ... + integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power) + integrator.add_transfer_rate(...) + +Defining transfer rates +----------------------- + +The :meth:`~openmc.deplete.abc.Integrator.add_transfer_rate()` method requires a +:class:`~openmc.Material` instance (alternatively, a material id or +the name) as the depletable material from which nuclides are processed, +a list of elements that share the same transfer rate, and a transfer rate itself. + +.. caution:: + + Make sure you set the transfer rate value with the right sign. + A positive transfer rate assumes removal, while a negative one assumes feed. + +The ``transfer_rate_units`` argument specifies the units for the transfer rate. +The default is `1/s`, but '1/min', '1/h', '1/d' and '1/a' are also valid +options. + +For example, to define continuous removal of xenon from one material with a +removal rate value of 0.1 s\ :sup:`-1` (or a cycle time of 10 s), you'd use:: + + mat1 = openmc.Material(material_id=1, name='fuel') + + ... + + integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power) + # by openmc.Material object + integrator.add_transfer_rate(mat1, ['Xe'], 0.1) + # or by material id + integrator.add_transfer_rate(1, ['Xe'], 0.1) + # or by material name + integrator.add_transfer_rate('fuel', ['Xe'], 0.1) + +Note that in this case the xenon isotopes that are removed will not be tracked. + +Defining a destination material +------------------------------- + +To transfer elements from one depletable material to another, the +``destination_material`` parameter needs to be passed to the +:meth:`~openmc.deplete.abc.Integrator.add_transfer_rate()` method. For example, +to transfer xenon from one material to another, you'd use:: + + ... + mat2 = openmc.Material(name='storage') + + ... + + integrator.add_transfer_rate(mat1, ['Xe'], 0.1, destination_material=mat2) diff --git a/openmc/deplete/__init__.py b/openmc/deplete/__init__.py index 329ee8b52..8a9509e90 100644 --- a/openmc/deplete/__init__.py +++ b/openmc/deplete/__init__.py @@ -16,6 +16,7 @@ from .atom_number import * from .stepresult import * from .results import * from .integrators import * +from .transfer_rates import * from . import abc from . import cram from . import helpers diff --git a/openmc/deplete/abc.py b/openmc/deplete/abc.py index 5b405b1d5..e871eccfe 100644 --- a/openmc/deplete/abc.py +++ b/openmc/deplete/abc.py @@ -24,6 +24,7 @@ from .stepresult import StepResult from .chain import Chain from .results import Results from .pool import deplete +from .transfer_rates import TransferRates __all__ = [ @@ -547,7 +548,6 @@ class Integrator(ABC): :attr:`solver`. .. versionadded:: 0.12 - Attributes ---------- operator : openmc.deplete.abc.TransportOperator @@ -574,7 +574,10 @@ class Integrator(ABC): * ``n1`` is a :class:`numpy.ndarray` of compositions at the next time step. Expected to be of the same shape as ``n0`` - .. versionadded:: 0.12 + transfer_rates : openmc.deplete.TransferRates + Instance of TransferRates class to perform continuous transfer during depletion + + .. versionadded:: 0.13.4 """ @@ -654,6 +657,8 @@ class Integrator(ABC): self.timesteps = asarray(seconds) self.source_rates = asarray(source_rates) + self.transfer_rates = None + if isinstance(solver, str): # Delay importing of cram module, which requires this file if solver == "cram48": @@ -704,7 +709,8 @@ class Integrator(ABC): def _timed_deplete(self, concs, rates, dt, matrix_func=None): start = time.time() results = deplete( - self._solver, self.chain, concs, rates, dt, matrix_func) + self._solver, self.chain, concs, rates, dt, matrix_func, + self.transfer_rates) return time.time() - start, results @abstractmethod @@ -835,6 +841,31 @@ class Integrator(ABC): self.operator.finalize() + def add_transfer_rate(self, material, elements, transfer_rate, + transfer_rate_units='1/s', destination_material=None): + """Add transfer rates to depletable material. + + Parameters + ---------- + material : openmc.Material or str or int + Depletable material + elements : list of str + List of strings of elements that share transfer rate + transfer_rate : float + Rate at which elements are transferred. A positive or negative values + set removal of feed rates, respectively. + destination_material : openmc.Material or str or int, Optional + Destination material to where nuclides get fed. + transfer_rate_units : {'1/s', '1/min', '1/h', '1/d', '1/a'} + Units for values specified in the transfer_rate argument. 's' means + seconds, 'min' means minutes, 'h' means hours, 'a' means Julian years. + + """ + if self.transfer_rates is None: + self.transfer_rates = TransferRates(self.operator, self.operator.model) + + self.transfer_rates.set_transfer_rate(material, elements, transfer_rate, + transfer_rate_units, destination_material) @add_params class SIIntegrator(Integrator): diff --git a/openmc/deplete/chain.py b/openmc/deplete/chain.py index 372ed3510..5e3157f90 100644 --- a/openmc/deplete/chain.py +++ b/openmc/deplete/chain.py @@ -698,6 +698,64 @@ class Chain: dict.update(matrix_dok, matrix) return matrix_dok.tocsr() + def form_rr_term(self, transfer_rates, materials): + """Function to form the transfer rate term matrices. + + .. versionadded:: 0.13.4 + + Parameters + ---------- + transfer_rates : openmc.deplete.TransferRates + Instance of openmc.deplete.TransferRates + materials : string or two-tuple of strings + Two cases are possible: + + 1) Material ID as string: + Nuclide transfer only. In this case the transfer rate terms will be + subtracted from the respective depletion matrix + + 2) Two-tuple of material IDs as strings: + Nuclide transfer from one material into another. + The pair is assumed to be + ``(destination_material, source_material)``, where + ``destination_material`` and ``source_material`` are the nuclide + receiving and losing materials, respectively. + The transfer rate terms get placed in the final matrix with indexing + position corresponding to the ID of the materials set. + + Returns + ------- + scipy.sparse.csr_matrix + Sparse matrix representing transfer term. + + """ + matrix = defaultdict(float) + + for i, nuclide in enumerate(self.nuclides): + element = re.split(r'\d+', nuclide.name)[0] + # Build transfer terms matrices + if isinstance(materials, str): + material = materials + if element in transfer_rates.get_elements(material): + matrix[i, i] = transfer_rates.get_transfer_rate(material, element) + else: + matrix[i, i] = 0.0 + #Build transfer terms matrices + elif isinstance(materials, tuple): + destination_material, material = materials + if transfer_rates.get_destination_material(material, element) == destination_material: + matrix[i, i] = transfer_rates.get_transfer_rate(material, element) + else: + warn(f'Material {destination_material} is not defined ' + f'as a destination material for Material {material}. ' + 'Setting transfer rate to 0.0') + matrix[i, i] = 0.0 + #Nothing else is allowed + n = len(self) + matrix_dok = sp.dok_matrix((n, n)) + dict.update(matrix_dok, matrix) + return matrix_dok.tocsr() + def get_branch_ratios(self, reaction="(n,gamma)"): """Return a dictionary with reaction branching ratios diff --git a/openmc/deplete/pool.py b/openmc/deplete/pool.py index aba60b558..fc006a884 100644 --- a/openmc/deplete/pool.py +++ b/openmc/deplete/pool.py @@ -4,6 +4,8 @@ Provided to avoid some circular imports """ from itertools import repeat, starmap from multiprocessing import Pool +from scipy.sparse import bmat +import numpy as np # Configurable switch that enables / disables the use of @@ -15,7 +17,8 @@ USE_MULTIPROCESSING = True NUM_PROCESSES = None -def deplete(func, chain, x, rates, dt, matrix_func=None, *matrix_args): +def deplete(func, chain, x, rates, dt, matrix_func=None, transfer_rates=None, + *matrix_args): """Deplete materials using given reaction rates for a specified time Parameters @@ -37,7 +40,11 @@ def deplete(func, chain, x, rates, dt, matrix_func=None, *matrix_args): ``fission_yields = {parent: {product: yield_frac}}`` Expected to return the depletion matrix required by ``func`` - matrix_args : Any, optional + transfer_rates : openmc.deplete.TransferRates, Optional + Object to perform continuous reprocessing. + + .. versionadded:: 0.13.4 + matrix_args: Any, optional Additional arguments passed to matrix_func Returns @@ -62,6 +69,51 @@ def deplete(func, chain, x, rates, dt, matrix_func=None, *matrix_args): matrices = map(matrix_func, repeat(chain), rates, fission_yields, *matrix_args) + if transfer_rates is not None: + # Calculate transfer rate terms as diagonal matrices + transfers = map(chain.form_rr_term, repeat(transfer_rates), + transfer_rates.burnable_mats) + # Subtract transfer rate terms from Bateman matrices + matrices = [matrix - transfer for (matrix, transfer) in zip(matrices, + transfers)] + + if len(transfer_rates.index_transfer) > 0: + # Calculate transfer rate terms as diagonal matrices + transfer_pair = { + mat_pair: chain.form_rr_term(transfer_rates, mat_pair) + for mat_pair in transfer_rates.index_transfer + } + + # Combine all matrices together in a single matrix of matrices + # to be solved in one go + n_rows = n_cols = len(transfer_rates.burnable_mats) + rows = [] + for row in range(n_rows): + cols = [] + for col in range(n_cols): + mat_pair = (transfer_rates.burnable_mats[row], + transfer_rates.burnable_mats[col]) + if row == col: + # Fill the diagonals with the Bateman matrices + cols.append(matrices[row]) + elif mat_pair in transfer_rates.index_transfer: + # Fill the off-diagonals with the transfer pair matrices + cols.append(transfer_pair[mat_pair]) + else: + cols.append(None) + + rows.append(cols) + matrix = bmat(rows) + + # Concatenate vectors of nuclides in one + x_multi = np.concatenate([xx for xx in x]) + x_result = func(matrix, x_multi, dt) + + # Split back the nuclide vector result into the original form + x_result = np.split(x_result, np.cumsum([len(i) for i in x])[:-1]) + + return x_result + inputs = zip(matrices, x, repeat(dt)) if USE_MULTIPROCESSING: diff --git a/openmc/deplete/transfer_rates.py b/openmc/deplete/transfer_rates.py new file mode 100644 index 000000000..c6f0cd94d --- /dev/null +++ b/openmc/deplete/transfer_rates.py @@ -0,0 +1,188 @@ +from numbers import Real + +from openmc.checkvalue import check_type, check_value +from openmc import Material +from openmc.data import ELEMENT_SYMBOL + + +class TransferRates: + """Class for defining continuous removals and feeds. + + Molten Salt Reactors (MSRs) benefit from continuous reprocessing, + which removes fission products and feeds fresh fuel into the system. MSRs + inspired the development of this class. + + An instance of this class can be passed directly to an instance of one of + the :class:`openmc.deplete.Integrator` classes. + + .. versionadded:: 0.13.4 + + Parameters + ---------- + operator : openmc.TransportOperator + Depletion operator + model : openmc.Model + OpenMC model containing materials and geometry. If using + :class:`openmc.deplete.CoupledOperator`, the model must also contain + a :class:`opnemc.Settings` object. + + Attributes + ---------- + burnable_mats : list of str + All burnable material IDs. + transfer_rates : dict of str to dict + Container of transfer rates, elements and destination material + index_transfer : Set of pair of str + Pair of strings needed to build final matrix (destination_material, mat) + """ + + def __init__(self, operator, model): + + self.materials = model.materials + self.burnable_mats = operator.burnable_mats + + #initialize transfer rates container dict + self.transfer_rates = {mat: {} for mat in self.burnable_mats} + self.index_transfer = set() + + def _get_material_id(self, val): + """Helper method for getting material id from Material obj or name. + + Parameters + ---------- + val : openmc.Material or str or int representing material name/id + + Returns + ------- + material_id : str + + """ + if isinstance(val, Material): + check_value('Depeletable Material', str(val.id), self.burnable_mats) + val = val.id + + elif isinstance(val, str): + if val.isnumeric(): + check_value('Material ID', str(val), self.burnable_mats) + else: + check_value('Material name', val, + [mat.name for mat in self.materials if mat.depletable]) + val = [mat.id for mat in self.materials if mat.name == val][0] + + elif isinstance(val, int): + check_value('Material ID', str(val), self.burnable_mats) + + return str(val) + + def get_transfer_rate(self, material, element): + """Return transfer rate for given material and element. + + Parameters + ---------- + material : openmc.Material or str or int + Depletable material + element : str + Element to get transfer rate value + + Returns + ------- + transfer_rate : float + Transfer rate value + + """ + material_id = self._get_material_id(material) + check_value('element', element, ELEMENT_SYMBOL.values()) + return self.transfer_rates[material_id][element][0] + + def get_destination_material(self, material, element): + """Return destination (or transfer) material for given material and + element, if defined. + + Parameters + ---------- + material : openmc.Material or str or int + Depletable material + element : str + Element that gets transferred to another material. + + Returns + ------- + destination_material_id : str + Depletable material ID to where the element gets transferred + + """ + material_id = self._get_material_id(material) + check_value('element', element, ELEMENT_SYMBOL.values()) + if element in self.transfer_rates[material_id]: + return self.transfer_rates[material_id][element][1] + + def get_elements(self, material): + """Extract removing elements for a given material + + Parameters + ---------- + material : openmc.Material or str or int + Depletable material + + Returns + ------- + elements : list + List of elements where transfer rates exist + + """ + material_id = self._get_material_id(material) + if material_id in self.transfer_rates: + return self.transfer_rates[material_id].keys() + + def set_transfer_rate(self, material, elements, transfer_rate, transfer_rate_units='1/s', + destination_material=None): + """Set element transfer rates in a depletable material. + + Parameters + ---------- + material : openmc.Material or str or int + Depletable material + elements : list of str + List of strings of elements that share transfer rate + transfer_rate : float + Rate at which elements are transferred. A positive or negative values + set removal of feed rates, respectively. + destination_material : openmc.Material or str or int, Optional + Destination material to where nuclides get fed. + transfer_rate_units : {'1/s', '1/min', '1/h', '1/d', '1/a'} + Units for values specified in the transfer_rate argument. 's' means + seconds, 'min' means minutes, 'h' means hours, 'a' means Julian years. + + """ + material_id = self._get_material_id(material) + check_type('transfer_rate', transfer_rate, Real) + + if destination_material is not None: + destination_material_id = self._get_material_id(destination_material) + if len(self.burnable_mats) > 1: + check_value('destination_material', str(destination_material_id), + self.burnable_mats) + else: + raise ValueError(f'Transfer to material {destination_material_id} '\ + 'is set, but there is only one depletable material') + else: + destination_material_id = None + + if transfer_rate_units in ('1/s', '1/sec'): + unit_conv = 1 + elif transfer_rate_units in ('1/min', '1/minute'): + unit_conv = 60 + elif transfer_rate_units in ('1/h', '1/hr', '1/hour'): + unit_conv = 60*60 + elif transfer_rate_units in ('1/d', '1/day'): + unit_conv = 24*60*60 + elif transfer_rate_units in ('1/a', '1/year'): + unit_conv = 365.25*24*60*60 + else: + raise ValueError("Invalid transfer rate unit '{}'".format(transfer_rate_units)) + + for element in elements: + check_value('element', element, ELEMENT_SYMBOL.values()) + self.transfer_rates[material_id][element] = transfer_rate / unit_conv, destination_material_id + if destination_material_id is not None: + self.index_transfer.add((destination_material_id, material_id)) diff --git a/tests/regression_tests/deplete_with_transfer_rates/__init__.py b/tests/regression_tests/deplete_with_transfer_rates/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_feed.h5 b/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_feed.h5 new file mode 100644 index 0000000000000000000000000000000000000000..cb8d49bfc65eb7dce6c4ab104320603d7ae4f3d6 GIT binary patch literal 37328 zcmeHPeN+_J6@SA5S(Pf1s7Q@n6RPW1vI;a<&JF_$5n7LdsBvv}xH0{CDCP7ba;8e6T^X{EtM`o8@ zz)*Adk9~9ByYIgHe((P7yqTTZ9XVMKe`RFCNTz(ORu;)()kFNc37@vzyar`#FMtC9 z%6cf9BKUzITU6TOp*|Mu+bq1I8R{25nw`y91Vj5+LNjZK=0~Yektb^_p91xpzV37c zalsRI&V4M|5Lp$@a5U)5M#&rzKQ3n>7mBq->z0O-%Ra{=iVy!Zb)E4fTAaG^? z9pnh-ALI*`F%0zC2)VzNtu8L{R1 ziVpG({qF+4#RHEx-vm>>HS&yL$Oh;Qm3KM<$Tv4k6yqrs;s@tH2M*>!iSa3+LWt#} zg^y$jEQ-}aKPfA)0W~SeDzO=X^*ARPYXcsjpGIEbJbh!hJ=E&}*&4^|Q(-ts3_R`^ z)oOO$Z6C!wY5>Q09$=<<;EVE~#vL9A81n$~NP_bNdWk$iK4N{iCC>v9a4u!+!d=`F zj)M(0!+EMA#$T@yc*QV>Ix&Go%lj8AD=u5*RdVhv`xnFXV;IkT z0=tXdBllZZUQ$yk?^4Fv(qRaRC68jsHw}5VaKp;nGbb*ko8=`3y((64a$bhdu{@|SylA@iy4{l0O(672d~8pmO_0ot{&q08wL zeX6|E5I{aP!^Ds;r(lC`7H37tyhHsrcquT8-}NWn?`i)gZ%4Z@jI7)D{h|8tT%!pT zX`hOSQP~#)7{6|4#C>YVJ^tF5ca3rAuhoBv?^7djd&TPD`Vhxy0{Im95eudr6Ziq@ zxk0X-$m`=luQ9=2&j)h*BwlX=JvP~2&kk~qrhtmWv&8eN5ik9|Se&{~^;DOa+RAE5 zS9&XLxMx*nWVUaGZK&K=QYS9RyH;htTj})_^4n*gYHu$~!R(mjt4>%q%2z7P^c4gF zvwYPEafID5KMYpBLY^7(jn0qocoqc?>vqE98QL-9nfqQo5b+Engvvk&nDLB39AV#~ z;+gS$t@9&1o<)Mgx}ETNhIY(&=CkpEh-VlfR0cx8jAunQ|6W7IGvj_&=SO%vL$7r^ z;qeUZnDNY>$Tf&(h!mB95HRCe+x`B%hKgqmQ}{T|FrLnj@OXwlC+l{?;~Cm9<5|Nr z?l?J|uih`re1?ANcEaNs+A-sqFNyn0Ji`c~G7ti0JgdM5WuY_1-C*62 zE}Fs*a9rbcrUUfaU|)=S#(gq>I*qYZFpPS}^SlK1_Ca1kJ!78A0X^%3`~dYFAm7gc+~{Y6`_;% zvRP|9rKO(fP_eqWN`3DJd5Qc-{a__MR}>7PFL%2!F0BBE1o;Zl!T1}j^SSZ)_dLUM zMNB!V%O9H|Z#vR=J^JGYxz@?+zii)lel+40>hrpYo7wY?oq#W{KjjrfjLNQe05paOWQ` z=`!Mh?ZI$Cz0mRD7eby{%FGG}X$LB9t0f?34S@(Ur{aa8UOS(Vcj}J+N03kpK5CVh% zAwUQa0))VDLLl7hhqiuRH#D1EFVJNqSCi_Oz$Nt8CBohLJ03kpK5CVh%AwUQWF9br}H+BTsFP<{lC({0Kc&$5m zLkJK8ga9Ex2!w*b`q{5|*5+@P`dda4xmeB$FT9s9Z_FYtbL?8@>Cam%u<>GyB%aen-w@1d!4?PA^D6W5C;Wr|Y{J}~Bq zt3A&1Ex%m0>CF$FKlxzMwGS3d(fY56FZ!;hFjD)Mr-M<_vF-)?HGk5cow0D|AZ>$fa8@0|3`*Dmg_xRO5qvBsy0p5LbVq2d*PO;SsrS1nJ9@1^hm z*cV?Zu@oIj#yuNR&{i-zL?gVM`{C|Gru}qYrj%WY(?3UG%J@feC zDIXU|M^X}m%`@|*WvRQ^iyIe9?it5#Y>?yO`gTj@%imutITDucd`*^ry~N>YnW(RC z+*`Aw7R~yXhDSu+?zdO`xIx47vHDA|9sOv6);`hd+P~!9LzbeDc_gZsh_=T9$>TKIa^ht4Bg4@B2( z)Yn%!n)KVf?^!flg|~lRz2iUn`YxF>*ZxATzP?qV5uOi^OLRJ$YvApOeMc-xdGv9ev(;`5&vF Y&M&=^Ufwh6YH`(e?gy?bmc7>Ze>f4d8UO$Q literal 0 HcmV?d00001 diff --git a/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_removal.h5 b/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_removal.h5 new file mode 100644 index 0000000000000000000000000000000000000000..b20bbc645ef97c997d19515a82745e3cac3c1559 GIT binary patch literal 37328 zcmeHPeQ;FO6~8YZu+bC(LI@G-2BJo&vLQ(grn~#d21qQ0#Kb0mlq|_YBKa~O7*sw5 z8mA@(beLG6grRDzc8ntxkZBvVY4M9DV37`XCOVy(nGP0g(PRj*@4b7@d)Y_cZg!K! zh-d%U`_8@Ro_p@^p5J-*y?t+=Oix?zy&Ix$VDiUgVj(O-ImEB8;nT5`*Px8;`EVdW zSqEiH5I+zki*j2O)JK4Qb1<)Hh5DTP=g(&>h@pKfp_w&~;YZ0(ktu1*pM2$;zOHu! z(&x|1l;?=?f^n8F26pBXAYO4ijO%cgt_&!x`waFsZ+Ztu8KcRj;pbfjYZMQk) z90wcBy7N>;xVv5x$VbNc>&dsfAD;j9c!e{IGBJUTk@hbpR$Nx-l5_4=`xo8x!x_(f z0vpRFN&VK9m(-L>yOarLI}9PR#2Z=S*Sb7g@c62X&SRh0jq;KO@+;*f9pEGmihxmG zI;-P3=8KWaOUNgEUh$Pr=$F5II#n0>;7eByNqic7zpvk~c)Pw9Cv%w1fOaix=+Xy8 zpDM361dvayFfruI3E1Ep#aU4z?@;$cUJ5MO?fMh%_jKLM+tDrzBk2x)f2b~!YqWqO z?NbpkDhERV<2Mf)ai7{f$z2=su09UkwT3S7eQF49uh?d|KE!caKt2I}M1X0_6n=nu zc_7zL<@HgZ*A(rp*9LOeG+u88J@$Qfy*QB5RRvTWoh6=E^>}Ib#U?8IRA+U0skyAC zbd{^pjC)pjMn?Nq*oI1dMP=fGv}=|3yOl0y0l$6btac5sU&+~b6#!4WXxSu(^C_6iiwdf>be1>%a$4A%<;|-T*CU97@2zF@yrBqgna|W zvlFxU2;Fz_Qsakz=QH$Ev*RDn(2fz$+7h_G#50T#Dm@`!#Ip)~P!>94+>O+Iee}su z6`znK2@4D3ljNV|%;44RfjJ4vFZNO?tI4wRKv$lNE6AnX5B+d==X~S9 zc2NvJz;X50nHJFNfPFFQ>G#R}=`_ZY!7%FS&+{VKYlggpdip$-4tl29`~dYVAZO~< z71Xa;fxqYI{n7HGfB6;PQ_$@AmtWD2*Y^$Pa>(#^?`sAi4Kykc(C1Avc+~>W6`_-H zWj|QsEG>232^Fi0tCaU{keA4R)c01>b4A_|`fj)D))k9gM$`I-l#Gf6vrC zSH$F#y7V#Y@}?z)*P}mqAlF)X{desf&yRY%e0^RQaWi_p@dDtB>rZ(F5uzpMu(6zc)e*1{As`8r30@q4eOBq^wKbmVG z-hB3Y?=<=9_2cswf<3RF(R%s5K4!(4;Ghr-pNU)i4p^6b&FKIhU5y7vPzZ9z1LF~K zK-u^CN5TAqD<0e{Z&e*h2ayf4}{)WG)vL>M#E?s@PyXpYnE@AdE16oH<_;n`k3G@DaB>y7qSdH*a7 z@_#UmJ~|%Jz&l_+HrwELWz@&(exvTqJKbAQA4|GF;EoSTga9Ex2oM5< z03kpK5CVk2XhOi>>xYgZUN^KFTrbdd!f1*D@`Vr}1PB2_fDj-A2m$X180{OqJr~vm z`$p8m`^Nea+YI)NsE_^i{YWMR2mwNX5Fi8y0YZQf7+nbXx^L|EvR^!5uur7@;pkd- z@`ex~1PB2_fDrHnf%SL&#`#dz7ICQMn5r~K{T-c-MJMae_4QfPPGv1ls{Pxxw0*+o z1&24~h{G+nyuPcxJ#xGa=MXQROwQQZ=gt{Bu2@g+ZM*08IdQg4`;PS$PfNAMJb(B2 zwO?GZ{^RiPmTr9GLu=>h>`SK?#i;$WqKZ~H3qsU?d1?q1kMu4&p!&1viJ6P*m&dC5 zsn4H(e~10JDtEuq`}3yW1XX_TtBQAXk7e2Zaw%h1QHEg4e)Hh;Gvnsk3ce`6&>5L+ z``?7z?uqRQHb?uj>$fiY$Qu6VThAS+xM<(Dz3I`S-P=_^6uhD?iAve#Qp$~c-nBQM z|Dd<^pvAGLsr?r}I~wbFZdc;U-plto-hAt)*$3Z`bxdtfiQW3tJE?JJ9Swh~IG6fD z%k;VeX?>qf$$h}}l{M}6{V!*gOY8gkLpO_IH${tk7XJHFXYf>!9i4yq)YjGFtkk2g z#QZZ~+?NzBY?+lME=_J=dpG2Wc{7jpKQ6^X-;UtQmmbLxEz!C4zm??Jr+Z%;m!_?6 z*ySIuJ%2!3-{Zgja`)`UMitLJhjQbWZHQFyvQNKnY;*OCs%+c!SYuew?W(+@aPJFe z|7Eo;KJ>^hr+nOJ6;p2hvVV7+ZB;|(j)m4#Tk$i?pQv~BTTi9G-Thg_ht`AxKUy=p z*R1y6``@f73tZazK6+*Aa~-V*Rex6gq<`9)su`;OTPy7qS^a-hW$~8#UVZq3L{;AS z*^5t>&s}2sv?}S$NxNVxoOC??53aekb6djqSGV17OZw}tUX2c$ZCkeDwsjqUKVw~% b{c%>_s*85hd-1QmleUeoFRm+=1J?Kd4VSmH literal 0 HcmV?d00001 diff --git a/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_transfer.h5 b/tests/regression_tests/deplete_with_transfer_rates/ref_depletion_with_transfer.h5 new file mode 100644 index 0000000000000000000000000000000000000000..d428cdd7cc9a856b83f4b8517e1192c1336fcd06 GIT binary patch literal 37328 zcmeHPeN+_J6@SBrt_f8HB@}J47-CSXqF}^$b{W_;LJv_0P&d^!Agp1-XZcWQ)ruOW zttgnZZc@|u(X=EU&+)W1YRw02w56@#pr(G*e&9bXrza-b98cgxv@`SWonc32mtDb7 zbN7#ZbKkq~zWaXf{_eb)o!RfF+n@N-eKGej`C~G%5EiZ+;@2(sbhhvsl(9V@4g@G0 zp==G}2ZCf#Zi|NcaIkL<<`svaKKH3bix>-HXdg>xX3eAcQ8H9yN!s!!U-_o5dmVxF zMX6cx95G%n&ho|Ko%wi(R~!%HI*g?&0}IT=!ldY4?X0SHRg|-ElSCt>g*zq)oEbm| zDZ=>&`GRE_1AR6|>ThDJT&2#M`bsCLvq(vCydDH-Ai=RXGy=0b__DWZa(WzLeN zc7cVmI_M|by&F)Hf}|3j>Rpd>lCe(U0s5)u1hKV>}Nq!#wbL`A_4H30)Zn|>JMnNMJ& z*#xQI>Wb3ZGHI7G-kbtMh%D)TmUK&(XA9S^%s74ga*9!2N`U-Ic}WL2i31{Fl$S2( zc#ipEsPYo>NuO8zTRW;38&2OGNd0nw+* zdkq2P(;=7`@}(U%_(pM7oWwgcKE_Ld1$$h7;{Bek$9X&2g<&M!f$tABMsbZ+P^5h- zB1Yvv2w?m=pb__}-4i^uG4JZ*&{J#h65pqW;P#4bfa^mXrxoOO@FN^dCp^FpP|pE! z-IsZNH0ZfwJoVZ@?wZW&&7j8~^wf(7IbBsi#gSR!c~y^>c3&()*{2rNRFs*^Ys*$T ztIW7(m1ktMZ-s5B)K^p{E=appdB0obEGXo+&kAas{Ve&iW0bG1!n#quQedVp9|#!b zD>uXucE|iMRQU>drq4GTKLX=fC^)Ry35;iG$B1W+ukeA0XBZ(=dPBg7XAI&9`wkb+ z^yh1hAA#{K1RU1v1jaM8W5lyIGara}h7m%gHw27$R&4g{HC#N??{_tR1jaM;TC)=v z&(Mw$&*I~_2JsA$qS6}zMm+0$*t6Gg@vJ$Pk8=pd)A$h>&+zAD%}!uELpw%1Ykrjb zO+53A6Osl-z=&tb5J%W+xOmnB=Y?n(Pvb{mJcFN4@Ay6l+A-o;&vcFq@eCt`N^b}l z@oW#o5%wJ}p5@2$5jMbhcgiypIIP(T?0klHjChtmiw{IR!w8|$8v;f=GeI0--{Im} z`wTup_e@@D{0QuPhJI>x0^=FlG2&TUBKMbgh7m%gHw27$R*4VFLT8M-p}Mb+KG~}0 zB+iqB#YJ=G$v?@<<$S5k;MFW|Zldywy;RC-vaCGN<#2KZxs>~%fBCY~FnsVi-*~WH z9Lo=IT>W)s0_b(Zz8LlN`(*xf8e_>|81?k$c@gaGfxLuz`aG2mdZyX@0QC|;&eE+b zsNb~$f6vkTqvgfG@+-clpxFs5zoH$V?;9-Okm2v%cMU)qXha~O&zolOsui9qLMP+P zzxr%JSy{m}s95EyR^GcoULyZd-&aY`75PHwi`}k|OB28$LcT(DF#d+>e6D}~Jxljo z5tC2q(#NdJn+XefJ^JGSxz57tzi8ihe$?aT@AJBdo6+-)R{>vKf66O}7?lGdV3b$f zI-VkrF|G!(ZAD!U-R;m$u? z(qzN~+kN4Jdj8|XFN}hNMs?+CuJp%8HjSq_>CoRD6*sE?E|VU!ta25U;(N;FwS}dw zA}5OK^O#t{AMs>ry}3#H&&pOLj7nIgQ*BY1n_7P*%6}44`&SzvTWoYgFXs&^H z^V#dY)8yONkI!ES_P%~b>*fFYm=$M&gF+JcOg!SZ-@4>$P7mYJ?Ranmg&MK7G|WG^^4ok~7@IKnx?!BHZ_|Ut1LFUv(Ae?J~H`d+lVsDbYTh%jbWiub|ep*dc+zt_K4Qw(|uHt+WGq1l|~U2mkf%KK+w zfd7MO_S5n3?Ysl_Lmv7vb9INO2KK`e^BbPG`vs(DN0JBuLVyq;1PB2_fDj-A2mwNX z5FiBZ3Iaygk!HM*^V5qxa6O6r{WY@zFGOKn?2o1Xx)No&zVvgUNSzQM1PB2_fDj-A z2mwNX5Fi8y0YZQfxB~=??k9u~(qSIaw2R!j1i4Y(J2mwNX z5Fi8y0YZQf7)b~Odi~Hj$m@nf2GEr#5RL{BkE&+eLs>30YZQfAOr{jLVyq;1V$DD{_Y#Qee4(84fct&KO9-> zPTmj#ga9Ex2oM7PAW%Q^je_U08^pnuqpIz;rbkcPGCP_s-MC?~pUGY_ukJ6K?XL(Q z6(0Fsu6U>A#b54hIvF+2igP%a(vh6;>PMy@zHB5QmrXZ_e^_uOP}RuyW8eZof~g$dFA*G*W@&-`L$W& zO8(ntDL-=f>Gf}ZU}-p;bNy^)tlIyD=;Ch{6o#n(@^mUxJldPNSM_Jprs+$XmdC03 zj@Pc7+LH2HRZcqCyVc#BsLHc@D*v?Xc(%3qddAM;48i*DZ}v|)H+F$l_@v_M>8Kp* z+vAsYN1RNwo_O!Y`i+_AEf4(up-X!!ucaKby1!HW+ZA$Bv z_j?cRPq4k?KKa9M9gDN=-kJ1F?`MzO{{FkK=j=ZfXFGRdVcf=Te@y$&dE53Am6y_f z-a4hRP+H#?7cR?lermD5ee*zeg|xnhpPM9(iHs3{vaTfc(my^;izu0P@~6==ZGVe8 zbN-g&lQj3d7@=WCwzwjBJKMc3S9DB2c5|&14}DvLt6pD|E6#{n)^t#k(|5hKZPi=a z`u?b6X3xZWZGBJAk8iTyY*z7fKDlUg-ZgD~PaMdLXju9S6~E14bN6-}o~p{xJ0>mO zyg1cbGBsnbYi6ItwrkC%3(4`;3DJ>}?bFh%k$+e`d+Cmw7UwL-`ug=BSYGb?_dn0Z zPEq^UY+YFr^W-QM*Ri+$y0QJw`&556efGXR>O)Q6eZ2GXtPi#Iy_i_D=A}tVsy)}G zi<$3Eddg~g`W=sCeA1VG- zu^5qR3O+_1RYo;he?h{)PW@} zx^2}|Uq$xqsH&Kt`jaOP9~L4aU>}m2#YmYts-=o1sqKCe-Zy{U?g%s72kqUA|yEoN0(9z#Vik8%Hb-*|Hf1UVNO*~S3 zvy6NjRT)9c2JnXE?T!HW)8Zmrvrt6@`dBiH&~)n3{gL)A-4f9 zDM*#7gSqt-Cxw_M9>7nX7ntYdJ<1;RW~m%ouj*@QIOjm_ac@AX8rp*$KP$5b8= z#yoJV{O56t0|7k`fJY9UAK)c;0zRTXtWI-LzrpT~flj$gS!dVL z5RTYWA@=#&hJ0mC8%tlb%d%M@+!O>x$~ zM|Bu}Qb|ihGp;{zzh~xY)egHfjMUA3e`t7<(ikU2-lqaFme~+M{I*ad_NlWQGqsU- z^*GGbT6#(CQ%kVD66fjq5aWzfd5Zj~BGd8v)B*Hbs62GPs;?%!(ao89lT@DBqUvqZ z6AxzU)lj){UIB{5DRExa@$&DBt@rk+$^Pz6yKA8HRI1m;p4FX^Y2QlQP}$ei_2Wv-$L;$QKLpUPvb3@XEAcvw-fGs z20JF6B_3A;G0zYoEOQ}X;#rL12z?91v#ITBgxMWR`uqrYK7*gWop3yZ9TU$ccPoFH zXNVA%xeze%tOpOuQfI{7O5N9oPxXB}ckh*jqis9)x}pyy4Syc(zHir{2T*JEdrot?=?sN!_G&wK9%yafNDpHK30MfniAv)lE! zj1e3Tkl7IqAEy5FghxPIK}w=#EMm)xSIC zW0t+?WCz|;?iy(ANVla>oPQn@73vXBF<)=@N&~Vnl|GZ|8ccQ!P}lyn`}UC#ecc1S zt*KXBEpKT5{b;2Dyw&V=@BH%O^;7e=qPf@4v|hp2$M&3wMoQvpCK>U&U|kZv(*-fPW#BRXnvaKTu1N7v z>!<(k0dNRE|H#Y*fBp-2Fr0tUR#5}*12{BhQ(f-CYpA(K+u!T&)wGjdygs+R1U1`- za_bd)tGj=;viu*+NI=J{K2#mh4?N6gX6~a*4fKOVe#3ctK|p497y(9r5hx}E!d^d2FY&rz!f?I7*9pZG1?&qWzz8q`i~u9R2rvTq5iso= z^F0?fL;FVP;l8mx;-q2U2z~U|`*E2OU<4QeMt~7u1Q-EEptukSx^JA#vtOJt>=Su^ zSX}GQ-Y^1;03*N%FakjkIJe`=>N6knB<>)WiYul4!YFaNxLvopA9^77fQ4mqFFdjI)e`IbNG)AME2 Uo8&_O_;`s3U|k`po4daM0hw;}0RR91 literal 0 HcmV?d00001 diff --git a/tests/regression_tests/deplete_with_transfer_rates/ref_no_depletion_only_removal.h5 b/tests/regression_tests/deplete_with_transfer_rates/ref_no_depletion_only_removal.h5 new file mode 100644 index 0000000000000000000000000000000000000000..81c97c755abf982fe95517e5ba48f6fd296a4897 GIT binary patch literal 37328 zcmeHPOKcoP5bfFZXOkFkJ|rdzXd;keK7q{-98mBic7h|J;6RXk0I?f;5^v-m{>c#_ zf*`^i0uCoj5y(dqP~-y=NV)ii;eeb-h#)S81M&&^f|Ns0=$WqP@wT^jJ!7x4B(;aF z>F(<4>euyZdS++#lcvTa_pRN!R(PL8LX?SW{}5ld=`(#%)u4>_6dhP3he?jd)PW@} zdTljSUrqMyxT=_-`V+4nJ}g8`z&<23i;)U-R7VxfQrr8a{BQoc+Yx9we6ZP@BjN>d zmU4C7S$!Es+paGP(!d-o8wCm#9ui8u=~CG0RfSAUenhS0C^NvPzJ? zSR?x<#Hmb&+kdvlC3Ue;Dz29W=(o8&QbCSr_qu%p9sPZzXh{uE2YiG7H;HdG#3RKw z%gDD;l@YXT0B=~{?FfKxEi_TYQys+*=D&*$c9TSWI#38HUpw|nqb1745cP8!^BeG! zf>fzKm|stEQiy5d0sPc?fq70opzJ|!mdc3@s=khfb6(6p?wx8iop)!}DUU|TF_i~| zF%R4+|9RZfKtRs};E_Y;2Y3mdfRCt;Dpekc(YaKJ8+O0P z!V!DdiM_YAJlp#2$)ne<{!wqrOG(PFoR>7fxtJ3HQ(n5E@f`VLrScN^r012ee1cz5 z`E-4F^BZ6N`>o_t?)|>LUnz8b?L4Spwh8Sa+R!!SM4v0~HUz+@37Q!AGDRDFQ=GN$ zRUL+(Rnijitm{wQ@0odCwZkqABXx7%9~$1QG{#Ah_o+aPWiA8|zb({=ed_F{Y;ELS zJr1+AmR?f()G}=S-lyV_N}xHm3L+n>%dMlp23cZXSI(h4dxk;Vwn#C6VIj}&F)n! zo{el*<4n+a0Y9Sg3_m9acB1hNc1%1Q*`fSqo@K`gmvSRu;#nQV5qcGiXE*7*P($Md z{D{Ug`uTL(_d#IC#Iu`EC}fyth!B?f5HRsb^dFYUtawXP+#*)V6D%_mjLQpOf;asur*B-sAsbFQupE zmHh|0T3n^zNv|LLOL#AdegKohMvApR!^r1QAdWMr=RB? zvNuL~33_^-Y9hVFGwJ|(Nh&vM>q_dYCh&WX?oW^xqvcn;rx4hQmS16~;QI#q6=e9` zd({Bsfl2}aJ#X6N)i^y@1Se~|o;scG>`XsS6{j+N{(CpzCHN2hLXw{=DumFz-LA)F zg5Yo{UjZG&-%6d&_2=K4wdab2cT$%hTg#it7gRm`X`%Aa0absmedGM7;}!OKUEpSV zzHyH5#rktz0b(q3Az;cYqZ&`aW5iW1TVd+J2R$!@Cac_Y!F$79eN-tZ0}r$p!UcNa@lh8>@mR&oiNKYG_;{vqnv;iNcT~!#{@p1b zv+T{JJMf-z*FbAWrp-lh{&`GPs7E}-e7)T(4ai0}bK311Om_@W*Zz$6_K^^M-2=U? z?pvOgKXmYZw9)|HYWBK!L3#1|srg&+{Of00ukh<*`_9B;Wl1%YtoU89E~&ujA|5Tq zLm_CzvhhGX0tb}C&%cJ|AMbfv%?q(fa~~?=Yzv!SDjtY8;DPco@ECt5z{4|FqIjtF z)BpDXI7FX+Z03A0{{=i4F1%=~sDbwZ92&E^KL6mg)Lg6W@Ada;+DR|jkl$X4n(agR z^-8_f+do@5{tsp(q~q0}st)J}9u_h)_gS_E`avSU;k>;dAhSDii4kA~7y(9r5nu!u z0Y-okU<4QeMqpJCFkMI5xR49!#TZ>rqJLO3=WroP3%}>5^YY>@5<1}b-$_m*v0G?=tJ`Tfoy!Z#0W3~ zi~u9R2rvSS03*N%loA3_uOFtDc-=5zxL)Avgi?wE_Jt8(1Q-EEfDvE>7=gkFnD&i@ zo(r3yeIxX6-&h}U(y(uYKKkqZxXcJJ0*nA7zz8q`i~u80S_p*QH_jH=FHRZuiM&57 zt#xN_7y(9r5nu!ufiMW1-FZ2E=J*HBQe{P7W5b0V*BV~_`NHqFZXIa6e*DCYTf%WKr=6cZ{PpXDyE~oh zSGHVA_I&QVN9zw+PV0MITVLzTD{ua|VXHH^dGhL+Zw@(^XubdXxN_T{4Vn2e>P>Q? Re|)?|1hB4<)XiJp{{Z#A^V_6dhP32T6{`)PW@} zdTrHIUq$xqxT+YV`s1%3IwV9)z&<23i=hg2R7(|2Qrr8a{BQnx&=F`nbfC$bBjN>d zmU4DH{EOEs+qFGP=*WJ-wN(PEnPRH1b=NW0s|uL3EH2u0G%kWThZ| zv0U~~h?ALix9@DXOX^~+R9r0!&~J6SrGgyM?s0qj+xvP+(UKaT4)_NDZxi3DiARcW zmXUA6DkEsw0N${C&=CONnrWhlr&@|1%zq~x>>`Qybf6GYzH#)m21}HS0qW;8top+~JDvyT9F_i~| zF%R4?|9RZfKtRs};E_Y;2Y3mdfRCt;Dpekc(YaKJTg#LsjDrSSJ5P02X6p@8`MXD> zdh-A7hx5OVSEWe$6SG8x+`lA5rnA-ca_*e{i#GjAp)#K(mWg$;-)-n7-)|(Ja_{%`{Ys(hYv(ZqvrT9Z(1xxtC;D9ZpdkQ0jnTxwmkHY7o8qi( zkLob^ypooPXI+2de$Uj4svUM|7^$24{?Oohr7=p1yiWyUEOQ}%_-&>}>{F-LWoskv z>T#H@HUE;@rAtQGyR*OJgxh0d&+5&{v~Q(tsO;G3GG@ikZiD%7^tAUtjh!B?f5HRsfP#mFev3RDRuLFKW<5?Lw9N3A*GuSclY}{4@ zG0zYoEb}2?;#r%W-K$tU)AzdpKcevrUI%uf@eFoMJga$1X)w=#6w7=Fn0Pk%WOlD& z@oZ>|8fT2g3-}R@XZSfeuoI1Euw&xc&{pL)^DH|~xRe_K6VGZXj?k-EJiATjg=!iv z;72r`(a)y~z7GOBCZ65iu8?7#AwpQ@L%_td5sD-9Ef&vGHEM(lG~PmamLP`%JJHT( zuw&v`>RB}q^9&KfG9Lmao+T)b(6?ATn|NA{Fug;`fFIG$XYezy6OCuEW8&HPZsjlY z3=zUI9|9(xb>l%<>WsKsto!=#slIpT?!B_`O6$(O-cRzLd_l@1s#?6cYq$T4y_BAs zSN0$1YIc=^C%t~~FXbr};e+RVYshxn7IlDe_3O+e=}pqU7<&3XSv{R5L@gPHo_?Np z$leI$CFtpSs*&^(+GC8^w`tt+W7nZWNkx<5f)jFw;Vo$R0yGmyIqgV z1i|4@z5+Ukzr{MA>(9S8Y0ni2@1!n2ww5=OFR6O?(@f=o{i^<9`^NcE$1Cjfy1>o! zeB%t^i}mNc0>oJ6Lco+)hBcmo$B3(3w!+kb4|-k*%Ln-7y{qRRQ-y{1V#`zo_J7bv zc{w?a_q&`BjA$f7qarSQcVe{kuax zX4#WTx8ptK&iyiqbF5=N# zJQRXfEE^BRBXB@D{QRqF{_*a2)w~d!H1~lb&bF}W`Qm|i10EoGz#;nlV^bG{`Oo6PVBtkuLJhnR;Lw;&b@>Ocq2?NGf3Lq+(?)v9`uz4%)NCKj zuUG1=-u~Ij@qaKwAsw&!M0G$v@UW1XnNPDd&<_&%4d?Ay0h!&AON;;`zz8q`i~u9R z2rvSS03*N%Fak@0fayBY#)Vu+FGlEk68*!Pxqu5%8W;T`_3KKM`T8>CLa{m{zz8q` zi~u9R2rvSS03*N%FanGKBd`DjO!pJ2=4o@BepiM*uKP{hN3LYIKp&Ft4`k!RB}RY| zU<4QeMt~7u1Q-EEpp+1Zdi^jt&+CRU!}S7RCzMhYurG`NBftnS0*nA7zz7sZz_f2H z^jz2!?Hi$o`^Nf+c)5kt^<|`|D8|p7_y;1+_ub2P0b7z0U&0|OR4*d3U z!{^ptE#H1{+*zm$<2&nW2vq;=%&tkTI$zg+aqYg9%MLp;n-h+EE$#gB(eGa$*wx|O zyuRUTvimdVeOiCWa$4VG+WJ~wTz~uLHJhAE^{31JxqQNTZ)59EUv4|=w4VF%`|Sq@ YoU63{h2(7i_;{WOU|k`po43CI0qyDW3jhEB literal 0 HcmV?d00001 diff --git a/tests/regression_tests/deplete_with_transfer_rates/test.py b/tests/regression_tests/deplete_with_transfer_rates/test.py new file mode 100644 index 000000000..7e20a7607 --- /dev/null +++ b/tests/regression_tests/deplete_with_transfer_rates/test.py @@ -0,0 +1,80 @@ +""" TransferRates depletion test suite """ + +from pathlib import Path + +import numpy as np +import pytest +import openmc +import openmc.deplete +from openmc.deplete import CoupledOperator + +@pytest.fixture +def model(): + f = openmc.Material(name="f") + f.add_element("U", 1, percent_type="ao", enrichment=4.25) + f.add_element("O", 2) + f.set_density("g/cc", 10.4) + + w = openmc.Material(name="w") + w.add_element("O", 1) + w.add_element("H", 2) + w.set_density("g/cc", 1.0) + w.depletable = True + + radii = [0.42, 0.45] + f.volume = np.pi * radii[0] ** 2 + w.volume = np.pi * (radii[1]**2 - radii[0]**2) + + materials = openmc.Materials([f, w]) + + surf_f = openmc.Sphere(r=radii[0]) + surf_w = openmc.Sphere(r=radii[1], boundary_type='reflective') + cell_f = openmc.Cell(fill=f, region=-surf_f) + cell_w = openmc.Cell(fill=w, region=+surf_f & -surf_w) + geometry = openmc.Geometry([cell_f, cell_w]) + + settings = openmc.Settings() + settings.particles = 100 + settings.inactive = 0 + settings.batches = 10 + + return openmc.Model(geometry, materials, settings) + +@pytest.mark.parametrize("rate, dest_mat, power, ref_result", [ + (1e-5, None, 0.0, 'no_depletion_only_removal'), + (-1e-5, None, 0.0, 'no_depletion_only_feed'), + (1e-5, None, 174.0, 'depletion_with_removal'), + (-1e-5, None, 174.0, 'depletion_with_feed'), + (-1e-5, 'w', 0.0, 'no_depletion_with_transfer'), + (1e-5, 'w', 174.0, 'depletion_with_transfer'), + ]) +def test_transfer_rates(run_in_tmpdir, model, rate, dest_mat, power, ref_result): + """Tests transfer_rates depletion class with transfer rates""" + + chain_file = Path(__file__).parents[2] / 'chain_simple.xml' + + transfer_elements = ['Xe'] + + op = CoupledOperator(model, chain_file) + integrator = openmc.deplete.PredictorIntegrator( + op, [1], power, timestep_units = 'd') + integrator.add_transfer_rate('f', transfer_elements, rate, + destination_material=dest_mat) + integrator.integrate() + + # Get path to test and reference results + path_test = op.output_dir / 'depletion_results.h5' + path_reference = Path(__file__).with_name(f'ref_{ref_result}.h5') + + # Load the reference/test results + res_ref = openmc.deplete.Results(path_reference) + res_test = openmc.deplete.Results(path_test) + + for mat in res_test[0].rates[0].index_mat: + for nuc in res_test[0].rates[0].index_nuc: + for rx in res_test[0].rates[0].index_rx: + y_test = res_test.get_reaction_rate(mat, nuc, rx)[1] / \ + res_test.get_atoms(mat, nuc)[1] + y_ref = res_ref.get_reaction_rate(mat, nuc, rx)[1] / \ + res_ref.get_atoms(mat, nuc)[1] + assert y_test == pytest.approx(y_ref, abs=1e-5) diff --git a/tests/unit_tests/test_deplete_transfer_rates.py b/tests/unit_tests/test_deplete_transfer_rates.py new file mode 100644 index 000000000..1572755c2 --- /dev/null +++ b/tests/unit_tests/test_deplete_transfer_rates.py @@ -0,0 +1,118 @@ +""" Tests for TransferRates class """ + +from pathlib import Path +from math import exp + +import pytest +import numpy as np + +import openmc +from openmc.deplete import CoupledOperator +from openmc.deplete.transfer_rates import TransferRates +from openmc.deplete.abc import (_SECONDS_PER_MINUTE, _SECONDS_PER_HOUR, + _SECONDS_PER_DAY, _SECONDS_PER_JULIAN_YEAR) + +CHAIN_PATH = Path(__file__).parents[1] / "chain_simple.xml" + +@pytest.fixture +def model(): + f = openmc.Material(name="f") + f.add_element("U", 1, percent_type="ao", enrichment=4.25) + f.add_element("O", 2) + f.set_density("g/cc", 10.4) + + w = openmc.Material(name="w") + w.add_element("O", 1) + w.add_element("H", 2) + w.set_density("g/cc", 1.0) + w.depletable = True + + radii = [0.42, 0.45] + f.volume = np.pi * radii[0] ** 2 + w.volume = np.pi * (radii[1]**2 - radii[0]**2) + materials = openmc.Materials([f, w]) + + surf_f = openmc.Sphere(r=radii[0]) + surf_w = openmc.Sphere(r=radii[1], boundary_type='vacuum') + cell_f = openmc.Cell(fill=f, region=-surf_f) + cell_w = openmc.Cell(fill=w, region=+surf_f & -surf_w) + geometry = openmc.Geometry([cell_f, cell_w]) + + settings = openmc.Settings() + settings.particles = 1000 + settings.inactive = 10 + settings.batches = 50 + + return openmc.Model(geometry, materials, settings) + + +def test_get_set(model): + """Tests the get/set methods""" + + # create transfer rates for U and Xe + transfer_rates = {'U': 0.01, 'Xe': 0.1} + op = CoupledOperator(model, CHAIN_PATH) + transfer = TransferRates(op, model) + + # Test by Openmc material, material name and material id + material, dest_material = [m for m in model.materials if m.depletable] + + for material_input in [material, material.name, material.id]: + for dest_material_input in [dest_material, dest_material.name, + dest_material.id]: + for element, transfer_rate in transfer_rates.items(): + transfer.set_transfer_rate(material_input, [element], transfer_rate, + destination_material=dest_material_input) + assert transfer.get_transfer_rate(material_input, + element) == transfer_rate + assert transfer.get_destination_material(material_input, + element) == str(dest_material.id) + assert transfer.get_elements(material_input) == transfer_rates.keys() + + +@pytest.mark.parametrize("transfer_rate_units, unit_conv", [ + ('1/s', 1), + ('1/sec', 1), + ('1/min', _SECONDS_PER_MINUTE), + ('1/minute', _SECONDS_PER_MINUTE), + ('1/h', _SECONDS_PER_HOUR), + ('1/hr', _SECONDS_PER_HOUR), + ('1/hour', _SECONDS_PER_HOUR), + ('1/d', _SECONDS_PER_DAY), + ('1/day', _SECONDS_PER_DAY), + ('1/a', _SECONDS_PER_JULIAN_YEAR), + ('1/year', _SECONDS_PER_JULIAN_YEAR), + ]) +def test_units(transfer_rate_units, unit_conv, model): + """ Units testing""" + # create transfer rate Xe + element = 'Xe' + transfer_rate = 1e-5 + op = CoupledOperator(model, CHAIN_PATH) + transfer = TransferRates(op, model) + + transfer.set_transfer_rate('f', [element], transfer_rate * unit_conv, + transfer_rate_units=transfer_rate_units) + assert transfer.get_transfer_rate('f', element) == transfer_rate + + +def test_transfer(run_in_tmpdir, model): + """Tests transfer depletion class without neither reaction rates nor decay + but only transfer rates""" + + # create transfer rate for U + element = ['U'] + transfer_rate = 1e-5 + op = CoupledOperator(model, CHAIN_PATH) + integrator = openmc.deplete.PredictorIntegrator( + op, [1,1], 0.0, timestep_units = 'd') + integrator.add_transfer_rate('f', element, transfer_rate) + integrator.integrate() + + # Get number of U238 atoms from results + results = openmc.deplete.Results('depletion_results.h5') + _, atoms = results.get_atoms(model.materials[0], "U238") + + # Ensure number of atoms equal transfer decay + assert atoms[1] == pytest.approx(atoms[0]*exp(-transfer_rate*3600*24)) + assert atoms[2] == pytest.approx(atoms[1]*exp(-transfer_rate*3600*24))