diff --git a/docs/source/pythonapi/base.rst b/docs/source/pythonapi/base.rst
index 3ffcac69cda..f248f3ed3c7 100644
--- a/docs/source/pythonapi/base.rst
+++ b/docs/source/pythonapi/base.rst
@@ -268,6 +268,7 @@ Variance Reduction
openmc.WeightWindows
openmc.WeightWindowsList
openmc.WeightWindowGenerator
+ openmc.WeightWindowsExodus
.. autosummary::
:toctree: generated
@@ -278,6 +279,7 @@ Variance Reduction
openmc.wwinp_to_wws
+
Coarse Mesh Finite Difference Acceleration
------------------------------------------
diff --git a/docs/source/usersguide/variance_reduction.rst b/docs/source/usersguide/variance_reduction.rst
index d551195f5cc..0f8d6ba0f63 100644
--- a/docs/source/usersguide/variance_reduction.rst
+++ b/docs/source/usersguide/variance_reduction.rst
@@ -227,6 +227,83 @@ mesh does not need to be redefined. Monte Carlo solves that load a weight window
file as above will utilize weight windows to reduce the variance of the
simulation.
+------------------------------------------------------------
+Generating Weight Windows from an Exodus II Adjoint Solution
+------------------------------------------------------------
+
+Weight windows can be directly generated from a multigroup adjoint flux
+solution stored as elemental data (as CONSTANT MONOMIAL) in an Exodus II file.
+This provides a convenient workflow for generating CADIS or FW-CADIS weight windows when the
+adjoint problem is solved using an external deterministic transport solver
+that produces Exodus II output, such as Griffin or any other MOOSE-based solver.
+
+In this workflow, the deterministic solver is used to solve the adjoint
+problem and write the resulting multigroup adjoint flux to an Exodus II
+file. OpenMC then reads the mesh and adjoint flux directly from the Exodus
+II file and constructs the corresponding weight windows during simulation
+initialization. No intermediate mesh or flux conversion is required.
+
+The Exodus II file must contain one elemental variable for each energy group,
+representing the adjoint flux for that group. The variables are specified
+using the ``adjoint_flux_variables`` input and must be listed in the same
+order as the energy groups specified by ``energy_bounds``. The energy bounds
+must be given in ascending order. Therefore, if the deterministic solver
+numbers its energy groups from the highest energy to the lowest energy
+(e.g., group 0 is the fastest group), the adjoint flux variables should be
+listed in the opposite order, from the lowest-energy group to the
+highest-energy group.
+
+The following example reads a two-group adjoint solution from
+``adjoint.e`` and generates neutron weight windows using the FW-CADIS
+normalization:
+
+.. code-block:: xml
+
+
+ adjoint.e
+ flux_g1 flux_g2
+ 0.0 1.0e5 2.0e7
+
+
+The optional ``timestep``, ``particle_type``, ``survival_ratio``,
+``upper_bound_ratio``, and ``max_split`` elements can be used to control
+the resulting weight windows. If these elements are omitted, their default
+values are used.
+
+The same functionality is available through the Python API using
+:class:`openmc.WeightWindowsExodus`:
+
+.. code-block:: python
+
+ import openmc
+
+ settings = openmc.Settings()
+ settings.weight_windows_exodus = openmc.WeightWindowsExodus(
+ file='adjoint.e',
+ adjoint_flux_variables=['flux_g1', 'flux_g2'],
+ energy_bounds=[0.0, 1.0e5, 2.0e7],
+ )
+
+The ``energy_bounds`` argument may also be specified using an
+:class:`openmc.mgxs.EnergyGroups` object.
+
+This capability can be used as part of both CADIS and FW-CADIS variance
+reduction workflows. For CADIS, the adjoint solution represents the
+importance of particles with respect to a specified tally response. For
+FW-CADIS, the adjoint solution can be used to construct weight windows that
+account for the desired spatial and energy-dependent importance distribution.
+In either case, the adjoint calculation can be performed by an external
+deterministic solver, while OpenMC uses the resulting adjoint flux to
+construct the Monte Carlo weight windows.
+
+.. warning::
+
+ When ```` is specified, it cannot be used together
+ with other weight-window specifications. OpenMC will use the weight
+ windows generated from the Exodus II adjoint solution by default.
+ Other weight-window inputs should not be provided when using
+ ````.
+
.. _source_biasing:
--------------
diff --git a/include/openmc/mesh.h b/include/openmc/mesh.h
index 1c6044514bd..88d451f7ded 100644
--- a/include/openmc/mesh.h
+++ b/include/openmc/mesh.h
@@ -1005,6 +1005,15 @@ class LibMesh : public UnstructuredMesh {
LibMesh(const std::string& filename, double length_multiplier = 1.0);
LibMesh(libMesh::MeshBase& input_mesh, double length_multiplier = 1.0);
+ //! Create a mesh from an externally constructed libMesh mesh, transferring
+ //! ownership of the mesh to OpenMC
+ //
+ //! \param[in] input_mesh Externally built mesh (must be replicated)
+ //! \param[in] length_multiplier Multiplier applied to mesh coordinates
+ //! \param[in] filename Name of the file the mesh was read from, if any
+ LibMesh(unique_ptr input_mesh,
+ double length_multiplier = 1.0, const std::string& filename = "");
+
static const std::string mesh_lib_type;
// Overridden Methods
diff --git a/include/openmc/weight_windows.h b/include/openmc/weight_windows.h
index d0b385d169d..71f413ec1c5 100644
--- a/include/openmc/weight_windows.h
+++ b/include/openmc/weight_windows.h
@@ -248,6 +248,11 @@ void apply_weight_window(Particle& p, WeightWindow weight_window);
//! Free memory associated with weight windows
void free_memory_weight_windows();
+//! Build a WeightWindows object from multigroup adjoint flux stored as
+//! elemental data in an Exodus II file (requires libMesh support)
+//! \param[in] node XML node for in settings.xml
+void read_weight_windows_exodus(pugi::xml_node node);
+
//! Search weight window that apply to a particle
//! \param[in] p Particle to search weight window for
std::pair search_weight_window(const Particle& p);
diff --git a/openmc/settings.py b/openmc/settings.py
index 3d75ee80dc2..4536a4f761e 100644
--- a/openmc/settings.py
+++ b/openmc/settings.py
@@ -17,18 +17,23 @@
from .source import SourceBase, MeshSource, IndependentSource
from .utility_funcs import input_path, set_xml_input_path
from .volume import VolumeCalculation
-from .weight_windows import WeightWindows, WeightWindowGenerator, WeightWindowsList
+from .weight_windows import (
+ WeightWindows,
+ WeightWindowGenerator,
+ WeightWindowsList,
+ WeightWindowsExodus,
+)
class RunMode(Enum):
- EIGENVALUE = 'eigenvalue'
- FIXED_SOURCE = 'fixed source'
- PLOT = 'plot'
- VOLUME = 'volume'
- PARTICLE_RESTART = 'particle restart'
+ EIGENVALUE = "eigenvalue"
+ FIXED_SOURCE = "fixed source"
+ PLOT = "plot"
+ VOLUME = "volume"
+ PARTICLE_RESTART = "particle restart"
-_RES_SCAT_METHODS = {'dbrc', 'rvs'}
+_RES_SCAT_METHODS = {"dbrc", "rvs"}
class Settings:
@@ -395,6 +400,12 @@ class Settings:
Path to a weight window file to load during simulation initialization
.. versionadded::0.14.0
+ weight_windows_exodus : openmc.WeightWindowsExodus
+ Specification for building weight windows from multigroup adjoint flux
+ stored as elemental data in an Exodus file. Requires OpenMC to be
+ built with libMesh support.
+
+ .. versionadded:: 0.15.4
write_initial_source : bool
Indicate whether to write the initial source distribution to file
"""
@@ -415,7 +426,7 @@ def __init__(self, **kwargs):
self._max_order = None
# Source subelement
- self._source = cv.CheckedList(SourceBase, 'source distributions')
+ self._source = cv.CheckedList(SourceBase, "source distributions")
self._source_rejection_fraction = None
self._confidence_intervals = None
@@ -475,7 +486,8 @@ def __init__(self, **kwargs):
self._resonance_scattering = {}
self._volume_calculations = cv.CheckedList(
- VolumeCalculation, 'volume calculations')
+ VolumeCalculation, "volume calculations"
+ )
self._create_fission_neutrons = None
self._create_delayed_neutrons = None
@@ -489,10 +501,12 @@ def __init__(self, **kwargs):
self._write_initial_source = None
self._weight_windows = WeightWindowsList()
self._weight_window_generators = cv.CheckedList(
- WeightWindowGenerator, 'weight window generators')
+ WeightWindowGenerator, "weight window generators"
+ )
self._weight_windows_on = None
self._shared_secondary_bank = None
self._weight_windows_file = None
+ self._weight_windows_exodus = None
self._weight_window_checkpoints = {}
self._max_history_splits = None
self._max_tracks = None
@@ -505,11 +519,11 @@ def __init__(self, **kwargs):
setattr(self, key, value)
def __setattr__(self, name: str, value):
- if not name.startswith('_'):
+ if not name.startswith("_"):
try:
getattr(self, name)
except AttributeError as e:
- msg, = traceback.format_exception_only(e)
+ (msg,) = traceback.format_exception_only(e)
msg = msg.strip().split(maxsplit=1)[-1]
warnings.warn(msg, stacklevel=2)
super().__setattr__(name, value)
@@ -520,7 +534,7 @@ def run_mode(self) -> str:
@run_mode.setter
def run_mode(self, run_mode: str):
- cv.check_value('run mode', run_mode, {x.value for x in RunMode})
+ cv.check_value("run mode", run_mode, {x.value for x in RunMode})
for mode in RunMode:
if mode.value == run_mode:
self._run_mode = mode
@@ -531,8 +545,8 @@ def batches(self) -> int:
@batches.setter
def batches(self, batches: int):
- cv.check_type('batches', batches, Integral)
- cv.check_greater_than('batches', batches, 0)
+ cv.check_type("batches", batches, Integral)
+ cv.check_greater_than("batches", batches, 0)
self._batches = batches
@property
@@ -541,9 +555,8 @@ def generations_per_batch(self) -> int:
@generations_per_batch.setter
def generations_per_batch(self, generations_per_batch: int):
- cv.check_type('generations per batch', generations_per_batch, Integral)
- cv.check_greater_than('generations per batch',
- generations_per_batch, 0)
+ cv.check_type("generations per batch", generations_per_batch, Integral)
+ cv.check_greater_than("generations per batch", generations_per_batch, 0)
self._generations_per_batch = generations_per_batch
@property
@@ -552,8 +565,8 @@ def inactive(self) -> int:
@inactive.setter
def inactive(self, inactive: int):
- cv.check_type('inactive batches', inactive, Integral)
- cv.check_greater_than('inactive batches', inactive, 0, True)
+ cv.check_type("inactive batches", inactive, Integral)
+ cv.check_greater_than("inactive batches", inactive, 0, True)
self._inactive = inactive
@property
@@ -562,8 +575,8 @@ def max_lost_particles(self) -> int:
@max_lost_particles.setter
def max_lost_particles(self, max_lost_particles: int):
- cv.check_type('max_lost_particles', max_lost_particles, Integral)
- cv.check_greater_than('max_lost_particles', max_lost_particles, 0)
+ cv.check_type("max_lost_particles", max_lost_particles, Integral)
+ cv.check_greater_than("max_lost_particles", max_lost_particles, 0)
self._max_lost_particles = max_lost_particles
@property
@@ -572,10 +585,9 @@ def rel_max_lost_particles(self) -> float:
@rel_max_lost_particles.setter
def rel_max_lost_particles(self, rel_max_lost_particles: float):
- cv.check_type('rel_max_lost_particles', rel_max_lost_particles, Real)
- cv.check_greater_than('rel_max_lost_particles',
- rel_max_lost_particles, 0)
- cv.check_less_than('rel_max_lost_particles', rel_max_lost_particles, 1)
+ cv.check_type("rel_max_lost_particles", rel_max_lost_particles, Real)
+ cv.check_greater_than("rel_max_lost_particles", rel_max_lost_particles, 0)
+ cv.check_less_than("rel_max_lost_particles", rel_max_lost_particles, 1)
self._rel_max_lost_particles = rel_max_lost_particles
@property
@@ -584,10 +596,8 @@ def max_write_lost_particles(self) -> int:
@max_write_lost_particles.setter
def max_write_lost_particles(self, max_write_lost_particles: int):
- cv.check_type('max_write_lost_particles',
- max_write_lost_particles, Integral)
- cv.check_greater_than('max_write_lost_particles',
- max_write_lost_particles, 0)
+ cv.check_type("max_write_lost_particles", max_write_lost_particles, Integral)
+ cv.check_greater_than("max_write_lost_particles", max_write_lost_particles, 0)
self._max_write_lost_particles = max_write_lost_particles
@property
@@ -596,8 +606,8 @@ def particles(self) -> int:
@particles.setter
def particles(self, particles: int):
- cv.check_type('particles', particles, Integral)
- cv.check_greater_than('particles', particles, 0)
+ cv.check_type("particles", particles, Integral)
+ cv.check_greater_than("particles", particles, 0)
self._particles = particles
@property
@@ -607,28 +617,36 @@ def keff_trigger(self) -> dict:
@keff_trigger.setter
def keff_trigger(self, keff_trigger: dict):
if not isinstance(keff_trigger, dict):
- msg = f'Unable to set a trigger on keff from "{keff_trigger}" ' \
- 'which is not a Python dictionary'
+ msg = (
+ f'Unable to set a trigger on keff from "{keff_trigger}" '
+ "which is not a Python dictionary"
+ )
raise ValueError(msg)
- elif 'type' not in keff_trigger:
- msg = f'Unable to set a trigger on keff from "{keff_trigger}" ' \
+ elif "type" not in keff_trigger:
+ msg = (
+ f'Unable to set a trigger on keff from "{keff_trigger}" '
'which does not have a "type" key'
+ )
raise ValueError(msg)
- elif keff_trigger['type'] not in ['variance', 'std_dev', 'rel_err']:
- msg = 'Unable to set a trigger on keff with ' \
- 'type "{0}"'.format(keff_trigger['type'])
+ elif keff_trigger["type"] not in ["variance", "std_dev", "rel_err"]:
+ msg = "Unable to set a trigger on keff with " 'type "{0}"'.format(
+ keff_trigger["type"]
+ )
raise ValueError(msg)
- elif 'threshold' not in keff_trigger:
- msg = f'Unable to set a trigger on keff from "{keff_trigger}" ' \
+ elif "threshold" not in keff_trigger:
+ msg = (
+ f'Unable to set a trigger on keff from "{keff_trigger}" '
'which does not have a "threshold" key'
+ )
raise ValueError(msg)
- elif not isinstance(keff_trigger['threshold'], Real):
- msg = 'Unable to set a trigger on keff with ' \
- 'threshold "{0}"'.format(keff_trigger['threshold'])
+ elif not isinstance(keff_trigger["threshold"], Real):
+ msg = "Unable to set a trigger on keff with " 'threshold "{0}"'.format(
+ keff_trigger["threshold"]
+ )
raise ValueError(msg)
self._keff_trigger = keff_trigger
@@ -639,8 +657,7 @@ def energy_mode(self) -> str:
@energy_mode.setter
def energy_mode(self, energy_mode: str):
- cv.check_value('energy mode', energy_mode,
- ['continuous-energy', 'multi-group'])
+ cv.check_value("energy mode", energy_mode, ["continuous-energy", "multi-group"])
self._energy_mode = energy_mode
@property
@@ -650,9 +667,8 @@ def max_order(self) -> int:
@max_order.setter
def max_order(self, max_order: int | None):
if max_order is not None:
- cv.check_type('maximum scattering order', max_order, Integral)
- cv.check_greater_than('maximum scattering order', max_order, 0,
- True)
+ cv.check_type("maximum scattering order", max_order, Integral)
+ cv.check_greater_than("maximum scattering order", max_order, 0, True)
self._max_order = max_order
@property
@@ -663,8 +679,7 @@ def source(self) -> list[SourceBase]:
def source(self, source: SourceBase | Iterable[SourceBase]):
if not isinstance(source, MutableSequence):
source = [source]
- self._source = cv.CheckedList(
- SourceBase, 'source distributions', source)
+ self._source = cv.CheckedList(SourceBase, "source distributions", source)
@property
def confidence_intervals(self) -> bool:
@@ -672,7 +687,7 @@ def confidence_intervals(self) -> bool:
@confidence_intervals.setter
def confidence_intervals(self, confidence_intervals: bool):
- cv.check_type('confidence interval', confidence_intervals, bool)
+ cv.check_type("confidence interval", confidence_intervals, bool)
self._confidence_intervals = confidence_intervals
@property
@@ -681,8 +696,7 @@ def electron_treatment(self) -> str:
@electron_treatment.setter
def electron_treatment(self, electron_treatment: str):
- cv.check_value('electron treatment',
- electron_treatment, ['led', 'ttb'])
+ cv.check_value("electron treatment", electron_treatment, ["led", "ttb"])
self._electron_treatment = electron_treatment
@property
@@ -691,7 +705,7 @@ def atomic_relaxation(self) -> bool:
@atomic_relaxation.setter
def atomic_relaxation(self, atomic_relaxation: bool):
- cv.check_type('atomic relaxation', atomic_relaxation, bool)
+ cv.check_type("atomic relaxation", atomic_relaxation, bool)
self._atomic_relaxation = atomic_relaxation
@property
@@ -700,7 +714,7 @@ def ptables(self) -> bool:
@ptables.setter
def ptables(self, ptables: bool):
- cv.check_type('probability tables', ptables, bool)
+ cv.check_type("probability tables", ptables, bool)
self._ptables = ptables
@property
@@ -709,7 +723,7 @@ def photon_transport(self) -> bool:
@photon_transport.setter
def photon_transport(self, photon_transport: bool):
- cv.check_type('photon transport', photon_transport, bool)
+ cv.check_type("photon transport", photon_transport, bool)
self._photon_transport = photon_transport
@property
@@ -718,7 +732,7 @@ def uniform_source_sampling(self) -> bool:
@uniform_source_sampling.setter
def uniform_source_sampling(self, uniform_source_sampling: bool):
- cv.check_type('strength as weights', uniform_source_sampling, bool)
+ cv.check_type("strength as weights", uniform_source_sampling, bool)
self._uniform_source_sampling = uniform_source_sampling
@property
@@ -727,8 +741,8 @@ def plot_seed(self):
@plot_seed.setter
def plot_seed(self, seed):
- cv.check_type('random plot color seed', seed, Integral)
- cv.check_greater_than('random plot color seed', seed, 0)
+ cv.check_type("random plot color seed", seed, Integral)
+ cv.check_greater_than("random plot color seed", seed, 0)
self._plot_seed = seed
@property
@@ -737,8 +751,8 @@ def seed(self) -> int:
@seed.setter
def seed(self, seed: int):
- cv.check_type('random number generator seed', seed, Integral)
- cv.check_greater_than('random number generator seed', seed, 0)
+ cv.check_type("random number generator seed", seed, Integral)
+ cv.check_greater_than("random number generator seed", seed, 0)
self._seed = seed
@property
@@ -747,8 +761,8 @@ def stride(self) -> int:
@stride.setter
def stride(self, stride: int):
- cv.check_type('random number generator stride', stride, Integral)
- cv.check_greater_than('random number generator stride', stride, 0)
+ cv.check_type("random number generator stride", stride, Integral)
+ cv.check_greater_than("random number generator stride", stride, 0)
self._stride = stride
@property
@@ -757,9 +771,9 @@ def surface_grazing_cutoff(self) -> float:
@surface_grazing_cutoff.setter
def surface_grazing_cutoff(self, surface_grazing_cutoff: float):
- cv.check_type('surface grazing cutoff', surface_grazing_cutoff, float)
- cv.check_greater_than('surface grazing cutoff', surface_grazing_cutoff, 0.0)
- cv.check_less_than('surface grazing cutoff', surface_grazing_cutoff, 1.0)
+ cv.check_type("surface grazing cutoff", surface_grazing_cutoff, float)
+ cv.check_greater_than("surface grazing cutoff", surface_grazing_cutoff, 0.0)
+ cv.check_less_than("surface grazing cutoff", surface_grazing_cutoff, 1.0)
self._surface_grazing_cutoff = surface_grazing_cutoff
@property
@@ -768,8 +782,8 @@ def surface_grazing_ratio(self) -> float:
@surface_grazing_ratio.setter
def surface_grazing_ratio(self, surface_grazing_ratio: float):
- cv.check_type('surface grazing ratio', surface_grazing_ratio, float)
- cv.check_greater_than('surface grazing ratio', surface_grazing_ratio, 0.0)
+ cv.check_type("surface grazing ratio", surface_grazing_ratio, float)
+ cv.check_greater_than("surface grazing ratio", surface_grazing_ratio, 0.0)
self._surface_grazing_ratio = surface_grazing_ratio
@property
@@ -778,7 +792,7 @@ def survival_biasing(self) -> bool:
@survival_biasing.setter
def survival_biasing(self, survival_biasing: bool):
- cv.check_type('survival biasing', survival_biasing, bool)
+ cv.check_type("survival biasing", survival_biasing, bool)
self._survival_biasing = survival_biasing
@property
@@ -787,7 +801,7 @@ def entropy_mesh(self) -> RegularMesh:
@entropy_mesh.setter
def entropy_mesh(self, entropy: RegularMesh):
- cv.check_type('entropy mesh', entropy, RegularMesh)
+ cv.check_type("entropy mesh", entropy, RegularMesh)
self._entropy_mesh = entropy
@property
@@ -796,7 +810,7 @@ def trigger_active(self) -> bool:
@trigger_active.setter
def trigger_active(self, trigger_active: bool):
- cv.check_type('trigger active', trigger_active, bool)
+ cv.check_type("trigger active", trigger_active, bool)
self._trigger_active = trigger_active
@property
@@ -805,9 +819,8 @@ def trigger_max_batches(self) -> int:
@trigger_max_batches.setter
def trigger_max_batches(self, trigger_max_batches: int):
- cv.check_type('trigger maximum batches', trigger_max_batches, Integral)
- cv.check_greater_than('trigger maximum batches',
- trigger_max_batches, 0)
+ cv.check_type("trigger maximum batches", trigger_max_batches, Integral)
+ cv.check_greater_than("trigger maximum batches", trigger_max_batches, 0)
self._trigger_max_batches = trigger_max_batches
@property
@@ -816,10 +829,8 @@ def trigger_batch_interval(self) -> int:
@trigger_batch_interval.setter
def trigger_batch_interval(self, trigger_batch_interval: int):
- cv.check_type('trigger batch interval',
- trigger_batch_interval, Integral)
- cv.check_greater_than('trigger batch interval',
- trigger_batch_interval, 0)
+ cv.check_type("trigger batch interval", trigger_batch_interval, Integral)
+ cv.check_greater_than("trigger batch interval", trigger_batch_interval, 0)
self._trigger_batch_interval = trigger_batch_interval
@property
@@ -828,10 +839,10 @@ def output(self) -> dict:
@output.setter
def output(self, output: dict):
- cv.check_type('output', output, Mapping)
+ cv.check_type("output", output, Mapping)
for key, value in output.items():
- cv.check_value('output key', key, ('summary', 'tallies', 'path'))
- if key in ('summary', 'tallies'):
+ cv.check_value("output key", key, ("summary", "tallies", "path"))
+ if key in ("summary", "tallies"):
cv.check_type(f"output['{key}']", value, bool)
else:
cv.check_type("output['path']", value, str)
@@ -843,23 +854,25 @@ def sourcepoint(self) -> dict:
@sourcepoint.setter
def sourcepoint(self, sourcepoint: dict):
- cv.check_type('sourcepoint options', sourcepoint, Mapping)
+ cv.check_type("sourcepoint options", sourcepoint, Mapping)
for key, value in sourcepoint.items():
- if key == 'batches':
- cv.check_type('sourcepoint batches', value, Iterable, Integral)
+ if key == "batches":
+ cv.check_type("sourcepoint batches", value, Iterable, Integral)
for batch in value:
- cv.check_greater_than('sourcepoint batch', batch, 0)
- elif key == 'separate':
- cv.check_type('sourcepoint separate', value, bool)
- elif key == 'write':
- cv.check_type('sourcepoint write', value, bool)
- elif key == 'overwrite':
- cv.check_type('sourcepoint overwrite', value, bool)
- elif key == 'mcpl':
- cv.check_type('sourcepoint mcpl', value, bool)
+ cv.check_greater_than("sourcepoint batch", batch, 0)
+ elif key == "separate":
+ cv.check_type("sourcepoint separate", value, bool)
+ elif key == "write":
+ cv.check_type("sourcepoint write", value, bool)
+ elif key == "overwrite":
+ cv.check_type("sourcepoint overwrite", value, bool)
+ elif key == "mcpl":
+ cv.check_type("sourcepoint mcpl", value, bool)
else:
- raise ValueError(f"Unknown key '{key}' encountered when "
- "setting sourcepoint options.")
+ raise ValueError(
+ f"Unknown key '{key}' encountered when "
+ "setting sourcepoint options."
+ )
self._sourcepoint = sourcepoint
@property
@@ -868,15 +881,17 @@ def statepoint(self) -> dict:
@statepoint.setter
def statepoint(self, statepoint: dict):
- cv.check_type('statepoint options', statepoint, Mapping)
+ cv.check_type("statepoint options", statepoint, Mapping)
for key, value in statepoint.items():
- if key == 'batches':
- cv.check_type('statepoint batches', value, Iterable, Integral)
+ if key == "batches":
+ cv.check_type("statepoint batches", value, Iterable, Integral)
for batch in value:
- cv.check_greater_than('statepoint batch', batch, 0)
+ cv.check_greater_than("statepoint batch", batch, 0)
else:
- raise ValueError(f"Unknown key '{key}' encountered when "
- "setting statepoint options.")
+ raise ValueError(
+ f"Unknown key '{key}' encountered when "
+ "setting statepoint options."
+ )
self._statepoint = statepoint
@property
@@ -885,17 +900,16 @@ def surf_source_read(self) -> dict:
@surf_source_read.setter
def surf_source_read(self, ssr: dict):
- cv.check_type('surface source reading options', ssr, Mapping)
+ cv.check_type("surface source reading options", ssr, Mapping)
for key, value in ssr.items():
- cv.check_value('surface source reading key', key,
- ('path'))
- if key == 'path':
- cv.check_type('path to surface source file', value, PathLike)
+ cv.check_value("surface source reading key", key, ("path"))
+ if key == "path":
+ cv.check_type("path to surface source file", value, PathLike)
self._surf_source_read = dict(ssr)
# Resolve path to surface source file
- if 'path' in ssr:
- self._surf_source_read['path'] = input_path(ssr['path'])
+ if "path" in ssr:
+ self._surf_source_read["path"] = input_path(ssr["path"])
@property
def surf_source_write(self) -> dict:
@@ -903,26 +917,37 @@ def surf_source_write(self) -> dict:
@surf_source_write.setter
def surf_source_write(self, surf_source_write: dict):
- cv.check_type("surface source writing options",
- surf_source_write, Mapping)
+ cv.check_type("surface source writing options", surf_source_write, Mapping)
for key, value in surf_source_write.items():
cv.check_value(
"surface source writing key",
key,
- ("surface_ids", "max_particles", "max_source_files",
- "mcpl", "cell", "cellfrom", "cellto"),
+ (
+ "surface_ids",
+ "max_particles",
+ "max_source_files",
+ "mcpl",
+ "cell",
+ "cellfrom",
+ "cellto",
+ ),
)
if key == "surface_ids":
cv.check_type(
"surface ids for source banking", value, Iterable, Integral
)
for surf_id in value:
- cv.check_greater_than(
- "surface id for source banking", surf_id, 0)
+ cv.check_greater_than("surface id for source banking", surf_id, 0)
elif key == "mcpl":
cv.check_type("write to an MCPL-format file", value, bool)
- elif key in ("max_particles", "max_source_files", "cell", "cellfrom", "cellto"):
+ elif key in (
+ "max_particles",
+ "max_source_files",
+ "cell",
+ "cellfrom",
+ "cellto",
+ ):
name = {
"max_particles": "maximum particle banks on surfaces per process",
"max_source_files": "maximun surface source files to be written",
@@ -941,24 +966,42 @@ def collision_track(self) -> dict:
@collision_track.setter
def collision_track(self, collision_track: dict):
- cv.check_type('Collision tracking options', collision_track, Mapping)
+ cv.check_type("Collision tracking options", collision_track, Mapping)
for key, value in collision_track.items():
- cv.check_value('collision_track key', key,
- ('cell_ids', 'reactions', 'universe_ids', 'material_ids', 'nuclides',
- 'deposited_E_threshold', 'max_collisions', 'max_collision_track_files', 'mcpl'))
- if key == 'cell_ids':
- cv.check_type('cell ids for collision tracking data banking', value,
- Iterable, Integral)
+ cv.check_value(
+ "collision_track key",
+ key,
+ (
+ "cell_ids",
+ "reactions",
+ "universe_ids",
+ "material_ids",
+ "nuclides",
+ "deposited_E_threshold",
+ "max_collisions",
+ "max_collision_track_files",
+ "mcpl",
+ ),
+ )
+ if key == "cell_ids":
+ cv.check_type(
+ "cell ids for collision tracking data banking",
+ value,
+ Iterable,
+ Integral,
+ )
for cell_id in value:
- cv.check_greater_than('cell id for collision tracking data banking',
- cell_id, 0)
- elif key == 'reactions':
- cv.check_type('MT numbers for collision tracking data banking', value,
- Iterable)
+ cv.check_greater_than(
+ "cell id for collision tracking data banking", cell_id, 0
+ )
+ elif key == "reactions":
+ cv.check_type(
+ "MT numbers for collision tracking data banking", value, Iterable
+ )
for reaction in value:
if isinstance(reaction, int):
cv.check_greater_than(
- 'MT number for collision tracking data banking', reaction, 0
+ "MT number for collision tracking data banking", reaction, 0
)
elif isinstance(reaction, str):
# check against allowed strings? so far let C++ code handle it
@@ -966,45 +1009,67 @@ def collision_track(self, collision_track: dict):
else:
raise TypeError(
f"MT number for collision tracking data banking must be a positive int or string, "
- f"got {type(reaction).__name__}")
- elif key == 'universe_ids':
- cv.check_type('universe ids for collision tracking data banking', value,
- Iterable, Integral)
+ f"got {type(reaction).__name__}"
+ )
+ elif key == "universe_ids":
+ cv.check_type(
+ "universe ids for collision tracking data banking",
+ value,
+ Iterable,
+ Integral,
+ )
for universe_id in value:
- cv.check_greater_than('universe id for collision tracking data banking',
- universe_id, 0)
- elif key == 'material_ids':
- cv.check_type('material ids for collision tracking data banking', value,
- Iterable, Integral)
+ cv.check_greater_than(
+ "universe id for collision tracking data banking",
+ universe_id,
+ 0,
+ )
+ elif key == "material_ids":
+ cv.check_type(
+ "material ids for collision tracking data banking",
+ value,
+ Iterable,
+ Integral,
+ )
for material_id in value:
- cv.check_greater_than('material id for collision tracking data banking',
- material_id, 0)
- elif key == 'nuclides':
- cv.check_type('nuclides for collision tracking data banking', value,
- Iterable, str)
+ cv.check_greater_than(
+ "material id for collision tracking data banking",
+ material_id,
+ 0,
+ )
+ elif key == "nuclides":
+ cv.check_type(
+ "nuclides for collision tracking data banking", value, Iterable, str
+ )
for nuclide in value:
# If nuclide name doesn't look valid, give a warning
try:
openmc.data.zam(nuclide)
except ValueError:
warnings.warn(f"Nuclide {nuclide} is not valid")
- elif key == 'deposited_E_threshold':
- cv.check_type('Deposited Energy Threshold for collision tracking data banking',
- value, Real)
- cv.check_greater_than('Deposited Energy Threshold for collision tracking data banking',
- value, 0)
- elif key == 'max_collisions':
- cv.check_type('maximum collisions banks per file',
- value, Integral)
- cv.check_greater_than('maximum collisions banks in collision tracking',
- value, 0)
- elif key == 'max_collision_track_files':
- cv.check_type('maximum collisions banks',
- value, Integral)
- cv.check_greater_than('maximum number of collision_track files ',
- value, 0)
- elif key == 'mcpl':
- cv.check_type('write to an MCPL-format file', value, bool)
+ elif key == "deposited_E_threshold":
+ cv.check_type(
+ "Deposited Energy Threshold for collision tracking data banking",
+ value,
+ Real,
+ )
+ cv.check_greater_than(
+ "Deposited Energy Threshold for collision tracking data banking",
+ value,
+ 0,
+ )
+ elif key == "max_collisions":
+ cv.check_type("maximum collisions banks per file", value, Integral)
+ cv.check_greater_than(
+ "maximum collisions banks in collision tracking", value, 0
+ )
+ elif key == "max_collision_track_files":
+ cv.check_type("maximum collisions banks", value, Integral)
+ cv.check_greater_than(
+ "maximum number of collision_track files ", value, 0
+ )
+ elif key == "mcpl":
+ cv.check_type("write to an MCPL-format file", value, bool)
self._collision_track = collision_track
@@ -1014,7 +1079,7 @@ def no_reduce(self) -> bool:
@no_reduce.setter
def no_reduce(self, no_reduce: bool):
- cv.check_type('no reduction option', no_reduce, bool)
+ cv.check_type("no reduction option", no_reduce, bool)
self._no_reduce = no_reduce
@property
@@ -1023,9 +1088,9 @@ def verbosity(self) -> int:
@verbosity.setter
def verbosity(self, verbosity: int):
- cv.check_type('verbosity', verbosity, Integral)
- cv.check_greater_than('verbosity', verbosity, 1, True)
- cv.check_less_than('verbosity', verbosity, 10, True)
+ cv.check_type("verbosity", verbosity, Integral)
+ cv.check_greater_than("verbosity", verbosity, 1, True)
+ cv.check_less_than("verbosity", verbosity, 10, True)
self._verbosity = verbosity
@property
@@ -1045,15 +1110,14 @@ def tabular_legendre(self) -> dict:
@tabular_legendre.setter
def tabular_legendre(self, tabular_legendre: dict):
- cv.check_type('tabular_legendre settings', tabular_legendre, Mapping)
+ cv.check_type("tabular_legendre settings", tabular_legendre, Mapping)
for key, value in tabular_legendre.items():
- cv.check_value('tabular_legendre key', key,
- ['enable', 'num_points'])
- if key == 'enable':
- cv.check_type('enable tabular_legendre', value, bool)
- elif key == 'num_points':
- cv.check_type('num_points tabular_legendre', value, Integral)
- cv.check_greater_than('num_points tabular_legendre', value, 0)
+ cv.check_value("tabular_legendre key", key, ["enable", "num_points"])
+ if key == "enable":
+ cv.check_type("enable tabular_legendre", value, bool)
+ elif key == "num_points":
+ cv.check_type("num_points tabular_legendre", value, Integral)
+ cv.check_greater_than("num_points tabular_legendre", value, 0)
self._tabular_legendre = tabular_legendre
@property
@@ -1063,24 +1127,27 @@ def temperature(self) -> dict:
@temperature.setter
def temperature(self, temperature: dict):
- cv.check_type('temperature settings', temperature, Mapping)
+ cv.check_type("temperature settings", temperature, Mapping)
for key, value in temperature.items():
- cv.check_value('temperature key', key,
- ['default', 'method', 'tolerance', 'multipole',
- 'range'])
- if key == 'default':
- cv.check_type('default temperature', value, Real)
- elif key == 'method':
- cv.check_value('temperature method', value,
- ['nearest', 'interpolation'])
- elif key == 'tolerance':
- cv.check_type('temperature tolerance', value, Real)
- elif key == 'multipole':
- cv.check_type('temperature multipole', value, bool)
- elif key == 'range':
- cv.check_length('temperature range', value, 2)
+ cv.check_value(
+ "temperature key",
+ key,
+ ["default", "method", "tolerance", "multipole", "range"],
+ )
+ if key == "default":
+ cv.check_type("default temperature", value, Real)
+ elif key == "method":
+ cv.check_value(
+ "temperature method", value, ["nearest", "interpolation"]
+ )
+ elif key == "tolerance":
+ cv.check_type("temperature tolerance", value, Real)
+ elif key == "multipole":
+ cv.check_type("temperature multipole", value, bool)
+ elif key == "range":
+ cv.check_length("temperature range", value, 2)
for T in value:
- cv.check_type('temperature', T, Real)
+ cv.check_type("temperature", T, Real)
self._temperature = temperature
@@ -1093,7 +1160,7 @@ def properties_file(self, value: PathLike | None):
if value is None:
self._properties_file = None
else:
- cv.check_type('properties file', value, PathLike)
+ cv.check_type("properties file", value, PathLike)
self._properties_file = input_path(value)
@property
@@ -1102,11 +1169,11 @@ def trace(self) -> Iterable:
@trace.setter
def trace(self, trace: Iterable):
- cv.check_type('trace', trace, Iterable, Integral)
- cv.check_length('trace', trace, 3)
- cv.check_greater_than('trace batch', trace[0], 0)
- cv.check_greater_than('trace generation', trace[1], 0)
- cv.check_greater_than('trace particle', trace[2], 0)
+ cv.check_type("trace", trace, Iterable, Integral)
+ cv.check_length("trace", trace, 3)
+ cv.check_greater_than("trace batch", trace[0], 0)
+ cv.check_greater_than("trace generation", trace[1], 0)
+ cv.check_greater_than("trace particle", trace[2], 0)
self._trace = trace
@property
@@ -1115,17 +1182,17 @@ def track(self) -> Iterable[Iterable[int]]:
@track.setter
def track(self, track: Iterable[Iterable[int]]):
- cv.check_type('track', track, Sequence)
+ cv.check_type("track", track, Sequence)
for t in track:
if len(t) != 3:
msg = f'Unable to set the track to "{t}" since its length is not 3'
raise ValueError(msg)
- cv.check_greater_than('track batch', t[0], 0)
- cv.check_greater_than('track generation', t[1], 0)
- cv.check_greater_than('track particle', t[2], 0)
- cv.check_type('track batch', t[0], Integral)
- cv.check_type('track generation', t[1], Integral)
- cv.check_type('track particle', t[2], Integral)
+ cv.check_greater_than("track batch", t[0], 0)
+ cv.check_greater_than("track generation", t[1], 0)
+ cv.check_greater_than("track particle", t[2], 0)
+ cv.check_type("track batch", t[0], Integral)
+ cv.check_type("track generation", t[1], Integral)
+ cv.check_type("track particle", t[2], Integral)
self._track = track
@property
@@ -1135,26 +1202,32 @@ def cutoff(self) -> dict:
@cutoff.setter
def cutoff(self, cutoff: dict):
if not isinstance(cutoff, Mapping):
- msg = f'Unable to set cutoff from "{cutoff}" which is not a '\
- 'Python dictionary'
+ msg = (
+ f'Unable to set cutoff from "{cutoff}" which is not a '
+ "Python dictionary"
+ )
raise ValueError(msg)
for key in cutoff:
- if key == 'weight':
- cv.check_type('weight cutoff', cutoff[key], Real)
- cv.check_greater_than('weight cutoff', cutoff[key], 0.0)
- elif key == 'weight_avg':
- cv.check_type('average survival weight', cutoff[key], Real)
- cv.check_greater_than('average survival weight',
- cutoff[key], 0.0)
- elif key == 'survival_normalization':
- cv.check_type('survival normalization', cutoff[key], bool)
- elif key in ['energy_neutron', 'energy_photon', 'energy_electron',
- 'energy_positron']:
- cv.check_type('energy cutoff', cutoff[key], Real)
- cv.check_greater_than('energy cutoff', cutoff[key], 0.0)
+ if key == "weight":
+ cv.check_type("weight cutoff", cutoff[key], Real)
+ cv.check_greater_than("weight cutoff", cutoff[key], 0.0)
+ elif key == "weight_avg":
+ cv.check_type("average survival weight", cutoff[key], Real)
+ cv.check_greater_than("average survival weight", cutoff[key], 0.0)
+ elif key == "survival_normalization":
+ cv.check_type("survival normalization", cutoff[key], bool)
+ elif key in [
+ "energy_neutron",
+ "energy_photon",
+ "energy_electron",
+ "energy_positron",
+ ]:
+ cv.check_type("energy cutoff", cutoff[key], Real)
+ cv.check_greater_than("energy cutoff", cutoff[key], 0.0)
else:
- msg = f'Unable to set cutoff to "{key}" which is unsupported ' \
- 'by OpenMC'
+ msg = (
+ f'Unable to set cutoff to "{key}" which is unsupported ' "by OpenMC"
+ )
self._cutoff = cutoff
@@ -1164,10 +1237,10 @@ def ufs_mesh(self) -> RegularMesh:
@ufs_mesh.setter
def ufs_mesh(self, ufs_mesh: RegularMesh):
- cv.check_type('UFS mesh', ufs_mesh, RegularMesh)
- cv.check_length('UFS mesh dimension', ufs_mesh.dimension, 3)
- cv.check_length('UFS mesh lower-left corner', ufs_mesh.lower_left, 3)
- cv.check_length('UFS mesh upper-right corner', ufs_mesh.upper_right, 3)
+ cv.check_type("UFS mesh", ufs_mesh, RegularMesh)
+ cv.check_length("UFS mesh dimension", ufs_mesh.dimension, 3)
+ cv.check_length("UFS mesh lower-left corner", ufs_mesh.lower_left, 3)
+ cv.check_length("UFS mesh upper-right corner", ufs_mesh.upper_right, 3)
self._ufs_mesh = ufs_mesh
@property
@@ -1176,26 +1249,24 @@ def resonance_scattering(self) -> dict:
@resonance_scattering.setter
def resonance_scattering(self, res: dict):
- cv.check_type('resonance scattering settings', res, Mapping)
- keys = ('enable', 'method', 'energy_min', 'energy_max', 'nuclides')
+ cv.check_type("resonance scattering settings", res, Mapping)
+ keys = ("enable", "method", "energy_min", "energy_max", "nuclides")
for key, value in res.items():
- cv.check_value('resonance scattering dictionary key', key, keys)
- if key == 'enable':
- cv.check_type('resonance scattering enable', value, bool)
- elif key == 'method':
- cv.check_value('resonance scattering method', value,
- _RES_SCAT_METHODS)
- elif key == 'energy_min':
- name = 'resonance scattering minimum energy'
+ cv.check_value("resonance scattering dictionary key", key, keys)
+ if key == "enable":
+ cv.check_type("resonance scattering enable", value, bool)
+ elif key == "method":
+ cv.check_value("resonance scattering method", value, _RES_SCAT_METHODS)
+ elif key == "energy_min":
+ name = "resonance scattering minimum energy"
cv.check_type(name, value, Real)
cv.check_greater_than(name, value, 0)
- elif key == 'energy_max':
- name = 'resonance scattering minimum energy'
+ elif key == "energy_max":
+ name = "resonance scattering minimum energy"
cv.check_type(name, value, Real)
cv.check_greater_than(name, value, 0)
- elif key == 'nuclides':
- cv.check_type('resonance scattering nuclides', value,
- Iterable, str)
+ elif key == "nuclides":
+ cv.check_type("resonance scattering nuclides", value, Iterable, str)
self._resonance_scattering = res
@property
@@ -1209,7 +1280,8 @@ def volume_calculations(
if not isinstance(vol_calcs, MutableSequence):
vol_calcs = [vol_calcs]
self._volume_calculations = cv.CheckedList(
- VolumeCalculation, 'stochastic volume calculations', vol_calcs)
+ VolumeCalculation, "stochastic volume calculations", vol_calcs
+ )
@property
def create_fission_neutrons(self) -> bool:
@@ -1217,8 +1289,7 @@ def create_fission_neutrons(self) -> bool:
@create_fission_neutrons.setter
def create_fission_neutrons(self, create_fission_neutrons: bool):
- cv.check_type('Whether create fission neutrons',
- create_fission_neutrons, bool)
+ cv.check_type("Whether create fission neutrons", create_fission_neutrons, bool)
self._create_fission_neutrons = create_fission_neutrons
@property
@@ -1227,8 +1298,9 @@ def create_delayed_neutrons(self) -> bool:
@create_delayed_neutrons.setter
def create_delayed_neutrons(self, create_delayed_neutrons: bool):
- cv.check_type('Whether create only prompt neutrons',
- create_delayed_neutrons, bool)
+ cv.check_type(
+ "Whether create only prompt neutrons", create_delayed_neutrons, bool
+ )
self._create_delayed_neutrons = create_delayed_neutrons
@property
@@ -1237,7 +1309,7 @@ def delayed_photon_scaling(self) -> bool:
@delayed_photon_scaling.setter
def delayed_photon_scaling(self, value: bool):
- cv.check_type('delayed photon scaling', value, bool)
+ cv.check_type("delayed photon scaling", value, bool)
self._delayed_photon_scaling = value
@property
@@ -1246,7 +1318,7 @@ def material_cell_offsets(self) -> bool:
@material_cell_offsets.setter
def material_cell_offsets(self, value: bool):
- cv.check_type('material cell offsets', value, bool)
+ cv.check_type("material cell offsets", value, bool)
self._material_cell_offsets = value
@property
@@ -1255,8 +1327,8 @@ def log_grid_bins(self) -> int:
@log_grid_bins.setter
def log_grid_bins(self, log_grid_bins: int):
- cv.check_type('log grid bins', log_grid_bins, Real)
- cv.check_greater_than('log grid bins', log_grid_bins, 0)
+ cv.check_type("log grid bins", log_grid_bins, Real)
+ cv.check_greater_than("log grid bins", log_grid_bins, 0)
self._log_grid_bins = log_grid_bins
@property
@@ -1265,7 +1337,7 @@ def event_based(self) -> bool:
@event_based.setter
def event_based(self, value: bool):
- cv.check_type('event based', value, bool)
+ cv.check_type("event based", value, bool)
self._event_based = value
@property
@@ -1274,8 +1346,8 @@ def max_particles_in_flight(self) -> int:
@max_particles_in_flight.setter
def max_particles_in_flight(self, value: int):
- cv.check_type('max particles in flight', value, Integral)
- cv.check_greater_than('max particles in flight', value, 0)
+ cv.check_type("max particles in flight", value, Integral)
+ cv.check_greater_than("max particles in flight", value, 0)
self._max_particles_in_flight = value
@property
@@ -1284,8 +1356,8 @@ def max_particle_events(self) -> int:
@max_particle_events.setter
def max_particle_events(self, value: int):
- cv.check_type('max particle events', value, Integral)
- cv.check_greater_than('max particle events', value, 0)
+ cv.check_type("max particle events", value, Integral)
+ cv.check_greater_than("max particle events", value, 0)
self._max_particle_events = value
@property
@@ -1294,7 +1366,7 @@ def write_initial_source(self) -> bool:
@write_initial_source.setter
def write_initial_source(self, value: bool):
- cv.check_type('write initial source', value, bool)
+ cv.check_type("write initial source", value, bool)
self._write_initial_source = value
@property
@@ -1313,7 +1385,7 @@ def weight_windows_on(self) -> bool:
@weight_windows_on.setter
def weight_windows_on(self, value: bool):
- cv.check_type('weight windows on', value, bool)
+ cv.check_type("weight windows on", value, bool)
self._weight_windows_on = value
@property
@@ -1322,7 +1394,7 @@ def shared_secondary_bank(self) -> bool:
@shared_secondary_bank.setter
def shared_secondary_bank(self, value: bool):
- cv.check_type('shared secondary bank', value, bool)
+ cv.check_type("shared secondary bank", value, bool)
self._shared_secondary_bank = value
@property
@@ -1332,14 +1404,14 @@ def weight_window_checkpoints(self) -> dict:
@weight_window_checkpoints.setter
def weight_window_checkpoints(self, weight_window_checkpoints: dict):
for key in weight_window_checkpoints.keys():
- cv.check_value('weight_window_checkpoints',
- key, ('collision', 'surface'))
+ cv.check_value("weight_window_checkpoints", key, ("collision", "surface"))
self._weight_window_checkpoints = weight_window_checkpoints
@property
def max_splits(self):
raise AttributeError(
- 'max_splits has been deprecated. Please use max_history_splits instead')
+ "max_splits has been deprecated. Please use max_history_splits instead"
+ )
@property
def max_history_splits(self) -> int:
@@ -1347,8 +1419,8 @@ def max_history_splits(self) -> int:
@max_history_splits.setter
def max_history_splits(self, value: int):
- cv.check_type('maximum particle splits', value, Integral)
- cv.check_greater_than('max particle splits', value, 0)
+ cv.check_type("maximum particle splits", value, Integral)
+ cv.check_greater_than("max particle splits", value, 0)
self._max_history_splits = value
@property
@@ -1357,8 +1429,8 @@ def max_secondaries(self) -> int:
@max_secondaries.setter
def max_secondaries(self, value: int):
- cv.check_type('maximum secondary bank size', value, Integral)
- cv.check_greater_than('max secondary bank size', value, 0)
+ cv.check_type("maximum secondary bank size", value, Integral)
+ cv.check_greater_than("max secondary bank size", value, 0)
self._max_secondaries = value
@property
@@ -1367,8 +1439,8 @@ def max_tracks(self) -> int:
@max_tracks.setter
def max_tracks(self, value: int):
- cv.check_type('maximum particle tracks', value, Integral)
- cv.check_greater_than('maximum particle tracks', value, 0, True)
+ cv.check_type("maximum particle tracks", value, Integral)
+ cv.check_greater_than("maximum particle tracks", value, 0, True)
self._max_tracks = value
@property
@@ -1380,9 +1452,19 @@ def weight_windows_file(self, value: PathLike | None):
if value is None:
self._weight_windows_file = None
else:
- cv.check_type('weight windows file', value, PathLike)
+ cv.check_type("weight windows file", value, PathLike)
self._weight_windows_file = input_path(value)
+ @property
+ def weight_windows_exodus(self) -> WeightWindowsExodus | None:
+ return self._weight_windows_exodus
+
+ @weight_windows_exodus.setter
+ def weight_windows_exodus(self, value: WeightWindowsExodus | None):
+ if value is not None:
+ cv.check_type("weight windows exodus", value, WeightWindowsExodus)
+ self._weight_windows_exodus = value
+
@property
def weight_window_generators(self) -> list[WeightWindowGenerator]:
return self._weight_window_generators
@@ -1392,7 +1474,8 @@ def weight_window_generators(self, wwgs):
if not isinstance(wwgs, MutableSequence):
wwgs = [wwgs]
self._weight_window_generators = cv.CheckedList(
- WeightWindowGenerator, 'weight window generators', wwgs)
+ WeightWindowGenerator, "weight window generators", wwgs
+ )
@property
def random_ray(self) -> dict:
@@ -1401,59 +1484,61 @@ def random_ray(self) -> dict:
@random_ray.setter
def random_ray(self, random_ray: dict):
if not isinstance(random_ray, Mapping):
- raise ValueError(f'Unable to set random_ray from "{random_ray}" '
- 'which is not a dict.')
+ raise ValueError(
+ f'Unable to set random_ray from "{random_ray}" ' "which is not a dict."
+ )
for key, value in random_ray.items():
- if key == 'distance_active':
- cv.check_type('active ray length', value, Real)
- cv.check_greater_than('active ray length', value, 0.0)
- elif key == 'distance_inactive':
- cv.check_type('inactive ray length', value, Real)
- cv.check_greater_than('inactive ray length',
- value, 0.0, True)
- elif key == 'ray_source':
- cv.check_type('random ray source', value, SourceBase)
- elif key == 'volume_estimator':
- cv.check_value('volume estimator', value,
- ('naive', 'simulation_averaged',
- 'hybrid'))
- elif key == 'source_shape':
- cv.check_value('source shape', value,
- ('flat', 'linear', 'linear_xy'))
- elif key == 'volume_normalized_flux_tallies':
- cv.check_type('volume normalized flux tallies', value, bool)
- elif key == 'adjoint':
- cv.check_type('adjoint', value, bool)
- elif key == 'source_region_meshes':
- cv.check_type('source region meshes', value, Iterable)
+ if key == "distance_active":
+ cv.check_type("active ray length", value, Real)
+ cv.check_greater_than("active ray length", value, 0.0)
+ elif key == "distance_inactive":
+ cv.check_type("inactive ray length", value, Real)
+ cv.check_greater_than("inactive ray length", value, 0.0, True)
+ elif key == "ray_source":
+ cv.check_type("random ray source", value, SourceBase)
+ elif key == "volume_estimator":
+ cv.check_value(
+ "volume estimator",
+ value,
+ ("naive", "simulation_averaged", "hybrid"),
+ )
+ elif key == "source_shape":
+ cv.check_value("source shape", value, ("flat", "linear", "linear_xy"))
+ elif key == "volume_normalized_flux_tallies":
+ cv.check_type("volume normalized flux tallies", value, bool)
+ elif key == "adjoint":
+ cv.check_type("adjoint", value, bool)
+ elif key == "source_region_meshes":
+ cv.check_type("source region meshes", value, Iterable)
for mesh, domains in value:
- cv.check_type('mesh', mesh, MeshBase)
- cv.check_type('domains', domains, Iterable)
- valid_types = (openmc.Material,
- openmc.Cell, openmc.Universe)
+ cv.check_type("mesh", mesh, MeshBase)
+ cv.check_type("domains", domains, Iterable)
+ valid_types = (openmc.Material, openmc.Cell, openmc.Universe)
for domain in domains:
if not isinstance(domain, valid_types):
raise ValueError(
- f'Invalid domain type: {type(domain)}. Expected '
- 'openmc.Material, openmc.Cell, or openmc.Universe.')
- elif key == 'sample_method':
- cv.check_value('sample method', value,
- ('prng', 'halton', 's2'))
- elif key == 'diagonal_stabilization_rho':
- cv.check_type('diagonal stabilization rho', value, Real)
- cv.check_greater_than('diagonal stabilization rho',
- value, 0.0, True)
- elif key == 'adjoint_source':
+ f"Invalid domain type: {type(domain)}. Expected "
+ "openmc.Material, openmc.Cell, or openmc.Universe."
+ )
+ elif key == "sample_method":
+ cv.check_value("sample method", value, ("prng", "halton", "s2"))
+ elif key == "diagonal_stabilization_rho":
+ cv.check_type("diagonal stabilization rho", value, Real)
+ cv.check_greater_than("diagonal stabilization rho", value, 0.0, True)
+ elif key == "adjoint_source":
if not isinstance(value, MutableSequence):
value = [value]
for source in value:
if not isinstance(source, SourceBase):
raise ValueError(
- f'Invalid adjoint source type: {type(source)}. '
- 'Expected openmc.SourceBase.')
+ f"Invalid adjoint source type: {type(source)}. "
+ "Expected openmc.SourceBase."
+ )
else:
- raise ValueError(f'Unable to set random ray to "{key}" which is '
- 'unsupported by OpenMC')
+ raise ValueError(
+ f'Unable to set random ray to "{key}" which is '
+ "unsupported by OpenMC"
+ )
self._random_ray = random_ray
@@ -1463,7 +1548,7 @@ def use_decay_photons(self) -> bool:
@use_decay_photons.setter
def use_decay_photons(self, value):
- cv.check_type('use decay photons', value, bool)
+ cv.check_type("use decay photons", value, bool)
self._use_decay_photons = value
@property
@@ -1472,12 +1557,9 @@ def source_rejection_fraction(self) -> float:
@source_rejection_fraction.setter
def source_rejection_fraction(self, source_rejection_fraction: float):
- cv.check_type('source_rejection_fraction',
- source_rejection_fraction, Real)
- cv.check_greater_than('source_rejection_fraction',
- source_rejection_fraction, 0)
- cv.check_less_than('source_rejection_fraction',
- source_rejection_fraction, 1)
+ cv.check_type("source_rejection_fraction", source_rejection_fraction, Real)
+ cv.check_greater_than("source_rejection_fraction", source_rejection_fraction, 0)
+ cv.check_less_than("source_rejection_fraction", source_rejection_fraction, 1)
self._source_rejection_fraction = source_rejection_fraction
@property
@@ -1487,8 +1569,8 @@ def free_gas_threshold(self) -> float | None:
@free_gas_threshold.setter
def free_gas_threshold(self, free_gas_threshold: float | None):
if free_gas_threshold is not None:
- cv.check_type('free gas threshold', free_gas_threshold, Real)
- cv.check_greater_than('free gas threshold', free_gas_threshold, 0.0)
+ cv.check_type("free gas threshold", free_gas_threshold, Real)
+ cv.check_greater_than("free gas threshold", free_gas_threshold, 0.0)
self._free_gas_threshold = free_gas_threshold
def _create_run_mode_subelement(self, root):
@@ -1550,7 +1632,9 @@ def _create_max_order_subelement(self, root):
def _create_source_subelement(self, root, mesh_memo=None):
for source in self.source:
root.append(source.to_xml_element())
- if isinstance(source, IndependentSource) and isinstance(source.space, MeshSpatial):
+ if isinstance(source, IndependentSource) and isinstance(
+ source.space, MeshSpatial
+ ):
path = f"./mesh[@id='{source.space.mesh.id}']"
if root.find(path) is None:
root.append(source.space.mesh.to_xml_element())
@@ -1570,7 +1654,7 @@ def _create_output_subelement(self, root):
element = ET.SubElement(root, "output")
for key, value in sorted(self._output.items()):
subelement = ET.SubElement(element, key)
- if key in ('summary', 'tallies'):
+ if key in ("summary", "tallies"):
subelement.text = str(value).lower()
else:
subelement.text = value
@@ -1583,10 +1667,9 @@ def _create_verbosity_subelement(self, root):
def _create_statepoint_subelement(self, root):
if self._statepoint:
element = ET.SubElement(root, "state_point")
- if 'batches' in self._statepoint:
+ if "batches" in self._statepoint:
subelement = ET.SubElement(element, "batches")
- subelement.text = ' '.join(
- str(x) for x in self._statepoint['batches'])
+ subelement.text = " ".join(str(x) for x in self._statepoint["batches"])
def _create_uniform_source_sampling_subelement(self, root):
if self._uniform_source_sampling is not None:
@@ -1597,34 +1680,33 @@ def _create_sourcepoint_subelement(self, root):
if self._sourcepoint:
element = ET.SubElement(root, "source_point")
- if 'batches' in self._sourcepoint:
+ if "batches" in self._sourcepoint:
subelement = ET.SubElement(element, "batches")
- subelement.text = ' '.join(
- str(x) for x in self._sourcepoint['batches'])
+ subelement.text = " ".join(str(x) for x in self._sourcepoint["batches"])
- if 'separate' in self._sourcepoint:
+ if "separate" in self._sourcepoint:
subelement = ET.SubElement(element, "separate")
- subelement.text = str(self._sourcepoint['separate']).lower()
+ subelement.text = str(self._sourcepoint["separate"]).lower()
- if 'write' in self._sourcepoint:
+ if "write" in self._sourcepoint:
subelement = ET.SubElement(element, "write")
- subelement.text = str(self._sourcepoint['write']).lower()
+ subelement.text = str(self._sourcepoint["write"]).lower()
# Overwrite latest subelement
- if 'overwrite' in self._sourcepoint:
+ if "overwrite" in self._sourcepoint:
subelement = ET.SubElement(element, "overwrite_latest")
- subelement.text = str(self._sourcepoint['overwrite']).lower()
+ subelement.text = str(self._sourcepoint["overwrite"]).lower()
- if 'mcpl' in self._sourcepoint:
+ if "mcpl" in self._sourcepoint:
subelement = ET.SubElement(element, "mcpl")
- subelement.text = str(self._sourcepoint['mcpl']).lower()
+ subelement.text = str(self._sourcepoint["mcpl"]).lower()
def _create_surf_source_read_subelement(self, root):
if self._surf_source_read:
element = ET.SubElement(root, "surf_source_read")
- if 'path' in self._surf_source_read:
+ if "path" in self._surf_source_read:
subelement = ET.SubElement(element, "path")
- subelement.text = str(self._surf_source_read['path'])
+ subelement.text = str(self._surf_source_read["path"])
def _create_surf_source_write_subelement(self, root):
if self._surf_source_write:
@@ -1637,7 +1719,13 @@ def _create_surf_source_write_subelement(self, root):
if "mcpl" in self._surf_source_write:
subelement = ET.SubElement(element, "mcpl")
subelement.text = str(self._surf_source_write["mcpl"]).lower()
- for key in ("max_particles", "max_source_files", "cell", "cellfrom", "cellto"):
+ for key in (
+ "max_particles",
+ "max_source_files",
+ "cell",
+ "cellfrom",
+ "cellto",
+ ):
if key in self._surf_source_write:
subelement = ET.SubElement(element, key)
subelement.text = str(self._surf_source_write[key])
@@ -1645,41 +1733,45 @@ def _create_surf_source_write_subelement(self, root):
def _create_collision_track_subelement(self, root):
if self._collision_track:
element = ET.SubElement(root, "collision_track")
- if 'cell_ids' in self._collision_track:
+ if "cell_ids" in self._collision_track:
subelement = ET.SubElement(element, "cell_ids")
- subelement.text = ' '.join(
- str(x) for x in self._collision_track['cell_ids'])
- if 'reactions' in self._collision_track:
+ subelement.text = " ".join(
+ str(x) for x in self._collision_track["cell_ids"]
+ )
+ if "reactions" in self._collision_track:
subelement = ET.SubElement(element, "reactions")
- subelement.text = ' '.join(
- str(x) for x in self._collision_track['reactions'])
- if 'universe_ids' in self._collision_track:
+ subelement.text = " ".join(
+ str(x) for x in self._collision_track["reactions"]
+ )
+ if "universe_ids" in self._collision_track:
subelement = ET.SubElement(element, "universe_ids")
- subelement.text = ' '.join(
- str(x) for x in self._collision_track['universe_ids'])
- if 'material_ids' in self._collision_track:
+ subelement.text = " ".join(
+ str(x) for x in self._collision_track["universe_ids"]
+ )
+ if "material_ids" in self._collision_track:
subelement = ET.SubElement(element, "material_ids")
- subelement.text = ' '.join(
- str(x) for x in self._collision_track['material_ids'])
- if 'nuclides' in self._collision_track:
+ subelement.text = " ".join(
+ str(x) for x in self._collision_track["material_ids"]
+ )
+ if "nuclides" in self._collision_track:
subelement = ET.SubElement(element, "nuclides")
- subelement.text = ' '.join(
- str(x) for x in self._collision_track['nuclides'])
- if 'deposited_E_threshold' in self._collision_track:
+ subelement.text = " ".join(
+ str(x) for x in self._collision_track["nuclides"]
+ )
+ if "deposited_E_threshold" in self._collision_track:
subelement = ET.SubElement(element, "deposited_E_threshold")
- subelement.text = str(
- self._collision_track['deposited_E_threshold'])
- if 'max_collisions' in self._collision_track:
+ subelement.text = str(self._collision_track["deposited_E_threshold"])
+ if "max_collisions" in self._collision_track:
subelement = ET.SubElement(element, "max_collisions")
- subelement.text = str(self._collision_track['max_collisions'])
- if 'max_collision_track_files' in self._collision_track:
- subelement = ET.SubElement(
- element, "max_collision_track_files")
+ subelement.text = str(self._collision_track["max_collisions"])
+ if "max_collision_track_files" in self._collision_track:
+ subelement = ET.SubElement(element, "max_collision_track_files")
subelement.text = str(
- self._collision_track['max_collision_track_files'])
- if 'mcpl' in self._collision_track:
+ self._collision_track["max_collision_track_files"]
+ )
+ if "mcpl" in self._collision_track:
subelement = ET.SubElement(element, "mcpl")
- subelement.text = str(self._collision_track['mcpl']).lower()
+ subelement.text = str(self._collision_track["mcpl"]).lower()
def _create_confidence_intervals(self, root):
if self._confidence_intervals is not None:
@@ -1741,8 +1833,11 @@ def _create_cutoff_subelement(self, root):
element = ET.SubElement(root, "cutoff")
for key, value in self._cutoff.items():
subelement = ET.SubElement(element, key)
- subelement.text = str(value) if key != 'survival_normalization' \
+ subelement.text = (
+ str(value)
+ if key != "survival_normalization"
else str(value).lower()
+ )
def _create_entropy_mesh_subelement(self, root, mesh_memo=None):
if self.entropy_mesh is None:
@@ -1751,12 +1846,14 @@ def _create_entropy_mesh_subelement(self, root, mesh_memo=None):
# use default heuristic for entropy mesh if not set by user
if self.entropy_mesh.dimension is None:
if self.particles is None:
- raise RuntimeError("Number of particles must be set in order to "
- "use entropy mesh dimension heuristic")
+ raise RuntimeError(
+ "Number of particles must be set in order to "
+ "use entropy mesh dimension heuristic"
+ )
else:
- n = ceil((self.particles / 20.0)**(1.0 / 3.0))
+ n = ceil((self.particles / 20.0) ** (1.0 / 3.0))
d = len(self.entropy_mesh.lower_left)
- self.entropy_mesh.dimension = (n,)*d
+ self.entropy_mesh.dimension = (n,) * d
# add mesh ID to this element
subelement = ET.SubElement(root, "entropy_mesh")
@@ -1802,10 +1899,10 @@ def _create_tabular_legendre_subelements(self, root):
if self.tabular_legendre:
element = ET.SubElement(root, "tabular_legendre")
subelement = ET.SubElement(element, "enable")
- subelement.text = str(self._tabular_legendre['enable']).lower()
- if 'num_points' in self._tabular_legendre:
+ subelement.text = str(self._tabular_legendre["enable"]).lower()
+ if "num_points" in self._tabular_legendre:
subelement = ET.SubElement(element, "num_points")
- subelement.text = str(self._tabular_legendre['num_points'])
+ subelement.text = str(self._tabular_legendre["num_points"])
def _create_temperature_subelements(self, root):
if self.temperature:
@@ -1813,8 +1910,8 @@ def _create_temperature_subelements(self, root):
element = ET.SubElement(root, f"temperature_{key}")
if isinstance(value, bool):
element.text = str(value).lower()
- elif key == 'range':
- element.text = ' '.join(str(T) for T in value)
+ elif key == "range":
+ element.text = " ".join(str(T) for T in value)
else:
element.text = str(value)
@@ -1827,12 +1924,12 @@ def _create_properties_file_element(self, root):
def _create_trace_subelement(self, root):
if self._trace is not None:
element = ET.SubElement(root, "trace")
- element.text = ' '.join(map(str, self._trace))
+ element.text = " ".join(map(str, self._trace))
def _create_track_subelement(self, root):
if self._track is not None:
element = ET.SubElement(root, "track")
- element.text = ' '.join(map(str, itertools.chain(*self._track)))
+ element.text = " ".join(map(str, itertools.chain(*self._track)))
def _create_ufs_mesh_subelement(self, root, mesh_memo=None):
if self.ufs_mesh is None:
@@ -1859,22 +1956,22 @@ def _create_use_decay_photons_subelement(self, root):
def _create_resonance_scattering_subelement(self, root):
res = self.resonance_scattering
if res:
- elem = ET.SubElement(root, 'resonance_scattering')
- if 'enable' in res:
- subelem = ET.SubElement(elem, 'enable')
- subelem.text = str(res['enable']).lower()
- if 'method' in res:
- subelem = ET.SubElement(elem, 'method')
- subelem.text = res['method']
- if 'energy_min' in res:
- subelem = ET.SubElement(elem, 'energy_min')
- subelem.text = str(res['energy_min'])
- if 'energy_max' in res:
- subelem = ET.SubElement(elem, 'energy_max')
- subelem.text = str(res['energy_max'])
- if 'nuclides' in res:
- subelem = ET.SubElement(elem, 'nuclides')
- subelem.text = ' '.join(res['nuclides'])
+ elem = ET.SubElement(root, "resonance_scattering")
+ if "enable" in res:
+ subelem = ET.SubElement(elem, "enable")
+ subelem.text = str(res["enable"]).lower()
+ if "method" in res:
+ subelem = ET.SubElement(elem, "method")
+ subelem.text = res["method"]
+ if "energy_min" in res:
+ subelem = ET.SubElement(elem, "energy_min")
+ subelem.text = str(res["energy_min"])
+ if "energy_max" in res:
+ subelem = ET.SubElement(elem, "energy_max")
+ subelem.text = str(res["energy_max"])
+ if "nuclides" in res:
+ subelem = ET.SubElement(elem, "nuclides")
+ subelem.text = " ".join(res["nuclides"])
def _create_create_fission_neutrons_subelement(self, root):
if self._create_fission_neutrons is not None:
@@ -1951,7 +2048,7 @@ def _create_shared_secondary_bank_subelement(self, root):
def _create_weight_window_generators_subelement(self, root, mesh_memo=None):
if not self.weight_window_generators:
return
- elem = ET.SubElement(root, 'weight_window_generators')
+ elem = ET.SubElement(root, "weight_window_generators")
for wwg in self.weight_window_generators:
elem.append(wwg.to_xml_element())
@@ -1973,20 +2070,22 @@ def _create_weight_windows_file_element(self, root):
element.text = str(self.weight_windows_file)
root.append(element)
+ def _create_weight_windows_exodus_subelement(self, root):
+ if self._weight_windows_exodus is not None:
+ root.append(self._weight_windows_exodus.to_xml_element())
+
def _create_weight_window_checkpoints_subelement(self, root):
if not self._weight_window_checkpoints:
return
element = ET.SubElement(root, "weight_window_checkpoints")
- if 'collision' in self._weight_window_checkpoints:
+ if "collision" in self._weight_window_checkpoints:
subelement = ET.SubElement(element, "collision")
- subelement.text = str(
- self._weight_window_checkpoints['collision']).lower()
+ subelement.text = str(self._weight_window_checkpoints["collision"]).lower()
- if 'surface' in self._weight_window_checkpoints:
+ if "surface" in self._weight_window_checkpoints:
subelement = ET.SubElement(element, "surface")
- subelement.text = str(
- self._weight_window_checkpoints['surface']).lower()
+ subelement.text = str(self._weight_window_checkpoints["surface"]).lower()
def _create_max_history_splits_subelement(self, root):
if self._max_history_splits is not None:
@@ -2007,34 +2106,34 @@ def _create_random_ray_subelement(self, root, mesh_memo=None):
if self._random_ray:
element = ET.SubElement(root, "random_ray")
for key, value in self._random_ray.items():
- if key == 'ray_source' and isinstance(value, SourceBase):
- subelement = ET.SubElement(element, 'ray_source')
+ if key == "ray_source" and isinstance(value, SourceBase):
+ subelement = ET.SubElement(element, "ray_source")
source_element = value.to_xml_element()
- if source_element.find('bias') is not None:
+ if source_element.find("bias") is not None:
raise RuntimeError(
- "Ray source distributions should not be biased.")
+ "Ray source distributions should not be biased."
+ )
subelement.append(source_element)
- elif key == 'source_region_meshes':
- subelement = ET.SubElement(element, 'source_region_meshes')
+ elif key == "source_region_meshes":
+ subelement = ET.SubElement(element, "source_region_meshes")
for mesh, domains in value:
- mesh_elem = ET.SubElement(subelement, 'mesh')
- mesh_elem.set('id', str(mesh.id))
+ mesh_elem = ET.SubElement(subelement, "mesh")
+ mesh_elem.set("id", str(mesh.id))
for domain in domains:
- domain_elem = ET.SubElement(mesh_elem, 'domain')
- domain_elem.set('id', str(domain.id))
- domain_elem.set(
- 'type', domain.__class__.__name__.lower())
+ domain_elem = ET.SubElement(mesh_elem, "domain")
+ domain_elem.set("id", str(domain.id))
+ domain_elem.set("type", domain.__class__.__name__.lower())
if mesh_memo is not None and mesh.id not in mesh_memo:
- domain_elem.set('type', domain.__class__.__name__.lower())
+ domain_elem.set("type", domain.__class__.__name__.lower())
# See if a element already exists -- if not, add it
path = f"./mesh[@id='{mesh.id}']"
if root.find(path) is None:
root.append(mesh.to_xml_element())
if mesh_memo is not None:
mesh_memo.add(mesh.id)
- elif key == 'adjoint_source':
- subelement = ET.SubElement(element, 'adjoint_source')
+ elif key == "adjoint_source":
+ subelement = ET.SubElement(element, "adjoint_source")
# Check that all entries are valid SourceBase instances, in case
# the random_ray setter was not used to populate dict entries.
if not isinstance(value, MutableSequence):
@@ -2042,8 +2141,9 @@ def _create_random_ray_subelement(self, root, mesh_memo=None):
for source in value:
if not isinstance(source, SourceBase):
raise ValueError(
- f'Invalid adjoint source type: {type(source)}. '
- 'Expected openmc.SourceBase.')
+ f"Invalid adjoint source type: {type(source)}. "
+ "Expected openmc.SourceBase."
+ )
subelement.append(source.to_xml_element())
elif isinstance(value, bool):
subelement = ET.SubElement(element, key)
@@ -2063,7 +2163,7 @@ def _create_free_gas_threshold_subelement(self, root):
element.text = str(self._free_gas_threshold)
def _eigenvalue_from_xml_element(self, root):
- elem = root.find('eigenvalue')
+ elem = root.find("eigenvalue")
if elem is not None:
self._run_mode_from_xml_element(elem)
self._particles_from_xml_element(elem)
@@ -2075,54 +2175,54 @@ def _eigenvalue_from_xml_element(self, root):
self._generations_per_batch_from_xml_element(elem)
def _run_mode_from_xml_element(self, root):
- text = get_text(root, 'run_mode')
+ text = get_text(root, "run_mode")
if text is not None:
self.run_mode = text
def _particles_from_xml_element(self, root):
- text = get_text(root, 'particles')
+ text = get_text(root, "particles")
if text is not None:
self.particles = int(text)
def _batches_from_xml_element(self, root):
- text = get_text(root, 'batches')
+ text = get_text(root, "batches")
if text is not None:
self.batches = int(text)
def _inactive_from_xml_element(self, root):
- text = get_text(root, 'inactive')
+ text = get_text(root, "inactive")
if text is not None:
self.inactive = int(text)
def _max_lost_particles_from_xml_element(self, root):
- text = get_text(root, 'max_lost_particles')
+ text = get_text(root, "max_lost_particles")
if text is not None:
self.max_lost_particles = int(text)
def _rel_max_lost_particles_from_xml_element(self, root):
- text = get_text(root, 'rel_max_lost_particles')
+ text = get_text(root, "rel_max_lost_particles")
if text is not None:
self.rel_max_lost_particles = float(text)
def _max_write_lost_particles_from_xml_element(self, root):
- text = get_text(root, 'max_write_lost_particles')
+ text = get_text(root, "max_write_lost_particles")
if text is not None:
self.max_write_lost_particles = int(text)
def _generations_per_batch_from_xml_element(self, root):
- text = get_text(root, 'generations_per_batch')
+ text = get_text(root, "generations_per_batch")
if text is not None:
self.generations_per_batch = int(text)
def _keff_trigger_from_xml_element(self, root):
- elem = root.find('keff_trigger')
+ elem = root.find("keff_trigger")
if elem is not None:
- trigger = get_text(elem, 'type')
- threshold = float(get_text(elem, 'threshold'))
- self.keff_trigger = {'type': trigger, 'threshold': threshold}
+ trigger = get_text(elem, "type")
+ threshold = float(get_text(elem, "threshold"))
+ self.keff_trigger = {"type": trigger, "threshold": threshold}
def _source_from_xml_element(self, root, meshes=None):
- for elem in root.findall('source'):
+ for elem in root.findall("source"):
src = SourceBase.from_xml_element(elem, meshes)
# add newly constructed source object to the list
self.source.append(src)
@@ -2130,171 +2230,204 @@ def _source_from_xml_element(self, root, meshes=None):
def _volume_calcs_from_xml_element(self, root):
volume_elems = root.findall("volume_calc")
if volume_elems:
- self.volume_calculations = [VolumeCalculation.from_xml_element(elem)
- for elem in volume_elems]
+ self.volume_calculations = [
+ VolumeCalculation.from_xml_element(elem) for elem in volume_elems
+ ]
def _output_from_xml_element(self, root):
- elem = root.find('output')
+ elem = root.find("output")
if elem is not None:
self.output = {}
- for key in ('summary', 'tallies', 'path'):
+ for key in ("summary", "tallies", "path"):
value = get_text(elem, key)
if value is not None:
- if key in ('summary', 'tallies'):
- value = value in ('true', '1')
+ if key in ("summary", "tallies"):
+ value = value in ("true", "1")
self.output[key] = value
def _statepoint_from_xml_element(self, root):
- elem = root.find('state_point')
+ elem = root.find("state_point")
if elem is not None:
batches = get_elem_list(elem, "batches", int)
if batches is not None:
- self.statepoint['batches'] = batches
+ self.statepoint["batches"] = batches
def _sourcepoint_from_xml_element(self, root):
- elem = root.find('source_point')
+ elem = root.find("source_point")
if elem is not None:
- for key in ('separate', 'write', 'overwrite_latest', 'batches', 'mcpl'):
- if key in ('separate', 'write', 'mcpl', 'overwrite_latest'):
- value = get_text(elem, key) in ('true', '1')
- if key == 'overwrite_latest':
- key = 'overwrite'
+ for key in ("separate", "write", "overwrite_latest", "batches", "mcpl"):
+ if key in ("separate", "write", "mcpl", "overwrite_latest"):
+ value = get_text(elem, key) in ("true", "1")
+ if key == "overwrite_latest":
+ key = "overwrite"
else:
value = get_elem_list(elem, key, int)
if value is not None:
self.sourcepoint[key] = value
def _surf_source_read_from_xml_element(self, root):
- elem = root.find('surf_source_read')
+ elem = root.find("surf_source_read")
if elem is not None:
ssr = {}
- value = get_text(elem, 'path')
+ value = get_text(elem, "path")
if value is not None:
- ssr['path'] = value
+ ssr["path"] = value
self.surf_source_read = ssr
def _surf_source_write_from_xml_element(self, root):
- elem = root.find('surf_source_write')
+ elem = root.find("surf_source_write")
if elem is None:
return
- for key in ('surface_ids', 'max_particles', 'max_source_files', 'mcpl', 'cell', 'cellto', 'cellfrom'):
- if key == 'surface_ids':
+ for key in (
+ "surface_ids",
+ "max_particles",
+ "max_source_files",
+ "mcpl",
+ "cell",
+ "cellto",
+ "cellfrom",
+ ):
+ if key == "surface_ids":
value = get_elem_list(elem, key, int)
else:
value = get_text(elem, key)
if value is not None:
- if key == 'mcpl':
- value = value in ('true', '1')
- elif key in ('max_particles', 'max_source_files', 'cell', 'cellfrom', 'cellto'):
+ if key == "mcpl":
+ value = value in ("true", "1")
+ elif key in (
+ "max_particles",
+ "max_source_files",
+ "cell",
+ "cellfrom",
+ "cellto",
+ ):
value = int(value)
self.surf_source_write[key] = value
def _collision_track_from_xml_element(self, root):
- elem = root.find('collision_track')
+ elem = root.find("collision_track")
if elem is not None:
- for key in ('cell_ids', 'reactions', 'universe_ids', 'material_ids', 'nuclides',
- 'deposited_E_threshold', 'max_collisions', "max_collision_track_files", 'mcpl'):
+ for key in (
+ "cell_ids",
+ "reactions",
+ "universe_ids",
+ "material_ids",
+ "nuclides",
+ "deposited_E_threshold",
+ "max_collisions",
+ "max_collision_track_files",
+ "mcpl",
+ ):
value = get_text(elem, key)
if value is not None:
- if key in ('cell_ids', 'universe_ids', 'material_ids'):
+ if key in ("cell_ids", "universe_ids", "material_ids"):
value = [int(x) for x in value.split()]
- elif key in ('reactions', 'nuclides'):
+ elif key in ("reactions", "nuclides"):
value = value.split()
- elif key in ('max_collisions', 'max_collision_track_files'):
+ elif key in ("max_collisions", "max_collision_track_files"):
value = int(value)
- elif key == 'deposited_E_threshold':
+ elif key == "deposited_E_threshold":
value = float(value)
- elif key == 'mcpl':
- value = value in ('true', '1')
+ elif key == "mcpl":
+ value = value in ("true", "1")
self.collision_track[key] = value
def _confidence_intervals_from_xml_element(self, root):
- text = get_text(root, 'confidence_intervals')
+ text = get_text(root, "confidence_intervals")
if text is not None:
- self.confidence_intervals = text in ('true', '1')
+ self.confidence_intervals = text in ("true", "1")
def _electron_treatment_from_xml_element(self, root):
- text = get_text(root, 'electron_treatment')
+ text = get_text(root, "electron_treatment")
if text is not None:
self.electron_treatment = text
def _atomic_relaxation_from_xml_element(self, root):
- text = get_text(root, 'atomic_relaxation')
+ text = get_text(root, "atomic_relaxation")
if text is not None:
- self.atomic_relaxation = text in ('true', '1')
+ self.atomic_relaxation = text in ("true", "1")
def _energy_mode_from_xml_element(self, root):
- text = get_text(root, 'energy_mode')
+ text = get_text(root, "energy_mode")
if text is not None:
self.energy_mode = text
def _max_order_from_xml_element(self, root):
- text = get_text(root, 'max_order')
+ text = get_text(root, "max_order")
if text is not None:
self.max_order = int(text)
def _photon_transport_from_xml_element(self, root):
- text = get_text(root, 'photon_transport')
+ text = get_text(root, "photon_transport")
if text is not None:
- self.photon_transport = text in ('true', '1')
+ self.photon_transport = text in ("true", "1")
def _uniform_source_sampling_from_xml_element(self, root):
- text = get_text(root, 'uniform_source_sampling')
+ text = get_text(root, "uniform_source_sampling")
if text is not None:
- self.uniform_source_sampling = text in ('true', '1')
+ self.uniform_source_sampling = text in ("true", "1")
def _plot_seed_from_xml_element(self, root):
- text = get_text(root, 'plot_seed')
+ text = get_text(root, "plot_seed")
if text is not None:
self.plot_seed = int(text)
def _ptables_from_xml_element(self, root):
- text = get_text(root, 'ptables')
+ text = get_text(root, "ptables")
if text is not None:
- self.ptables = text in ('true', '1')
+ self.ptables = text in ("true", "1")
def _seed_from_xml_element(self, root):
- text = get_text(root, 'seed')
+ text = get_text(root, "seed")
if text is not None:
self.seed = int(text)
def _stride_from_xml_element(self, root):
- text = get_text(root, 'stride')
+ text = get_text(root, "stride")
if text is not None:
self.stride = int(text)
def _surface_grazing_cutoff_from_xml_element(self, root):
- text = get_text(root, 'surface_grazing_cutoff')
+ text = get_text(root, "surface_grazing_cutoff")
if text is not None:
self.surface_grazing_cutoff = float(text)
def _surface_grazing_ratio_from_xml_element(self, root):
- text = get_text(root, 'surface_grazing_ratio')
+ text = get_text(root, "surface_grazing_ratio")
if text is not None:
self.surface_grazing_ratio = float(text)
def _survival_biasing_from_xml_element(self, root):
- text = get_text(root, 'survival_biasing')
+ text = get_text(root, "survival_biasing")
if text is not None:
- self.survival_biasing = text in ('true', '1')
+ self.survival_biasing = text in ("true", "1")
def _cutoff_from_xml_element(self, root):
- elem = root.find('cutoff')
+ elem = root.find("cutoff")
if elem is not None:
self.cutoff = {}
- for key in ('energy_neutron', 'energy_photon', 'energy_electron',
- 'energy_positron', 'weight', 'weight_avg', 'time_neutron',
- 'time_photon', 'time_electron', 'time_positron',
- 'survival_normalization'):
+ for key in (
+ "energy_neutron",
+ "energy_photon",
+ "energy_electron",
+ "energy_positron",
+ "weight",
+ "weight_avg",
+ "time_neutron",
+ "time_photon",
+ "time_electron",
+ "time_positron",
+ "survival_normalization",
+ ):
value = get_text(elem, key)
if value is not None:
- if key == 'survival_normalization':
- self.cutoff[key] = value in ('true', '1')
+ if key == "survival_normalization":
+ self.cutoff[key] = value in ("true", "1")
else:
self.cutoff[key] = float(value)
def _entropy_mesh_from_xml_element(self, root, meshes):
- text = get_text(root, 'entropy_mesh')
+ text = get_text(root, "entropy_mesh")
if text is None:
return
mesh_id = int(text)
@@ -2303,59 +2436,59 @@ def _entropy_mesh_from_xml_element(self, root, meshes):
self.entropy_mesh = meshes[mesh_id]
def _trigger_from_xml_element(self, root):
- elem = root.find('trigger')
+ elem = root.find("trigger")
if elem is not None:
- self.trigger_active = get_text(elem, 'active') in ('true', '1')
- text = get_text(elem, 'max_batches')
+ self.trigger_active = get_text(elem, "active") in ("true", "1")
+ text = get_text(elem, "max_batches")
if text is not None:
self.trigger_max_batches = int(text)
- text = get_text(elem, 'batch_interval')
+ text = get_text(elem, "batch_interval")
if text is not None:
self.trigger_batch_interval = int(text)
def _no_reduce_from_xml_element(self, root):
- text = get_text(root, 'no_reduce')
+ text = get_text(root, "no_reduce")
if text is not None:
- self.no_reduce = text in ('true', '1')
+ self.no_reduce = text in ("true", "1")
def _verbosity_from_xml_element(self, root):
- text = get_text(root, 'verbosity')
+ text = get_text(root, "verbosity")
if text is not None:
self.verbosity = int(text)
def _ifp_n_generation_from_xml_element(self, root):
- text = get_text(root, 'ifp_n_generation')
+ text = get_text(root, "ifp_n_generation")
if text is not None:
self.ifp_n_generation = int(text)
def _tabular_legendre_from_xml_element(self, root):
- elem = root.find('tabular_legendre')
+ elem = root.find("tabular_legendre")
if elem is not None:
- text = get_text(elem, 'enable')
- self.tabular_legendre['enable'] = text in ('true', '1')
- text = get_text(elem, 'num_points')
+ text = get_text(elem, "enable")
+ self.tabular_legendre["enable"] = text in ("true", "1")
+ text = get_text(elem, "num_points")
if text is not None:
- self.tabular_legendre['num_points'] = int(text)
+ self.tabular_legendre["num_points"] = int(text)
def _temperature_from_xml_element(self, root):
- text = get_text(root, 'temperature_default')
+ text = get_text(root, "temperature_default")
if text is not None:
- self.temperature['default'] = float(text)
- text = get_text(root, 'temperature_tolerance')
+ self.temperature["default"] = float(text)
+ text = get_text(root, "temperature_tolerance")
if text is not None:
- self.temperature['tolerance'] = float(text)
- text = get_text(root, 'temperature_method')
+ self.temperature["tolerance"] = float(text)
+ text = get_text(root, "temperature_method")
if text is not None:
- self.temperature['method'] = text
+ self.temperature["method"] = text
text = get_elem_list(root, "temperature_range", float)
if text is not None:
- self.temperature['range'] = text
- text = get_text(root, 'temperature_multipole')
+ self.temperature["range"] = text
+ text = get_text(root, "temperature_multipole")
if text is not None:
- self.temperature['multipole'] = text in ('true', '1')
+ self.temperature["multipole"] = text in ("true", "1")
def _properties_file_from_xml_element(self, root):
- text = get_text(root, 'properties_file')
+ text = get_text(root, "properties_file")
if text is not None:
self.properties_file = text
@@ -2370,7 +2503,7 @@ def _track_from_xml_element(self, root):
self.track = list(zip(values[::3], values[1::3], values[2::3]))
def _ufs_mesh_from_xml_element(self, root, meshes):
- text = get_text(root, 'ufs_mesh')
+ text = get_text(root, "ufs_mesh")
if text is None:
return
mesh_id = int(text)
@@ -2379,184 +2512,192 @@ def _ufs_mesh_from_xml_element(self, root, meshes):
self.ufs_mesh = meshes[mesh_id]
def _resonance_scattering_from_xml_element(self, root):
- elem = root.find('resonance_scattering')
+ elem = root.find("resonance_scattering")
if elem is not None:
- keys = ('enable', 'method', 'energy_min', 'energy_max', 'nuclides')
+ keys = ("enable", "method", "energy_min", "energy_max", "nuclides")
for key in keys:
- if key == 'nuclides':
+ if key == "nuclides":
value = get_elem_list(elem, key, str)
else:
value = get_text(elem, key)
if value is not None:
- if key == 'enable':
- value = value in ('true', '1')
- elif key in ('energy_min', 'energy_max'):
+ if key == "enable":
+ value = value in ("true", "1")
+ elif key in ("energy_min", "energy_max"):
value = float(value)
self.resonance_scattering[key] = value
def _create_fission_neutrons_from_xml_element(self, root):
- text = get_text(root, 'create_fission_neutrons')
+ text = get_text(root, "create_fission_neutrons")
if text is not None:
- self.create_fission_neutrons = text in ('true', '1')
+ self.create_fission_neutrons = text in ("true", "1")
def _create_delayed_neutrons_from_xml_element(self, root):
- text = get_text(root, 'create_delayed_neutrons')
+ text = get_text(root, "create_delayed_neutrons")
if text is not None:
- self.create_delayed_neutrons = text in ('true', '1')
+ self.create_delayed_neutrons = text in ("true", "1")
def _delayed_photon_scaling_from_xml_element(self, root):
- text = get_text(root, 'delayed_photon_scaling')
+ text = get_text(root, "delayed_photon_scaling")
if text is not None:
- self.delayed_photon_scaling = text in ('true', '1')
+ self.delayed_photon_scaling = text in ("true", "1")
def _event_based_from_xml_element(self, root):
- text = get_text(root, 'event_based')
+ text = get_text(root, "event_based")
if text is not None:
- self.event_based = text in ('true', '1')
+ self.event_based = text in ("true", "1")
def _max_particles_in_flight_from_xml_element(self, root):
- text = get_text(root, 'max_particles_in_flight')
+ text = get_text(root, "max_particles_in_flight")
if text is not None:
self.max_particles_in_flight = int(text)
def _max_particle_events_from_xml_element(self, root):
- text = get_text(root, 'max_particle_events')
+ text = get_text(root, "max_particle_events")
if text is not None:
self.max_particle_events = int(text)
def _material_cell_offsets_from_xml_element(self, root):
- text = get_text(root, 'material_cell_offsets')
+ text = get_text(root, "material_cell_offsets")
if text is not None:
- self.material_cell_offsets = text in ('true', '1')
+ self.material_cell_offsets = text in ("true", "1")
def _log_grid_bins_from_xml_element(self, root):
- text = get_text(root, 'log_grid_bins')
+ text = get_text(root, "log_grid_bins")
if text is not None:
self.log_grid_bins = int(text)
def _write_initial_source_from_xml_element(self, root):
- text = get_text(root, 'write_initial_source')
+ text = get_text(root, "write_initial_source")
if text is not None:
- self.write_initial_source = text in ('true', '1')
+ self.write_initial_source = text in ("true", "1")
def _weight_window_generators_from_xml_element(self, root, meshes=None):
- for elem in root.iter('weight_windows_generator'):
+ for elem in root.iter("weight_windows_generator"):
wwg = WeightWindowGenerator.from_xml_element(elem, meshes)
self.weight_window_generators.append(wwg)
def _weight_windows_from_xml_element(self, root, meshes=None):
- for elem in root.findall('weight_windows'):
+ for elem in root.findall("weight_windows"):
ww = WeightWindows.from_xml_element(elem, meshes)
self.weight_windows.append(ww)
def _weight_windows_on_from_xml_element(self, root):
- text = get_text(root, 'weight_windows_on')
+ text = get_text(root, "weight_windows_on")
if text is not None:
- self.weight_windows_on = text in ('true', '1')
+ self.weight_windows_on = text in ("true", "1")
def _shared_secondary_bank_from_xml_element(self, root):
- text = get_text(root, 'shared_secondary_bank')
+ text = get_text(root, "shared_secondary_bank")
if text is not None:
- self.shared_secondary_bank = text in ('true', '1')
+ self.shared_secondary_bank = text in ("true", "1")
def _weight_windows_file_from_xml_element(self, root):
- text = get_text(root, 'weight_windows_file')
+ text = get_text(root, "weight_windows_file")
if text is not None:
self.weight_windows_file = text
+ def _weight_windows_exodus_from_xml_element(self, root):
+ elem = root.find("weight_windows_exodus")
+ if elem is not None:
+ self.weight_windows_exodus = WeightWindowsExodus.from_xml_element(elem)
+
def _weight_window_checkpoints_from_xml_element(self, root):
- elem = root.find('weight_window_checkpoints')
+ elem = root.find("weight_window_checkpoints")
if elem is None:
return
- for key in ('collision', 'surface'):
+ for key in ("collision", "surface"):
value = get_text(elem, key)
if value is not None:
- value = value in ('true', '1')
+ value = value in ("true", "1")
self.weight_window_checkpoints[key] = value
def _max_history_splits_from_xml_element(self, root):
- text = get_text(root, 'max_history_splits')
+ text = get_text(root, "max_history_splits")
if text is not None:
self.max_history_splits = int(text)
def _max_secondaries_from_xml_element(self, root):
- text = get_text(root, 'max_secondaries')
+ text = get_text(root, "max_secondaries")
if text is not None:
self.max_secondaries = int(text)
def _max_tracks_from_xml_element(self, root):
- text = get_text(root, 'max_tracks')
+ text = get_text(root, "max_tracks")
if text is not None:
self.max_tracks = int(text)
def _random_ray_from_xml_element(self, root, meshes=None):
- elem = root.find('random_ray')
+ elem = root.find("random_ray")
if elem is not None:
self.random_ray = {}
for child in elem:
- if child.tag in ('distance_inactive', 'distance_active', 'diagonal_stabilization_rho'):
+ if child.tag in (
+ "distance_inactive",
+ "distance_active",
+ "diagonal_stabilization_rho",
+ ):
self.random_ray[child.tag] = float(child.text)
- elif child.tag == 'ray_source':
- source_element = child.find('source')
+ elif child.tag == "ray_source":
+ source_element = child.find("source")
source = SourceBase.from_xml_element(source_element)
- if child.find('bias') is not None:
+ if child.find("bias") is not None:
raise RuntimeError(
- "Ray source distributions should not be biased.")
- self.random_ray['ray_source'] = source
- elif child.tag == 'volume_estimator':
- self.random_ray['volume_estimator'] = child.text
- elif child.tag == 'source_shape':
- self.random_ray['source_shape'] = child.text
- elif child.tag == 'volume_normalized_flux_tallies':
- self.random_ray['volume_normalized_flux_tallies'] = (
- child.text in ('true', '1')
- )
- elif child.tag == 'adjoint':
- self.random_ray['adjoint'] = (
- child.text in ('true', '1')
+ "Ray source distributions should not be biased."
+ )
+ self.random_ray["ray_source"] = source
+ elif child.tag == "volume_estimator":
+ self.random_ray["volume_estimator"] = child.text
+ elif child.tag == "source_shape":
+ self.random_ray["source_shape"] = child.text
+ elif child.tag == "volume_normalized_flux_tallies":
+ self.random_ray["volume_normalized_flux_tallies"] = child.text in (
+ "true",
+ "1",
)
- elif child.tag == 'adjoint_source':
- self.random_ray['adjoint_source'] = []
- for subelem in child.findall('source'):
+ elif child.tag == "adjoint":
+ self.random_ray["adjoint"] = child.text in ("true", "1")
+ elif child.tag == "adjoint_source":
+ self.random_ray["adjoint_source"] = []
+ for subelem in child.findall("source"):
src = SourceBase.from_xml_element(subelem)
# add newly constructed source object to the list
- self.random_ray['adjoint_source'].append(src)
- elif child.tag == 'sample_method':
- self.random_ray['sample_method'] = child.text
- elif child.tag == 'source_region_meshes':
- self.random_ray['source_region_meshes'] = []
- for mesh_elem in child.findall('mesh'):
- mesh_id = int(get_text(mesh_elem, 'id'))
+ self.random_ray["adjoint_source"].append(src)
+ elif child.tag == "sample_method":
+ self.random_ray["sample_method"] = child.text
+ elif child.tag == "source_region_meshes":
+ self.random_ray["source_region_meshes"] = []
+ for mesh_elem in child.findall("mesh"):
+ mesh_id = int(get_text(mesh_elem, "id"))
if meshes and mesh_id in meshes:
mesh = meshes[mesh_id]
else:
mesh = MeshBase.from_xml_element(mesh_elem)
domains = []
- for domain_elem in mesh_elem.findall('domain'):
+ for domain_elem in mesh_elem.findall("domain"):
domain_id = int(get_text(domain_elem, "id"))
domain_type = get_text(domain_elem, "type")
- if domain_type == 'material':
+ if domain_type == "material":
domain = openmc.Material(domain_id)
- elif domain_type == 'cell':
+ elif domain_type == "cell":
domain = openmc.Cell(domain_id)
- elif domain_type == 'universe':
+ elif domain_type == "universe":
domain = openmc.Universe(domain_id)
domains.append(domain)
- self.random_ray['source_region_meshes'].append(
- (mesh, domains))
+ self.random_ray["source_region_meshes"].append((mesh, domains))
def _use_decay_photons_from_xml_element(self, root):
- text = get_text(root, 'use_decay_photons')
+ text = get_text(root, "use_decay_photons")
if text is not None:
- self.use_decay_photons = text in ('true', '1')
+ self.use_decay_photons = text in ("true", "1")
def _source_rejection_fraction_from_xml_element(self, root):
- text = get_text(root, 'source_rejection_fraction')
+ text = get_text(root, "source_rejection_fraction")
if text is not None:
self.source_rejection_fraction = float(text)
def _free_gas_threshold_from_xml_element(self, root):
- text = get_text(root, 'free_gas_threshold')
+ text = get_text(root, "free_gas_threshold")
if text is not None:
self.free_gas_threshold = float(text)
@@ -2629,6 +2770,7 @@ def to_xml_element(self, mesh_memo=None):
self._create_shared_secondary_bank_subelement(element)
self._create_weight_window_generators_subelement(element, mesh_memo)
self._create_weight_windows_file_element(element)
+ self._create_weight_windows_exodus_subelement(element)
self._create_weight_window_checkpoints_subelement(element)
self._create_max_history_splits_subelement(element)
self._create_max_tracks_subelement(element)
@@ -2643,7 +2785,7 @@ def to_xml_element(self, mesh_memo=None):
return element
- def export_to_xml(self, path: PathLike = 'settings.xml'):
+ def export_to_xml(self, path: PathLike = "settings.xml"):
"""Export simulation settings to an XML file.
Parameters
@@ -2657,11 +2799,11 @@ def export_to_xml(self, path: PathLike = 'settings.xml'):
# Check if path is a directory
p = Path(path)
if p.is_dir():
- p /= 'settings.xml'
+ p /= "settings.xml"
# Write the XML Tree to the settings.xml file
tree = ET.ElementTree(root_element)
- tree.write(str(p), xml_declaration=True, encoding='utf-8')
+ tree.write(str(p), xml_declaration=True, encoding="utf-8")
@classmethod
def from_xml_element(cls, elem, meshes=None):
@@ -2746,6 +2888,7 @@ def from_xml_element(cls, elem, meshes=None):
settings._weight_windows_on_from_xml_element(elem)
settings._shared_secondary_bank_from_xml_element(elem)
settings._weight_windows_file_from_xml_element(elem)
+ settings._weight_windows_exodus_from_xml_element(elem)
settings._weight_window_generators_from_xml_element(elem, meshes)
settings._weight_window_checkpoints_from_xml_element(elem)
settings._max_history_splits_from_xml_element(elem)
@@ -2759,7 +2902,7 @@ def from_xml_element(cls, elem, meshes=None):
return settings
@classmethod
- def from_xml(cls, path: PathLike = 'settings.xml'):
+ def from_xml(cls, path: PathLike = "settings.xml"):
"""Generate settings from XML file
.. versionadded:: 0.13.0
diff --git a/openmc/weight_windows.py b/openmc/weight_windows.py
index 63af2596efc..cebe3f44ec1 100644
--- a/openmc/weight_windows.py
+++ b/openmc/weight_windows.py
@@ -10,13 +10,20 @@
import h5py
import openmc
-from openmc.mesh import MeshBase, RectilinearMesh, CylindricalMesh, SphericalMesh, UnstructuredMesh
+from openmc.mesh import (
+ MeshBase,
+ RectilinearMesh,
+ CylindricalMesh,
+ SphericalMesh,
+ UnstructuredMesh,
+)
from openmc.tallies import Tallies
import openmc.checkvalue as cv
from openmc.checkvalue import PathLike
from ._xml import get_elem_list, get_text, clean_indentation
from .mixin import IDManagerMixin
from .particle_type import ParticleType
+from .utility_funcs import input_path
class WeightWindows(IDManagerMixin):
@@ -106,6 +113,7 @@ class WeightWindows(IDManagerMixin):
openmc.Settings
"""
+
next_id = 1
used_ids = set()
@@ -116,12 +124,12 @@ def __init__(
upper_ww_bounds: Iterable[float] | None = None,
upper_bound_ratio: float | None = None,
energy_bounds: Iterable[Real] | None = None,
- particle_type: str | int | openmc.ParticleType = 'neutron',
+ particle_type: str | int | openmc.ParticleType = "neutron",
survival_ratio: float = 3.0,
max_lower_bound_ratio: float | None = None,
max_split: int = 10,
- weight_cutoff: float = 1.e-38,
- id: int | None = None
+ weight_cutoff: float = 1.0e-38,
+ id: int | None = None,
):
self.mesh = mesh
self.id = id
@@ -132,12 +140,16 @@ def __init__(
self.lower_ww_bounds = lower_ww_bounds
if upper_ww_bounds is not None and upper_bound_ratio:
- raise ValueError("Exactly one of upper_ww_bounds and "
- "upper_bound_ratio must be present.")
+ raise ValueError(
+ "Exactly one of upper_ww_bounds and "
+ "upper_bound_ratio must be present."
+ )
if upper_ww_bounds is None and upper_bound_ratio is None:
- raise ValueError("Exactly one of upper_ww_bounds and "
- "upper_bound_ratio must be present.")
+ raise ValueError(
+ "Exactly one of upper_ww_bounds and "
+ "upper_bound_ratio must be present."
+ )
if upper_bound_ratio:
self.upper_ww_bounds = [
@@ -148,8 +160,9 @@ def __init__(
self.upper_ww_bounds = upper_ww_bounds
if len(self.lower_ww_bounds) != len(self.upper_ww_bounds):
- raise ValueError('Size of the lower and upper weight '
- 'window bounds do not match')
+ raise ValueError(
+ "Size of the lower and upper weight " "window bounds do not match"
+ )
self.survival_ratio = survival_ratio
@@ -161,17 +174,19 @@ def __init__(
self.weight_cutoff = weight_cutoff
def __repr__(self) -> str:
- string = type(self).__name__ + '\n'
- string += '{: <16}=\t{}\n'.format('\tID', self._id)
- string += '{: <16}=\t{}\n'.format('\tMesh', self.mesh)
- string += '{: <16}=\t{}\n'.format('\tParticle Type', self._particle_type)
- string += '{: <16}=\t{}\n'.format('\tEnergy Bounds', self._energy_bounds)
- string += '{: <16}=\t{}\n'.format('\tMax lower bound ratio', self.max_lower_bound_ratio)
- string += '{: <16}=\t{}\n'.format('\tLower WW Bounds', self._lower_ww_bounds)
- string += '{: <16}=\t{}\n'.format('\tUpper WW Bounds', self._upper_ww_bounds)
- string += '{: <16}=\t{}\n'.format('\tSurvival Ratio', self._survival_ratio)
- string += '{: <16}=\t{}\n'.format('\tMax Split', self._max_split)
- string += '{: <16}=\t{}\n'.format('\tWeight Cutoff', self._weight_cutoff)
+ string = type(self).__name__ + "\n"
+ string += "{: <16}=\t{}\n".format("\tID", self._id)
+ string += "{: <16}=\t{}\n".format("\tMesh", self.mesh)
+ string += "{: <16}=\t{}\n".format("\tParticle Type", self._particle_type)
+ string += "{: <16}=\t{}\n".format("\tEnergy Bounds", self._energy_bounds)
+ string += "{: <16}=\t{}\n".format(
+ "\tMax lower bound ratio", self.max_lower_bound_ratio
+ )
+ string += "{: <16}=\t{}\n".format("\tLower WW Bounds", self._lower_ww_bounds)
+ string += "{: <16}=\t{}\n".format("\tUpper WW Bounds", self._upper_ww_bounds)
+ string += "{: <16}=\t{}\n".format("\tSurvival Ratio", self._survival_ratio)
+ string += "{: <16}=\t{}\n".format("\tMax Split", self._max_split)
+ string += "{: <16}=\t{}\n".format("\tWeight Cutoff", self._weight_cutoff)
return string
def __eq__(self, other: WeightWindows) -> bool:
@@ -182,11 +197,13 @@ def __eq__(self, other: WeightWindows) -> bool:
# TODO: add ability to check mesh equality
# check several attributes directly
- attrs = ('particle_type',
- 'survival_ratio',
- 'max_lower_bound_ratio',
- 'max_split',
- 'weight_cutoff')
+ attrs = (
+ "particle_type",
+ "survival_ratio",
+ "max_lower_bound_ratio",
+ "max_split",
+ "weight_cutoff",
+ )
for attr in attrs:
if getattr(self, attr) != getattr(other, attr):
return False
@@ -209,7 +226,7 @@ def mesh(self) -> MeshBase:
@mesh.setter
def mesh(self, mesh: MeshBase):
- cv.check_type('Weight window mesh', mesh, MeshBase)
+ cv.check_type("Weight window mesh", mesh, MeshBase)
self._mesh = mesh
@property
@@ -220,7 +237,9 @@ def particle_type(self) -> ParticleType:
def particle_type(self, pt):
ptype = ParticleType(pt)
if ptype not in {ParticleType.NEUTRON, ParticleType.PHOTON}:
- raise ValueError("Weight windows can only be applied for neutrons or photons")
+ raise ValueError(
+ "Weight windows can only be applied for neutrons or photons"
+ )
self._particle_type = ptype
@property
@@ -229,7 +248,7 @@ def energy_bounds(self) -> Iterable[Real]:
@energy_bounds.setter
def energy_bounds(self, bounds: Iterable[float]):
- cv.check_type('Energy bounds', bounds, Iterable, Real)
+ cv.check_type("Energy bounds", bounds, Iterable, Real)
self._energy_bounds = np.asarray(bounds)
@property
@@ -244,11 +263,9 @@ def lower_ww_bounds(self) -> np.ndarray:
@lower_ww_bounds.setter
def lower_ww_bounds(self, bounds: Iterable[float]):
- cv.check_iterable_type('Lower WW bounds',
- bounds,
- Real,
- min_depth=1,
- max_depth=4)
+ cv.check_iterable_type(
+ "Lower WW bounds", bounds, Real, min_depth=1, max_depth=4
+ )
# reshape data according to mesh and energy bins
bounds = np.asarray(bounds)
if isinstance(self.mesh, UnstructuredMesh):
@@ -263,11 +280,9 @@ def upper_ww_bounds(self) -> np.ndarray:
@upper_ww_bounds.setter
def upper_ww_bounds(self, bounds: Iterable[float]):
- cv.check_iterable_type('Upper WW bounds',
- bounds,
- Real,
- min_depth=1,
- max_depth=4)
+ cv.check_iterable_type(
+ "Upper WW bounds", bounds, Real, min_depth=1, max_depth=4
+ )
# reshape data according to mesh and energy bins
bounds = np.asarray(bounds)
if isinstance(self.mesh, UnstructuredMesh):
@@ -282,8 +297,8 @@ def survival_ratio(self) -> float:
@survival_ratio.setter
def survival_ratio(self, val: float):
- cv.check_type('Survival ratio', val, Real)
- cv.check_greater_than('Survival ratio', val, 1.0, True)
+ cv.check_type("Survival ratio", val, Real)
+ cv.check_greater_than("Survival ratio", val, 1.0, True)
self._survival_ratio = val
@property
@@ -292,8 +307,8 @@ def max_lower_bound_ratio(self) -> float:
@max_lower_bound_ratio.setter
def max_lower_bound_ratio(self, val: float):
- cv.check_type('Maximum lower bound ratio', val, Real)
- cv.check_greater_than('Maximum lower bound ratio', val, 1.0, equality=True)
+ cv.check_type("Maximum lower bound ratio", val, Real)
+ cv.check_greater_than("Maximum lower bound ratio", val, 1.0, equality=True)
self._max_lower_bound_ratio = val
@property
@@ -302,7 +317,7 @@ def max_split(self) -> int:
@max_split.setter
def max_split(self, val: int):
- cv.check_type('Max split', val, Integral)
+ cv.check_type("Max split", val, Integral)
self._max_split = val
@property
@@ -311,8 +326,8 @@ def weight_cutoff(self) -> float:
@weight_cutoff.setter
def weight_cutoff(self, cutoff: float):
- cv.check_type('Weight cutoff', cutoff, Real)
- cv.check_greater_than('Weight cutoff', cutoff, 0.0, True)
+ cv.check_type("Weight cutoff", cutoff, Real)
+ cv.check_greater_than("Weight cutoff", cutoff, 0.0, True)
self._weight_cutoff = cutoff
def to_xml_element(self) -> ET.Element:
@@ -323,37 +338,37 @@ def to_xml_element(self) -> ET.Element:
element : lxml.etree._Element
XML element containing the weight window information
"""
- element = ET.Element('weight_windows')
+ element = ET.Element("weight_windows")
- element.set('id', str(self._id))
+ element.set("id", str(self._id))
- subelement = ET.SubElement(element, 'mesh')
+ subelement = ET.SubElement(element, "mesh")
subelement.text = str(self.mesh.id)
- subelement = ET.SubElement(element, 'particle_type')
+ subelement = ET.SubElement(element, "particle_type")
subelement.text = str(self.particle_type)
if self.energy_bounds is not None:
- subelement = ET.SubElement(element, 'energy_bounds')
- subelement.text = ' '.join(str(e) for e in self.energy_bounds)
+ subelement = ET.SubElement(element, "energy_bounds")
+ subelement.text = " ".join(str(e) for e in self.energy_bounds)
- subelement = ET.SubElement(element, 'lower_ww_bounds')
- subelement.text = ' '.join(str(b) for b in self.lower_ww_bounds.ravel('F'))
+ subelement = ET.SubElement(element, "lower_ww_bounds")
+ subelement.text = " ".join(str(b) for b in self.lower_ww_bounds.ravel("F"))
- subelement = ET.SubElement(element, 'upper_ww_bounds')
- subelement.text = ' '.join(str(b) for b in self.upper_ww_bounds.ravel('F'))
+ subelement = ET.SubElement(element, "upper_ww_bounds")
+ subelement.text = " ".join(str(b) for b in self.upper_ww_bounds.ravel("F"))
- subelement = ET.SubElement(element, 'survival_ratio')
+ subelement = ET.SubElement(element, "survival_ratio")
subelement.text = str(self.survival_ratio)
if self.max_lower_bound_ratio is not None:
- subelement = ET.SubElement(element, 'max_lower_bound_ratio')
+ subelement = ET.SubElement(element, "max_lower_bound_ratio")
subelement.text = str(self.max_lower_bound_ratio)
- subelement = ET.SubElement(element, 'max_split')
+ subelement = ET.SubElement(element, "max_split")
subelement.text = str(self.max_split)
- subelement = ET.SubElement(element, 'weight_cutoff')
+ subelement = ET.SubElement(element, "weight_cutoff")
subelement.text = str(self.weight_cutoff)
return element
@@ -375,7 +390,7 @@ def from_xml_element(cls, elem: ET.Element, meshes: dict[int, MeshBase]) -> Self
Weight windows object
"""
# Get mesh for weight windows
- mesh_id = int(get_text(elem, 'mesh'))
+ mesh_id = int(get_text(elem, "mesh"))
if mesh_id not in meshes:
raise ValueError(f'Could not locate mesh with ID "{mesh_id}"')
mesh = meshes[mesh_id]
@@ -384,20 +399,20 @@ def from_xml_element(cls, elem: ET.Element, meshes: dict[int, MeshBase]) -> Self
lower_ww_bounds = get_elem_list(elem, "lower_ww_bounds", float)
upper_ww_bounds = get_elem_list(elem, "upper_ww_bounds", float)
e_bounds = get_elem_list(elem, "energy_bounds", float)
- particle_type = get_text(elem, 'particle_type')
- survival_ratio = float(get_text(elem, 'survival_ratio'))
+ particle_type = get_text(elem, "particle_type")
+ survival_ratio = float(get_text(elem, "survival_ratio"))
ww_shape = (len(e_bounds) - 1,) + mesh.dimension[::-1]
lower_ww_bounds = np.array(lower_ww_bounds).reshape(ww_shape).T
upper_ww_bounds = np.array(upper_ww_bounds).reshape(ww_shape).T
max_lower_bound_ratio = None
- if get_text(elem, 'max_lower_bound_ratio'):
- max_lower_bound_ratio = float(get_text(elem, 'max_lower_bound_ratio'))
+ if get_text(elem, "max_lower_bound_ratio"):
+ max_lower_bound_ratio = float(get_text(elem, "max_lower_bound_ratio"))
- max_split = int(get_text(elem, 'max_split'))
- weight_cutoff = float(get_text(elem, 'weight_cutoff'))
- id = int(get_text(elem, 'id'))
+ max_split = int(get_text(elem, "max_split"))
+ weight_cutoff = float(get_text(elem, "weight_cutoff"))
+ id = int(get_text(elem, "id"))
return cls(
mesh=mesh,
@@ -409,7 +424,7 @@ def from_xml_element(cls, elem: ET.Element, meshes: dict[int, MeshBase]) -> Self
max_lower_bound_ratio=max_lower_bound_ratio,
max_split=max_split,
weight_cutoff=weight_cutoff,
- id=id
+ id=id,
)
@classmethod
@@ -429,25 +444,25 @@ def from_hdf5(cls, group: h5py.Group, meshes: dict[int, MeshBase]) -> Self:
A weight window object
"""
- id = int(group.name.split('/')[-1].lstrip('weight_windows'))
- mesh_id = group['mesh'][()]
+ id = int(group.name.split("/")[-1].lstrip("weight_windows"))
+ mesh_id = group["mesh"][()]
mesh = meshes[mesh_id]
- ptype = group['particle_type'][()].decode()
- e_bounds = group['energy_bounds'][()]
+ ptype = group["particle_type"][()].decode()
+ e_bounds = group["energy_bounds"][()]
# weight window bounds are stored with the shape (e, k, j, i)
# in C++ and HDF5 -- the opposite of how they are stored here
- shape = (e_bounds.size - 1, *mesh.dimension[::-1])
- lower_ww_bounds = group['lower_ww_bounds'][()].reshape(shape).T
- upper_ww_bounds = group['upper_ww_bounds'][()].reshape(shape).T
- survival_ratio = group['survival_ratio'][()]
+ shape = (e_bounds.size - 1, *mesh.dimension[::-1])
+ lower_ww_bounds = group["lower_ww_bounds"][()].reshape(shape).T
+ upper_ww_bounds = group["upper_ww_bounds"][()].reshape(shape).T
+ survival_ratio = group["survival_ratio"][()]
max_lower_bound_ratio = None
- if group.get('max_lower_bound_ratio') is not None:
- max_lower_bound_ratio = group['max_lower_bound_ratio'][()]
+ if group.get("max_lower_bound_ratio") is not None:
+ max_lower_bound_ratio = group["max_lower_bound_ratio"][()]
- max_split = group['max_split'][()]
- weight_cutoff = group['weight_cutoff'][()]
+ max_split = group["max_split"][()]
+ weight_cutoff = group["weight_cutoff"][()]
return cls(
mesh=mesh,
@@ -459,7 +474,7 @@ def from_hdf5(cls, group: h5py.Group, meshes: dict[int, MeshBase]) -> Self:
max_lower_bound_ratio=max_lower_bound_ratio,
max_split=max_split,
weight_cutoff=weight_cutoff,
- id=id
+ id=id,
)
@@ -479,7 +494,7 @@ def wwinp_to_wws(path: PathLike) -> WeightWindowsList:
"""
warnings.warn(
"This function is deprecated in favor of 'WeightWindowsList.from_wwinp'",
- FutureWarning
+ FutureWarning,
)
return WeightWindowsList.from_wwinp(path)
@@ -534,18 +549,18 @@ class WeightWindowGenerator:
Whether or not to apply weight windows on the fly.
"""
- _WWG_PARAMS = {'value': str, 'threshold': float, 'ratio': float}
+ _WWG_PARAMS = {"value": str, "threshold": float, "ratio": float}
def __init__(
self,
mesh: openmc.MeshBase,
energy_bounds: Sequence[float] | None = None,
- particle_type: str | int | openmc.ParticleType = 'neutron',
- method: str = 'magic',
+ particle_type: str | int | openmc.ParticleType = "neutron",
+ method: str = "magic",
targets: openmc.Tallies | Iterable[int] | None = None,
max_realizations: int = 1,
update_interval: int = 1,
- on_the_fly: bool = True
+ on_the_fly: bool = True,
):
self._update_parameters = None
@@ -561,7 +576,7 @@ def __init__(
self.on_the_fly = on_the_fly
def __repr__(self):
- string = type(self).__name__ + '\n'
+ string = type(self).__name__ + "\n"
string += f'\t{"Mesh":<20}=\t{self.mesh.id}\n'
string += f'\t{"Particle:":<20}=\t{str(self.particle_type)}\n'
string += f'\t{"Energy Bounds:":<20}=\t{self.energy_bounds}\n'
@@ -581,7 +596,7 @@ def mesh(self) -> openmc.MeshBase:
@mesh.setter
def mesh(self, m: openmc.MeshBase):
- cv.check_type('mesh', m, openmc.MeshBase)
+ cv.check_type("mesh", m, openmc.MeshBase)
self._mesh = m
@property
@@ -590,7 +605,7 @@ def energy_bounds(self) -> Iterable[Real]:
@energy_bounds.setter
def energy_bounds(self, eb: Iterable[float]):
- cv.check_type('energy bounds', eb, Iterable, Real)
+ cv.check_type("energy bounds", eb, Iterable, Real)
self._energy_bounds = eb
@property
@@ -601,7 +616,9 @@ def particle_type(self) -> ParticleType:
def particle_type(self, pt):
ptype = ParticleType(pt)
if ptype not in {ParticleType.NEUTRON, ParticleType.PHOTON}:
- raise ValueError("Weight windows can only be applied for neutrons or photons")
+ raise ValueError(
+ "Weight windows can only be applied for neutrons or photons"
+ )
self._particle_type = ptype
@property
@@ -610,15 +627,15 @@ def method(self) -> str:
@method.setter
def method(self, m: str):
- cv.check_type('generation method', m, str)
- cv.check_value('generation method', m, ('magic', 'fw_cadis'))
+ cv.check_type("generation method", m, str)
+ cv.check_value("generation method", m, ("magic", "fw_cadis"))
self._method = m
if self._update_parameters is not None:
try:
self._check_update_parameters()
except (TypeError, KeyError):
warnings.warn(f'Update parameters are invalid for the "{m}" method.')
-
+
@property
def targets(self) -> openmc.Tallies:
return self._targets
@@ -628,10 +645,10 @@ def targets(self, t):
if t is None:
self._targets = t
else:
- cv.check_type('Local FW-CADIS target tallies', t, Iterable)
- cv.check_greater_than('Local FW-CADIS target tallies', len(t), 0)
+ cv.check_type("Local FW-CADIS target tallies", t, Iterable)
+ cv.check_greater_than("Local FW-CADIS target tallies", len(t), 0)
if not isinstance(t, openmc.Tallies):
- cv.check_iterable_type('Local FW-CADIS target tallies', t, int)
+ cv.check_iterable_type("Local FW-CADIS target tallies", t, int)
t = np.asarray(list(t), dtype=int)
self._targets = t
@@ -641,8 +658,8 @@ def max_realizations(self) -> int:
@max_realizations.setter
def max_realizations(self, m: int):
- cv.check_type('max tally realizations', m, Integral)
- cv.check_greater_than('max tally realizations', m, 0)
+ cv.check_type("max tally realizations", m, Integral)
+ cv.check_greater_than("max tally realizations", m, 0)
self._max_realizations = m
@property
@@ -651,8 +668,8 @@ def update_interval(self) -> int:
@update_interval.setter
def update_interval(self, ui: int):
- cv.check_type('update interval', ui, Integral)
- cv.check_greater_than('update interval', ui , 0)
+ cv.check_type("update interval", ui, Integral)
+ cv.check_greater_than("update interval", ui, 0)
self._update_interval = ui
@property
@@ -660,14 +677,18 @@ def update_parameters(self) -> dict:
return self._update_parameters
def _check_update_parameters(self, params: dict):
- if self.method == 'magic' or self.method == 'fw_cadis':
+ if self.method == "magic" or self.method == "fw_cadis":
check_params = self._WWG_PARAMS
for key, val in params.items():
if key not in check_params:
- raise ValueError(f'Invalid param "{key}" for {self.method} '
- 'weight window generation')
- cv.check_type(f'weight window generation param: "{key}"', val, self._WWG_PARAMS[key])
+ raise ValueError(
+ f'Invalid param "{key}" for {self.method} '
+ "weight window generation"
+ )
+ cv.check_type(
+ f'weight window generation param: "{key}"', val, self._WWG_PARAMS[key]
+ )
@update_parameters.setter
def update_parameters(self, params: dict):
@@ -680,13 +701,13 @@ def on_the_fly(self) -> bool:
@on_the_fly.setter
def on_the_fly(self, otf: bool):
- cv.check_type('on the fly generation', otf, bool)
+ cv.check_type("on the fly generation", otf, bool)
self._on_the_fly = otf
def _update_parameters_subelement(self, element: ET.Element):
if not self.update_parameters:
return
- params_element = ET.SubElement(element, 'update_parameters')
+ params_element = ET.SubElement(element, "update_parameters")
for pname, value in self.update_parameters.items():
param_element = ET.SubElement(params_element, pname)
param_element.text = str(value)
@@ -703,7 +724,7 @@ def _sanitize_update_parameters(cls, method: str, update_parameters: dict):
update_parameters : dict
The update parameters as-read from the XML node (keys: str, values: str)
"""
- if method == 'magic' or method == 'fw_cadis':
+ if method == "magic" or method == "fw_cadis":
check_params = cls._WWG_PARAMS
for param, param_type in check_params.items():
@@ -711,38 +732,39 @@ def _sanitize_update_parameters(cls, method: str, update_parameters: dict):
update_parameters[param] = param_type(update_parameters[param])
def to_xml_element(self):
- """Creates a 'weight_window_generator' element to be written to an XML file.
- """
- element = ET.Element('weight_windows_generator')
+ """Creates a 'weight_window_generator' element to be written to an XML file."""
+ element = ET.Element("weight_windows_generator")
- mesh_elem = ET.SubElement(element, 'mesh')
+ mesh_elem = ET.SubElement(element, "mesh")
mesh_elem.text = str(self.mesh.id)
if self.energy_bounds is not None:
- subelement = ET.SubElement(element, 'energy_bounds')
- subelement.text = ' '.join(str(e) for e in self.energy_bounds)
- particle_elem = ET.SubElement(element, 'particle_type')
+ subelement = ET.SubElement(element, "energy_bounds")
+ subelement.text = " ".join(str(e) for e in self.energy_bounds)
+ particle_elem = ET.SubElement(element, "particle_type")
particle_elem.text = str(self.particle_type)
- realizations_elem = ET.SubElement(element, 'max_realizations')
+ realizations_elem = ET.SubElement(element, "max_realizations")
realizations_elem.text = str(self.max_realizations)
- update_interval_elem = ET.SubElement(element, 'update_interval')
+ update_interval_elem = ET.SubElement(element, "update_interval")
update_interval_elem.text = str(self.update_interval)
- otf_elem = ET.SubElement(element, 'on_the_fly')
+ otf_elem = ET.SubElement(element, "on_the_fly")
otf_elem.text = str(self.on_the_fly).lower()
- method_elem = ET.SubElement(element, 'method')
+ method_elem = ET.SubElement(element, "method")
method_elem.text = self.method
if self.targets is not None:
- if self.method != 'fw_cadis':
+ if self.method != "fw_cadis":
raise ValueError(
- "FW-CADIS update method is required in order to use " \
- "target tallies for WeightWindowGenerator.")
+ "FW-CADIS update method is required in order to use "
+ "target tallies for WeightWindowGenerator."
+ )
elif isinstance(self.targets, openmc.Tallies):
raise RuntimeError(
- "FW-CADIS target tallies must be checked to ensure they are " \
- "present on model.tallies. Use model.export_to_xml() or " \
- "model.export_to_model_xml() to link FW-CADIS target tallies.")
+ "FW-CADIS target tallies must be checked to ensure they are "
+ "present on model.tallies. Use model.export_to_xml() or "
+ "model.export_to_model_xml() to link FW-CADIS target tallies."
+ )
else:
- targets_elem = ET.SubElement(element, 'targets')
- targets_elem.text = ' '.join(str(tally_id) for tally_id in self.targets)
+ targets_elem = ET.SubElement(element, "targets")
+ targets_elem.text = " ".join(str(tally_id) for tally_id in self.targets)
if self.update_parameters is not None:
self._update_parameters_subelement(element)
@@ -768,30 +790,31 @@ def from_xml_element(cls, elem: ET.Element, meshes: dict) -> Self:
openmc.WeightWindowGenerator
"""
- mesh_id = int(get_text(elem, 'mesh'))
+ mesh_id = int(get_text(elem, "mesh"))
mesh = meshes[mesh_id]
-
+
energy_bounds = get_elem_list(elem, "energy_bounds", float)
- particle_type = get_text(elem, 'particle_type')
+ particle_type = get_text(elem, "particle_type")
wwg = cls(mesh, energy_bounds, particle_type)
- wwg.max_realizations = int(get_text(elem, 'max_realizations'))
- wwg.update_interval = int(get_text(elem, 'update_interval'))
- wwg.on_the_fly = bool(get_text(elem, 'on_the_fly'))
- wwg.method = get_text(elem, 'method')
- targets_elem = elem.find('targets')
+ wwg.max_realizations = int(get_text(elem, "max_realizations"))
+ wwg.update_interval = int(get_text(elem, "update_interval"))
+ wwg.on_the_fly = bool(get_text(elem, "on_the_fly"))
+ wwg.method = get_text(elem, "method")
+ targets_elem = elem.find("targets")
if targets_elem is not None:
- if wwg.method != 'fw_cadis':
+ if wwg.method != "fw_cadis":
raise ValueError(
- "FW-CADIS update method is required in order to use " \
- "target tallies for WeightWindowGenerator.")
+ "FW-CADIS update method is required in order to use "
+ "target tallies for WeightWindowGenerator."
+ )
else:
wwg.targets = get_elem_list(elem, "targets")
- if elem.find('update_parameters') is not None:
+ if elem.find("update_parameters") is not None:
update_parameters = {}
- params_elem = elem.find('update_parameters')
+ params_elem = elem.find("update_parameters")
for entry in params_elem:
update_parameters[entry.tag] = entry.text
@@ -800,7 +823,323 @@ def from_xml_element(cls, elem: ET.Element, meshes: dict) -> Self:
return wwg
-def hdf5_to_wws(path='weight_windows.h5') -> WeightWindowsList:
+
+class WeightWindowsExodus:
+ """Specification for building weight windows from an Exodus file.
+
+ The Exodus file is expected to contain multigroup adjoint flux stored as
+ CONSTANT MONOMIAL elemental variables (one variable per energy group).
+ At simulation initialization, OpenMC reads the mesh and flux, registers the
+ mesh as an unstructured (libMesh) mesh, applies FW-CADIS-style normalization
+ and creates the corresponding weight windows. An instance of this class can
+ be assigned to the :attr:`openmc.Settings.weight_windows_exodus` attribute.
+
+ Requires OpenMC to be built with libMesh support.
+
+ .. versionadded:: 0.16.1
+
+ Parameters
+ ----------
+ file : path-like
+ Path to the Exodus file containing the mesh and adjoint flux
+ adjoint_flux_variables : iterable of str
+ Names of the elemental variables containing the adjoint flux, one per
+ energy group, ordered consistently with `energy_bounds` (ascending
+ energy). Solvers that write group 0 as the fastest group (e.g.
+ Griffin) require the variables to be listed thermal-first.
+ energy_bounds : iterable of float or openmc.mgxs.EnergyGroups
+ Monotonically increasing energy group boundaries in [eV]. The number
+ of boundaries must be one more than the number of flux variables. An
+ :class:`openmc.mgxs.EnergyGroups` instance may be passed directly.
+ timestep : int, optional
+ Zero-based index of the Exodus time step to read the flux from. If
+ not given, the last time step in the file is used.
+ particle_type : str or int or openmc.ParticleType
+ Particle type the weight windows apply to (default: 'neutron')
+ survival_ratio : float, optional
+ Ratio of the survival weight to the lower weight window bound for
+ rouletting. If not given, the default of the transport code (3.0)
+ applies.
+ upper_bound_ratio : float, optional
+ Ratio of the upper to lower weight window bounds. If not given, the
+ default of the transport code (5.0) applies.
+ max_split : int, optional
+ Maximum allowable number of particles when splitting. If not given,
+ the default of the transport code (10) applies.
+
+ Attributes
+ ----------
+ file : pathlib.Path
+ Path to the Exodus file containing the mesh and adjoint flux
+ adjoint_flux_variables : list of str
+ Names of the elemental variables containing the adjoint flux
+ energy_bounds : numpy.ndarray of float
+ Monotonically increasing energy group boundaries in [eV]
+ timestep : int or None
+ Zero-based index of the Exodus time step to read the flux from
+ particle_type : openmc.ParticleType
+ Particle type the weight windows apply to
+ survival_ratio : float or None
+ Ratio of the survival weight to the lower weight window bound
+ upper_bound_ratio : float or None
+ Ratio of the upper to lower weight window bounds
+ max_split : int or None
+ Maximum allowable number of particles when splitting
+
+ See Also
+ --------
+ openmc.Settings.weight_windows_exodus
+
+ """
+
+ def __init__(
+ self,
+ file: PathLike,
+ adjoint_flux_variables: Iterable[str],
+ energy_bounds,
+ timestep: int | None = None,
+ particle_type: str | int | openmc.ParticleType = "neutron",
+ survival_ratio: float | None = None,
+ upper_bound_ratio: float | None = None,
+ max_split: int | None = None,
+ ):
+ self.file = file
+ self.adjoint_flux_variables = adjoint_flux_variables
+ self.energy_bounds = energy_bounds
+ self.timestep = timestep
+ self.particle_type = particle_type
+ self.survival_ratio = survival_ratio
+ self.upper_bound_ratio = upper_bound_ratio
+ self.max_split = max_split
+ self._check_consistency()
+
+ def _check_consistency(self):
+ """Cross-attribute checks mirroring those performed by the C++ layer"""
+ n_groups = len(self.adjoint_flux_variables)
+ if self.energy_bounds.size != n_groups + 1:
+ raise ValueError(
+ f"Number of energy bounds ({self.energy_bounds.size}) must be "
+ f"one more than the number of adjoint flux variables "
+ f"({n_groups})."
+ )
+ # compare using the transport code defaults when a value is unset
+ survival = 3.0 if self.survival_ratio is None else self.survival_ratio
+ upper = 5.0 if self.upper_bound_ratio is None else self.upper_bound_ratio
+ if upper <= survival:
+ raise ValueError(
+ f"Upper bound ratio ({upper}) must be larger than the "
+ f"survival ratio ({survival})."
+ )
+
+ def __repr__(self) -> str:
+ string = type(self).__name__ + "\n"
+ string += f'\t{"File":<20}=\t{self.file}\n'
+ string += f'\t{"Flux variables":<20}=\t{self.adjoint_flux_variables}\n'
+ string += f'\t{"Energy bounds":<20}=\t{self.energy_bounds}\n'
+ string += f'\t{"Timestep":<20}=\t{self.timestep}\n'
+ string += f'\t{"Particle":<20}=\t{str(self.particle_type)}\n'
+ string += f'\t{"Survival ratio":<20}=\t{self.survival_ratio}\n'
+ string += f'\t{"Upper bound ratio":<20}=\t{self.upper_bound_ratio}\n'
+ string += f'\t{"Max split":<20}=\t{self.max_split}\n'
+ return string
+
+ def __eq__(self, other) -> bool:
+ if not isinstance(other, WeightWindowsExodus):
+ return False
+ attrs = (
+ "file",
+ "adjoint_flux_variables",
+ "timestep",
+ "particle_type",
+ "survival_ratio",
+ "upper_bound_ratio",
+ "max_split",
+ )
+ for attr in attrs:
+ if getattr(self, attr) != getattr(other, attr):
+ return False
+ return np.array_equal(self.energy_bounds, other.energy_bounds)
+
+ @property
+ def file(self) -> Path:
+ return self._file
+
+ @file.setter
+ def file(self, value: PathLike):
+ cv.check_type("Exodus weight windows file", value, PathLike)
+ self._file = input_path(value)
+
+ @property
+ def adjoint_flux_variables(self) -> list[str]:
+ return self._adjoint_flux_variables
+
+ @adjoint_flux_variables.setter
+ def adjoint_flux_variables(self, variables: Iterable[str]):
+ cv.check_type("adjoint flux variables", variables, Iterable, str)
+ variables = list(variables)
+ cv.check_greater_than("number of adjoint flux variables", len(variables), 0)
+ self._adjoint_flux_variables = variables
+
+ @property
+ def energy_bounds(self) -> np.ndarray:
+ return self._energy_bounds
+
+ @energy_bounds.setter
+ def energy_bounds(self, bounds):
+ # accept an openmc.mgxs.EnergyGroups directly; local import avoids a
+ # circular import between openmc.weight_windows and openmc.mgxs
+ from openmc.mgxs import EnergyGroups
+
+ if isinstance(bounds, EnergyGroups):
+ bounds = bounds.group_edges
+ cv.check_type("energy bounds", bounds, Iterable, Real)
+ bounds = np.asarray(bounds, dtype=float)
+ if bounds.ndim != 1 or bounds.size < 2:
+ raise ValueError("At least two energy bounds must be provided.")
+ if np.any(np.diff(bounds) <= 0.0):
+ raise ValueError("Energy bounds must be strictly increasing.")
+ self._energy_bounds = bounds
+
+ @property
+ def timestep(self) -> int | None:
+ return self._timestep
+
+ @timestep.setter
+ def timestep(self, value: int | None):
+ if value is not None:
+ cv.check_type("timestep", value, Integral)
+ cv.check_greater_than("timestep", value, 0, equality=True)
+ self._timestep = value
+
+ @property
+ def particle_type(self) -> ParticleType:
+ return self._particle_type
+
+ @particle_type.setter
+ def particle_type(self, pt):
+ ptype = ParticleType(pt)
+ if ptype not in {ParticleType.NEUTRON, ParticleType.PHOTON}:
+ raise ValueError(
+ "Weight windows can only be applied for neutrons or photons"
+ )
+ self._particle_type = ptype
+
+ @property
+ def survival_ratio(self) -> float | None:
+ return self._survival_ratio
+
+ @survival_ratio.setter
+ def survival_ratio(self, value: float | None):
+ if value is not None:
+ cv.check_type("survival ratio", value, Real)
+ cv.check_greater_than("survival ratio", value, 1.0)
+ self._survival_ratio = value
+
+ @property
+ def upper_bound_ratio(self) -> float | None:
+ return self._upper_bound_ratio
+
+ @upper_bound_ratio.setter
+ def upper_bound_ratio(self, value: float | None):
+ if value is not None:
+ cv.check_type("upper bound ratio", value, Real)
+ cv.check_greater_than("upper bound ratio", value, 1.0)
+ self._upper_bound_ratio = value
+
+ @property
+ def max_split(self) -> int | None:
+ return self._max_split
+
+ @max_split.setter
+ def max_split(self, value: int | None):
+ if value is not None:
+ cv.check_type("max split", value, Integral)
+ cv.check_greater_than("max split", value, 1)
+ self._max_split = value
+
+ def to_xml_element(self) -> ET.Element:
+ """Create a 'weight_windows_exodus' element to be written to an XML file."""
+ self._check_consistency()
+
+ element = ET.Element("weight_windows_exodus")
+
+ subelement = ET.SubElement(element, "file")
+ subelement.text = str(self.file)
+
+ subelement = ET.SubElement(element, "adjoint_flux_variables")
+ subelement.text = " ".join(self.adjoint_flux_variables)
+
+ subelement = ET.SubElement(element, "energy_bounds")
+ subelement.text = " ".join(str(e) for e in self.energy_bounds)
+
+ if self.timestep is not None:
+ subelement = ET.SubElement(element, "timestep")
+ subelement.text = str(self.timestep)
+
+ subelement = ET.SubElement(element, "particle_type")
+ subelement.text = str(self.particle_type)
+
+ # optional values are omitted so that the transport code defaults apply
+ if self.survival_ratio is not None:
+ subelement = ET.SubElement(element, "survival_ratio")
+ subelement.text = str(self.survival_ratio)
+
+ if self.upper_bound_ratio is not None:
+ subelement = ET.SubElement(element, "upper_bound_ratio")
+ subelement.text = str(self.upper_bound_ratio)
+
+ if self.max_split is not None:
+ subelement = ET.SubElement(element, "max_split")
+ subelement.text = str(self.max_split)
+
+ clean_indentation(element)
+
+ return element
+
+ @classmethod
+ def from_xml_element(cls, elem: ET.Element) -> Self:
+ """Create a WeightWindowsExodus object from an XML element
+
+ Parameters
+ ----------
+ elem : lxml.etree._Element
+ XML element
+
+ Returns
+ -------
+ openmc.WeightWindowsExodus
+ """
+ file = get_text(elem, "file")
+ variables = get_elem_list(elem, "adjoint_flux_variables", str)
+ energy_bounds = get_elem_list(elem, "energy_bounds", float)
+
+ wwe = cls(file, variables, energy_bounds)
+
+ timestep = get_text(elem, "timestep")
+ if timestep is not None:
+ wwe.timestep = int(timestep)
+
+ particle_type = get_text(elem, "particle_type")
+ if particle_type is not None:
+ wwe.particle_type = particle_type
+
+ survival_ratio = get_text(elem, "survival_ratio")
+ if survival_ratio is not None:
+ wwe.survival_ratio = float(survival_ratio)
+
+ upper_bound_ratio = get_text(elem, "upper_bound_ratio")
+ if upper_bound_ratio is not None:
+ wwe.upper_bound_ratio = float(upper_bound_ratio)
+
+ max_split = get_text(elem, "max_split")
+ if max_split is not None:
+ wwe.max_split = int(max_split)
+
+ wwe._check_consistency()
+ return wwe
+
+
+def hdf5_to_wws(path="weight_windows.h5") -> WeightWindowsList:
"""Create a WeightWindowsList from a weight windows HDF5 file
.. versionadded:: 0.14.0
@@ -816,7 +1155,7 @@ def hdf5_to_wws(path='weight_windows.h5') -> WeightWindowsList:
"""
warnings.warn(
"This function is deprecated in favor of 'WeightWindowsList.from_hdf5'",
- FutureWarning
+ FutureWarning,
)
return WeightWindowsList.from_hdf5(path)
@@ -832,11 +1171,12 @@ class WeightWindowsList(list):
An iterable of WeightWindows objects to initialize the list with
"""
+
def __init__(self, iterable: Iterable[WeightWindows] = ()):
super().__init__(iterable)
@classmethod
- def from_hdf5(cls, path: PathLike = 'weight_windows.h5') -> Self:
+ def from_hdf5(cls, path: PathLike = "weight_windows.h5") -> Self:
"""Create WeightWindowsList from a weight windows HDF5 file.
Parameters
@@ -853,12 +1193,12 @@ def from_hdf5(cls, path: PathLike = 'weight_windows.h5') -> Self:
with h5py.File(path) as h5_file:
# read in all of the meshes in the mesh node
meshes = {}
- for mesh_group in h5_file['meshes']:
- mesh = MeshBase.from_hdf5(h5_file['meshes'][mesh_group])
+ for mesh_group in h5_file["meshes"]:
+ mesh = MeshBase.from_hdf5(h5_file["meshes"][mesh_group])
meshes[mesh.id] = mesh
wws = [
WeightWindows.from_hdf5(ww, meshes)
- for ww in h5_file['weight_windows'].values()
+ for ww in h5_file["weight_windows"].values()
]
return cls(wws)
@@ -887,35 +1227,36 @@ def from_wwinp(cls, path: PathLike) -> Self:
# header value checks
if _if != 1:
- raise ValueError(f'Found incorrect file type, if: {_if}')
+ raise ValueError(f"Found incorrect file type, if: {_if}")
if iv > 1:
# read number of time bins for each particle, 'nt(1...ni)'
- nt = np.fromstring(wwinp.readline(), sep=' ', dtype=int)
+ nt = np.fromstring(wwinp.readline(), sep=" ", dtype=int)
# raise error if time bins are present for now
- raise ValueError('Time-dependent weight windows '
- 'are not yet supported')
+ raise ValueError(
+ "Time-dependent weight windows " "are not yet supported"
+ )
else:
nt = ni * [1]
# read number of energy bins for each particle, 'ne(1...ni)'
- ne = np.fromstring(wwinp.readline(), sep=' ', dtype=int)
+ ne = np.fromstring(wwinp.readline(), sep=" ", dtype=int)
# read coarse mesh dimensions and lower left corner
- mesh_description = np.fromstring(wwinp.readline(), sep=' ')
+ mesh_description = np.fromstring(wwinp.readline(), sep=" ")
nfx, nfy, nfz = mesh_description[:3].astype(int)
xyz0 = mesh_description[3:]
# read cylindrical and spherical mesh vectors if present
if nr == 16:
# read number of coarse bins
- line_arr = np.fromstring(wwinp.readline(), sep=' ')
+ line_arr = np.fromstring(wwinp.readline(), sep=" ")
ncx, ncy, ncz = line_arr[:3].astype(int)
# read polar vector (x1, y1, z1)
xyz1 = line_arr[3:]
# read azimuthal vector (x2, y2, z2)
- line_arr = np.fromstring(wwinp.readline(), sep=' ')
+ line_arr = np.fromstring(wwinp.readline(), sep=" ")
xyz2 = line_arr[:3]
# Get polar and azimuthal axes
@@ -924,13 +1265,17 @@ def from_wwinp(cls, path: PathLike) -> Self:
# Check for polar axis other than (0, 0, 1)
norm = np.linalg.norm(polar_axis)
- if not np.isclose(polar_axis[2]/norm, 1.0):
- raise NotImplementedError('Polar axis not aligned to z-axis not supported')
+ if not np.isclose(polar_axis[2] / norm, 1.0):
+ raise NotImplementedError(
+ "Polar axis not aligned to z-axis not supported"
+ )
# Check for azimuthal axis other than (1, 0, 0)
norm = np.linalg.norm(azimuthal_axis)
- if not np.isclose(azimuthal_axis[0]/norm, 1.0):
- raise NotImplementedError('Azimuthal axis not aligned to x-axis not supported')
+ if not np.isclose(azimuthal_axis[0] / norm, 1.0):
+ raise NotImplementedError(
+ "Azimuthal axis not aligned to x-axis not supported"
+ )
# read geometry type
nwg = int(line_arr[-1])
@@ -938,13 +1283,14 @@ def from_wwinp(cls, path: PathLike) -> Self:
elif nr == 10:
# read rectilinear data:
# number of coarse mesh bins and mesh type
- ncx, ncy, ncz, nwg = \
- np.fromstring(wwinp.readline(), sep=' ').astype(int)
+ ncx, ncy, ncz, nwg = np.fromstring(wwinp.readline(), sep=" ").astype(
+ int
+ )
else:
- raise RuntimeError(f'Invalid mesh description (nr) found: {nr}')
+ raise RuntimeError(f"Invalid mesh description (nr) found: {nr}")
# read BLOCK 2 and BLOCK 3 data into a single array
- ww_data = np.fromstring(wwinp.read(), sep=' ')
+ ww_data = np.fromstring(wwinp.read(), sep=" ")
# extract mesh data from the ww_data array
start_idx = 0
@@ -952,26 +1298,30 @@ def from_wwinp(cls, path: PathLike) -> Self:
# first values in the mesh definition arrays are the first
# coordinate of the grid
end_idx = start_idx + 1 + 3 * ncx
- i0, i_vals = ww_data[start_idx], ww_data[start_idx+1:end_idx]
+ i0, i_vals = ww_data[start_idx], ww_data[start_idx + 1 : end_idx]
start_idx = end_idx
end_idx = start_idx + 1 + 3 * ncy
- j0, j_vals = ww_data[start_idx], ww_data[start_idx+1:end_idx]
+ j0, j_vals = ww_data[start_idx], ww_data[start_idx + 1 : end_idx]
start_idx = end_idx
end_idx = start_idx + 1 + 3 * ncz
- k0, k_vals = ww_data[start_idx], ww_data[start_idx+1:end_idx]
+ k0, k_vals = ww_data[start_idx], ww_data[start_idx + 1 : end_idx]
start_idx = end_idx
# mesh consistency checks
if nr == 16 and nwg == 1 or nr == 10 and nwg != 1:
- raise ValueError(f'Mesh description in header ({nr}) '
- f'does not match the mesh type ({nwg})')
+ raise ValueError(
+ f"Mesh description in header ({nr}) "
+ f"does not match the mesh type ({nwg})"
+ )
if nr == 10 and (xyz0 != (i0, j0, k0)).any():
- raise ValueError(f'Mesh origin in the header ({xyz0}) '
- f' does not match the origin in the mesh '
- f' description ({i0, j0, k0})')
+ raise ValueError(
+ f"Mesh origin in the header ({xyz0}) "
+ f" does not match the origin in the mesh "
+ f" description ({i0, j0, k0})"
+ )
# create openmc mesh object
grids = []
@@ -979,13 +1329,16 @@ def from_wwinp(cls, path: PathLike) -> Self:
for grid0, grid_vals, n_pnts in mesh_definition:
# file spec checks for the mesh definition
if (grid_vals[2::3] != 1.0).any():
- raise ValueError('One or more mesh ratio value, qx, '
- 'is not equal to one')
+ raise ValueError(
+ "One or more mesh ratio value, qx, " "is not equal to one"
+ )
s = int(grid_vals[::3].sum())
if s != n_pnts:
- raise ValueError(f'Sum of the fine bin entries, {s}, does '
- f'not match the number of fine bins, {n_pnts}')
+ raise ValueError(
+ f"Sum of the fine bin entries, {s}, does "
+ f"not match the number of fine bins, {n_pnts}"
+ )
# extend the grid based on the next coarse bin endpoint, px
# and the number of fine bins in the coarse bin, sx
@@ -1004,19 +1357,16 @@ def from_wwinp(cls, path: PathLike) -> Self:
r_grid=grids[0],
z_grid=grids[1],
phi_grid=grids[2],
- origin = xyz0,
+ origin=xyz0,
)
elif nwg == 3:
mesh = SphericalMesh(
- r_grid=grids[0],
- theta_grid=grids[1],
- phi_grid=grids[2],
- origin = xyz0
+ r_grid=grids[0], theta_grid=grids[1], phi_grid=grids[2], origin=xyz0
)
# extract weight window values from array
wws = cls()
- for ne_i, nt_i, particle_type in zip(ne, nt, ('neutron', 'photon')):
+ for ne_i, nt_i, particle_type in zip(ne, nt, ("neutron", "photon")):
# no information to read for this particle if
# either the energy bins or time bins are empty
if ne_i == 0 or nt_i == 0:
@@ -1050,17 +1400,19 @@ def from_wwinp(cls, path: PathLike) -> Self:
start_idx = end_idx
# create a weight window object
- ww = WeightWindows(id=None,
- mesh=mesh,
- lower_ww_bounds=ww_values,
- upper_bound_ratio=5.0,
- energy_bounds=energy_bounds,
- particle_type=particle_type)
+ ww = WeightWindows(
+ id=None,
+ mesh=mesh,
+ lower_ww_bounds=ww_values,
+ upper_bound_ratio=5.0,
+ energy_bounds=energy_bounds,
+ particle_type=particle_type,
+ )
wws.append(ww)
return wws
- def export_to_hdf5(self, path: PathLike = 'weight_windows.h5', **init_kwargs):
+ def export_to_hdf5(self, path: PathLike = "weight_windows.h5", **init_kwargs):
"""Write weight windows to an HDF5 file.
Parameters
@@ -1072,11 +1424,12 @@ def export_to_hdf5(self, path: PathLike = 'weight_windows.h5', **init_kwargs):
"""
import openmc.lib
- cv.check_type('path', path, PathLike)
+
+ cv.check_type("path", path, PathLike)
# Create a temporary model with the weight windows
model = openmc.Model()
- sph = openmc.Sphere(boundary_type='vacuum')
+ sph = openmc.Sphere(boundary_type="vacuum")
cell = openmc.Cell(region=-sph)
model.geometry = openmc.Geometry([cell])
model.settings.weight_windows = self
diff --git a/src/mesh.cpp b/src/mesh.cpp
index 46698a6e390..2e3aa0136fc 100644
--- a/src/mesh.cpp
+++ b/src/mesh.cpp
@@ -3626,6 +3626,23 @@ LibMesh::LibMesh(libMesh::MeshBase& input_mesh, double length_multiplier)
initialize();
}
+// create the mesh from an externally constructed libMesh mesh, transferring
+// ownership to OpenMC
+LibMesh::LibMesh(unique_ptr input_mesh,
+ double length_multiplier, const std::string& filename)
+{
+ if (!input_mesh->is_replicated()) {
+ fatal_error("At present LibMesh tallies require a replicated mesh. Please "
+ "ensure 'input_mesh' is a libMesh::ReplicatedMesh.");
+ }
+
+ unique_m_ = std::move(input_mesh);
+ m_ = unique_m_.get();
+ filename_ = filename;
+ set_length_multiplier(length_multiplier);
+ initialize();
+}
+
// create the mesh from an input file
LibMesh::LibMesh(const std::string& filename, double length_multiplier)
{
@@ -3842,7 +3859,7 @@ void LibMesh::set_score_data(const std::string& var_name,
unsigned int std_dev_num = variable_map_.at(std_dev_name);
for (auto it = m_->local_elements_begin(); it != m_->local_elements_end();
- it++) {
+ it++) {
if (!(*it)->active()) {
continue;
}
diff --git a/src/settings.cpp b/src/settings.cpp
index 15a7b9c27c9..008d7934702 100644
--- a/src/settings.cpp
+++ b/src/settings.cpp
@@ -1240,6 +1240,11 @@ void read_settings_xml(pugi::xml_node root)
std::make_unique(node_ww));
}
+ // Weight windows built from adjoint flux in an Exodus file (libMesh)
+ if (check_for_node(root, "weight_windows_exodus")) {
+ read_weight_windows_exodus(root.child("weight_windows_exodus"));
+ }
+
// Enable weight windows by default if one or more are present
if (variance_reduction::weight_windows.size() > 0)
settings::weight_windows_on = true;
diff --git a/src/weight_windows.cpp b/src/weight_windows.cpp
index d3565eaaf63..92a33534069 100644
--- a/src/weight_windows.cpp
+++ b/src/weight_windows.cpp
@@ -30,6 +30,17 @@
#include
+#ifdef OPENMC_LIBMESH_ENABLED
+#include "libmesh/dof_map.h"
+#include "libmesh/elem.h"
+#include "libmesh/equation_systems.h"
+#include "libmesh/exodusII_io.h"
+#include "libmesh/explicit_system.h"
+#include "libmesh/mesh_communication.h"
+#include "libmesh/numeric_vector.h"
+#include "libmesh/replicated_mesh.h"
+#endif
+
namespace openmc {
//==============================================================================
@@ -930,6 +941,275 @@ void WeightWindowsGenerator::update() const
// Non-member functions
//==============================================================================
+//! Compute FW-CADIS weight window bounds from multigroup adjoint flux
+//
+//! Mirrors the FW_CADIS branch of WeightWindows::update_weights(): positive
+//! flux values are inverted and normalized by twice the global maximum of the
+//! inverted values. Elements with non-positive flux keep the sentinel -1.0.
+//! \param[in] flux flux[g][e]: adjoint flux for group g, element e
+//! \param[in] upper_bound_ratio ratio of upper to lower ww bounds
+//! \param[out] flat_lower lower bounds, flat layout [g * n_elem + e]
+//! \param[out] flat_upper upper bounds, same layout
+//! \return false if no positive flux value exists anywhere
+static bool fw_cadis_bounds(const vector>& flux,
+ double upper_bound_ratio, vector& flat_lower,
+ vector& flat_upper)
+{
+ const size_t n_groups = flux.size();
+ const size_t n_elem = n_groups ? flux[0].size() : 0;
+
+ flat_lower.assign(n_groups * n_elem, -1.0);
+ flat_upper.assign(n_groups * n_elem, -1.0);
+
+ // Invert positive flux values and track the global maximum
+ double inv_max = 0.0;
+ for (size_t g = 0; g < n_groups; ++g) {
+ for (size_t e = 0; e < n_elem; ++e) {
+ if (flux[g][e] > 0.0) {
+ double inv = 1.0 / flux[g][e];
+ flat_lower[g * n_elem + e] = inv;
+ inv_max = std::max(inv_max, inv);
+ }
+ }
+ }
+
+ if (inv_max <= 0.0)
+ return false;
+
+ const double norm_factor = 1.0 / (2.0 * inv_max);
+ for (size_t i = 0; i < n_groups * n_elem; ++i) {
+ if (flat_lower[i] >= 0.0) {
+ flat_lower[i] *= norm_factor;
+ flat_upper[i] = flat_lower[i] * upper_bound_ratio;
+ }
+ }
+ return true;
+}
+
+void read_weight_windows_exodus(pugi::xml_node node)
+{
+#ifndef OPENMC_LIBMESH_ENABLED
+ (void)node;
+ fatal_error(" requires OpenMC to be compiled "
+ "with libMesh support (-DOPENMC_USE_LIBMESH=on).");
+#else
+ // Make sure required elements are present
+ const vector required_elems {
+ "file", "adjoint_flux_variables", "energy_bounds"};
+ for (const auto& elem : required_elems) {
+ if (!check_for_node(node, elem.c_str())) {
+ fatal_error(
+ fmt::format("Must specify <{}> for .", elem));
+ }
+ }
+
+ const std::string file = get_node_value(node, "file", true);
+ if (!file_exists(file))
+ fatal_error(fmt::format(
+ ": mesh file '{}' does not exist.", file));
+
+ // One elemental variable per energy group, ordered by ascending energy
+ // consistently with
+ const vector flux_vars =
+ get_node_array(node, "adjoint_flux_variables");
+ if (flux_vars.empty())
+ fatal_error(": must "
+ "list at least one variable.");
+ const int n_groups = static_cast(flux_vars.size());
+
+ const vector e_bounds = get_node_array(node, "energy_bounds");
+ if (static_cast(e_bounds.size()) != n_groups + 1)
+ fatal_error(fmt::format(
+ ": must have exactly {} values "
+ "for {} group(s), but {} were provided.",
+ n_groups + 1, n_groups, e_bounds.size()));
+ for (int g = 0; g < n_groups; ++g) {
+ if (e_bounds[g] >= e_bounds[g + 1])
+ fatal_error(
+ fmt::format(": must be strictly "
+ "increasing; bounds[{}] = {} >= bounds[{}] = {}.",
+ g, e_bounds[g], g + 1, e_bounds[g + 1]));
+ }
+
+ // is 0-based; default -1 selects the last step in the file
+ const int ts_user = check_for_node(node, "timestep")
+ ? std::stoi(get_node_value(node, "timestep", true))
+ : -1;
+
+ const std::string p_type_str = check_for_node(node, "particle_type")
+ ? get_node_value(node, "particle_type", true)
+ : "neutron";
+
+ const double survival_ratio =
+ check_for_node(node, "survival_ratio")
+ ? std::stod(get_node_value(node, "survival_ratio", true))
+ : 3.0;
+ if (survival_ratio <= 1)
+ fatal_error("Survival to lower weight window ratio must bigger than 1 "
+ "and less than the upper to lower weight window ratio.");
+
+ const double upper_bound_ratio =
+ check_for_node(node, "upper_bound_ratio")
+ ? std::stod(get_node_value(node, "upper_bound_ratio", true))
+ : 5.0;
+ if (upper_bound_ratio <= survival_ratio)
+ fatal_error(fmt::format(
+ ": ({}) must be larger "
+ "than ({}).",
+ upper_bound_ratio, survival_ratio));
+
+ const int max_split = check_for_node(node, "max_split")
+ ? std::stoi(get_node_value(node, "max_split", true))
+ : 10;
+ if (max_split <= 1)
+ fatal_error("max split must be larger than 1");
+
+ // Read the mesh and all group flux variables in a single pass. Note that
+ // copy_elemental_solution() must be called on the same ExodusII_IO object
+ // that performed read(), and allow_renumbering(false) must be set before
+ // read() so that element IDs match the Exodus element block entries.
+ if (!settings::libmesh_comm)
+ fatal_error(": no libMesh communicator is "
+ "initialized.");
+
+ auto mesh = make_unique(*settings::libmesh_comm, 3);
+ mesh->allow_renumbering(false);
+
+ libMesh::ExodusII_IO exo(*mesh);
+ exo.read(file);
+
+ // The reader only populates rank 0, so replicate the mesh to the other MPI
+ // ranks before use (no-op in serial)
+ libMesh::MeshCommunication().broadcast(*mesh);
+ mesh->prepare_for_use();
+
+ const int n_elem = static_cast(mesh->n_active_elem());
+ if (n_elem == 0)
+ fatal_error(fmt::format(
+ ": mesh file '{}' has no elements.", file));
+
+ // Resolve the requested time step (Exodus steps are 1-based internally).
+ // The file is only open on rank 0, so query metadata there and broadcast.
+ const auto& comm = mesh->comm();
+ int n_steps = 0;
+ if (comm.rank() == 0)
+ n_steps = static_cast(exo.get_time_steps().size());
+ comm.broadcast(n_steps);
+
+ const int ts_1based = (ts_user < 0) ? n_steps : (ts_user + 1);
+ if (ts_1based < 1 || ts_1based > n_steps)
+ fatal_error(fmt::format(
+ ": requested timestep {} is out of range "
+ "[0, {}) for file '{}'.",
+ (ts_user < 0 ? n_steps - 1 : ts_user), n_steps, file));
+
+ // Verify every requested variable exists before reading any of them
+ if (comm.rank() == 0) {
+ const auto& exo_elem_vars = exo.get_elem_var_names();
+ for (const auto& vname : flux_vars) {
+ if (std::find(exo_elem_vars.begin(), exo_elem_vars.end(), vname) ==
+ exo_elem_vars.end()) {
+ std::string available;
+ for (size_t vi = 0; vi < exo_elem_vars.size(); ++vi) {
+ if (vi)
+ available += ", ";
+ available += exo_elem_vars[vi];
+ }
+ fatal_error(fmt::format(
+ ": variable '{}' not found in '{}'.\n"
+ " Available element variables: [{}]",
+ vname, file, available));
+ }
+ }
+ }
+
+ // Index flux arrays by elem->id() - first_id so that the flux index matches
+ // the mesh bin computed by LibMesh::get_bin_from_element()
+ const auto first_id = (*mesh->elements_begin())->id();
+
+ // flux[g][e] matches the (n_energy_bins, n_mesh_bins) layout of lower_ww_
+ vector> flux(n_groups);
+
+ for (int g = 0; g < n_groups; ++g) {
+ // Use a fresh EquationSystems per group to avoid DOF conflicts from
+ // multiple active variables
+ libMesh::EquationSystems eq_sys(*mesh);
+ auto& sys = eq_sys.add_system("adjoint_ww");
+ sys.add_variable(flux_vars[g], libMesh::CONSTANT, libMesh::MONOMIAL);
+ eq_sys.init();
+
+ exo.copy_elemental_solution(sys, flux_vars[g], flux_vars[g], ts_1based);
+
+ // Under MPI the solution vector is distributed; gather the full vector
+ // onto every rank (collective, no-op in serial)
+ std::vector soln_local;
+ sys.solution->localize(soln_local);
+
+ const libMesh::DofMap& dof_map = sys.get_dof_map();
+ flux[g].assign(n_elem, 0.0);
+ for (const auto* elem : mesh->active_element_ptr_range()) {
+ std::vector dofs;
+ dof_map.dof_indices(elem, dofs);
+ if (dofs.size() != 1)
+ fatal_error(fmt::format(
+ ": expected one DOF per element but found "
+ "{} for element {}.",
+ dofs.size(), elem->id()));
+ const auto bin = elem->id() - first_id;
+ if (bin >= static_cast(n_elem))
+ fatal_error(fmt::format(
+ ": element IDs in '{}' are not contiguous "
+ "(element {} with first ID {}).",
+ file, elem->id(), first_id));
+ flux[g][bin] = soln_local[dofs[0]];
+ }
+ }
+
+ // Register the mesh with OpenMC, transferring ownership
+ int32_t mesh_id = 1;
+ for (const auto& m : model::meshes)
+ mesh_id = std::max(mesh_id, m->id_ + 1);
+
+ model::meshes.push_back(make_unique(std::move(mesh), 1.0, file));
+ model::meshes.back()->set_id(mesh_id);
+
+ // Normalize (FW-CADIS) and build the WeightWindows object
+ vector flat_lower;
+ vector flat_upper;
+ if (!fw_cadis_bounds(flux, upper_bound_ratio, flat_lower, flat_upper))
+ fatal_error(fmt::format(
+ ": all adjoint flux values across all {} "
+ "group(s) in '{}' are zero or negative -- cannot compute FW-CADIS "
+ "weight windows.",
+ n_groups, file));
+
+ // set_mesh() and set_energy_bounds() must precede set_bounds() since both
+ // trigger allocate_ww_bounds()
+ WeightWindows* wws = WeightWindows::create();
+ wws->set_mesh(model::mesh_map.at(mesh_id));
+ wws->set_particle_type(ParticleType {p_type_str});
+ wws->set_energy_bounds(span(e_bounds.data(), e_bounds.size()));
+ wws->survival_ratio() = survival_ratio;
+ wws->max_split() = max_split;
+ wws->set_bounds(span(flat_lower.data(), flat_lower.size()),
+ span(flat_upper.data(), flat_upper.size()));
+
+ std::string varlist;
+ for (int g = 0; g < n_groups; ++g) {
+ if (g)
+ varlist += ", ";
+ varlist += flux_vars[g];
+ }
+ write_message(
+ fmt::format("Loaded {}-group adjoint weight windows from '{}':\n"
+ " {} elements, variables [{}], timestep {}, "
+ "upper_bound_ratio {:g}.",
+ n_groups, file, n_elem, varlist, (ts_user < 0 ? n_steps - 1 : ts_user),
+ upper_bound_ratio),
+ 5);
+#endif // OPENMC_LIBMESH_ENABLED
+}
+
std::pair search_weight_window(const Particle& p)
{
// TODO: this is a linear search - should do something more clever
@@ -1419,4 +1699,4 @@ extern "C" int openmc_weight_windows_import(const char* filename)
return 0;
}
-} // namespace openmc
+} // namespace openmc
\ No newline at end of file
diff --git a/tests/regression_tests/weightwindows_exodus/inputs_true.dat b/tests/regression_tests/weightwindows_exodus/inputs_true.dat
new file mode 100644
index 00000000000..70067349c7e
--- /dev/null
+++ b/tests/regression_tests/weightwindows_exodus/inputs_true.dat
@@ -0,0 +1,36 @@
+
+
+
+
+
+
+
+
+
+ |
+
+
+
+
+
+
+
+
+ fixed source
+ 100
+ 3
+ 0
+
+
+ 0.5 0.8 0.2
+
+
+ 12345
+
+ test_out.e
+ adjoint_flux_g1 adjoint_flux_g0
+ 0.0 0.625 20000000.0
+ neutron
+
+
+
diff --git a/tests/regression_tests/weightwindows_exodus/results_true.dat b/tests/regression_tests/weightwindows_exodus/results_true.dat
new file mode 100644
index 00000000000..e69de29bb2d
diff --git a/tests/regression_tests/weightwindows_exodus/test.py b/tests/regression_tests/weightwindows_exodus/test.py
new file mode 100644
index 00000000000..2997a90b638
--- /dev/null
+++ b/tests/regression_tests/weightwindows_exodus/test.py
@@ -0,0 +1,42 @@
+import openmc
+from tests.testing_harness import PyAPITestHarness
+
+
+def test_weight_windows_exodus():
+ # Build a minimal basic model shell to trigger initial setup parsing
+ model = openmc.Model()
+
+ # Simple fuel material
+ uo2 = openmc.Material(name="uo2")
+ uo2.add_nuclide("U235", 1.0)
+ uo2.set_density("g/cm3", 10.0)
+ model.materials.append(uo2)
+
+ # Tiny box geometry
+ box = openmc.model.RectangularPrism(
+ width=0.6, height=0.9, origin=(0.3, 0.45), boundary_type="vacuum"
+ )
+ z_top = openmc.ZPlane(z0=0.3, boundary_type="vacuum")
+ z_bot = openmc.ZPlane(z0=-0.0, boundary_type="vacuum")
+ cell = openmc.Cell(fill=uo2, region=-box & +z_bot & -z_top)
+ model.geometry = openmc.Geometry([cell])
+
+ # Setup standard simulation run settings
+ model.settings.batches = 3
+ model.settings.inactive = 0
+ model.settings.particles = 100
+ model.settings.run_mode = "fixed source"
+ model.settings.seed = 12345
+ model.settings.source = openmc.IndependentSource(
+ space=openmc.stats.Point((0.5, 0.8, 0.2))
+ )
+ # Configure your new Exodus feature
+ model.settings.weight_windows_exodus = openmc.WeightWindowsExodus(
+ file="test_out.e", # Path to your small test mesh in this folder
+ adjoint_flux_variables=["adjoint_flux_g1", "adjoint_flux_g0"],
+ energy_bounds=[0.0, 0.625, 2.0e7],
+ )
+
+ # Use OpenMC's test harness to run the simulation
+ harness = PyAPITestHarness("statepoint.3.h5", model)
+ harness.main()
diff --git a/tests/regression_tests/weightwindows_exodus/test_out.e b/tests/regression_tests/weightwindows_exodus/test_out.e
new file mode 100755
index 00000000000..8773b4f67fd
Binary files /dev/null and b/tests/regression_tests/weightwindows_exodus/test_out.e differ
diff --git a/tests/unit_tests/test_settings.py b/tests/unit_tests/test_settings.py
index bdb3ea8fe9f..fb98f3d2304 100644
--- a/tests/unit_tests/test_settings.py
+++ b/tests/unit_tests/test_settings.py
@@ -10,85 +10,108 @@
def test_export_to_xml(run_in_tmpdir):
- tmp_properties_file = 'properties_test.h5'
+ tmp_properties_file = "properties_test.h5"
- s = openmc.Settings(run_mode='fixed source', batches=1000, seed=17)
+ s = openmc.Settings(run_mode="fixed source", batches=1000, seed=17)
s.generations_per_batch = 10
s.inactive = 100
s.particles = 1000000
s.max_lost_particles = 5
s.rel_max_lost_particles = 1e-4
- s.keff_trigger = {'type': 'std_dev', 'threshold': 0.001}
- s.energy_mode = 'continuous-energy'
+ s.keff_trigger = {"type": "std_dev", "threshold": 0.001}
+ s.energy_mode = "continuous-energy"
s.max_order = 5
s.max_tracks = 1234
s.source = openmc.IndependentSource(space=openmc.stats.Point())
- s.output = {'summary': True, 'tallies': False, 'path': 'here'}
+ s.output = {"summary": True, "tallies": False, "path": "here"}
s.verbosity = 7
- s.sourcepoint = {'batches': [50, 150, 500, 1000], 'separate': True,
- 'write': True, 'overwrite': True, 'mcpl': True}
- s.statepoint = {'batches': [50, 150, 500, 1000]}
- s.surf_source_read = {'path': 'surface_source_1.h5'}
- s.surf_source_write = {'surface_ids': [2], 'max_particles': 200}
+ s.sourcepoint = {
+ "batches": [50, 150, 500, 1000],
+ "separate": True,
+ "write": True,
+ "overwrite": True,
+ "mcpl": True,
+ }
+ s.statepoint = {"batches": [50, 150, 500, 1000]}
+ s.surf_source_read = {"path": "surface_source_1.h5"}
+ s.surf_source_write = {"surface_ids": [2], "max_particles": 200}
s.surface_grazing_ratio = 0.7
s.surface_grazing_cutoff = 0.1
s.confidence_intervals = True
s.ptables = True
s.plot_seed = 100
s.survival_biasing = True
- s.cutoff = {'weight': 0.25, 'weight_avg': 0.5, 'energy_neutron': 1.0e-5,
- 'survival_normalization': True,
- 'energy_photon': 1000.0, 'energy_electron': 1.0e-5,
- 'energy_positron': 1.0e-5, 'time_neutron': 1.0e-5,
- 'time_photon': 1.0e-5, 'time_electron': 1.0e-5,
- 'time_positron': 1.0e-5}
+ s.cutoff = {
+ "weight": 0.25,
+ "weight_avg": 0.5,
+ "energy_neutron": 1.0e-5,
+ "survival_normalization": True,
+ "energy_photon": 1000.0,
+ "energy_electron": 1.0e-5,
+ "energy_positron": 1.0e-5,
+ "time_neutron": 1.0e-5,
+ "time_photon": 1.0e-5,
+ "time_electron": 1.0e-5,
+ "time_positron": 1.0e-5,
+ }
mesh = openmc.RegularMesh()
- mesh.lower_left = (-10., -10., -10.)
- mesh.upper_right = (10., 10., 10.)
+ mesh.lower_left = (-10.0, -10.0, -10.0)
+ mesh.upper_right = (10.0, 10.0, 10.0)
mesh.dimension = (5, 5, 5)
s.entropy_mesh = mesh
s.trigger_active = True
s.trigger_max_batches = 10000
s.trigger_batch_interval = 50
s.no_reduce = False
- s.tabular_legendre = {'enable': True, 'num_points': 50}
- s.temperature = {'default': 293.6, 'method': 'interpolation',
- 'multipole': True, 'range': (200., 1000.)}
+ s.tabular_legendre = {"enable": True, "num_points": 50}
+ s.temperature = {
+ "default": 293.6,
+ "method": "interpolation",
+ "multipole": True,
+ "range": (200.0, 1000.0),
+ }
s.properties_file = tmp_properties_file
s.trace = (10, 1, 20)
s.track = [(1, 1, 1), (2, 1, 1)]
s.ufs_mesh = mesh
- s.resonance_scattering = {'enable': True, 'method': 'rvs',
- 'energy_min': 1.0, 'energy_max': 1000.0,
- 'nuclides': ['U235', 'U238', 'Pu239']}
+ s.resonance_scattering = {
+ "enable": True,
+ "method": "rvs",
+ "energy_min": 1.0,
+ "energy_max": 1000.0,
+ "nuclides": ["U235", "U238", "Pu239"],
+ }
s.volume_calculations = openmc.VolumeCalculation(
- domains=[openmc.Cell()], samples=1000, lower_left=(-10., -10., -10.),
- upper_right = (10., 10., 10.))
+ domains=[openmc.Cell()],
+ samples=1000,
+ lower_left=(-10.0, -10.0, -10.0),
+ upper_right=(10.0, 10.0, 10.0),
+ )
s.create_fission_neutrons = True
s.create_delayed_neutrons = False
s.log_grid_bins = 2000
s.photon_transport = False
- s.electron_treatment = 'led'
+ s.electron_treatment = "led"
s.atomic_relaxation = False
s.write_initial_source = True
- s.weight_window_checkpoints = {'surface': True, 'collision': False}
+ s.weight_window_checkpoints = {"surface": True, "collision": False}
source_region_mesh = openmc.RegularMesh()
source_region_mesh.dimension = [2, 2, 2]
source_region_mesh.lower_left = [-2, -2, -2]
source_region_mesh.upper_right = [2, 2, 2]
root_universe = openmc.Universe()
s.random_ray = {
- 'distance_inactive': 10.0,
- 'distance_active': 100.0,
- 'ray_source': openmc.IndependentSource(
- space=openmc.stats.Box((-1., -1., -1.), (1., 1., 1.))
+ "distance_inactive": 10.0,
+ "distance_active": 100.0,
+ "ray_source": openmc.IndependentSource(
+ space=openmc.stats.Box((-1.0, -1.0, -1.0), (1.0, 1.0, 1.0))
),
- 'source_region_meshes': [(source_region_mesh, [root_universe])],
- 'volume_estimator': 'hybrid',
- 'source_shape': 'linear',
- 'volume_normalized_flux_tallies': True,
- 'adjoint': False,
- 'sample_method': 'halton'
+ "source_region_meshes": [(source_region_mesh, [root_universe])],
+ "volume_estimator": "hybrid",
+ "source_shape": "linear",
+ "volume_normalized_flux_tallies": True,
+ "adjoint": False,
+ "sample_method": "halton",
}
s.max_particle_events = 100
s.max_secondaries = 1_000_000
@@ -98,29 +121,33 @@ def test_export_to_xml(run_in_tmpdir):
# Make sure exporting XML works
s.export_to_xml()
-
# Generate settings from XML
s = openmc.Settings.from_xml()
- assert s.run_mode == 'fixed source'
+ assert s.run_mode == "fixed source"
assert s.batches == 1000
assert s.generations_per_batch == 10
assert s.inactive == 100
assert s.particles == 1000000
assert s.max_lost_particles == 5
assert s.rel_max_lost_particles == 1e-4
- assert s.keff_trigger == {'type': 'std_dev', 'threshold': 0.001}
- assert s.energy_mode == 'continuous-energy'
+ assert s.keff_trigger == {"type": "std_dev", "threshold": 0.001}
+ assert s.energy_mode == "continuous-energy"
assert s.max_order == 5
assert s.max_tracks == 1234
assert isinstance(s.source[0], openmc.IndependentSource)
assert isinstance(s.source[0].space, openmc.stats.Point)
- assert s.output == {'summary': True, 'tallies': False, 'path': 'here'}
+ assert s.output == {"summary": True, "tallies": False, "path": "here"}
assert s.verbosity == 7
- assert s.sourcepoint == {'batches': [50, 150, 500, 1000], 'separate': True,
- 'write': True, 'overwrite': True, 'mcpl': True}
- assert s.statepoint == {'batches': [50, 150, 500, 1000]}
- assert s.surf_source_read['path'].name == 'surface_source_1.h5'
- assert s.surf_source_write == {'surface_ids': [2], 'max_particles': 200}
+ assert s.sourcepoint == {
+ "batches": [50, 150, 500, 1000],
+ "separate": True,
+ "write": True,
+ "overwrite": True,
+ "mcpl": True,
+ }
+ assert s.statepoint == {"batches": [50, 150, 500, 1000]}
+ assert s.surf_source_read["path"].name == "surface_source_1.h5"
+ assert s.surf_source_write == {"surface_ids": [2], "max_particles": 200}
assert s.surface_grazing_ratio == 0.7
assert s.surface_grazing_cutoff == 0.1
assert s.confidence_intervals
@@ -128,65 +155,80 @@ def test_export_to_xml(run_in_tmpdir):
assert s.plot_seed == 100
assert s.seed == 17
assert s.survival_biasing
- assert s.cutoff == {'weight': 0.25, 'weight_avg': 0.5,
- 'survival_normalization': True,
- 'energy_neutron': 1.0e-5, 'energy_photon': 1000.0,
- 'energy_electron': 1.0e-5, 'energy_positron': 1.0e-5,
- 'time_neutron': 1.0e-5, 'time_photon': 1.0e-5,
- 'time_electron': 1.0e-5, 'time_positron': 1.0e-5}
+ assert s.cutoff == {
+ "weight": 0.25,
+ "weight_avg": 0.5,
+ "survival_normalization": True,
+ "energy_neutron": 1.0e-5,
+ "energy_photon": 1000.0,
+ "energy_electron": 1.0e-5,
+ "energy_positron": 1.0e-5,
+ "time_neutron": 1.0e-5,
+ "time_photon": 1.0e-5,
+ "time_electron": 1.0e-5,
+ "time_positron": 1.0e-5,
+ }
assert isinstance(s.entropy_mesh, openmc.RegularMesh)
- assert s.entropy_mesh.lower_left == [-10., -10., -10.]
- assert s.entropy_mesh.upper_right == [10., 10., 10.]
+ assert s.entropy_mesh.lower_left == [-10.0, -10.0, -10.0]
+ assert s.entropy_mesh.upper_right == [10.0, 10.0, 10.0]
assert s.entropy_mesh.dimension == (5, 5, 5)
assert s.trigger_active
assert s.trigger_max_batches == 10000
assert s.trigger_batch_interval == 50
assert not s.no_reduce
- assert s.tabular_legendre == {'enable': True, 'num_points': 50}
- assert s.temperature == {'default': 293.6, 'method': 'interpolation',
- 'multipole': True, 'range': [200., 1000.]}
+ assert s.tabular_legendre == {"enable": True, "num_points": 50}
+ assert s.temperature == {
+ "default": 293.6,
+ "method": "interpolation",
+ "multipole": True,
+ "range": [200.0, 1000.0],
+ }
assert s.properties_file == Path(tmp_properties_file)
assert s.trace == [10, 1, 20]
assert s.track == [(1, 1, 1), (2, 1, 1)]
assert isinstance(s.ufs_mesh, openmc.RegularMesh)
- assert s.ufs_mesh.lower_left == [-10., -10., -10.]
- assert s.ufs_mesh.upper_right == [10., 10., 10.]
+ assert s.ufs_mesh.lower_left == [-10.0, -10.0, -10.0]
+ assert s.ufs_mesh.upper_right == [10.0, 10.0, 10.0]
assert s.ufs_mesh.dimension == (5, 5, 5)
- assert s.resonance_scattering == {'enable': True, 'method': 'rvs',
- 'energy_min': 1.0, 'energy_max': 1000.0,
- 'nuclides': ['U235', 'U238', 'Pu239']}
+ assert s.resonance_scattering == {
+ "enable": True,
+ "method": "rvs",
+ "energy_min": 1.0,
+ "energy_max": 1000.0,
+ "nuclides": ["U235", "U238", "Pu239"],
+ }
assert s.create_fission_neutrons
assert not s.create_delayed_neutrons
assert s.log_grid_bins == 2000
assert not s.photon_transport
- assert s.electron_treatment == 'led'
+ assert s.electron_treatment == "led"
assert not s.atomic_relaxation
assert s.write_initial_source
assert len(s.volume_calculations) == 1
vol = s.volume_calculations[0]
- assert vol.domain_type == 'cell'
+ assert vol.domain_type == "cell"
assert len(vol.ids) == 1
assert vol.samples == 1000
- assert vol.lower_left == (-10., -10., -10.)
- assert vol.upper_right == (10., 10., 10.)
- assert s.weight_window_checkpoints == {'surface': True, 'collision': False}
+ assert vol.lower_left == (-10.0, -10.0, -10.0)
+ assert vol.upper_right == (10.0, 10.0, 10.0)
+ assert s.weight_window_checkpoints == {"surface": True, "collision": False}
assert s.max_particle_events == 100
- assert s.random_ray['distance_inactive'] == 10.0
- assert s.random_ray['distance_active'] == 100.0
- assert s.random_ray['ray_source'].space.lower_left == [-1., -1., -1.]
- assert s.random_ray['ray_source'].space.upper_right == [1., 1., 1.]
- assert 'source_region_meshes' in s.random_ray
- assert len(s.random_ray['source_region_meshes']) == 1
- mesh_and_domains = s.random_ray['source_region_meshes'][0]
+ assert s.random_ray["distance_inactive"] == 10.0
+ assert s.random_ray["distance_active"] == 100.0
+ assert s.random_ray["ray_source"].space.lower_left == [-1.0, -1.0, -1.0]
+ assert s.random_ray["ray_source"].space.upper_right == [1.0, 1.0, 1.0]
+ assert "source_region_meshes" in s.random_ray
+ assert len(s.random_ray["source_region_meshes"]) == 1
+ mesh_and_domains = s.random_ray["source_region_meshes"][0]
recovered_mesh = mesh_and_domains[0]
assert recovered_mesh.dimension == (2, 2, 2)
- assert recovered_mesh.lower_left == [-2., -2., -2.]
- assert recovered_mesh.upper_right == [2., 2., 2.]
- assert s.random_ray['volume_estimator'] == 'hybrid'
- assert s.random_ray['source_shape'] == 'linear'
- assert s.random_ray['volume_normalized_flux_tallies']
- assert not s.random_ray['adjoint']
- assert s.random_ray['sample_method'] == 'halton'
+ assert recovered_mesh.lower_left == [-2.0, -2.0, -2.0]
+ assert recovered_mesh.upper_right == [2.0, 2.0, 2.0]
+ assert s.random_ray["volume_estimator"] == "hybrid"
+ assert s.random_ray["source_shape"] == "linear"
+ assert s.random_ray["volume_normalized_flux_tallies"]
+ assert not s.random_ray["adjoint"]
+ assert s.random_ray["sample_method"] == "halton"
assert s.max_secondaries == 1_000_000
assert s.source_rejection_fraction == 0.01
assert s.free_gas_threshold == 800.0
@@ -197,10 +239,10 @@ def test_properties_file_load(tmp_path, mpi_intracomm):
# Session 1: export a structurally valid properties file via the C++ API,
# then collect the cell/material structure so we can patch it with h5py.
- cell_instances = {} # {cell_id: n_instances} — material cells only
- mat_densities = {} # {mat_id: original atom/b-cm density}
+ cell_instances = {} # {cell_id: n_instances} — material cells only
+ mat_densities = {} # {mat_id: original atom/b-cm density}
- props_path = tmp_path / 'properties.h5'
+ props_path = tmp_path / "properties.h5"
with openmc.lib.TemporarySession(model, intracomm=mpi_intracomm):
openmc.lib.export_properties(str(props_path))
for cell_id, cell in openmc.lib.cells.items():
@@ -210,25 +252,26 @@ def test_properties_file_load(tmp_path, mpi_intracomm):
except NotImplementedError:
pass
for mat_id, mat in openmc.lib.materials.items():
- mat_densities[mat_id] = mat.get_density('atom/b-cm')
+ mat_densities[mat_id] = mat.get_density("atom/b-cm")
assert any(n > 1 for n in cell_instances.values())
# Patch the exported properties file overwriting temperatures
# with per-instance values and scale material atom densities.
density_factor = 0.75
- with h5py.File(props_path, 'r+') as f:
- cells_grp = f['geometry/cells']
+ with h5py.File(props_path, "r+") as f:
+ cells_grp = f["geometry/cells"]
for cell_id, n in cell_instances.items():
- cell_grp = cells_grp[f'cell {cell_id}']
- del cell_grp['temperature']
+ cell_grp = cells_grp[f"cell {cell_id}"]
+ del cell_grp["temperature"]
cell_grp.create_dataset(
- 'temperature', data=[500.0 + 5.0 * i for i in range(n)]
+ "temperature", data=[500.0 + 5.0 * i for i in range(n)]
)
for mat_id, orig_density in mat_densities.items():
- f['materials'][f'material {mat_id}'].attrs['atom_density'] = \
+ f["materials"][f"material {mat_id}"].attrs["atom_density"] = (
orig_density * density_factor
+ )
# now apply the newly patched properties file using the settings
# and load the model again, checking that the new temperature and
@@ -243,6 +286,38 @@ def test_properties_file_load(tmp_path, mpi_intracomm):
for mat_id, orig_density in mat_densities.items():
mat = openmc.lib.materials[mat_id]
- assert mat.get_density('atom/b-cm') == pytest.approx(
+ assert mat.get_density("atom/b-cm") == pytest.approx(
orig_density * density_factor, rel=1e-5
)
+
+
+def test_weight_windows_exodus():
+ # Instantiate the new class with dummy parameters
+ wwe = openmc.WeightWindowsExodus(
+ file="test_out.e",
+ adjoint_flux_variables=["adjoint_flux_g1", "adjoint_flux_g1"],
+ energy_bounds=[0.0, 1.0e6, 2.0e7],
+ timestep=2,
+ particle_type="neutron",
+ survival_ratio=3.5,
+ upper_bound_ratio=5.5,
+ max_split=12,
+ )
+
+ settings = openmc.Settings()
+ settings.weight_windows_exodus = wwe
+
+ # Convert settings to XML element
+ root = settings.to_xml_element()
+ node = root.find("weight_windows_exodus")
+ assert node is not None
+
+ # Assert all values map cleanly to the expected XML tags
+ assert node.find("file").text == "test_out.e"
+ assert node.find("adjoint_flux_variables").text == "adjoint_flux_g1 adjoint_flux_g1"
+ assert node.find("energy_bounds").text == "0.0 1000000.0 20000000.0"
+ assert node.find("timestep").text == "2"
+ assert node.find("particle_type").text == "neutron"
+ assert node.find("survival_ratio").text == "3.5"
+ assert node.find("upper_bound_ratio").text == "5.5"
+ assert node.find("max_split").text == "12"