mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 13:45:36 -04:00
Add Surface.evaluate() methods and Region.__contains__ methods
This commit is contained in:
Paul Romano
parent
e6d396a198
commit
07be6306ce
2 changed files with 298 additions and 0 deletions
|
|
@ -28,6 +28,10 @@ class Region(object):
|
|||
def __invert__(self):
|
||||
return Complement(self)
|
||||
|
||||
@abstractmethod
|
||||
def __contains__(self, point):
|
||||
return False
|
||||
|
||||
@abstractmethod
|
||||
def __str__(self):
|
||||
return ''
|
||||
|
|
@ -229,6 +233,22 @@ class Intersection(Region):
|
|||
def __init__(self, *nodes):
|
||||
self.nodes = list(nodes)
|
||||
|
||||
def __contains__(self, point):
|
||||
"""Check whether a point is contained in the region.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
Cartesian coordinates, :math:`(x',y',z')`, of the point
|
||||
|
||||
Returns
|
||||
-------
|
||||
bool
|
||||
Whether the point is in the region
|
||||
|
||||
"""
|
||||
return all(point in n for n in self.nodes)
|
||||
|
||||
def __str__(self):
|
||||
return '(' + ' '.join(map(str, self.nodes)) + ')'
|
||||
|
||||
|
|
@ -281,6 +301,22 @@ class Union(Region):
|
|||
def __init__(self, *nodes):
|
||||
self.nodes = list(nodes)
|
||||
|
||||
def __contains__(self, point):
|
||||
"""Check whether a point is contained in the region.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
Cartesian coordinates, :math:`(x',y',z')`, of the point
|
||||
|
||||
Returns
|
||||
-------
|
||||
bool
|
||||
Whether the point is in the region
|
||||
|
||||
"""
|
||||
return any(point in n for n in self.nodes)
|
||||
|
||||
def __str__(self):
|
||||
return '(' + ' | '.join(map(str, self.nodes)) + ')'
|
||||
|
||||
|
|
@ -336,6 +372,22 @@ class Complement(Region):
|
|||
def __init__(self, node):
|
||||
self.node = node
|
||||
|
||||
def __contains__(self, point):
|
||||
"""Check whether a point is contained in the region.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
Cartesian coordinates, :math:`(x',y',z')`, of the point
|
||||
|
||||
Returns
|
||||
-------
|
||||
bool
|
||||
Whether the point is in the region
|
||||
|
||||
"""
|
||||
return point not in self.node
|
||||
|
||||
def __str__(self):
|
||||
return '~' + str(self.node)
|
||||
|
||||
|
|
|
|||
|
|
@ -295,6 +295,25 @@ class Plane(Surface):
|
|||
self._periodic_surface = periodic_surface
|
||||
periodic_surface._periodic_surface = self
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`Ax' + By' + Cz' - d`
|
||||
|
||||
"""
|
||||
|
||||
x, y, z = point
|
||||
return self.a*x + self.b*y + self.c*z - self.d
|
||||
|
||||
def create_xml_subelement(self):
|
||||
element = super(Plane, self).create_xml_subelement()
|
||||
|
||||
|
|
@ -392,6 +411,23 @@ class XPlane(Plane):
|
|||
return (np.array([self.x0, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`x' - x_0`
|
||||
|
||||
"""
|
||||
return point[0] - self.x0
|
||||
|
||||
|
||||
class YPlane(Plane):
|
||||
"""A plane perpendicular to the y axis of the form :math:`y - y_0 = 0`
|
||||
|
|
@ -481,6 +517,23 @@ class YPlane(Plane):
|
|||
return (np.array([-np.inf, self.y0, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`y' - y_0`
|
||||
|
||||
"""
|
||||
return point[1] - self.y0
|
||||
|
||||
|
||||
class ZPlane(Plane):
|
||||
"""A plane perpendicular to the z axis of the form :math:`z - z_0 = 0`
|
||||
|
|
@ -570,6 +623,23 @@ class ZPlane(Plane):
|
|||
return (np.array([-np.inf, -np.inf, self.z0]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`z' - z_0`
|
||||
|
||||
"""
|
||||
return point[2] - self.z0
|
||||
|
||||
|
||||
class Cylinder(Surface):
|
||||
"""A cylinder whose length is parallel to the x-, y-, or z-axis.
|
||||
|
|
@ -728,6 +798,25 @@ class XCylinder(Cylinder):
|
|||
return (np.array([-np.inf, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2`
|
||||
|
||||
"""
|
||||
y = point[1] - self.y0
|
||||
z = point[2] - self.z0
|
||||
return y**2 + z**2 - self.r**2
|
||||
|
||||
|
||||
class YCylinder(Cylinder):
|
||||
"""An infinite cylinder whose length is parallel to the y-axis of the form
|
||||
|
|
@ -831,6 +920,25 @@ class YCylinder(Cylinder):
|
|||
return (np.array([-np.inf, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
z = point[2] - self.z0
|
||||
return x**2 + z**2 - self.r**2
|
||||
|
||||
|
||||
class ZCylinder(Cylinder):
|
||||
"""An infinite cylinder whose length is parallel to the z-axis of the form
|
||||
|
|
@ -934,6 +1042,25 @@ class ZCylinder(Cylinder):
|
|||
return (np.array([-np.inf, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
y = point[1] - self.y0
|
||||
return x**2 + y**2 - self.r**2
|
||||
|
||||
|
||||
class Sphere(Surface):
|
||||
"""A sphere of the form :math:`(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 = R^2`.
|
||||
|
|
@ -1062,6 +1189,26 @@ class Sphere(Surface):
|
|||
return (np.array([-np.inf, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(x' - x_0)^2 + (y' - y_0)^2 + (z' - z_0)^2 - R^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
y = point[1] - self.y0
|
||||
z = point[2] - self.z0
|
||||
return x**2 + y**2 + z**2 - self.r**2
|
||||
|
||||
|
||||
class Cone(Surface):
|
||||
"""A conical surface parallel to the x-, y-, or z-axis.
|
||||
|
|
@ -1214,6 +1361,26 @@ class XCone(Cone):
|
|||
|
||||
self._type = 'x-cone'
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2(x' - x_0)^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
y = point[1] - self.y0
|
||||
z = point[2] - self.z0
|
||||
return y**2 + z**2 - self.r2*x**2
|
||||
|
||||
|
||||
class YCone(Cone):
|
||||
"""A cone parallel to the y-axis of the form :math:`(x - x_0)^2 + (z - z_0)^2 =
|
||||
|
|
@ -1270,6 +1437,26 @@ class YCone(Cone):
|
|||
|
||||
self._type = 'y-cone'
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2(y' - y_0)^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
y = point[1] - self.y0
|
||||
z = point[2] - self.z0
|
||||
return x**2 + z**2 - self.r2*y**2
|
||||
|
||||
|
||||
class ZCone(Cone):
|
||||
"""A cone parallel to the x-axis of the form :math:`(x - x_0)^2 + (y - y_0)^2 =
|
||||
|
|
@ -1326,6 +1513,26 @@ class ZCone(Cone):
|
|||
|
||||
self._type = 'z-cone'
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2(z' - z_0)^2`
|
||||
|
||||
"""
|
||||
x = point[0] - self.x0
|
||||
y = point[1] - self.y0
|
||||
z = point[2] - self.z0
|
||||
return x**2 + y**2 - self.r2*z**2
|
||||
|
||||
|
||||
class Quadric(Surface):
|
||||
"""A surface of the form :math:`Ax^2 + By^2 + Cz^2 + Dxy + Eyz + Fxz + Gx + Hy +
|
||||
|
|
@ -1471,6 +1678,27 @@ class Quadric(Surface):
|
|||
check_type('k coefficient', k, Real)
|
||||
self._coefficients['k'] = k
|
||||
|
||||
def evaluate(self, point):
|
||||
"""Evaluate the surface equation at a given point.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
The Cartesian coordinates, :math:`(x',y',z')`, at which the surface
|
||||
equation should be evaluated.
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`Ax'^2 + By'^2 + Cz'^2 + Dx'y' + Ey'z' + Fx'z' + Gx' + Hy' +
|
||||
Jz' + K = 0`
|
||||
|
||||
"""
|
||||
x, y, z = point
|
||||
return x*(self.a*x + self.d*y + self.g) + \
|
||||
y*(self.b*y + self.e*z + self.h) + \
|
||||
z*(self.c*z + self.f*x + self.j) + self.k
|
||||
|
||||
|
||||
class Halfspace(Region):
|
||||
"""A positive or negative half-space region.
|
||||
|
|
@ -1516,6 +1744,24 @@ class Halfspace(Region):
|
|||
def __invert__(self):
|
||||
return -self.surface if self.side == '+' else +self.surface
|
||||
|
||||
def __contains__(self, point):
|
||||
"""Check whether a point is contained in the half-space.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
point : 3-tuple of float
|
||||
Cartesian coordinates, :math:`(x',y',z')`, of the point
|
||||
|
||||
Returns
|
||||
-------
|
||||
bool
|
||||
Whether the point is in the half-space
|
||||
|
||||
"""
|
||||
|
||||
val = self.surface.evaluate(point)
|
||||
return val >= 0. if self.side == '+' else val < 0.
|
||||
|
||||
@property
|
||||
def surface(self):
|
||||
return self._surface
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue