Generic filter for rejecting sources points when sampling

[shimwell]

As originally discussed on Discourse https://openmc.discourse.group/t/source-uniform-sampling-in-cell-of-given-material/1062

It would be great to add generic filter to allow the rejection of source points when being sampled during the creation of particles.

Source sampling should allow for rejection of the proposed sampled coordinates if the coordinates are not in a user specified material or not in a user specified cell.

While the generic source sampling filter would be separate to the distributions themselves it could be passed to the distributions and the distributions should make use of the filter to ascertain if the sampling point results in a neutron or not.

Perhaps the best distribution to demonstrate this on is the openmc.stats.Box distribution as it already has a similar feature.

This will allow us to make interesting shaped sources like the example picture below which shows source sampling randomly across a openmc.stats.Box rejecting sample points that don't appear within a certain material / cell

There are some related lines in OpenMC relating to a similar feature that allows sampling points to be rejected if they are not within a fissionable material with the only_fissionable keyword that might be useful for anyone considering this feature.

openmc/src/distribution_spatial.cpp

Lines 187 to 204 in efcc0c1

	SpatialBox::SpatialBox(pugi::xml_node node, bool fission)
	: only_fissionable_{fission}
	{
	// Read lower-right/upper-left coordinates
	auto params = get_node_array<double>(node, "parameters");
	if (params.size() != 6)
	openmc::fatal_error("Box/fission spatial source must have six "
	"parameters specified.");
	lower_left_ = Position{params[0], params[1], params[2]};
	upper_right_ = Position{params[3], params[4], params[5]};
	}
	Position SpatialBox::sample(uint64_t* seed) const
	{
	Position xi {prn(seed), prn(seed), prn(seed)};
	return lower_left_ + xi*(upper_right_ - lower_left_);
	}

openmc/src/source.cpp

Lines 171 to 173 in efcc0c1

	if (space_box) {
	if (space_box->only_fissionable()) {
	// Determine material

openmc/openmc/stats/multivariate.py

Lines 617 to 716 in efcc0c1

	class Box(Spatial):
	"""Uniform distribution of coordinates in a rectangular cuboid.
	Parameters
	----------
	lower_left : Iterable of float
	Lower-left coordinates of cuboid
	upper_right : Iterable of float
	Upper-right coordinates of cuboid
	only_fissionable : bool, optional
	Whether spatial sites should only be accepted if they occur in
	fissionable materials
	Attributes
	----------
	lower_left : Iterable of float
	Lower-left coordinates of cuboid
	upper_right : Iterable of float
	Upper-right coordinates of cuboid
	only_fissionable : bool, optional
	Whether spatial sites should only be accepted if they occur in
	fissionable materials
	"""
	def __init__(self, lower_left, upper_right, only_fissionable=False):
	self.lower_left = lower_left
	self.upper_right = upper_right
	self.only_fissionable = only_fissionable
	@property
	def lower_left(self):
	return self._lower_left
	@property
	def upper_right(self):
	return self._upper_right
	@property
	def only_fissionable(self):
	return self._only_fissionable
	@lower_left.setter
	def lower_left(self, lower_left):
	cv.check_type('lower left coordinate', lower_left, Iterable, Real)
	cv.check_length('lower left coordinate', lower_left, 3)
	self._lower_left = lower_left
	@upper_right.setter
	def upper_right(self, upper_right):
	cv.check_type('upper right coordinate', upper_right, Iterable, Real)
	cv.check_length('upper right coordinate', upper_right, 3)
	self._upper_right = upper_right
	@only_fissionable.setter
	def only_fissionable(self, only_fissionable):
	cv.check_type('only fissionable', only_fissionable, bool)
	self._only_fissionable = only_fissionable
	def to_xml_element(self):
	"""Return XML representation of the box distribution
	Returns
	-------
	element : xml.etree.ElementTree.Element
	XML element containing box distribution data
	"""
	element = ET.Element('space')
	if self.only_fissionable:
	element.set("type", "fission")
	else:
	element.set("type", "box")
	params = ET.SubElement(element, "parameters")
	params.text = ' '.join(map(str, self.lower_left)) + ' ' + \
	' '.join(map(str, self.upper_right))
	return element
	@classmethod
	def from_xml_element(cls, elem):
	"""Generate box distribution from an XML element
	Parameters
	----------
	elem : xml.etree.ElementTree.Element
	XML element
	Returns
	-------
	openmc.stats.Box
	Box distribution generated from XML element
	"""
	only_fissionable = get_text(elem, 'type') == 'fission'
	params = [float(x) for x in get_text(elem, 'parameters').split()]
	lower_left = params[:len(params)//2]
	upper_right = params[len(params)//2:]
	return cls(lower_left, upper_right, only_fissionable)

[shimwell]

I've issued a bounty reward for this issue https://www.bountysource.com/issues/98138076-generic-filter-for-rejecting-sources-points-when-sampling

[RemDelaporteMathurin]

Could the spawn probability be inhomogeneous and a function of x, y, z?

[paulromano]

I have a branch with this implemented -- just need to clean it up, make it a little more generic, and open a pull request. I should have it out there in the next few weeks.

[RemDelaporteMathurin]

@paulromano in some applications, the probability of birth is not homogeneous (see this example). Do you think this feature would allow users to specify a probability distribution or strengths distributions?

[shimwell]

Yes that would also work, the feature provides a cell used to exclude birth locations outside of the volume. You can throw any distribution of points and it will cut out the ones that are not within the cell.

[RemDelaporteMathurin]

Oh so the usage would "just" be 1) give a bunch of sources 2) it exclude the sources outside of this cell

I thought the usage would be 1) give a cell 2) openmc samples sources in random locations within that cell

[shimwell]

👍 user provides openmc.Source or list of openmc.Source as normal, openmc samples particle birth coordinates then birth coordinates outside of cell/material are excluded. So not quite excluding entire openmc.Sources but particles of sources which is slightly different

Oh so the usage would "just" be 1) give a bunch of sources 2) it exclude the sources outside of this cell

This approach is handy as it works for complex cell shapes. This issue also requests exclusion by material as well which will be handy for DAGMC use case where we don't always know the cell ID