mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-28 22:26:08 -04:00
Add translate methods to surfaces and regions
This commit is contained in:
Paul Romano
parent
2943e38a81
commit
394ac75560
4 changed files with 409 additions and 12 deletions
|
|
@ -229,14 +229,27 @@ class Region(metaclass=ABCMeta):
|
|||
clone : openmc.Region
|
||||
The clone of this region
|
||||
|
||||
Raises
|
||||
------
|
||||
NotImplementedError
|
||||
This method is not implemented for the abstract region class.
|
||||
"""
|
||||
pass
|
||||
|
||||
@abstractmethod
|
||||
def translate(self, vector, memo=None):
|
||||
"""Translate region in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which region should be translated
|
||||
memo : dict or None
|
||||
Dictionary used for memoization
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Region
|
||||
Translated region
|
||||
|
||||
"""
|
||||
raise NotImplementedError('The clone method is not implemented for '
|
||||
'the abstract region class.')
|
||||
pass
|
||||
|
||||
|
||||
class Intersection(Region, MutableSequence):
|
||||
|
|
@ -352,6 +365,26 @@ class Intersection(Region, MutableSequence):
|
|||
clone[:] = [n.clone(memo) for n in self]
|
||||
return clone
|
||||
|
||||
def translate(self, vector, memo=None):
|
||||
"""Translate region in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which region should be translated
|
||||
memo : dict or None
|
||||
Dictionary used for memoization
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Intersection
|
||||
Translated region
|
||||
|
||||
"""
|
||||
if memo is None:
|
||||
memo = {}
|
||||
return type(self)(n.translate(vector, memo) for n in self)
|
||||
|
||||
|
||||
class Union(Region, MutableSequence):
|
||||
r"""Union of two or more regions.
|
||||
|
|
@ -464,6 +497,26 @@ class Union(Region, MutableSequence):
|
|||
clone[:] = [n.clone(memo) for n in self]
|
||||
return clone
|
||||
|
||||
def translate(self, vector, memo=None):
|
||||
"""Translate region in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which region should be translated
|
||||
memo : dict or None
|
||||
Dictionary used for memoization
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Union
|
||||
Translated region
|
||||
|
||||
"""
|
||||
if memo is None:
|
||||
memo = {}
|
||||
return type(self)(n.translate(vector, memo) for n in self)
|
||||
|
||||
|
||||
class Complement(Region):
|
||||
"""Complement of a region.
|
||||
|
|
@ -584,3 +637,23 @@ class Complement(Region):
|
|||
clone = deepcopy(self)
|
||||
clone.node = self.node.clone(memo)
|
||||
return clone
|
||||
|
||||
def translate(self, vector, memo=None):
|
||||
"""Translate region in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which region should be translated
|
||||
memo : dict or None
|
||||
Dictionary used for memoization
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Complement
|
||||
Translated region
|
||||
|
||||
"""
|
||||
if memo is None:
|
||||
memo = {}
|
||||
return type(self)(self.node.translate(vector, memo))
|
||||
|
|
|
|||
|
|
@ -1,4 +1,4 @@
|
|||
from abc import ABCMeta
|
||||
from abc import ABCMeta, abstractmethod
|
||||
from collections import OrderedDict
|
||||
from copy import deepcopy
|
||||
from functools import partial
|
||||
|
|
@ -15,7 +15,7 @@ from openmc.mixin import IDManagerMixin
|
|||
_BOUNDARY_TYPES = ['transmission', 'vacuum', 'reflective', 'periodic']
|
||||
|
||||
|
||||
class Surface(IDManagerMixin):
|
||||
class Surface(IDManagerMixin, metaclass=ABCMeta):
|
||||
"""An implicit surface with an associated boundary condition.
|
||||
|
||||
An implicit surface is defined as the set of zeros of a function of the
|
||||
|
|
@ -142,7 +142,6 @@ class Surface(IDManagerMixin):
|
|||
desired half-space
|
||||
|
||||
"""
|
||||
|
||||
return (np.array([-np.inf, -np.inf, -np.inf]),
|
||||
np.array([np.inf, np.inf, np.inf]))
|
||||
|
||||
|
|
@ -175,6 +174,14 @@ class Surface(IDManagerMixin):
|
|||
|
||||
return memo[self]
|
||||
|
||||
@abstractmethod
|
||||
def evaluate(self, point):
|
||||
pass
|
||||
|
||||
@abstractmethod
|
||||
def translate(self, vector):
|
||||
pass
|
||||
|
||||
def to_xml_element(self):
|
||||
"""Return XML representation of the surface
|
||||
|
||||
|
|
@ -432,13 +439,34 @@ class Plane(Surface):
|
|||
Returns
|
||||
-------
|
||||
float
|
||||
:math:`Ax' + By' + Cz' - d`
|
||||
:math:`Ax' + By' + Cz' - D`
|
||||
|
||||
"""
|
||||
|
||||
x, y, z = point
|
||||
return self.a*x + self.b*y + self.c*z - self.d
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Plane
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
d = self.d + self.a*vx + self.b*vy + self.c*vz
|
||||
if d == self.d:
|
||||
return self
|
||||
else:
|
||||
return type(self)(A=self.a, B=self.b, C=self.c, D=d)
|
||||
|
||||
def to_xml_element(self):
|
||||
"""Return XML representation of the surface
|
||||
|
||||
|
|
@ -593,6 +621,26 @@ class XPlane(Plane):
|
|||
"""
|
||||
return point[0] - self.x0
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.XPlane
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx = vector[0]
|
||||
if vx == 0:
|
||||
return self
|
||||
else:
|
||||
return type(self)(x0=self.x0 + vx)
|
||||
|
||||
|
||||
class YPlane(Plane):
|
||||
"""A plane perpendicular to the y axis of the form :math:`y - y_0 = 0`
|
||||
|
|
@ -699,6 +747,26 @@ class YPlane(Plane):
|
|||
"""
|
||||
return point[1] - self.y0
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.YPlane
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vy = vector[1]
|
||||
if vy == 0.0:
|
||||
return self
|
||||
else:
|
||||
return type(self)(y0=self.y0 + vy)
|
||||
|
||||
|
||||
class ZPlane(Plane):
|
||||
"""A plane perpendicular to the z axis of the form :math:`z - z_0 = 0`
|
||||
|
|
@ -805,8 +873,28 @@ class ZPlane(Plane):
|
|||
"""
|
||||
return point[2] - self.z0
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
class Cylinder(Surface, metaclass=ABCMeta):
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.ZPlane
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vz = vector[2]
|
||||
if vz == 0.0:
|
||||
return self
|
||||
else:
|
||||
return type(self)(z0=self.z0 + vz)
|
||||
|
||||
|
||||
class Cylinder(Surface):
|
||||
"""A cylinder whose length is parallel to the x-, y-, or z-axis.
|
||||
|
||||
Parameters
|
||||
|
|
@ -977,6 +1065,28 @@ class XCylinder(Cylinder):
|
|||
z = point[2] - self.z0
|
||||
return y**2 + z**2 - self.r**2
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.XCylinder
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
if vy == 0.0 and vz == 0.0:
|
||||
return self
|
||||
else:
|
||||
y0 = self.y0 + vy
|
||||
z0 = self.z0 + vz
|
||||
return type(self)(y0=y0, z0=z0, R=self.r)
|
||||
|
||||
|
||||
class YCylinder(Cylinder):
|
||||
"""An infinite cylinder whose length is parallel to the y-axis of the form
|
||||
|
|
@ -1099,6 +1209,28 @@ class YCylinder(Cylinder):
|
|||
z = point[2] - self.z0
|
||||
return x**2 + z**2 - self.r**2
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.YCylinder
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
if vx == 0.0 and vz == 0.0:
|
||||
return self
|
||||
else:
|
||||
x0 = self.x0 + vx
|
||||
z0 = self.z0 + vz
|
||||
return type(self)(x0=x0, z0=z0, R=self.r)
|
||||
|
||||
|
||||
class ZCylinder(Cylinder):
|
||||
"""An infinite cylinder whose length is parallel to the z-axis of the form
|
||||
|
|
@ -1221,6 +1353,28 @@ class ZCylinder(Cylinder):
|
|||
y = point[1] - self.y0
|
||||
return x**2 + y**2 - self.r**2
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.ZCylinder
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
if vx == 0.0 and vy == 0.0:
|
||||
return self
|
||||
else:
|
||||
x0 = self.x0 + vx
|
||||
y0 = self.y0 + vy
|
||||
return type(self)(x0=x0, y0=y0, R=self.r)
|
||||
|
||||
|
||||
class Sphere(Surface):
|
||||
"""A sphere of the form :math:`(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 = R^2`.
|
||||
|
|
@ -1369,8 +1523,31 @@ class Sphere(Surface):
|
|||
z = point[2] - self.z0
|
||||
return x**2 + y**2 + z**2 - self.r**2
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
class Cone(Surface, metaclass=ABCMeta):
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Sphere
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
if vx == 0.0 and vy == 0.0 and vz == 0.0:
|
||||
return self
|
||||
else:
|
||||
x0 = self.x0 + vx
|
||||
y0 = self.y0 + vy
|
||||
z0 = self.z0 + vz
|
||||
return type(self)(x0=x0, y0=y0, z0=z0, R=self.r)
|
||||
|
||||
|
||||
class Cone(Surface):
|
||||
"""A conical surface parallel to the x-, y-, or z-axis.
|
||||
|
||||
Parameters
|
||||
|
|
@ -1463,6 +1640,29 @@ class Cone(Surface, metaclass=ABCMeta):
|
|||
check_type('R^2 coefficient', R2, Real)
|
||||
self._coefficients['R2'] = R2
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Cone
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
if vx == 0.0 and vy == 0.0 and vz == 0.0:
|
||||
return self
|
||||
else:
|
||||
x0 = self.x0 + vx
|
||||
y0 = self.y0 + vy
|
||||
z0 = self.z0 + vz
|
||||
return type(self)(x0=x0, y0=y0, z0=z0, R2=self.r2)
|
||||
|
||||
|
||||
class XCone(Cone):
|
||||
"""A cone parallel to the x-axis of the form :math:`(y - y_0)^2 + (z - z_0)^2 =
|
||||
|
|
@ -1842,6 +2042,30 @@ class Quadric(Surface):
|
|||
y*(self.b*y + self.e*z + self.h) + \
|
||||
z*(self.c*z + self.f*x + self.j) + self.k
|
||||
|
||||
def translate(self, vector):
|
||||
"""Translate surface in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which surface should be translated
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Quadric
|
||||
Translated surface
|
||||
|
||||
"""
|
||||
vx, vy, vz = vector
|
||||
a, b, c, d, e, f, g, h, j, k = (getattr(self, key) for key in
|
||||
self._coeff_keys)
|
||||
k = (k + vx*vx + vy*vy + vz*vz + d*vx*vy + e*vy*vz + f*vx*vz
|
||||
- g*vx - h*vy - j*vz)
|
||||
g = g - 2*a*vx - d*vy - f*vz
|
||||
h = h - 2*b*vy - d*vx - e*vz
|
||||
j = j - 2*c*vz - e*vy - f*vx
|
||||
return type(self)(a=a, b=b, c=c, d=d, e=e, f=f, g=g, h=h, j=j, k=k)
|
||||
|
||||
|
||||
class Halfspace(Region):
|
||||
"""A positive or negative half-space region.
|
||||
|
|
@ -1987,3 +2211,30 @@ class Halfspace(Region):
|
|||
clone = deepcopy(self)
|
||||
clone.surface = self.surface.clone(memo)
|
||||
return clone
|
||||
|
||||
def translate(self, vector, memo=None):
|
||||
"""Translate half-space in given direction
|
||||
|
||||
Parameters
|
||||
----------
|
||||
vector : iterable of float
|
||||
Direction in which region should be translated
|
||||
memo : dict or None
|
||||
Dictionary used for memoization
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Halfspace
|
||||
Translated half-space
|
||||
|
||||
"""
|
||||
if memo is None:
|
||||
memo = {}
|
||||
|
||||
# If translated surface not in memo, add it
|
||||
key = (self.surface, tuple(vector))
|
||||
if key not in memo:
|
||||
memo[key] = self.surface.translate(vector)
|
||||
|
||||
# Return translated surface
|
||||
return type(self)(memo[key], self.side)
|
||||
|
|
|
|||
|
|
@ -34,6 +34,12 @@ def test_union(reset):
|
|||
assert (6, -1, 0) not in reg2
|
||||
assert str(reg2) == '((1 | -2) 3)'
|
||||
|
||||
# translate method
|
||||
regt = region.translate((2.0, 0.0, 0.0))
|
||||
assert (-4, 0, 0) in regt
|
||||
assert (6, 0, 0) not in regt
|
||||
assert (8, 0, 0) in regt
|
||||
|
||||
|
||||
def test_intersection(reset):
|
||||
s1 = openmc.XPlane(surface_id=1, x0=5)
|
||||
|
|
@ -61,6 +67,12 @@ def test_intersection(reset):
|
|||
assert (-6, -2, 0) not in reg2
|
||||
assert str(reg2) == '((-1 2) | 3)'
|
||||
|
||||
# translate method
|
||||
regt = region.translate((2.0, 0.0, 0.0))
|
||||
assert (-4, 0, 0) not in regt
|
||||
assert (6, 0, 0) in regt
|
||||
assert (8, 0, 0) not in regt
|
||||
|
||||
|
||||
def test_complement(reset):
|
||||
zcyl = openmc.ZCylinder(surface_id=1, R=1.)
|
||||
|
|
@ -88,6 +100,12 @@ def test_complement(reset):
|
|||
assert (0, 0, 6) not in inside
|
||||
assert (0, 0, 6) in outside
|
||||
|
||||
# translate method
|
||||
inside_t = inside.translate((1.0, 1.0, 1.0))
|
||||
ll, ur = inside_t.bounding_box
|
||||
assert ll == pytest.approx((0., 0., -4.))
|
||||
assert ur == pytest.approx((2., 2., 6.))
|
||||
|
||||
|
||||
def test_get_surfaces():
|
||||
s1 = openmc.XPlane()
|
||||
|
|
|
|||
|
|
@ -29,6 +29,10 @@ def test_plane():
|
|||
x, y, z = (4, 3, 6)
|
||||
assert s.evaluate((x, y, z)) == pytest.approx(s.a*x + s.b*y + s.c*z - s.d)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 0.0, 0.0))
|
||||
assert (st.a, st.b, st.c, st.d) == (s.a, s.b, s.c, 4)
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -69,6 +73,10 @@ def test_xplane():
|
|||
# evaluate method
|
||||
assert s.evaluate((5., 0., 0.)) == pytest.approx(2.)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 0.0, 0.0))
|
||||
assert st.x0 == s.x0 + 1
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -94,6 +102,10 @@ def test_yplane():
|
|||
# evaluate method
|
||||
assert s.evaluate((0., 0., 0.)) == pytest.approx(-3.)
|
||||
|
||||
# translate method
|
||||
st = s.translate((0.0, 1.0, 0.0))
|
||||
assert st.y0 == s.y0 + 1
|
||||
|
||||
|
||||
def test_zplane():
|
||||
s = openmc.ZPlane(z0=3.)
|
||||
|
|
@ -116,6 +128,10 @@ def test_zplane():
|
|||
# evaluate method
|
||||
assert s.evaluate((0., 0., 10.)) == pytest.approx(7.)
|
||||
|
||||
# translate method
|
||||
st = s.translate((0.0, 0.0, 1.0))
|
||||
assert st.z0 == s.z0 + 1
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -138,6 +154,12 @@ def test_xcylinder():
|
|||
# evaluate method
|
||||
assert s.evaluate((0, y, z)) == pytest.approx(-r**2)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
assert st.y0 == s.y0 + 1
|
||||
assert st.z0 == s.z0 + 1
|
||||
assert st.r == s.r
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -169,6 +191,12 @@ def test_ycylinder():
|
|||
# evaluate method
|
||||
assert s.evaluate((x, 0, z)) == pytest.approx(-r**2)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
assert st.x0 == s.x0 + 1
|
||||
assert st.z0 == s.z0 + 1
|
||||
assert st.r == s.r
|
||||
|
||||
|
||||
def test_zcylinder():
|
||||
x, y, r = 3, 5, 2
|
||||
|
|
@ -188,6 +216,12 @@ def test_zcylinder():
|
|||
# evaluate method
|
||||
assert s.evaluate((x, y, 0)) == pytest.approx(-r**2)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
assert st.x0 == s.x0 + 1
|
||||
assert st.y0 == s.y0 + 1
|
||||
assert st.r == s.r
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -211,6 +245,13 @@ def test_sphere():
|
|||
# evaluate method
|
||||
assert s.evaluate((x, y, z)) == pytest.approx(-r**2)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
assert st.x0 == s.x0 + 1
|
||||
assert st.y0 == s.y0 + 1
|
||||
assert st.z0 == s.z0 + 1
|
||||
assert st.r == s.r
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -229,6 +270,13 @@ def cone_common(apex, r2, cls):
|
|||
# evaluate method -- should be zero at apex
|
||||
assert s.evaluate((x, y, z)) == pytest.approx(0.0)
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
assert st.x0 == s.x0 + 1
|
||||
assert st.y0 == s.y0 + 1
|
||||
assert st.z0 == s.z0 + 1
|
||||
assert st.r2 == s.r2
|
||||
|
||||
# Make sure repr works
|
||||
repr(s)
|
||||
|
||||
|
|
@ -272,6 +320,13 @@ def test_quadric():
|
|||
assert s.evaluate((0., 0., 0.)) == pytest.approx(coeffs['k'])
|
||||
assert s.evaluate((1., 1., 1.)) == pytest.approx(3 + coeffs['k'])
|
||||
|
||||
# translate method
|
||||
st = s.translate((1.0, 1.0, 1.0))
|
||||
for coeff in 'abcdef':
|
||||
assert getattr(s, coeff) == getattr(st, coeff)
|
||||
assert (st.g, st.h, st.j) == (-2, -2, -2)
|
||||
assert st.k == s.k + 3
|
||||
|
||||
|
||||
def test_cylinder_from_points():
|
||||
# Generate 45-degree rotated cylinder in x-y plane with radius 1
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue