Skip to content

Commit 7382b5d

Browse files
church89paulromano
church89
and
paulromano
authored
External transfer rates source term (#3088)
Co-authored-by: Paul Romano <paul.k.romano@gmail.com>
1 parent f615441 commit 7382b5d

17 files changed

Lines changed: 747 additions & 180 deletions

‎docs/source/pythonapi/deplete.rst‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -206,14 +206,15 @@ total system energy.
206206
The :class:`openmc.deplete.IndependentOperator` uses inner classes subclassed
207207
from those listed above to perform similar calculations.
208208

209-
The following classes are used to define transfer rates to model continuous
210-
removal or feed of nuclides during depletion.
209+
The following classes are used to define external source rates or transfer rates
210+
to model continuous removal or feed of nuclides during depletion.
211211

212212
.. autosummary::
213213
:toctree: generated
214214
:nosignatures:
215215
:template: myclass.rst
216216

217+
transfer_rates.ExternalSourceRates
217218
transfer_rates.TransferRates
218219

219220
Intermediate Classes

‎openmc/deplete/abc.py‎

Lines changed: 89 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -29,7 +29,7 @@
2929
_SECONDS_PER_DAY, _SECONDS_PER_JULIAN_YEAR
3030
from .pool import deplete
3131
from .reaction_rates import ReactionRates
32-
from .transfer_rates import TransferRates
32+
from .transfer_rates import TransferRates, ExternalSourceRates
3333

3434

3535
__all__ = [
@@ -607,9 +607,14 @@ class Integrator(ABC):
607607
next time step. Expected to be of the same shape as ``n0``
608608
609609
transfer_rates : openmc.deplete.TransferRates
610-
Instance of TransferRates class to perform continuous transfer during depletion
610+
Transfer rates for the depletion system used to model continuous
611+
removal/feed between materials.
611612
612613
.. versionadded:: 0.14.0
614+
external_source_rates : openmc.deplete.ExternalSourceRates
615+
External source rates for the depletion system.
616+
617+
.. versionadded:: 0.15.3
613618
614619
"""
615620

@@ -686,6 +691,7 @@ def __init__(
686691
self.source_rates = np.asarray(source_rates)
687692

688693
self.transfer_rates = None
694+
self.external_source_rates = None
689695

690696
if isinstance(solver, str):
691697
# Delay importing of cram module, which requires this file
@@ -731,11 +737,11 @@ def solver(self, func):
731737

732738
self._solver = func
733739

734-
def _timed_deplete(self, n, rates, dt, matrix_func=None):
740+
def _timed_deplete(self, n, rates, dt, i=None, matrix_func=None):
735741
start = time.time()
736742
results = deplete(
737-
self._solver, self.chain, n, rates, dt, matrix_func,
738-
self.transfer_rates)
743+
self._solver, self.chain, n, rates, dt, i, matrix_func,
744+
self.transfer_rates, self.external_source_rates)
739745
return time.time() - start, results
740746

741747
@abstractmethod
@@ -885,13 +891,14 @@ def integrate(
885891
self.operator.finalize()
886892

887893
def add_transfer_rate(
888-
self,
889-
material: Union[str, int, Material],
890-
components: Sequence[str],
891-
transfer_rate: float,
892-
transfer_rate_units: str = '1/s',
893-
destination_material: Optional[Union[str, int, Material]] = None
894-
):
894+
self,
895+
material: str | int | Material,
896+
components: Sequence[str],
897+
transfer_rate: float,
898+
transfer_rate_units: str = '1/s',
899+
timesteps: Sequence[int] | None = None,
900+
destination_material: str | int | Material | None = None
901+
):
895902
"""Add transfer rates to depletable material.
896903
897904
Parameters
@@ -905,18 +912,79 @@ def add_transfer_rate(
905912
transfer_rate : float
906913
Rate at which elements are transferred. A positive or negative values
907914
set removal of feed rates, respectively.
908-
destination_material : openmc.Material or str or int, Optional
909-
Destination material to where nuclides get fed.
910915
transfer_rate_units : {'1/s', '1/min', '1/h', '1/d', '1/a'}
911916
Units for values specified in the transfer_rate argument. 's' means
912917
seconds, 'min' means minutes, 'h' means hours, 'a' means Julian years.
918+
timesteps : list of int, optional
919+
List of timestep indices where to set external source rates.
920+
Defaults to None, which means the external source rate is set for
921+
all timesteps.
922+
destination_material : openmc.Material or str or int, Optional
923+
Destination material to where nuclides get fed.
913924
914925
"""
915926
if self.transfer_rates is None:
916-
self.transfer_rates = TransferRates(self.operator, self.operator.model)
927+
if hasattr(self.operator, 'model'):
928+
materials = self.operator.model.materials
929+
elif hasattr(self.operator, 'materials'):
930+
materials = self.operator.materials
931+
self.transfer_rates = TransferRates(
932+
self.operator, materials, len(self.timesteps))
933+
934+
if self.external_source_rates is not None and destination_material:
935+
raise ValueError('Currently is not possible to set a transfer rate '
936+
'with destination matrial in combination with '
937+
'external source rates.')
938+
939+
self.transfer_rates.set_transfer_rate(
940+
material, components, transfer_rate, transfer_rate_units,
941+
timesteps, destination_material)
942+
943+
def add_external_source_rate(
944+
self,
945+
material: str | int | Material,
946+
composition: dict[str, float],
947+
rate: float,
948+
rate_units: str = 'g/s',
949+
timesteps: Sequence[int] | None = None
950+
):
951+
"""Add external source rates to depletable material.
952+
953+
Parameters
954+
----------
955+
material : openmc.Material or str or int
956+
Depletable material
957+
composition : dict of str to float
958+
External source rate composition vector, where key can be an element
959+
or a nuclide and value the corresponding weight percent.
960+
rate : float
961+
External source rate in units of mass per time. A positive or
962+
negative value corresponds to a feed or removal rate, respectively.
963+
units : {'g/s', 'g/min', 'g/h', 'g/d', 'g/a'}
964+
Units for values specified in the `rate` argument. 's' for seconds,
965+
'min' for minutes, 'h' for hours, 'a' for Julian years.
966+
timesteps : list of int, optional
967+
List of timestep indices where to set external source rates.
968+
Defaults to None, which means the external source rate is set for
969+
all timesteps.
970+
971+
"""
972+
if self.external_source_rates is None:
973+
if hasattr(self.operator, 'model'):
974+
materials = self.operator.model.materials
975+
elif hasattr(self.operator, 'materials'):
976+
materials = self.operator.materials
977+
self.external_source_rates = ExternalSourceRates(
978+
self.operator, materials, len(self.timesteps))
979+
980+
if self.transfer_rates is not None and self.transfer_rates.index_transfer:
981+
raise ValueError('Currently is not possible to set an external '
982+
'source rate in combination with transfer rates '
983+
'with destination matrial.')
984+
985+
self.external_source_rates.set_external_source_rate(
986+
material, composition, rate, rate_units, timesteps)
917987

918-
self.transfer_rates.set_transfer_rate(material, components, transfer_rate,
919-
transfer_rate_units, destination_material)
920988

921989
@add_params
922990
class SIIntegrator(Integrator):
@@ -1047,10 +1115,10 @@ def _get_bos_data_from_operator(self, step_index, step_power, n_bos):
10471115
return inherited
10481116

10491117
def integrate(
1050-
self,
1051-
output: bool = True,
1052-
path: PathLike = "depletion_results.h5"
1053-
):
1118+
self,
1119+
output: bool = True,
1120+
path: PathLike = "depletion_results.h5"
1121+
):
10541122
"""Perform the entire depletion process across all steps
10551123
10561124
Parameters

‎openmc/deplete/chain.py‎

Lines changed: 57 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -703,7 +703,7 @@ def form_matrix(self, rates, fission_yields=None):
703703
# Return CSC representation instead of DOK
704704
return matrix.tocsc()
705705

706-
def form_rr_term(self, tr_rates, mats):
706+
def form_rr_term(self, tr_rates, current_timestep, mats):
707707
"""Function to form the transfer rate term matrices.
708708
709709
.. versionadded:: 0.14.0
@@ -712,6 +712,8 @@ def form_rr_term(self, tr_rates, mats):
712712
----------
713713
tr_rates : openmc.deplete.TransferRates
714714
Instance of openmc.deplete.TransferRates
715+
current_timestep : int
716+
Current timestep index
715717
mats : string or two-tuple of strings
716718
Two cases are possible:
717719
@@ -740,32 +742,74 @@ def form_rr_term(self, tr_rates, mats):
740742

741743
for i, nuc in enumerate(self.nuclides):
742744
elm = re.split(r'\d+', nuc.name)[0]
743-
# Build transfer terms matrices
745+
# Build transfer terms (nuclide transfer only)
744746
if isinstance(mats, str):
745747
mat = mats
746-
components = tr_rates.get_components(mat)
748+
components = tr_rates.get_components(mat, current_timestep)
749+
if not components:
750+
break
747751
if elm in components:
748-
matrix[i, i] = sum(tr_rates.get_transfer_rate(mat, elm))
752+
matrix[i, i] = sum(
753+
tr_rates.get_external_rate(mat, elm, current_timestep))
749754
elif nuc.name in components:
750-
matrix[i, i] = sum(tr_rates.get_transfer_rate(mat, nuc.name))
755+
matrix[i, i] = sum(
756+
tr_rates.get_external_rate(mat, nuc.name, current_timestep))
751757
else:
752758
matrix[i, i] = 0.0
753-
#Build transfer terms matrices
759+
760+
# Build transfer terms (transfer from one material into another)
754761
elif isinstance(mats, tuple):
755762
dest_mat, mat = mats
756-
if dest_mat in tr_rates.get_destination_material(mat, elm):
757-
dest_mat_idx = tr_rates.get_destination_material(mat, elm).index(dest_mat)
758-
matrix[i, i] = tr_rates.get_transfer_rate(mat, elm)[dest_mat_idx]
759-
elif dest_mat in tr_rates.get_destination_material(mat, nuc.name):
760-
dest_mat_idx = tr_rates.get_destination_material(mat, nuc.name).index(dest_mat)
761-
matrix[i, i] = tr_rates.get_transfer_rate(mat, nuc.name)[dest_mat_idx]
763+
components = tr_rates.get_components(mat, current_timestep, dest_mat)
764+
if elm in components:
765+
matrix[i, i] = tr_rates.get_external_rate(
766+
mat, elm, current_timestep, dest_mat)[0]
767+
elif nuc.name in components:
768+
matrix[i, i] = tr_rates.get_external_rate(
769+
mat, nuc.name, current_timestep, dest_mat)[0]
762770
else:
763771
matrix[i, i] = 0.0
764-
#Nothing else is allowed
765772

766773
# Return CSC instead of DOK
767774
return matrix.tocsc()
768775

776+
def form_ext_source_term(self, ext_source_rates, current_timestep, mat):
777+
"""Function to form the external source rate term vectors.
778+
779+
.. versionadded:: 0.15.3
780+
781+
Parameters
782+
----------
783+
ext_source_rates : openmc.deplete.ExternalSourceRates
784+
Instance of openmc.deplete.ExternalSourceRates
785+
current_timestep : int
786+
Current timestep index
787+
mat : string
788+
Material id
789+
790+
Returns
791+
-------
792+
scipy.sparse.csc_matrix
793+
Sparse vector representing external source term.
794+
795+
"""
796+
if not ext_source_rates.get_components(mat, current_timestep):
797+
return
798+
# Use DOK as intermediate representation
799+
n = len(self)
800+
vector = sp.dok_matrix((n, 1))
801+
802+
for i, nuc in enumerate(self.nuclides):
803+
# Build source term vector
804+
if nuc.name in ext_source_rates.get_components(mat, current_timestep):
805+
vector[i] = sum(ext_source_rates.get_external_rate(
806+
mat, nuc.name, current_timestep))
807+
else:
808+
vector[i] = 0.0
809+
810+
# Return CSC instead of DOK
811+
return vector.tocsc()
812+
769813
def get_branch_ratios(self, reaction="(n,gamma)"):
770814
"""Return a dictionary with reaction branching ratios
771815

0 commit comments

Comments
 (0)