crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-21 06:25:30 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Restricted file source (#2916)

Browse source 
Co-authored-by: church89 <l.chierici89@gmail.com>
Co-authored-by: Paul Romano <paul.k.romano@gmail.com>
This commit is contained in:
Erik B Knudsen 2024-05-15 22:49:04 +02:00 • committed by GitHub
parent a7d6939c11
commit 1702b4554b
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
32 changed files with 823 additions and 240 deletions
Download patch file Download diff file Expand all files Collapse all files

38
docs/source/io_formats/settings.rst
View file

@ -534,8 +534,9 @@ attributes/sub-elements:
*Default*: 1.0
:type:
Indicator of source type. One of ``independent``, ``file``, ``compiled``, or ``mesh``.
The type of the source will be determined by this attribute if it is present.
Indicator of source type. One of ``independent``, ``file``, ``compiled``, or
``mesh``. The type of the source will be determined by this attribute if it
is present.
:particle:
The source particle type, either ``neutron`` or ``photon``.
@ -716,6 +717,39 @@ attributes/sub-elements:
mesh element and follows the format for :ref:`source_element`. The number of
``<source>`` sub-elements should correspond to the number of mesh elements.
:constraints:
This sub-element indicates the presence of constraints on sampled source
sites (see :ref:`usersguide_source_constraints` for details). It may have
the following sub-elements:
:domain_ids:
The unique IDs of domains for which source sites must be within.
*Default*: None
:domain_type:
The type of each domain for source rejection ("cell", "material", or
"universe").
*Default*: None
:fissionable:
A boolean indicating whether source sites must be sampled within a
material that is fissionable in order to be accepted.
:time_bounds:
A pair of times in [s] indicating the lower and upper bound for a time
interval that source particles must be within.
:energy_bounds:
A pair of energies in [eV] indicating the lower and upper bound for an
energy interval that source particles must be within.
:rejection_strategy:
Either "resample", indicating that source sites should be resampled when
one is rejected, or "kill", indicating that a rejected source site is
assigned zero weight.
.. _univariate:
Univariate Probability Distributions

48
docs/source/usersguide/settings.rst
View file

@ -420,6 +420,54 @@ string, which gets passed down to the ``openmc_create_source()`` function::
settings.source = openmc.CompiledSource('libsource.so', '3.5e6')
.. _usersguide_source_constraints:
Source Constraints
------------------
All source classes in OpenMC have the ability to apply a set of "constraints"
that limit which sampled source sites are actually used for transport. The most
common use case is to sample source sites over some simple spatial distribution
(e.g., uniform over a box) and then only accept those that appear in a given
cell or material. This can be done with a domain constraint, which can be
specified as follows::
source_cell = openmc.Cell(...)
...
spatial_dist = openmc.stats.Box((-10., -10., -10.), (10., 10., 10.))
source = openmc.IndependentSource(
space=spatial_dist,
constraints={'domains': [source_cell]}
)
For k-eigenvalue problems, a convenient constraint is available that limits
source sites to those sampled in a fissionable material::
source = openmc.IndependentSource(
space=spatial_dist, constraints={'fissionable': True}
)
Constraints can also be placed on a range of energies or times::
# Only use source sites between 500 keV and 1 MeV and with times under 1 sec
source = openmc.FileSource(
'source.h5',
constraints={'energy_bounds': [500.0e3, 1.0e6], 'time_bounds': [0.0, 1.0]}
)
Normally, when a source site is rejected, a new one will be resampled until one
is found that meets the constraints. However, the rejection strategy can be
changed so that a rejected site will just not be simulated by specifying::
source = openmc.IndependentSource(
space=spatial_dist,
constraints={'domains': [cell], 'rejection_strategy': 'kill'}
)
In this case, the actual number of particles simulated may be less than what you
specified in :attr:`Settings.particles`.
.. _usersguide_entropy:
---------------

9
examples/assembly/assembly.py
View file

@ -112,11 +112,10 @@ def assembly_model():
model.settings.batches = 150
model.settings.inactive = 50
model.settings.particles = 1000
model.settings.source = openmc.IndependentSource(space=openmc.stats.Box(
(-pitch/2, -pitch/2, -1),
(pitch/2, pitch/2, 1),
only_fissionable=True
))
model.settings.source = openmc.IndependentSource(
space=openmc.stats.Box((-pitch/2, -pitch/2, -1), (pitch/2, pitch/2, 1)),
constraints={'fissionable': True}
)
# NOTE: We never actually created a Materials object. When you export/run
# using the Model object, if no materials were assigned it will look through

5
examples/pincell_depletion/restart_depletion.py
View file

@ -26,8 +26,9 @@ settings.particles = 10000
# Create an initial uniform spatial source distribution over fissionable zones
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:], only_fissionable=True)
settings.source = openmc.IndependentSource(space=uniform_dist)
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
settings.source = openmc.IndependentSource(
space=uniform_dist, constraints={'fissionable': True})
entropy_mesh = openmc.RegularMesh()
entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50]

5
examples/pincell_depletion/run_depletion.py
View file

@ -71,8 +71,9 @@ settings.particles = 1000
# Create an initial uniform spatial source distribution over fissionable zones
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:], only_fissionable=True)
settings.source = openmc.IndependentSource(space=uniform_dist)
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
settings.source = openmc.IndependentSource(
space=uniform_dist, constraints={'fissionable': True})
entropy_mesh = openmc.RegularMesh()
entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50]

5
include/openmc/distribution_spatial.h
View file

@ -116,6 +116,11 @@ public:
//! \return Sampled element index and position within that element
std::pair<int32_t, Position> sample_mesh(uint64_t* seed) const;
//! Sample a mesh element
//! \param seed Pseudorandom number seed pointer
//! \return Sampled element index
int32_t sample_element_index(uint64_t* seed) const;
//! For unstructured meshes, ensure that elements are all linear tetrahedra
void check_element_types() const;

90
include/openmc/source.h
View file

@ -4,6 +4,7 @@
#ifndef OPENMC_SOURCE_H
#define OPENMC_SOURCE_H
#include <limits>
#include <unordered_set>
#include "pugixml.hpp"
@ -39,19 +40,72 @@ extern vector<unique_ptr<Source>> external_sources;
//==============================================================================
//! Abstract source interface
//
//! The Source class provides the interface that must be implemented by derived
//! classes, namely the sample() method that returns a sampled source site. From
//! this base class, source rejection is handled within the
//! sample_with_constraints() method. However, note that some classes directly
//! check for constraints for efficiency reasons (like IndependentSource), in
//! which case the constraints_applied() method indicates that constraints
//! should not be checked a second time from the base class.
//==============================================================================
class Source {
public:
// Constructors, destructors
Source() = default;
explicit Source(pugi::xml_node node);
virtual ~Source() = default;
// Methods that must be implemented
// Methods that can be overridden
virtual double strength() const { return strength_; }
//! Sample a source site and apply constraints
//
//! \param[inout] seed Pseudorandom seed pointer
//! \return Sampled site
SourceSite sample_with_constraints(uint64_t* seed) const;
//! Sample a source site (without applying constraints)
//
//! Sample from the external source distribution
//! \param[inout] seed Pseudorandom seed pointer
//! \return Sampled site
virtual SourceSite sample(uint64_t* seed) const = 0;
// Methods that can be overridden
virtual double strength() const { return 1.0; }
static unique_ptr<Source> create(pugi::xml_node node);
protected:
// Domain types
enum class DomainType { UNIVERSE, MATERIAL, CELL };
// Strategy used for rejecting sites when constraints are applied. KILL means
// that sites are always accepted but if they don't satisfy constraints, they
// are given weight 0. RESAMPLE means that a new source site will be sampled
// until constraints are met.
enum class RejectionStrategy { KILL, RESAMPLE };
// Indicates whether derived class already handles constraints
virtual bool constraints_applied() const { return false; }
// Methods for constraints
void read_constraints(pugi::xml_node node);
bool satisfies_spatial_constraints(Position r) const;
bool satisfies_energy_constraints(double E) const;
bool satisfies_time_constraints(double time) const;
// Data members
double strength_ {1.0}; //!< Source strength
std::unordered_set<int32_t> domain_ids_; //!< Domains to reject from
DomainType domain_type_; //!< Domain type for rejection
std::pair<double, double> time_bounds_ {-std::numeric_limits<double>::max(),
std::numeric_limits<double>::max()}; //!< time limits
std::pair<double, double> energy_bounds_ {
0, std::numeric_limits<double>::max()}; //!< energy limits
bool only_fissionable_ {
false}; //!< Whether site must be in fissionable material
RejectionStrategy rejection_strategy_ {
RejectionStrategy::RESAMPLE}; //!< Procedure for rejecting
};
//==============================================================================
@ -73,7 +127,6 @@ public:
// Properties
ParticleType particle_type() const { return particle_; }
double strength() const override { return strength_; }
// Make observing pointers available
SpatialDistribution* space() const { return space_.get(); }
@ -81,19 +134,17 @@ public:
Distribution* energy() const { return energy_.get(); }
Distribution* time() const { return time_.get(); }
private:
// Domain types
enum class DomainType { UNIVERSE, MATERIAL, CELL };
protected:
// Indicates whether derived class already handles constraints
bool constraints_applied() const override { return true; }
private:
// Data members
ParticleType particle_ {ParticleType::neutron}; //!< Type of particle emitted
double strength_ {1.0}; //!< Source strength
UPtrSpace space_; //!< Spatial distribution
UPtrAngle angle_; //!< Angular distribution
UPtrDist energy_; //!< Energy distribution
UPtrDist time_; //!< Time distribution
DomainType domain_type_; //!< Domain type for rejection
std::unordered_set<int32_t> domain_ids_; //!< Domains to reject from
};
//==============================================================================
@ -107,9 +158,12 @@ public:
explicit FileSource(const std::string& path);
// Methods
SourceSite sample(uint64_t* seed) const override;
void load_sites_from_file(
const std::string& path); //!< Load source sites from file
protected:
SourceSite sample(uint64_t* seed) const override;
private:
vector<SourceSite> sites_; //!< Source sites from a file
};
@ -124,16 +178,17 @@ public:
CompiledSourceWrapper(pugi::xml_node node);
~CompiledSourceWrapper();
double strength() const override { return compiled_source_->strength(); }
void setup(const std::string& path, const std::string& parameters);
protected:
// Defer implementation to custom source library
SourceSite sample(uint64_t* seed) const override
{
return compiled_source_->sample(seed);
}
double strength() const override { return compiled_source_->strength(); }
void setup(const std::string& path, const std::string& parameters);
private:
void* shared_library_; //!< library from dlopen
unique_ptr<Source> compiled_source_;
@ -164,6 +219,9 @@ public:
return sources_.size() == 1 ? sources_[0] : sources_[i];
}
protected:
bool constraints_applied() const override { return true; }
private:
// Data members
unique_ptr<MeshSpatial> space_; //!< Mesh spatial

12
openmc/examples.py
View file

@ -76,8 +76,10 @@ def pwr_pin_cell():
model.settings.batches = 10
model.settings.inactive = 5
model.settings.particles = 100
model.settings.source = openmc.IndependentSource(space=openmc.stats.Box(
[-pitch/2, -pitch/2, -1], [pitch/2, pitch/2, 1], only_fissionable=True))
model.settings.source = openmc.IndependentSource(
space=openmc.stats.Box([-pitch/2, -pitch/2, -1], [pitch/2, pitch/2, 1]),
constraints={'fissionable': True}
)
plot = openmc.Plot.from_geometry(model.geometry)
plot.pixels = (300, 300)
@ -527,8 +529,10 @@ def pwr_assembly():
model.settings.batches = 10
model.settings.inactive = 5
model.settings.particles = 100
model.settings.source = openmc.IndependentSource(space=openmc.stats.Box(
[-pitch/2, -pitch/2, -1], [pitch/2, pitch/2, 1], only_fissionable=True))
model.settings.source = openmc.IndependentSource(
space=openmc.stats.Box([-pitch/2, -pitch/2, -1], [pitch/2, pitch/2, 1]),
constraints={'fissionable': True}
)
plot = openmc.Plot()
plot.origin = (0.0, 0.0, 0)

328
openmc/source.py
View file

@ -6,7 +6,7 @@ from numbers import Real
import warnings
import typing # imported separately as py3.8 requires typing.Iterable
# also required to prevent typing.Union namespace overwriting Union
from typing import Optional, Sequence
from typing import Optional, Sequence, Dict, Any
import lxml.etree as ET
import numpy as np
@ -28,6 +28,19 @@ class SourceBase(ABC):
----------
strength : float
Strength of the source
constraints : dict
Constraints on sampled source particles. Valid keys include 'domains',
'time_bounds', 'energy_bounds', 'fissionable', and 'rejection_strategy'.
For 'domains', the corresponding value is an iterable of
:class:`openmc.Cell`, :class:`openmc.Material`, or
:class:`openmc.Universe` for which sampled sites must be within. For
'time_bounds' and 'energy_bounds', the corresponding value is a sequence
of floats giving the lower and upper bounds on time in [s] or energy in
[eV] that the sampled particle must be within. For 'fissionable', the
value is a bool indicating that only sites in fissionable material
should be accepted. The 'rejection_strategy' indicates what should
happen when a source particle is rejected: either 'resample' (pick a new
particle) or 'kill' (accept and terminate).
Attributes
----------
@ -35,11 +48,20 @@ class SourceBase(ABC):
Indicator of source type.
strength : float
Strength of the source
constraints : dict
Constraints on sampled source particles. Valid keys include
'domain_type', 'domain_ids', 'time_bounds', 'energy_bounds',
'fissionable', and 'rejection_strategy'.
"""
def __init__(self, strength=1.0):
def __init__(
self,
strength: Optional[float] = 1.0,
constraints: Optional[Dict[str, Any]] = None
):
self.strength = strength
self.constraints = constraints
@property
def strength(self):
@ -47,10 +69,47 @@ class SourceBase(ABC):
@strength.setter
def strength(self, strength):
cv.check_type('source strength', strength, Real)
cv.check_greater_than('source strength', strength, 0.0, True)
cv.check_type('source strength', strength, Real, none_ok=True)
if strength is not None:
cv.check_greater_than('source strength', strength, 0.0, True)
self._strength = strength
@property
def constraints(self) -> Dict[str, Any]:
return self._constraints
@constraints.setter
def constraints(self, constraints: Optional[Dict[str, Any]]):
self._constraints = {}
if constraints is None:
return
for key, value in constraints.items():
if key == 'domains':
cv.check_type('domains', value, Iterable,
(openmc.Cell, openmc.Material, openmc.Universe))
if isinstance(value[0], openmc.Cell):
self._constraints['domain_type'] = 'cell'
elif isinstance(value[0], openmc.Material):
self._constraints['domain_type'] = 'material'
elif isinstance(value[0], openmc.Universe):
self._constraints['domain_type'] = 'universe'
self._constraints['domain_ids'] = [d.id for d in value]
elif key == 'time_bounds':
cv.check_type('time bounds', value, Iterable, Real)
self._constraints['time_bounds'] = tuple(value)
elif key == 'energy_bounds':
cv.check_type('energy bounds', value, Iterable, Real)
self._constraints['energy_bounds'] = tuple(value)
elif key == 'fissionable':
cv.check_type('fissionable', value, bool)
self._constraints['fissionable'] = value
elif key == 'rejection_strategy':
cv.check_value('rejection strategy', value, ('resample', 'kill'))
self._constraints['rejection_strategy'] = value
else:
raise ValueError('Unknown key in constraints dictionary: {key}')
@abstractmethod
def populate_xml_element(self, element):
"""Add necessary source information to an XML element
@ -73,8 +132,30 @@ class SourceBase(ABC):
"""
element = ET.Element("source")
element.set("type", self.type)
element.set("strength", str(self.strength))
if self.strength is not None:
element.set("strength", str(self.strength))
self.populate_xml_element(element)
constraints = self.constraints
if constraints:
constraints_elem = ET.SubElement(element, "constraints")
if "domain_ids" in constraints:
dt_elem = ET.SubElement(constraints_elem, "domain_type")
dt_elem.text = constraints["domain_type"]
id_elem = ET.SubElement(constraints_elem, "domain_ids")
id_elem.text = ' '.join(str(uid) for uid in constraints["domain_ids"])
if "time_bounds" in constraints:
dt_elem = ET.SubElement(constraints_elem, "time_bounds")
dt_elem.text = ' '.join(str(t) for t in constraints["time_bounds"])
if "energy_bounds" in constraints:
dt_elem = ET.SubElement(constraints_elem, "energy_bounds")
dt_elem.text = ' '.join(str(E) for E in constraints["energy_bounds"])
if "fissionable" in constraints:
dt_elem = ET.SubElement(constraints_elem, "fissionable")
dt_elem.text = str(constraints["fissionable"]).lower()
if "rejection_strategy" in constraints:
dt_elem = ET.SubElement(constraints_elem, "rejection_strategy")
dt_elem.text = constraints["rejection_strategy"]
return element
@classmethod
@ -118,6 +199,47 @@ class SourceBase(ABC):
else:
raise ValueError(f'Source type {source_type} is not recognized')
@staticmethod
def _get_constraints(elem: ET.Element) -> Dict[str, Any]:
# Find element containing constraints
constraints_elem = elem.find("constraints")
elem = constraints_elem if constraints_elem is not None else elem
constraints = {}
domain_type = get_text(elem, "domain_type")
if domain_type is not None:
domain_ids = [int(x) for x in get_text(elem, "domain_ids").split()]
# Instantiate some throw-away domains that are used by the
# constructor to assign IDs
with warnings.catch_warnings():
warnings.simplefilter('ignore', openmc.IDWarning)
if domain_type == 'cell':
domains = [openmc.Cell(uid) for uid in domain_ids]
elif domain_type == 'material':
domains = [openmc.Material(uid) for uid in domain_ids]
elif domain_type == 'universe':
domains = [openmc.Universe(uid) for uid in domain_ids]
constraints['domains'] = domains
time_bounds = get_text(elem, "time_bounds")
if time_bounds is not None:
constraints['time_bounds'] = [float(x) for x in time_bounds.split()]
energy_bounds = get_text(elem, "energy_bounds")
if energy_bounds is not None:
constraints['energy_bounds'] = [float(x) for x in energy_bounds.split()]
fissionable = get_text(elem, "fissionable")
if fissionable is not None:
constraints['fissionable'] = fissionable in ('true', '1')
rejection_strategy = get_text(elem, "rejection_strategy")
if rejection_strategy is not None:
constraints['rejection_strategy'] = rejection_strategy
return constraints
class IndependentSource(SourceBase):
"""Distribution of phase space coordinates for source sites.
@ -142,6 +264,22 @@ class IndependentSource(SourceBase):
Domains to reject based on, i.e., if a sampled spatial location is not
within one of these domains, it will be rejected.
.. deprecated:: 0.14.1
Use the `constraints` argument instead.
constraints : dict
Constraints on sampled source particles. Valid keys include 'domains',
'time_bounds', 'energy_bounds', 'fissionable', and 'rejection_strategy'.
For 'domains', the corresponding value is an iterable of
:class:`openmc.Cell`, :class:`openmc.Material`, or
:class:`openmc.Universe` for which sampled sites must be within. For
'time_bounds' and 'energy_bounds', the corresponding value is a sequence
of floats giving the lower and upper bounds on time in [s] or energy in
[eV] that the sampled particle must be within. For 'fissionable', the
value is a bool indicating that only sites in fissionable material
should be accepted. The 'rejection_strategy' indicates what should
happen when a source particle is rejected: either 'resample' (pick a new
particle) or 'kill' (accept and terminate).
Attributes
----------
space : openmc.stats.Spatial or None
@ -161,10 +299,10 @@ class IndependentSource(SourceBase):
particle : {'neutron', 'photon'}
Source particle type
ids : Iterable of int
IDs of domains to use for rejection
domain_type : {'cell', 'material', 'universe'}
Type of domain to use for rejection
constraints : dict
Constraints on sampled source particles. Valid keys include
'domain_type', 'domain_ids', 'time_bounds', 'energy_bounds',
'fissionable', and 'rejection_strategy'.
"""
@ -176,9 +314,15 @@ class IndependentSource(SourceBase):
time: Optional[openmc.stats.Univariate] = None,
strength: float = 1.0,
particle: str = 'neutron',
domains: Optional[Sequence[typing.Union[openmc.Cell, openmc.Material, openmc.Universe]]] = None
domains: Optional[Sequence[typing.Union[openmc.Cell, openmc.Material, openmc.Universe]]] = None,
constraints: Optional[Dict[str, Any]] = None
):
super().__init__(strength)
if domains is not None:
warnings.warn("The 'domains' arguments has been replaced by the "
"'constraints' argument.", FutureWarning)
constraints = {'domains': domains}
super().__init__(strength=strength, constraints=constraints)
self._space = None
self._angle = None
@ -193,20 +337,8 @@ class IndependentSource(SourceBase):
self.energy = energy
if time is not None:
self.time = time
self.strength = strength
self.particle = particle
self._domain_ids = []
self._domain_type = None
if domains is not None:
if isinstance(domains[0], openmc.Cell):
self.domain_type = 'cell'
elif isinstance(domains[0], openmc.Material):
self.domain_type = 'material'
elif isinstance(domains[0], openmc.Universe):
self.domain_type = 'universe'
self.domain_ids = [d.id for d in domains]
@property
def type(self) -> str:
return 'independent'
@ -273,24 +405,6 @@ class IndependentSource(SourceBase):
cv.check_value('source particle', particle, ['neutron', 'photon'])
self._particle = particle
@property
def domain_ids(self):
return self._domain_ids
@domain_ids.setter
def domain_ids(self, ids):
cv.check_type('domain IDs', ids, Iterable, Real)
self._domain_ids = ids
@property
def domain_type(self):
return self._domain_type
@domain_type.setter
def domain_type(self, domain_type):
cv.check_value('domain type', domain_type, ('cell', 'material', 'universe'))
self._domain_type = domain_type
def populate_xml_element(self, element):
"""Add necessary source information to an XML element
@ -309,11 +423,6 @@ class IndependentSource(SourceBase):
element.append(self.energy.to_xml_element('energy'))
if self.time is not None:
element.append(self.time.to_xml_element('time'))
if self.domain_ids:
dt_elem = ET.SubElement(element, "domain_type")
dt_elem.text = self.domain_type
id_elem = ET.SubElement(element, "domain_ids")
id_elem.text = ' '.join(str(uid) for uid in self.domain_ids)
@classmethod
def from_xml_element(cls, elem: ET.Element, meshes=None) -> SourceBase:
@ -333,24 +442,8 @@ class IndependentSource(SourceBase):
Source generated from XML element
"""
domain_type = get_text(elem, "domain_type")
if domain_type is not None:
domain_ids = [int(x) for x in get_text(elem, "domain_ids").split()]
# Instantiate some throw-away domains that are used by the
# constructor to assign IDs
with warnings.catch_warnings():
warnings.simplefilter('ignore', openmc.IDWarning)
if domain_type == 'cell':
domains = [openmc.Cell(uid) for uid in domain_ids]
elif domain_type == 'material':
domains = [openmc.Material(uid) for uid in domain_ids]
elif domain_type == 'universe':
domains = [openmc.Universe(uid) for uid in domain_ids]
else:
domains = None
source = cls(domains=domains)
constraints = cls._get_constraints(elem)
source = cls(constraints=constraints)
strength = get_text(elem, 'strength')
if strength is not None:
@ -360,10 +453,6 @@ class IndependentSource(SourceBase):
if particle is not None:
source.particle = particle
filename = get_text(elem, 'file')
if filename is not None:
source.file = filename
space = elem.find('space')
if space is not None:
source.space = Spatial.from_xml_element(space, meshes)
@ -405,6 +494,19 @@ class MeshSource(SourceBase):
multidimensional array whose shape matches the mesh shape. If spatial
distributions are set on any of the source objects, they will be ignored
during source site sampling.
constraints : dict
Constraints on sampled source particles. Valid keys include 'domains',
'time_bounds', 'energy_bounds', 'fissionable', and 'rejection_strategy'.
For 'domains', the corresponding value is an iterable of
:class:`openmc.Cell`, :class:`openmc.Material`, or
:class:`openmc.Universe` for which sampled sites must be within. For
'time_bounds' and 'energy_bounds', the corresponding value is a sequence
of floats giving the lower and upper bounds on time in [s] or energy in
[eV] that the sampled particle must be within. For 'fissionable', the
value is a bool indicating that only sites in fissionable material
should be accepted. The 'rejection_strategy' indicates what should
happen when a source particle is rejected: either 'resample' (pick a new
particle) or 'kill' (accept and terminate).
Attributes
----------
@ -416,9 +518,19 @@ class MeshSource(SourceBase):
Strength of the source
type : str
Indicator of source type: 'mesh'
constraints : dict
Constraints on sampled source particles. Valid keys include
'domain_type', 'domain_ids', 'time_bounds', 'energy_bounds',
'fissionable', and 'rejection_strategy'.
"""
def __init__(self, mesh: MeshBase, sources: Sequence[SourceBase]):
def __init__(
self,
mesh: MeshBase,
sources: Sequence[SourceBase],
constraints: Optional[Dict[str, Any]] = None,
):
super().__init__(strength=None, constraints=constraints)
self.mesh = mesh
self.sources = sources
@ -473,8 +585,9 @@ class MeshSource(SourceBase):
@strength.setter
def strength(self, val):
cv.check_type('mesh source strength', val, Real)
self.set_total_strength(val)
if val is not None:
cv.check_type('mesh source strength', val, Real)
self.set_total_strength(val)
def set_total_strength(self, strength: float):
"""Scales the element source strengths based on a desired total strength.
@ -531,7 +644,8 @@ class MeshSource(SourceBase):
mesh = meshes[mesh_id]
sources = [SourceBase.from_xml_element(e) for e in elem.iterchildren('source')]
return cls(mesh, sources)
constraints = cls._get_constraints(elem)
return cls(mesh, sources, constraints=constraints)
def Source(*args, **kwargs):
@ -556,6 +670,19 @@ class CompiledSource(SourceBase):
Parameters to be provided to the compiled shared library function
strength : float
Strength of the source
constraints : dict
Constraints on sampled source particles. Valid keys include 'domains',
'time_bounds', 'energy_bounds', 'fissionable', and 'rejection_strategy'.
For 'domains', the corresponding value is an iterable of
:class:`openmc.Cell`, :class:`openmc.Material`, or
:class:`openmc.Universe` for which sampled sites must be within. For
'time_bounds' and 'energy_bounds', the corresponding value is a sequence
of floats giving the lower and upper bounds on time in [s] or energy in
[eV] that the sampled particle must be within. For 'fissionable', the
value is a bool indicating that only sites in fissionable material
should be accepted. The 'rejection_strategy' indicates what should
happen when a source particle is rejected: either 'resample' (pick a new
particle) or 'kill' (accept and terminate).
Attributes
----------
@ -567,10 +694,20 @@ class CompiledSource(SourceBase):
Strength of the source
type : str
Indicator of source type: 'compiled'
constraints : dict
Constraints on sampled source particles. Valid keys include
'domain_type', 'domain_ids', 'time_bounds', 'energy_bounds',
'fissionable', and 'rejection_strategy'.
"""
def __init__(self, library: Optional[str] = None, parameters: Optional[str] = None, strength=1.0) -> None:
super().__init__(strength=strength)
def __init__(
self,
library: Optional[str] = None,
parameters: Optional[str] = None,
strength: float = 1.0,
constraints: Optional[Dict[str, Any]] = None
) -> None:
super().__init__(strength=strength, constraints=constraints)
self._library = None
if library is not None:
@ -634,9 +771,10 @@ class CompiledSource(SourceBase):
Source generated from XML element
"""
library = get_text(elem, 'library')
kwargs = {'constraints': cls._get_constraints(elem)}
kwargs['library'] = get_text(elem, 'library')
source = cls(library)
source = cls(**kwargs)
strength = get_text(elem, 'strength')
if strength is not None:
@ -660,6 +798,19 @@ class FileSource(SourceBase):
Path to the source file from which sites should be sampled
strength : float
Strength of the source (default is 1.0)
constraints : dict
Constraints on sampled source particles. Valid keys include 'domains',
'time_bounds', 'energy_bounds', 'fissionable', and 'rejection_strategy'.
For 'domains', the corresponding value is an iterable of
:class:`openmc.Cell`, :class:`openmc.Material`, or
:class:`openmc.Universe` for which sampled sites must be within. For
'time_bounds' and 'energy_bounds', the corresponding value is a sequence
of floats giving the lower and upper bounds on time in [s] or energy in
[eV] that the sampled particle must be within. For 'fissionable', the
value is a bool indicating that only sites in fissionable material
should be accepted. The 'rejection_strategy' indicates what should
happen when a source particle is rejected: either 'resample' (pick a new
particle) or 'kill' (accept and terminate).
Attributes
----------
@ -669,13 +820,21 @@ class FileSource(SourceBase):
Strength of the source
type : str
Indicator of source type: 'file'
constraints : dict
Constraints on sampled source particles. Valid keys include
'domain_type', 'domain_ids', 'time_bounds', 'energy_bounds',
'fissionable', and 'rejection_strategy'.
"""
def __init__(self, path: Optional[PathLike] = None, strength=1.0) -> None:
super().__init__(strength=strength)
def __init__(
self,
path: Optional[PathLike] = None,
strength: float = 1.0,
constraints: Optional[Dict[str, Any]] = None
):
super().__init__(strength=strength, constraints=constraints)
self._path = None
if path is not None:
self.path = path
@ -722,16 +881,13 @@ class FileSource(SourceBase):
Source generated from XML element
"""
filename = get_text(elem, 'file')
source = cls(filename)
kwargs = {'constraints': cls._get_constraints(elem)}
kwargs['path'] = get_text(elem, 'file')
strength = get_text(elem, 'strength')
if strength is not None:
source.strength = float(strength)
kwargs['strength'] = float(strength)
return source
return cls(**kwargs)
class ParticleType(IntEnum):

11
openmc/stats/multivariate.py
View file

@ -4,6 +4,7 @@ from abc import ABC, abstractmethod
from collections.abc import Iterable
from math import cos, pi
from numbers import Real
from warnings import warn
import lxml.etree as ET
import numpy as np
@ -768,6 +769,9 @@ class Box(Spatial):
Whether spatial sites should only be accepted if they occur in
fissionable materials
.. deprecated:: 0.14.1
Use the `constraints` argument when defining a source object instead.
Attributes
----------
lower_left : Iterable of float
@ -778,6 +782,9 @@ class Box(Spatial):
Whether spatial sites should only be accepted if they occur in
fissionable materials
.. deprecated:: 0.14.1
Use the `constraints` argument when defining a source object instead.
"""
def __init__(
@ -818,6 +825,10 @@ class Box(Spatial):
def only_fissionable(self, only_fissionable):
cv.check_type('only fissionable', only_fissionable, bool)
self._only_fissionable = only_fissionable
if only_fissionable:
warn("The 'only_fissionable' has been deprecated. Use the "
"'constraints' argument when defining a source instead.",
FutureWarning)
def to_xml_element(self):
"""Return XML representation of the box distribution

7
src/distribution_spatial.cpp
View file

@ -281,10 +281,15 @@ void MeshSpatial::check_element_types() const
}
}
int32_t MeshSpatial::sample_element_index(uint64_t* seed) const
{
return elem_idx_dist_.sample(seed);
}
std::pair<int32_t, Position> MeshSpatial::sample_mesh(uint64_t* seed) const
{
// Sample the CDF defined in initialization above
int32_t elem_idx = elem_idx_dist_.sample(seed);
int32_t elem_idx = this->sample_element_index(seed);
return {elem_idx, mesh()->sample_element(elem_idx, seed)};
}

4
src/simulation.cpp
View file

@ -750,7 +750,7 @@ void free_memory_simulation()
void transport_history_based_single_particle(Particle& p)
{
while (true) {
while (p.alive()) {
p.event_calculate_xs();
if (!p.alive())
break;
@ -761,8 +761,6 @@ void transport_history_based_single_particle(Particle& p)
p.event_collide();
}
p.event_revive_from_secondary();
if (!p.alive())
break;
}
p.event_death();
}

304
src/source.cpp
View file

@ -47,9 +47,23 @@ vector<unique_ptr<Source>> external_sources;
}
//==============================================================================
// Source create implementation
// Source implementation
//==============================================================================
Source::Source(pugi::xml_node node)
{
// Check for source strength
if (check_for_node(node, "strength")) {
strength_ = std::stod(get_node_value(node, "strength"));
if (strength_ < 0.0) {
fatal_error("Source strength is negative.");
}
}
// Check for additional defined constraints
read_constraints(node);
}
unique_ptr<Source> Source::create(pugi::xml_node node)
{
// if the source type is present, use it to determine the type
@ -79,6 +93,163 @@ unique_ptr<Source> Source::create(pugi::xml_node node)
}
}
void Source::read_constraints(pugi::xml_node node)
{
// Check for constraints node. For backwards compatibility, if no constraints
// node is given, still try searching for domain constraints from top-level
// node.
pugi::xml_node constraints_node = node.child("constraints");
if (constraints_node) {
node = constraints_node;
}
// Check for domains to reject from
if (check_for_node(node, "domain_type")) {
std::string domain_type = get_node_value(node, "domain_type");
if (domain_type == "cell") {
domain_type_ = DomainType::CELL;
} else if (domain_type == "material") {
domain_type_ = DomainType::MATERIAL;
} else if (domain_type == "universe") {
domain_type_ = DomainType::UNIVERSE;
} else {
fatal_error(
std::string("Unrecognized domain type for constraint: " + domain_type));
}
auto ids = get_node_array<int>(node, "domain_ids");
domain_ids_.insert(ids.begin(), ids.end());
}
if (check_for_node(node, "time_bounds")) {
auto ids = get_node_array<double>(node, "time_bounds");
if (ids.size() != 2) {
fatal_error("Time bounds must be represented by two numbers.");
}
time_bounds_ = std::make_pair(ids[0], ids[1]);
}
if (check_for_node(node, "energy_bounds")) {
auto ids = get_node_array<double>(node, "energy_bounds");
if (ids.size() != 2) {
fatal_error("Energy bounds must be represented by two numbers.");
}
energy_bounds_ = std::make_pair(ids[0], ids[1]);
}
if (check_for_node(node, "fissionable")) {
only_fissionable_ = get_node_value_bool(node, "fissionable");
}
// Check for how to handle rejected particles
if (check_for_node(node, "rejection_strategy")) {
std::string rejection_strategy = get_node_value(node, "rejection_strategy");
if (rejection_strategy == "kill") {
rejection_strategy_ = RejectionStrategy::KILL;
} else if (rejection_strategy == "resample") {
rejection_strategy_ = RejectionStrategy::RESAMPLE;
} else {
fatal_error(std::string(
"Unrecognized strategy source rejection: " + rejection_strategy));
}
}
}
SourceSite Source::sample_with_constraints(uint64_t* seed) const
{
bool accepted = false;
static int n_reject = 0;
static int n_accept = 0;
SourceSite site;
while (!accepted) {
// Sample a source site without considering constraints yet
site = this->sample(seed);
if (constraints_applied()) {
accepted = true;
} else {
// Check whether sampled site satisfies constraints
accepted = satisfies_spatial_constraints(site.r) &&
satisfies_energy_constraints(site.E) &&
satisfies_time_constraints(site.time);
if (!accepted) {
++n_reject;
if (n_reject >= EXTSRC_REJECT_THRESHOLD &&
static_cast<double>(n_accept) / n_reject <=
EXTSRC_REJECT_FRACTION) {
fatal_error("More than 95% of external source sites sampled were "
"rejected. Please check your source definition.");
}
// For the "kill" strategy, accept particle but set weight to 0 so that
// it is terminated immediately
if (rejection_strategy_ == RejectionStrategy::KILL) {
accepted = true;
site.wgt = 0.0;
}
}
}
}
// Increment number of accepted samples
++n_accept;
return site;
}
bool Source::satisfies_energy_constraints(double E) const
{
return E > energy_bounds_.first && E < energy_bounds_.second;
}
bool Source::satisfies_time_constraints(double time) const
{
return time > time_bounds_.first && time < time_bounds_.second;
}
bool Source::satisfies_spatial_constraints(Position r) const
{
GeometryState geom_state;
geom_state.r() = r;
// Reject particle if it's not in the geometry at all
bool found = exhaustive_find_cell(geom_state);
if (!found)
return false;
// Check the geometry state against specified domains
bool accepted = true;
if (!domain_ids_.empty()) {
if (domain_type_ == DomainType::MATERIAL) {
auto mat_index = geom_state.material();
if (mat_index != MATERIAL_VOID) {
accepted = contains(domain_ids_, model::materials[mat_index]->id());
}
} else {
for (int i = 0; i < geom_state.n_coord(); i++) {
auto id = (domain_type_ == DomainType::CELL)
? model::cells[geom_state.coord(i).cell]->id_
: model::universes[geom_state.coord(i).universe]->id_;
if ((accepted = contains(domain_ids_, id)))
break;
}
}
}
// Check if spatial site is in fissionable material
if (accepted && only_fissionable_) {
// Determine material
auto mat_index = geom_state.material();
if (mat_index == MATERIAL_VOID) {
accepted = false;
} else {
accepted = model::materials[mat_index]->fissionable();
}
}
return accepted;
}
//==============================================================================
// IndependentSource implementation
//==============================================================================
@ -89,7 +260,7 @@ IndependentSource::IndependentSource(
energy_ {std::move(energy)}, time_ {std::move(time)}
{}
IndependentSource::IndependentSource(pugi::xml_node node)
IndependentSource::IndependentSource(pugi::xml_node node) : Source(node)
{
// Check for particle type
if (check_for_node(node, "particle")) {
@ -104,11 +275,6 @@ IndependentSource::IndependentSource(pugi::xml_node node)
}
}
// Check for source strength
if (check_for_node(node, "strength")) {
strength_ = std::stod(get_node_value(node, "strength"));
}
// Check for external source file
if (check_for_node(node, "file")) {
@ -122,6 +288,15 @@ IndependentSource::IndependentSource(pugi::xml_node node)
space_ = UPtrSpace {new SpatialPoint()};
}
// For backwards compatibility, check for only fissionable setting on box
// source
auto space_box = dynamic_cast<SpatialBox*>(space_.get());
if (space_box) {
if (!only_fissionable_) {
only_fissionable_ = space_box->only_fissionable();
}
}
// Determine external source angular distribution
if (check_for_node(node, "angle")) {
angle_ = UnitSphereDistribution::create(node.child("angle"));
@ -148,24 +323,6 @@ IndependentSource::IndependentSource(pugi::xml_node node)
double p[] {1.0};
time_ = UPtrDist {new Discrete {T, p, 1}};
}
// Check for domains to reject from
if (check_for_node(node, "domain_type")) {
std::string domain_type = get_node_value(node, "domain_type");
if (domain_type == "cell") {
domain_type_ = DomainType::CELL;
} else if (domain_type == "material") {
domain_type_ = DomainType::MATERIAL;
} else if (domain_type == "universe") {
domain_type_ = DomainType::UNIVERSE;
} else {
fatal_error(std::string(
"Unrecognized domain type for source rejection: " + domain_type));
}
auto ids = get_node_array<int>(node, "domain_ids");
domain_ids_.insert(ids.begin(), ids.end());
}
}
}
@ -174,60 +331,21 @@ SourceSite IndependentSource::sample(uint64_t* seed) const
SourceSite site;
site.particle = particle_;
// Repeat sampling source location until a good site has been found
bool found = false;
int n_reject = 0;
// Repeat sampling source location until a good site has been accepted
bool accepted = false;
static int n_reject = 0;
static int n_accept = 0;
while (!found) {
// Set particle type
Particle p;
p.type() = particle_;
p.u() = {0.0, 0.0, 1.0};
while (!accepted) {
// Sample spatial distribution
p.r() = space_->sample(seed);
site.r = space_->sample(seed);
// Now search to see if location exists in geometry
found = exhaustive_find_cell(p);
// Check if spatial site is in fissionable material
if (found) {
auto space_box = dynamic_cast<SpatialBox*>(space_.get());
if (space_box) {
if (space_box->only_fissionable()) {
// Determine material
auto mat_index = p.material();
if (mat_index == MATERIAL_VOID) {
found = false;
} else {
found = model::materials[mat_index]->fissionable();
}
}
}
// Rejection based on cells/materials/universes
if (!domain_ids_.empty()) {
found = false;
if (domain_type_ == DomainType::MATERIAL) {
auto mat_index = p.material();
if (mat_index != MATERIAL_VOID) {
found = contains(domain_ids_, model::materials[mat_index]->id());
}
} else {
for (int i = 0; i < p.n_coord(); i++) {
auto id = (domain_type_ == DomainType::CELL)
? model::cells[p.coord(i).cell]->id_
: model::universes[p.coord(i).universe]->id_;
if ((found = contains(domain_ids_, id)))
break;
}
}
}
}
// Check if sampled position satisfies spatial constraints
accepted = satisfies_spatial_constraints(site.r);
// Check for rejection
if (!found) {
if (!accepted) {
++n_reject;
if (n_reject >= EXTSRC_REJECT_THRESHOLD &&
static_cast<double>(n_accept) / n_reject <= EXTSRC_REJECT_FRACTION) {
@ -236,8 +354,6 @@ SourceSite IndependentSource::sample(uint64_t* seed) const
"definition.");
}
}
site.r = p.r();
}
// Sample angle
@ -261,7 +377,8 @@ SourceSite IndependentSource::sample(uint64_t* seed) const
site.E = energy_->sample(seed);
// Resample if energy falls above maximum particle energy
if (site.E < data::energy_max[p])
if (site.E < data::energy_max[p] and
(satisfies_energy_constraints(site.E)))
break;
n_reject++;
@ -287,7 +404,8 @@ SourceSite IndependentSource::sample(uint64_t* seed) const
//==============================================================================
// FileSource implementation
//==============================================================================
FileSource::FileSource(pugi::xml_node node)
FileSource::FileSource(pugi::xml_node node) : Source(node)
{
auto path = get_node_value(node, "file", false, true);
if (ends_with(path, ".mcpl") || ends_with(path, ".mcpl.gz")) {
@ -332,6 +450,7 @@ void FileSource::load_sites_from_file(const std::string& path)
SourceSite FileSource::sample(uint64_t* seed) const
{
// Sample a particle randomly from list
size_t i_site = sites_.size() * prn(seed);
return sites_[i_site];
}
@ -339,7 +458,8 @@ SourceSite FileSource::sample(uint64_t* seed) const
//==============================================================================
// CompiledSourceWrapper implementation
//==============================================================================
CompiledSourceWrapper::CompiledSourceWrapper(pugi::xml_node node)
CompiledSourceWrapper::CompiledSourceWrapper(pugi::xml_node node) : Source(node)
{
// Get shared library path and parameters
auto path = get_node_value(node, "library", false, true);
@ -403,7 +523,7 @@ CompiledSourceWrapper::~CompiledSourceWrapper()
// MeshSource implementation
//==============================================================================
MeshSource::MeshSource(pugi::xml_node node)
MeshSource::MeshSource(pugi::xml_node node) : Source(node)
{
int32_t mesh_id = stoi(get_node_value(node, "mesh"));
int32_t mesh_idx = model::mesh_map.at(mesh_id);
@ -430,15 +550,27 @@ MeshSource::MeshSource(pugi::xml_node node)
SourceSite MeshSource::sample(uint64_t* seed) const
{
// sample location and element from mesh
auto mesh_location = space_->sample_mesh(seed);
// Sample the CDF defined in initialization above
int32_t element = space_->sample_element_index(seed);
// Sample source for the chosen element
int32_t element = mesh_location.first;
SourceSite site = source(element)->sample(seed);
// Sample position and apply rejection on spatial domains
Position r;
do {
r = space_->mesh()->sample_element(element, seed);
} while (!this->satisfies_spatial_constraints(r));
// Replace spatial position with the one already sampled
site.r = mesh_location.second;
SourceSite site;
while (true) {
// Sample source for the chosen element and replace the position
site = source(element)->sample_with_constraints(seed);
site.r = r;
// Apply other rejections
if (satisfies_energy_constraints(site.E) &&
satisfies_time_constraints(site.time)) {
break;
}
}
return site;
}
@ -493,7 +625,7 @@ SourceSite sample_external_source(uint64_t* seed)
}
// Sample source site from i-th source distribution
SourceSite site {model::external_sources[i]->sample(seed)};
SourceSite site {model::external_sources[i]->sample_with_constraints(seed)};
// If running in MG, convert site.E to group
if (!settings::run_CE) {

5
tests/regression_tests/asymmetric_lattice/inputs_true.dat
View file

@ -209,9 +209,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-32 -32 0 32 32 32</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

6
tests/regression_tests/asymmetric_lattice/test.py
View file

@ -54,8 +54,10 @@ class AsymmetricLatticeTestHarness(PyAPITestHarness):
self._model.tallies.append(tally)
# Specify summary output and correct source sampling box
self._model.settings.source = openmc.IndependentSource(space=openmc.stats.Box(
[-32, -32, 0], [32, 32, 32], only_fissionable = True))
self._model.settings.source = openmc.IndependentSource(
space=openmc.stats.Box([-32, -32, 0], [32, 32, 32]),
constraints={'fissionable': True}
)
def _get_results(self, hash_output=True):
"""Digest info in statepoint and summary and return as a string."""

5
tests/regression_tests/mg_temperature_multi/inputs_true.dat
View file

@ -28,9 +28,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
<energy_mode>multi-group</energy_mode>
</settings>

5
tests/regression_tests/mg_temperature_multi/test.py
View file

@ -133,8 +133,9 @@ def test_mg_temperature_multi():
# Create an initial uniform spatial source distribution over fissionable zones
lower_left = (-pitch/2, -pitch/2, -1)
upper_right = (pitch/2, pitch/2, 1)
uniform_dist = openmc.stats.Box(lower_left, upper_right, only_fissionable=True)
settings.source = openmc.IndependentSource(space=uniform_dist)
uniform_dist = openmc.stats.Box(lower_left, upper_right)
settings.source = openmc.IndependentSource(
space=uniform_dist, constraints={'fissionable': True})
###############################################################################
# Define tallies

5
tests/regression_tests/mgxs_library_ce_to_mg/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_ce_to_mg_nuclides/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_condense/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_correction/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_distribcell/inputs_true.dat
View file

@ -69,9 +69,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-10.71 -10.71 -1 10.71 10.71 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_hdf5/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_histogram/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_no_nuclides/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_nuclides/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/mgxs_library_specific_nuclides/inputs_true.dat
View file

@ -42,9 +42,12 @@
<batches>10</batches>
<inactive>5</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.63 -0.63 -1 0.63 0.63 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

5
tests/regression_tests/surface_tally/inputs_true.dat
View file

@ -30,9 +30,12 @@
<batches>10</batches>
<inactive>0</inactive>
<source particle="neutron" strength="1.0" type="independent">
<space type="fission">
<space type="box">
<parameters>-0.62992 -0.62992 -1 0.62992 0.62992 1</parameters>
</space>
<constraints>
<fissionable>true</fissionable>
</constraints>
</source>
</settings>
<tallies>

6
tests/regression_tests/surface_tally/test.py
View file

@ -72,9 +72,9 @@ class SurfaceTallyTestHarness(PyAPITestHarness):
# Create an initial uniform spatial source distribution
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:],\
only_fissionable=True)
settings_file.source = openmc.IndependentSource(space=uniform_dist)
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:],)
settings_file.source = openmc.IndependentSource(
space=uniform_dist, constraints={'fissionable': True})
self._model.settings = settings_file
# Tallies file

6
tests/unit_tests/test_lib.py
View file

@ -67,8 +67,10 @@ def uo2_trigger_model():
model.settings.batches = 10
model.settings.inactive = 5
model.settings.particles = 100
model.settings.source = openmc.IndependentSource(space=openmc.stats.Box(
[-0.5, -0.5, -1], [0.5, 0.5, 1], only_fissionable=True))
model.settings.source = openmc.IndependentSource(
space=openmc.stats.Box([-0.5, -0.5, -1], [0.5, 0.5, 1]),
constraints={'fissionable': True},
)
model.settings.verbosity = 1
model.settings.keff_trigger = {'type': 'std_dev', 'threshold': 0.001}
model.settings.trigger_active = True

6
tests/unit_tests/test_model.py
View file

@ -57,9 +57,9 @@ def pin_model_attributes():
# Create a uniform spatial source distribution over fissionable zones
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
uniform_dist = openmc.stats.Box(
bounds[:3], bounds[3:], only_fissionable=True)
settings.source = openmc.IndependentSource(space=uniform_dist)
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
settings.source = openmc.IndependentSource(
space=uniform_dist, constraints={'fissionable': True})
entropy_mesh = openmc.RegularMesh()
entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50]

103
tests/unit_tests/test_source.py
View file

@ -7,6 +7,8 @@ import numpy as np
import pytest
from pytest import approx
from tests.regression_tests import config
def test_source():
space = openmc.stats.Point()
@ -100,7 +102,8 @@ def test_source_xml_roundtrip():
assert new_src.strength == approx(src.strength)
def test_rejection(run_in_tmpdir):
@pytest.fixture
def sphere_box_model():
# Model with two spheres inside a box
mat = openmc.Material()
mat.add_nuclide('H1', 1.0)
@ -119,24 +122,108 @@ def test_rejection(run_in_tmpdir):
model.settings.batches = 10
model.settings.run_mode = 'fixed source'
return model, cell1, cell2, cell3
def test_constraints_independent(sphere_box_model, run_in_tmpdir):
model, cell1, cell2, cell3 = sphere_box_model
# Set up a box source with rejection on the spherical cell
space = openmc.stats.Box(*cell3.bounding_box)
model.settings.source = openmc.IndependentSource(space=space, domains=[cell1, cell2])
space = openmc.stats.Box((-4., -1., -1.), (4., 1., 1.))
model.settings.source = openmc.IndependentSource(
space=space, constraints={'domains': [cell1, cell2]}
)
# Load up model via openmc.lib and sample source
model.export_to_xml()
model.export_to_model_xml()
openmc.lib.init()
particles = openmc.lib.sample_external_source(1000)
# Make sure that all sampled sources are within one of the spheres
joint_region = cell1.region | cell2.region
for p in particles:
assert p.r in joint_region
assert p.r not in non_source_region
assert p.r in (cell1.region | cell2.region)
assert p.r not in cell3.region
openmc.lib.finalize()
def test_constraints_mesh(sphere_box_model, run_in_tmpdir):
model, cell1, cell2, cell3 = sphere_box_model
bbox = cell3.bounding_box
mesh = openmc.RegularMesh()
mesh.lower_left = bbox.lower_left
mesh.upper_right = bbox.upper_right
mesh.dimension = (2, 1, 1)
left_source = openmc.IndependentSource()
right_source = openmc.IndependentSource()
model.settings.source = openmc.MeshSource(
mesh, [left_source, right_source], constraints={'domains': [cell1, cell2]}
)
# Load up model via openmc.lib and sample source
model.export_to_model_xml()
openmc.lib.init()
particles = openmc.lib.sample_external_source(1000)
# Make sure that all sampled sources are within one of the spheres
for p in particles:
assert p.r in (cell1.region | cell2.region)
assert p.r not in cell3.region
openmc.lib.finalize()
def test_constraints_file(sphere_box_model, run_in_tmpdir):
model = sphere_box_model[0]
# Create source file with randomly sampled source sites
rng = np.random.default_rng()
energy = rng.uniform(0., 1e6, 10_000)
time = rng.uniform(0., 1., 10_000)
particles = [openmc.SourceParticle(E=e, time=t) for e, t in zip(energy, time)]
openmc.write_source_file(particles, 'uniform_source.h5')
# Use source file
model.settings.source = openmc.FileSource(
'uniform_source.h5',
constraints={
'time_bounds': [0.25, 0.75],
'energy_bounds': [500.e3, 1.0e6],
}
)
# Load up model via openmc.lib and sample source
model.export_to_model_xml()
openmc.lib.init()
particles = openmc.lib.sample_external_source(1000)
# Make sure that all sampled sources are within energy/time bounds
for p in particles:
assert 0.25 <= p.time <= 0.75
assert 500.e3 <= p.E <= 1.0e6
openmc.lib.finalize()
@pytest.mark.skipif(config['mpi'], reason='Not compatible with MPI')
def test_rejection_limit(sphere_box_model, run_in_tmpdir):
model, cell1 = sphere_box_model[:2]
# Define a point source that will get rejected 100% of the time
model.settings.source = openmc.IndependentSource(
space=openmc.stats.Point((-3., 0., 0.)),
constraints={'domains': [cell1]}
)
# Confirm that OpenMC doesn't run in an infinite loop. Note that this may
# work when running with MPI since it won't necessarily capture the error
# message correctly
with pytest.raises(RuntimeError, match="rejected"):
model.run(openmc_exec=config['exe'])
def test_exceptions():
with pytest.raises(AttributeError, match=r'Please use the FileSource class'):
@ -156,4 +243,4 @@ def test_exceptions():
with pytest.raises(AttributeError, match=r'has no attribute \'frisbee\''):
s = openmc.IndependentSource()
s.frisbee
s.frisbee