moved corner radius into get_rectangular_prism and get_hexagonal_prism functions

This commit is contained in:
Sam Shaner 2017-05-17 10:59:55 -04:00
parent 0fa5c63efb
commit 455c655c41

View file

@ -1960,7 +1960,7 @@ class Halfspace(Region):
return clone
def get_rectangular_prism(width, height, axis='z', origin=(0., 0.),
def get_rectangular_prism(width, height, corner_radius=0., axis='z', origin=(0., 0.),
boundary_type='transmission'):
"""Get an infinite rectangular prism from four planar surfaces.
@ -1972,6 +1972,8 @@ def get_rectangular_prism(width, height, axis='z', origin=(0., 0.),
height: float
Prism height in units of cm. The height is aligned with the z, z,
or y axes for prisms parallel to the x, y, or z axis, respectively.
corner_radius: float
Prism corner radius in units of cm. Defaults to 0.
axis : {'x', 'y', 'z'}
Axis with which the infinite length of the prism should be aligned.
Defaults to 'z'.
@ -1992,6 +1994,7 @@ def get_rectangular_prism(width, height, axis='z', origin=(0., 0.),
check_type('width', width, Real)
check_type('height', height, Real)
check_type('corner_radius', corner_radius, Real)
check_value('axis', axis, ['x','y','z'])
check_type('origin', origin, Iterable, Real)
@ -2026,10 +2029,120 @@ def get_rectangular_prism(width, height, axis='z', origin=(0., 0.),
boundary_type=boundary_type)
prism = +min_x & -max_x & +min_y & -max_y
if corner_radius > 0.:
if axis == 'x':
y_min_z_min = XCylinder(name='y min z min',
y0=origin[0] - width /2. + corner_radius,
z0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
y_min_z_max = XCylinder(name='y min z max',
y0=origin[0] - width /2. + corner_radius,
z0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
y_max_z_min = XCylinder(name='y max z min',
y0=origin[0] + width /2. - corner_radius,
z0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
y_max_z_max = XCylinder(name='y max z max',
y0=origin[0] + width /2. - corner_radius,
z0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
y_min = YPlane(name='min y', y0=-width/2.+origin[0]+corner_radius,
boundary_type=boundary_type)
y_max = YPlane(name='max y', y0=+width/2.+origin[0]-corner_radius,
boundary_type=boundary_type)
z_min = ZPlane(name='min z', z0=-height/2.+origin[1]+corner_radius,
boundary_type=boundary_type)
z_max = ZPlane(name='max z', z0=+height/2.+origin[1]-corner_radius,
boundary_type=boundary_type)
corners = (+y_min_z_min & -y_min & -z_min) | \
(+y_min_z_max & -y_min & +z_max) | \
(+y_max_z_min & +y_max & -z_min) | \
(+y_max_z_max & +y_max & +z_max)
elif axis == 'y':
x_min_z_min = YCylinder(name='x min z min',
x0=origin[0] - width /2. + corner_radius,
z0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_min_z_max = YCylinder(name='x min z max',
x0=origin[0] - width /2. + corner_radius,
z0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_max_z_min = YCylinder(name='x max z min',
x0=origin[0] + width /2. - corner_radius,
z0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_max_z_max = YCylinder(name='x max z max',
x0=origin[0] + width /2. - corner_radius,
z0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_min = XPlane(name='min x', x0=-width/2.+origin[0]+corner_radius,
boundary_type=boundary_type)
x_max = XPlane(name='max x', x0=+width/2.+origin[0]-corner_radius,
boundary_type=boundary_type)
z_min = ZPlane(name='min z', z0=-height/2.+origin[1]+corner_radius,
boundary_type=boundary_type)
z_max = ZPlane(name='max z', z0=+height/2.+origin[1]-corner_radius,
boundary_type=boundary_type)
corners = (+x_min_z_min & -x_min & -z_min) | \
(+x_min_z_max & -x_min & +z_max) | \
(+x_max_z_min & +x_max & -z_min) | \
(+x_max_z_max & +x_max & +z_max)
else:
x_min_y_min = ZCylinder(name='x min y min',
x0=origin[0] - width /2. + corner_radius,
y0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_min_y_max = ZCylinder(name='x min y max',
x0=origin[0] - width /2. + corner_radius,
y0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_max_y_min = ZCylinder(name='x max y min',
x0=origin[0] + width /2. - corner_radius,
y0=origin[1] - height/2. + corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_max_y_max = ZCylinder(name='x max y max',
x0=origin[0] + width /2. - corner_radius,
y0=origin[1] + height/2. - corner_radius,
R=corner_radius,
boundary_type=boundary_type)
x_min = XPlane(name='min x', x0=-width/2.+origin[0]+corner_radius,
boundary_type=boundary_type)
x_max = XPlane(name='max x', x0=+width/2.+origin[0]-corner_radius,
boundary_type=boundary_type)
y_min = YPlane(name='min y', y0=-height/2.+origin[1]+corner_radius,
boundary_type=boundary_type)
y_max = YPlane(name='max y', y0=+height/2.+origin[1]-corner_radius,
boundary_type=boundary_type)
corners = (+x_min_y_min & -x_min & -y_min) | \
(+x_min_y_max & -x_min & +y_max) | \
(+x_max_y_min & +x_max & -y_min) | \
(+x_max_y_max & +x_max & +y_max)
prism = prism & ~corners
return prism
def get_hexagonal_prism(edge_length=1., orientation='y',
def get_hexagonal_prism(edge_length=1., corner_radius=0., orientation='y',
boundary_type='transmission'):
"""Create a hexagon region from six surface planes.
@ -2037,6 +2150,8 @@ def get_hexagonal_prism(edge_length=1., orientation='y',
----------
edge_length : float
Length of a side of the hexagon in cm
corner_radius: float
Prism corner radius in units of cm. Defaults to 0.
orientation : {'x', 'y'}
An 'x' orientation means that two sides of the hexagon are parallel to
the x-axis and a 'y' orientation means that two sides of the hexagon are
@ -2070,8 +2185,8 @@ def get_hexagonal_prism(edge_length=1., orientation='y',
# y = -x/sqrt(3) - a
lower_left = Plane(A=c, B=1., D=-l, boundary_type=boundary_type)
return Intersection(-right, +left, -upper_right, -upper_left,
+lower_right, +lower_left)
prism = -right & +left & -upper_right & -upper_left & \
+lower_right & +lower_left
elif orientation == 'x':
top = YPlane(y0=sqrt(3.)/2.*l)
@ -2089,109 +2204,152 @@ def get_hexagonal_prism(edge_length=1., orientation='y',
# y = sqrt(3)*(x + a)
upper_left = Plane(A=-c, B=1., D=c*l, boundary_type=boundary_type)
return Intersection(-top, +bottom, -upper_right, +lower_right,
+lower_left, -upper_left)
prism = -top & +bottom & -upper_right & +lower_right & \
+lower_left & -upper_left
if corner_radius > 0.:
if orientation == 'x':
c = sqrt(3.)
x_min_y_min_in = ZCylinder(name='x min y min in',
x0=-0.5 * (l - 2./c*corner_radius),
y0=-c/2. * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_min_y_max_in = ZCylinder(name='x min y max in',
x0= 0.5 * (l - 2./c*corner_radius),
y0=-c/2. * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_max_y_min_in = ZCylinder(name='x max y min in',
x0=-0.5 * (l - 2./c*corner_radius),
y0= c/2. * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_max_y_max_in = ZCylinder(name='x max y max in',
x0= 0.5 * (l - 2./c*corner_radius),
y0= c/2. * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_min_in = ZCylinder(name='x min in',
x0=-(l - 2./c*corner_radius),
y0=0.,
R=corner_radius,
boundary_type=boundary_type)
x_max_in = ZCylinder(name='x max in',
x0= (l - 2./c*corner_radius),
y0=0.,
R=corner_radius,
boundary_type=boundary_type)
def get_square_cylinder(r_cylinder, r_corner, axis='z', origin=(0., 0.),
boundary_type='transmission'):
"""Get an infinite square cylinder using a collection of planes and cylinders.
x_min_y_min_out = ZCylinder(name='x min y min out',
x0=-0.5 * l,
y0=-c/2. * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_min_y_max_out = ZCylinder(name='x min y max out',
x0=+0.5 * l,
y0=-c/2. * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_max_y_min_out = ZCylinder(name='x max y min out',
x0=-0.5 * l,
y0= c/2. * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_max_y_max_out = ZCylinder(name='x max y max out',
x0= 0.5 * l,
y0= c/2. * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_min_out = ZCylinder(name='x min out',
x0=-l,
y0=0.,
R=corner_radius/c,
boundary_type=boundary_type)
x_max_out = ZCylinder(name='x max out',
x0= l,
y0=0.,
R=corner_radius/c,
boundary_type=boundary_type)
Parameters
----------
r_cylinder: float
Square cylinder radius in units of cm. The square cylinder is aligned
with the x, y, or z axis.
r_corner: float
Square cylinder corner radius in units of cm.
axis : {'x', 'y', 'z'}
Axis with which the infinite length of the square cylinder should be
aligned.
Defaults to 'z'.
origin: Iterable of two floats
Origin of the square cylinder. The two floats correspond to (y,z), (x,z)
or (x,y) for square cylinders parallel to the x, y or z axis,
respectively.
Defaults to (0., 0.).
boundary_type : {'transmission, 'vacuum', 'reflective', 'periodic'}
Boundary condition that defines the behavior for particles hitting the
surfaces comprising the square cylinder (default is 'transmission').
corners = (+x_min_y_min_in & -x_min_y_min_out) | \
(+x_min_y_max_in & -x_min_y_max_out) | \
(+x_max_y_min_in & -x_max_y_min_out) | \
(+x_max_y_max_in & -x_max_y_max_out) | \
(+x_min_in & -x_min_out ) | \
(+x_max_in & -x_max_out )
Returns
-------
openmc.Region
The inside of a square cylinder
elif orientation == 'y':
c = sqrt(3.)
x_min_y_min_in = ZCylinder(name='x min y min in',
x0=-c/2. * (l - 2./c*corner_radius),
y0=-0.5 * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_min_y_max_in = ZCylinder(name='x min y max in',
x0=-c/2. * (l - 2./c*corner_radius),
y0= 0.5 * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_max_y_min_in = ZCylinder(name='x max y min in',
x0= c/2. * (l - 2./c*corner_radius),
y0=-0.5 * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
x_max_y_max_in = ZCylinder(name='x max y max in',
x0= c/2. * (l - 2./c*corner_radius),
y0= 0.5 * (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
y_min_in = ZCylinder(name='y min in',
x0=0.,
y0=-(l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
y_max_in = ZCylinder(name='y max in',
x0=0.,
y0= (l - 2./c*corner_radius),
R=corner_radius,
boundary_type=boundary_type)
"""
x_min_y_min_out = ZCylinder(name='x min y min out',
x0=-c/2. * l,
y0=-0.5 * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_min_y_max_out = ZCylinder(name='x min y max out',
x0=-c/2. * l,
y0=+0.5 * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_max_y_min_out = ZCylinder(name='x max y min out',
x0= c/2. * l,
y0=-0.5 * l,
R=corner_radius/c,
boundary_type=boundary_type)
x_max_y_max_out = ZCylinder(name='x max y max out',
x0= c/2. * l,
y0= 0.5 * l,
R=corner_radius/c,
boundary_type=boundary_type)
y_min_out = ZCylinder(name='y min out',
x0=0.,
y0=-l,
R=corner_radius/c,
boundary_type=boundary_type)
y_max_out = ZCylinder(name='y max out',
x0=0.,
y0= l,
R=corner_radius/c,
boundary_type=boundary_type)
check_type('r_cylinder', r_cylinder, Real)
check_type('r_corner', r_corner, Real)
check_value('axis', axis, ['x','y','z'])
check_type('origin', origin, Iterable, Real)
corners = (+x_min_y_min_in & -x_min_y_min_out) | \
(+x_min_y_max_in & -x_min_y_max_out) | \
(+x_max_y_min_in & -x_max_y_min_out) | \
(+x_max_y_max_in & -x_max_y_max_out) | \
(+y_min_in & -y_min_out ) | \
(+y_max_in & -y_max_out )
if axis == 'x':
y_min_z_min = XCylinder(name='y min z min', y0=-r_cylinder+r_corner,
z0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
y_min_z_max = XCylinder(name='y min z max', y0=-r_cylinder+r_corner,
z0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
y_max_z_min = XCylinder(name='y max z min', y0= r_cylinder-r_corner,
z0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
y_max_z_max = XCylinder(name='y max z max', y0= r_cylinder-r_corner,
z0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
corners = -y_min_z_min | -y_min_z_max | -y_max_z_min | -y_max_z_max
elif axis == 'y':
x_min_z_min = YCylinder(name='x min z min', x0=-r_cylinder+r_corner,
z0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
x_min_z_max = YCylinder(name='x min z max', x0=-r_cylinder+r_corner,
z0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
x_max_z_min = YCylinder(name='x max z min', x0= r_cylinder-r_corner,
z0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
x_max_z_max = YCylinder(name='x max z max', x0= r_cylinder-r_corner,
z0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
corners = -x_min_z_min | -x_min_z_max | -x_max_z_min | -x_max_z_max
elif axis == 'z':
x_min_y_min = ZCylinder(name='x min y min', x0=-r_cylinder+r_corner,
y0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
x_min_y_max = ZCylinder(name='x min y max', x0=-r_cylinder+r_corner,
y0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
x_max_y_min = ZCylinder(name='x max y min', x0= r_cylinder-r_corner,
y0=-r_cylinder+r_corner, R=r_corner,
boundary_type=boundary_type)
x_max_y_max = ZCylinder(name='x max y max', x0= r_cylinder-r_corner,
y0= r_cylinder-r_corner, R=r_corner,
boundary_type=boundary_type)
corners = -x_min_y_min | -x_min_y_max | -x_max_y_min | -x_max_y_max
prism = prism & ~corners
prism_1 = get_rectangular_prism(width=2*r_cylinder,
height=2*(r_cylinder - r_corner),
axis=axis,
origin=origin,
boundary_type=boundary_type)
origin_2 = list(origin)
origin_2[1] += r_cylinder - r_corner/2.
prism_2 = get_rectangular_prism(width=2*(r_cylinder - r_corner),
height=r_corner,
axis=axis,
origin=origin_2,
boundary_type=boundary_type)
origin_3 = list(origin)
origin_3[1] -= r_cylinder - r_corner/2.
prism_3 = get_rectangular_prism(width=2*(r_cylinder - r_corner),
height=r_corner,
axis=axis,
origin=origin_3,
boundary_type=boundary_type)
sqc = corners | prism_1 | prism_2 | prism_3
return sqc
return prism