Allow transfer rates with nuclides.

openmc/deplete/abc.py

@@ -840,16 +840,18 @@ def integrate(self, final_step=True, output=True):
 
         self.operator.finalize()
 
-    def add_transfer_rate(self, material, elements, transfer_rate,
+    def add_transfer_rate(self, material, components, 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
+        components : list of str
+            List of strings of elements and/or nuclides that share transfer rate.
+            A transfer rate for a nuclide cannot be added to a material
+            alongside a transfer rate for its element and vice versa.
         transfer_rate : float
             Rate at which elements are transferred. A positive or negative values
             set removal of feed rates, respectively.
@@ -863,7 +865,7 @@ def add_transfer_rate(self, material, elements, transfer_rate,
         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,
+        self.transfer_rates.set_transfer_rate(material, components, transfer_rate,
                                       transfer_rate_units, destination_material)
 
 @add_params

openmc/deplete/chain.py

@@ -724,15 +724,20 @@ def form_rr_term(self, transfer_rates, materials):
             # Build transfer terms matrices
             if isinstance(materials, str):
                 material = materials
-                if element in transfer_rates.get_elements(material):
+                components = transfer_rates.get_components(material)
+                if element in components:
                     matrix[i, i] = transfer_rates.get_transfer_rate(material, element)
+                elif nuclide.name in components:
+                    matrix[i, i] = transfer_rates.get_transfer_rate(material, nuclide.name)
                 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)
+                elif transfer_rates.get_destination_material(material, nuclide.name) == destination_material:
+                    matrix[i, i] = transfer_rates.get_transfer_rate(material, nuclide.name)
                 else:
                     matrix[i, i] = 0.0
             #Nothing else is allowed

openmc/deplete/transfer_rates.py

@@ -1,4 +1,5 @@
 from numbers import Real
+import re
 
 from openmc.checkvalue import check_type, check_value
 from openmc import Material
@@ -31,7 +32,8 @@ class TransferRates:
     burnable_mats : list of str
         All burnable material IDs.
     transfer_rates : dict of str to dict
-        Container of transfer rates, elements and destination material
+        Container of transfer rates, components (elements and/or nuclides) and
+        destination material
     index_transfer : Set of pair of str
         Pair of strings needed to build final matrix (destination_material, mat)
     """
@@ -74,15 +76,15 @@ def _get_material_id(self, val):
 
         return str(val)
 
-    def get_transfer_rate(self, material, element):
+    def get_transfer_rate(self, material, component):
         """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
+        component : str
+            Element or nuclide to get transfer rate value
 
         Returns
         -------
@@ -91,33 +93,34 @@ def get_transfer_rate(self, material, element):
 
         """
         material_id = self._get_material_id(material)
-        check_value('element', element, ELEMENT_SYMBOL.values())
-        return self.transfer_rates[material_id][element][0]
+        check_type('component', component, str)
+        return self.transfer_rates[material_id][component][0]
 
-    def get_destination_material(self, material, element):
-        """Return destination (or transfer) material for given material and
-        element, if defined.
+    def get_destination_material(self, material, component):
+        """Return destination material for given material and
+        component, if defined.
 
         Parameters
         ----------
         material : openmc.Material or str or int
             Depletable material
-        element : str
-            Element that gets transferred to another material.
+        component : str
+            Element or nuclide that gets transferred to another material.
 
         Returns
         -------
         destination_material_id : str
-            Depletable material ID to where the element gets transferred
+            Depletable material ID to where the element or nuclide 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]
+        check_type('component', component, str)
+        if component in self.transfer_rates[material_id]:
+            return self.transfer_rates[material_id][component][1]
 
-    def get_elements(self, material):
-        """Extract removing elements for a given material
+    def get_components(self, material):
+        """Extract removing elements and/or nuclides for a given material
 
         Parameters
         ----------
@@ -127,44 +130,48 @@ def get_elements(self, material):
         Returns
         -------
         elements : list
-            List of elements where transfer rates exist
+            List of elements and nuclides 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.
+    def set_transfer_rate(self, material, components, transfer_rate,
+                          transfer_rate_units='1/s', destination_material=None):
+        """Set element and/or nuclide 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
+        components : list of str
+            List of strings of elements and/or nuclides that share transfer rate.
+            Cannot add transfer rates for nuclides to a material where a
+            transfer rate for its element is specified and vice versa.
         transfer_rate : float
-            Rate at which elements are transferred. A positive or negative values
-            set removal of feed rates, respectively.
+            Rate at which elements and/or nuclides are transferred. A positive or
+            negative value corresponds to a removal or feed rate, 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.
+            Units for values specified in the transfer_rate argument. 's' for
+            seconds, 'min' for minutes, 'h' for hours, 'a' for Julian years.
 
         """
         material_id = self._get_material_id(material)
         check_type('transfer_rate', transfer_rate, Real)
+        check_type('components', components, list, expected_iter_type=str)
 
         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')
+                raise ValueError('Transfer to material '
+                                 f'{destination_material_id} is set, but there '
+                                 'is only one depletable material')
         else:
             destination_material_id = None
 
@@ -179,10 +186,35 @@ def set_transfer_rate(self, material, elements, transfer_rate, transfer_rate_uni
         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
+            raise ValueError('Invalid transfer rate unit '
+                             f'"{transfer_rate_units}"')
+
+        for component in components:
+            current_components = self.transfer_rates[material_id].keys()
+            split_component = re.split(r'\d+', component)
+            element = split_component[0]
+            if element not in ELEMENT_SYMBOL.values():
+                raise ValueError(f'{component} is not a valid nuclide or '
+                                 'element.')
+            else:
+                if len(split_component) == 1:
+                    element_nucs = [c for c in current_components
+                                    if re.match(component + r'\d', c)]
+                    if len(element_nucs) > 0:
+                        nuc_str = ", ".join(element_nucs)
+                        raise ValueError('Cannot add transfer rate for element '
+                                         f'{component} to material {material_id} '
+                                         f'with transfer rate(s) for nuclide(s) '
+                                         f'{nuc_str}.')
+
+                else:
+                    if element in current_components:
+                        raise ValueError('Cannot add transfer rate for nuclide '
+                                         f'{component} to material {material_id} '
+                                         f'where element {element} already has '
+                                         'a transfer rate.')
+
+            self.transfer_rates[material_id][component] = \
+                transfer_rate / unit_conv, destination_material_id
             if destination_material_id is not None:
                 self.index_transfer.add((destination_material_id, material_id))

tests/unit_tests/test_deplete_transfer_rates.py

@@ -45,29 +45,73 @@ def model():
 
     return openmc.Model(geometry, materials, settings)
 
-
-def test_get_set(model):
+@pytest.mark.parametrize("case_name, transfer_rates", [
+    ('elements', {'U': 0.01, 'Xe': 0.1}),
+    ('nuclides', {'I135': 0.01, 'Gd156': 0.1, 'Gd157': 0.01}),
+    ('nuclides_elements', {'I135': 0.01, 'Gd156': 0.1, 'Gd157': 0.01, 'U': 0.01,
+                           'Xe': 0.1}),
+    ('elements_nuclides', {'U': 0.01, 'Xe': 0.1, 'I135': 0.01, 'Gd156': 0.1,
+                           'Gd157': 0.01}),
+    ('rates_invalid_1', {'Gd': 0.01, 'Gd157': 0.01, 'Gd156': 0.01}),
+    ('rates_invalid_2', {'Gd156': 0.01, 'Gd157': 0.01, 'Gd': 0.01}),
+    ('rates_invalid_3', {'Gb156': 0.01}),
+    ('rates_invalid_4', {'Gb': 0.01})
+    ])
+def test_get_set(model, case_name, transfer_rates):
     """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()
+            if case_name == 'rates_invalid_1':
+                with pytest.raises(ValueError, match='Cannot add transfer '
+                                    'rate for nuclide Gd157 to material 1 '
+                                    'where element Gd already has a '
+                                    'transfer rate.'):
+                    for component, transfer_rate in transfer_rates.items():
+                        transfer.set_transfer_rate(material_input,
+                                                   [component],
+                                                   transfer_rate)
+            elif case_name == 'rates_invalid_2':
+                with pytest.raises(ValueError, match='Cannot add transfer '
+                                    f'rate for element Gd to material 1 with '
+                                   'transfer rate\(s\) for nuclide\(s\) '
+                                   'Gd156, Gd157.'):
+                    for component, transfer_rate in transfer_rates.items():
+                        transfer.set_transfer_rate(material_input,
+                                                   [component],
+                                                   transfer_rate)
+            elif case_name == 'rates_invalid_3':
+                with pytest.raises(ValueError, match='Gb156 is not a valid '
+                                   'nuclide or element.'):
+                    for component, transfer_rate in transfer_rates.items():
+                        transfer.set_transfer_rate(material_input,
+                                                   [component],
+                                                   transfer_rate)
+            elif case_name == 'rates_invalid_4':
+                with pytest.raises(ValueError, match='Gb is not a valid '
+                                   'nuclide or element.'):
+                    for component, transfer_rate in transfer_rates.items():
+                        transfer.set_transfer_rate(material_input,
+                                                   [component],
+                                                   transfer_rate)
+            else:
+                for component, transfer_rate in transfer_rates.items():
+                    transfer.set_transfer_rate(material_input, [component],
+                                               transfer_rate,
+                                               destination_material=\
+                                               dest_material_input)
+                    assert transfer.get_transfer_rate(
+                        material_input, component) == transfer_rate
+                    assert transfer.get_destination_material(
+                        material_input, component) == str(dest_material.id)
+                assert transfer.get_components(material_input) == \
+                    transfer_rates.keys()
 
 
 @pytest.mark.parametrize("transfer_rate_units, unit_conv", [
@@ -86,14 +130,15 @@ def test_get_set(model):
 def test_units(transfer_rate_units, unit_conv, model):
     """ Units testing"""
     # create transfer rate Xe
-    element = 'Xe'
+    components = ['Xe', 'U235']
     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
+    for component in components:
+        transfer.set_transfer_rate('f', [component], transfer_rate * unit_conv,
+                                   transfer_rate_units=transfer_rate_units)
+        assert transfer.get_transfer_rate('f', component) == transfer_rate
 
 
 def test_transfer(run_in_tmpdir, model):