Skip to content

Commit

Permalink
Allow transfer rates with nuclides. (#2564)
Browse files Browse the repository at this point in the history
Co-authored-by: Paul Romano <[email protected]>
Co-authored-by: Lorenzo Chierici <[email protected]>
  • Loading branch information
3 people committed Jun 26, 2023
1 parent 19e4c37 commit f200e44
Show file tree
Hide file tree
Showing 4 changed files with 140 additions and 56 deletions.
10 changes: 6 additions & 4 deletions openmc/deplete/abc.py
Original file line number Diff line number Diff line change
Expand Up @@ -818,16 +818,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.
Expand All @@ -841,7 +843,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
Expand Down
7 changes: 6 additions & 1 deletion openmc/deplete/chain.py
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down
100 changes: 66 additions & 34 deletions openmc/deplete/transfer_rates.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
from numbers import Real
import re

from openmc.checkvalue import check_type, check_value
from openmc import Material
Expand Down Expand Up @@ -33,7 +34,8 @@ class TransferRates:
local_mats : list of str
All burnable material IDs being managed by a single process
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)
"""
Expand Down Expand Up @@ -77,15 +79,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
-------
Expand All @@ -94,33 +96,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
----------
Expand All @@ -130,44 +133,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

Expand All @@ -182,10 +189,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))
79 changes: 62 additions & 17 deletions tests/unit_tests/test_deplete_transfer_rates.py
Original file line number Diff line number Diff line change
Expand Up @@ -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", [
Expand All @@ -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):
Expand Down

0 comments on commit f200e44

Please sign in to comment.