Merge pull request #1188 from paulromano/translate-and-points

Add translate method for surfaces/regions and two convenience functions
This commit is contained in:
Adam Nelson 2019-03-08 14:52:54 -06:00 committed by GitHub
commit b2d6a06284
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
36 changed files with 749 additions and 212 deletions

View file

@ -53,9 +53,10 @@ extensions = ['sphinx.ext.autodoc',
'sphinx.ext.autosummary',
'sphinx.ext.intersphinx',
'sphinx.ext.viewcode',
'sphinx.ext.imgconverter',
'sphinx_numfig',
'notebook_sphinxext']
if not on_rtd:
extensions.append('sphinx.ext.imgconverter')
# Add any paths that contain templates here, relative to this directory.
templates_path = ['_templates']

View file

@ -11,6 +11,7 @@ Convenience Functions
:template: myfunction.rst
openmc.model.borated_water
openmc.model.cylinder_from_points
openmc.model.get_hexagonal_prism
openmc.model.get_rectangular_prism
openmc.model.subdivide

View file

@ -88,7 +88,7 @@ parameters for a sphere are the :math:`x,y,z` coordinates of the center of the
sphere and the radius of the sphere. All of these parameters can be set either
as optional keyword arguments to the class constructor or via attributes::
sphere = openmc.Sphere(R=10.0)
sphere = openmc.Sphere(r=10.0)
# This is equivalent
sphere = openmc.Sphere()
@ -98,7 +98,7 @@ Once a surface has been created, half-spaces can be obtained by applying the
unary ``-`` or ``+`` operators, corresponding to the negative and positive
half-spaces, respectively. For example::
>>> sphere = openmc.Sphere(R=10.0)
>>> sphere = openmc.Sphere(r=10.0)
>>> inside_sphere = -sphere
>>> outside_sphere = +sphere
>>> type(inside_sphere)
@ -140,10 +140,10 @@ may want to specify different behavior for particles passing through a
surface. To specify a vacuum boundary condition, simply change the
:attr:`Surface.boundary_type` attribute to 'vacuum'::
outer_surface = openmc.Sphere(R=100.0, boundary_type='vacuum')
outer_surface = openmc.Sphere(r=100.0, boundary_type='vacuum')
# This is equivalent
outer_surface = openmc.Sphere(R=100.0)
outer_surface = openmc.Sphere(r=100.0)
outer_surface.boundary_type = 'vacuum'
Reflective and periodic boundary conditions can be set with the strings
@ -154,8 +154,8 @@ can be determined automatically. For non-axis-aligned planes, it is necessary to
specify pairs explicitly using the :attr:`Surface.periodic_surface` attribute as
in the following example::
p1 = openmc.Plane(A=0.3, B=5.0, D=1.0, boundary_type='periodic')
p2 = openmc.Plane(A=0.3, B=5.0, D=-1.0, boundary_type='periodic')
p1 = openmc.Plane(a=0.3, b=5.0, d=1.0, boundary_type='periodic')
p2 = openmc.Plane(a=0.3, b=5.0, d=-1.0, boundary_type='periodic')
p1.periodic_surface = p2
Rotationally-periodic boundary conditions can be specified for a pair of

View file

@ -30,7 +30,7 @@ arguments are not necessary. For example,
::
sphere = openmc.Sphere(R=10.0)
sphere = openmc.Sphere(r=10.0)
cell = openm.Cell(region=-sphere)
vol_calc = openmc.VolumeCalculation([cell], 1000000)

View file

@ -36,9 +36,9 @@ materials_file.export_to_xml()
###############################################################################
# Instantiate ZCylinder surfaces
surf1 = openmc.ZCylinder(surface_id=1, x0=0, y0=0, R=7, name='surf 1')
surf2 = openmc.ZCylinder(surface_id=2, x0=0, y0=0, R=9, name='surf 2')
surf3 = openmc.ZCylinder(surface_id=3, x0=0, y0=0, R=11, name='surf 3')
surf1 = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=7, name='surf 1')
surf2 = openmc.ZCylinder(surface_id=2, x0=0, y0=0, r=9, name='surf 2')
surf3 = openmc.ZCylinder(surface_id=3, x0=0, y0=0, r=11, name='surf 3')
surf3.boundary_type = 'vacuum'
# Instantiate Cells

View file

@ -43,7 +43,7 @@ left = openmc.XPlane(surface_id=1, x0=-3, name='left')
right = openmc.XPlane(surface_id=2, x0=3, name='right')
bottom = openmc.YPlane(surface_id=3, y0=-4, name='bottom')
top = openmc.YPlane(surface_id=4, y0=4, name='top')
fuel_surf = openmc.ZCylinder(surface_id=5, x0=0, y0=0, R=0.4)
fuel_surf = openmc.ZCylinder(surface_id=5, x0=0, y0=0, r=0.4)
left.boundary_type = 'vacuum'
right.boundary_type = 'vacuum'

View file

@ -39,9 +39,9 @@ left = openmc.XPlane(surface_id=1, x0=-2, name='left')
right = openmc.XPlane(surface_id=2, x0=2, name='right')
bottom = openmc.YPlane(surface_id=3, y0=-2, name='bottom')
top = openmc.YPlane(surface_id=4, y0=2, name='top')
fuel1 = openmc.ZCylinder(surface_id=5, x0=0, y0=0, R=0.4)
fuel2 = openmc.ZCylinder(surface_id=6, x0=0, y0=0, R=0.3)
fuel3 = openmc.ZCylinder(surface_id=7, x0=0, y0=0, R=0.2)
fuel1 = openmc.ZCylinder(surface_id=5, x0=0, y0=0, r=0.4)
fuel2 = openmc.ZCylinder(surface_id=6, x0=0, y0=0, r=0.3)
fuel3 = openmc.ZCylinder(surface_id=7, x0=0, y0=0, r=0.2)
left.boundary_type = 'vacuum'
right.boundary_type = 'vacuum'

View file

@ -39,9 +39,9 @@ left = openmc.XPlane(surface_id=1, x0=-2, name='left')
right = openmc.XPlane(surface_id=2, x0=2, name='right')
bottom = openmc.YPlane(surface_id=3, y0=-2, name='bottom')
top = openmc.YPlane(surface_id=4, y0=2, name='top')
fuel1 = openmc.ZCylinder(surface_id=5, x0=0, y0=0, R=0.4)
fuel2 = openmc.ZCylinder(surface_id=6, x0=0, y0=0, R=0.3)
fuel3 = openmc.ZCylinder(surface_id=7, x0=0, y0=0, R=0.2)
fuel1 = openmc.ZCylinder(surface_id=5, x0=0, y0=0, r=0.4)
fuel2 = openmc.ZCylinder(surface_id=6, x0=0, y0=0, r=0.3)
fuel3 = openmc.ZCylinder(surface_id=7, x0=0, y0=0, r=0.2)
left.boundary_type = 'vacuum'
right.boundary_type = 'vacuum'

View file

@ -49,9 +49,9 @@ materials_file.export_to_xml()
###############################################################################
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, R=0.39218, name='Fuel OR')
clad_ir = openmc.ZCylinder(surface_id=2, x0=0, y0=0, R=0.40005, name='Clad IR')
clad_or = openmc.ZCylinder(surface_id=3, x0=0, y0=0, R=0.45720, name='Clad OR')
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=0.39218, name='Fuel OR')
clad_ir = openmc.ZCylinder(surface_id=2, x0=0, y0=0, r=0.40005, name='Clad IR')
clad_or = openmc.ZCylinder(surface_id=3, x0=0, y0=0, r=0.45720, name='Clad OR')
left = openmc.XPlane(surface_id=4, x0=-0.62992, name='left')
right = openmc.XPlane(surface_id=5, x0=0.62992, name='right')
bottom = openmc.YPlane(surface_id=6, y0=-0.62992, name='bottom')

View file

@ -54,9 +54,9 @@ borated_water.add_s_alpha_beta('c_H_in_H2O')
###############################################################################
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, R=0.39218, name='Fuel OR')
clad_ir = openmc.ZCylinder(surface_id=2, x0=0, y0=0, R=0.40005, name='Clad IR')
clad_or = openmc.ZCylinder(surface_id=3, x0=0, y0=0, R=0.45720, name='Clad OR')
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=0.39218, name='Fuel OR')
clad_ir = openmc.ZCylinder(surface_id=2, x0=0, y0=0, r=0.40005, name='Clad IR')
clad_or = openmc.ZCylinder(surface_id=3, x0=0, y0=0, r=0.45720, name='Clad OR')
left = openmc.XPlane(surface_id=4, x0=-0.62992, name='left')
right = openmc.XPlane(surface_id=5, x0=0.62992, name='right')
bottom = openmc.YPlane(surface_id=6, y0=-0.62992, name='bottom')

View file

@ -98,7 +98,7 @@ materials_file.export_to_xml()
###############################################################################
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, R=0.54, name='Fuel OR')
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=0.54, name='Fuel OR')
left = openmc.XPlane(surface_id=4, x0=-0.63, name='left')
right = openmc.XPlane(surface_id=5, x0=0.63, name='right')
bottom = openmc.YPlane(surface_id=6, y0=-0.63, name='bottom')

View file

@ -51,8 +51,8 @@ def pwr_pin_cell():
# Instantiate ZCylinder surfaces
pitch = 1.26
fuel_or = openmc.ZCylinder(x0=0, y0=0, R=0.39218, name='Fuel OR')
clad_or = openmc.ZCylinder(x0=0, y0=0, R=0.45720, name='Clad OR')
fuel_or = openmc.ZCylinder(x0=0, y0=0, r=0.39218, name='Fuel OR')
clad_or = openmc.ZCylinder(x0=0, y0=0, r=0.45720, name='Clad OR')
left = openmc.XPlane(x0=-pitch/2, name='left', boundary_type='reflective')
right = openmc.XPlane(x0=pitch/2, name='right', boundary_type='reflective')
bottom = openmc.YPlane(y0=-pitch/2, name='bottom',
@ -256,14 +256,14 @@ def pwr_core():
bot_nozzle, top_nozzle, top_fa, bot_fa)
# Define surfaces.
s1 = openmc.ZCylinder(R=0.41, surface_id=1)
s2 = openmc.ZCylinder(R=0.475, surface_id=2)
s3 = openmc.ZCylinder(R=0.56, surface_id=3)
s4 = openmc.ZCylinder(R=0.62, surface_id=4)
s5 = openmc.ZCylinder(R=187.6, surface_id=5)
s6 = openmc.ZCylinder(R=209.0, surface_id=6)
s7 = openmc.ZCylinder(R=229.0, surface_id=7)
s8 = openmc.ZCylinder(R=249.0, surface_id=8, boundary_type='vacuum')
s1 = openmc.ZCylinder(r=0.41, surface_id=1)
s2 = openmc.ZCylinder(r=0.475, surface_id=2)
s3 = openmc.ZCylinder(r=0.56, surface_id=3)
s4 = openmc.ZCylinder(r=0.62, surface_id=4)
s5 = openmc.ZCylinder(r=187.6, surface_id=5)
s6 = openmc.ZCylinder(r=209.0, surface_id=6)
s7 = openmc.ZCylinder(r=229.0, surface_id=7)
s8 = openmc.ZCylinder(r=249.0, surface_id=8, boundary_type='vacuum')
s31 = openmc.ZPlane(z0=-229.0, surface_id=31, boundary_type='vacuum')
s32 = openmc.ZPlane(z0=-199.0, surface_id=32)
@ -473,8 +473,8 @@ def pwr_assembly():
model.materials = (fuel, clad, hot_water)
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(x0=0, y0=0, R=0.39218, name='Fuel OR')
clad_or = openmc.ZCylinder(x0=0, y0=0, R=0.45720, name='Clad OR')
fuel_or = openmc.ZCylinder(x0=0, y0=0, r=0.39218, name='Fuel OR')
clad_or = openmc.ZCylinder(x0=0, y0=0, r=0.45720, name='Clad OR')
# Create boundary planes to surround the geometry
pitch = 21.42

View file

@ -3,7 +3,7 @@ from collections.abc import Iterable
from math import sqrt
from numbers import Real
from openmc import XPlane, YPlane, Plane, ZCylinder
from openmc import XPlane, YPlane, Plane, ZCylinder, Quadric
from openmc.checkvalue import check_type, check_value
import openmc.data
@ -237,16 +237,16 @@ def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.),
c = sqrt(3.)/3.
# y = -x/sqrt(3) + a
upper_right = Plane(A=c, B=1., D=l+x*c+y, boundary_type=boundary_type)
upper_right = Plane(a=c, b=1., d=l+x*c+y, boundary_type=boundary_type)
# y = x/sqrt(3) + a
upper_left = Plane(A=-c, B=1., D=l-x*c+y, boundary_type=boundary_type)
upper_left = Plane(a=-c, b=1., d=l-x*c+y, boundary_type=boundary_type)
# y = x/sqrt(3) - a
lower_right = Plane(A=-c, B=1., D=-l-x*c+y, boundary_type=boundary_type)
lower_right = Plane(a=-c, b=1., d=-l-x*c+y, boundary_type=boundary_type)
# y = -x/sqrt(3) - a
lower_left = Plane(A=c, B=1., D=-l+x*c+y, boundary_type=boundary_type)
lower_left = Plane(a=c, b=1., d=-l+x*c+y, boundary_type=boundary_type)
prism = -right & +left & -upper_right & -upper_left & \
+lower_right & +lower_left
@ -262,17 +262,17 @@ def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.),
c = sqrt(3.)
# y = -sqrt(3)*(x - a)
upper_right = Plane(A=c, B=1., D=c*l+x*c+y, boundary_type=boundary_type)
upper_right = Plane(a=c, b=1., d=c*l+x*c+y, boundary_type=boundary_type)
# y = sqrt(3)*(x + a)
lower_right = Plane(A=-c, B=1., D=-c*l-x*c+y,
lower_right = Plane(a=-c, b=1., d=-c*l-x*c+y,
boundary_type=boundary_type)
# y = -sqrt(3)*(x + a)
lower_left = Plane(A=c, B=1., D=-c*l+x*c+y, boundary_type=boundary_type)
lower_left = Plane(a=c, b=1., d=-c*l+x*c+y, boundary_type=boundary_type)
# y = sqrt(3)*(x + a)
upper_left = Plane(A=-c, B=1., D=c*l-x*c+y, boundary_type=boundary_type)
upper_left = Plane(a=-c, b=1., d=c*l-x*c+y, boundary_type=boundary_type)
prism = -top & +bottom & -upper_right & +lower_right & \
+lower_left & -upper_left
@ -292,8 +292,8 @@ def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.),
t = l - corner_radius/c
# Cylinder with corner radius and boundary type pre-applied
cyl1 = partial(ZCylinder, R=corner_radius, boundary_type=boundary_type)
cyl2 = partial(ZCylinder, R=corner_radius/(2*c),
cyl1 = partial(ZCylinder, r=corner_radius, boundary_type=boundary_type)
cyl2 = partial(ZCylinder, r=corner_radius/(2*c),
boundary_type=boundary_type)
if orientation == 'x':
@ -345,6 +345,54 @@ def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.),
return prism
def cylinder_from_points(p1, p2, r, **kwargs):
"""Return cylinder defined by two points passing through its center.
Parameters
----------
p1, p2 : 3-tuples
Coordinates of two points that pass through the center of the cylinder
r : float
Radius of the cylinder
kwargs : dict
Keyword arguments passed to the :class:`openmc.Quadric` constructor
Returns
-------
openmc.Quadric
Quadric surface representing the cylinder.
"""
# Get x, y, z coordinates of two points
x1, y1, z1 = p1
x2, y2, z2 = p2
# Define intermediate terms
dx = x2 - x1
dy = y2 - y1
dz = z2 - z1
cx = y1*z2 + y2*z1
cy = -(x1*z2 + x2*z1)
cz = x1*y2 + x2*y1
# Given p=(x,y,z), p1=(x1, y1, z1), p2=(x2, y2, z2), the equation for the
# cylinder can be derived as r = |(p - p1) (p - p2)| / |p2 - p1|.
# Expanding out all terms and grouping according to what Quadric expects
# gives the following coefficients.
kwargs['a'] = dy*dy + dz*dz
kwargs['b'] = dx*dx + dz*dz
kwargs['c'] = dx*dx + dy*dy
kwargs['d'] = -2*dx*dy
kwargs['e'] = -2*dy*dz
kwargs['f'] = -2*dx*dz
kwargs['g'] = cy*dz - cz*dy
kwargs['h'] = cz*dx - cx*dz
kwargs['j'] = cx*dy - cy*dx
kwargs['k'] = -(dx*dx + dy*dy + dz*dz)*r*r
return openmc.Quadric(**kwargs)
def subdivide(surfaces):
"""Create regions separated by a series of surfaces.

View file

@ -49,7 +49,7 @@ class TRISO(openmc.Cell):
"""
def __init__(self, outer_radius, fill, center=(0., 0., 0.)):
self._surface = openmc.Sphere(R=outer_radius)
self._surface = openmc.Sphere(r=outer_radius)
super().__init__(fill=fill, region=-self._surface)
self.center = np.asarray(center)

View file

@ -229,14 +229,28 @@ 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. This parameter is used internally
and should not be specified by the user.
Returns
-------
openmc.Region
Translated region
"""
raise NotImplementedError('The clone method is not implemented for '
'the abstract region class.')
pass
class Intersection(Region, MutableSequence):
@ -247,7 +261,7 @@ class Intersection(Region, MutableSequence):
following example:
>>> equator = openmc.ZPlane(z0=0.0)
>>> earth = openmc.Sphere(R=637.1e6)
>>> earth = openmc.Sphere(r=637.1e6)
>>> northern_hemisphere = -earth & +equator
>>> southern_hemisphere = -earth & -equator
>>> type(northern_hemisphere)
@ -352,6 +366,27 @@ 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. This parameter is used internally
and should not be specified by the user.
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.
@ -361,7 +396,7 @@ class Union(Region, MutableSequence):
example:
>>> s1 = openmc.ZPlane(z0=0.0)
>>> s2 = openmc.Sphere(R=637.1e6)
>>> s2 = openmc.Sphere(r=637.1e6)
>>> type(-s2 | +s1)
<class 'openmc.region.Union'>
@ -464,6 +499,27 @@ 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. This parameter is used internally
and should not be specified by the user.
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 +640,24 @@ 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. This parameter is used internally
and should not be specified by the user.
Returns
-------
openmc.Complement
Translated region
"""
if memo is None:
memo = {}
return type(self)(self.node.translate(vector, memo))

View file

@ -1,9 +1,10 @@
from abc import ABCMeta
from abc import ABCMeta, abstractmethod
from collections import OrderedDict
from copy import deepcopy
from functools import partial
from numbers import Real, Integral
from xml.etree import ElementTree as ET
from warnings import warn
import numpy as np
@ -14,8 +15,13 @@ from openmc.mixin import IDManagerMixin
_BOUNDARY_TYPES = ['transmission', 'vacuum', 'reflective', 'periodic']
_WARNING_UPPER = """\
"{}(...) accepts an argument named '{}', not '{}'. Future versions of OpenMC \
will not accept the capitalized version.\
"""
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 +148,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 +180,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
@ -280,29 +293,29 @@ class Surface(IDManagerMixin):
surface = Plane(surface_id, bc, A, B, C, D, name)
elif surf_type == 'x-cylinder':
y0, z0, R = coeffs
surface = XCylinder(surface_id, bc, y0, z0, R, name)
y0, z0, r = coeffs
surface = XCylinder(surface_id, bc, y0, z0, r, name)
elif surf_type == 'y-cylinder':
x0, z0, R = coeffs
surface = YCylinder(surface_id, bc, x0, z0, R, name)
x0, z0, r = coeffs
surface = YCylinder(surface_id, bc, x0, z0, r, name)
elif surf_type == 'z-cylinder':
x0, y0, R = coeffs
surface = ZCylinder(surface_id, bc, x0, y0, R, name)
x0, y0, r = coeffs
surface = ZCylinder(surface_id, bc, x0, y0, r, name)
elif surf_type == 'sphere':
x0, y0, z0, R = coeffs
surface = Sphere(surface_id, bc, x0, y0, z0, R, name)
x0, y0, z0, r = coeffs
surface = Sphere(surface_id, bc, x0, y0, z0, r, name)
elif surf_type in ['x-cone', 'y-cone', 'z-cone']:
x0, y0, z0, R2 = coeffs
x0, y0, z0, r2 = coeffs
if surf_type == 'x-cone':
surface = XCone(surface_id, bc, x0, y0, z0, R2, name)
surface = XCone(surface_id, bc, x0, y0, z0, r2, name)
elif surf_type == 'y-cone':
surface = YCone(surface_id, bc, x0, y0, z0, R2, name)
surface = YCone(surface_id, bc, x0, y0, z0, r2, name)
elif surf_type == 'z-cone':
surface = ZCone(surface_id, bc, x0, y0, z0, R2, name)
surface = ZCone(surface_id, bc, x0, y0, z0, r2, name)
elif surf_type == 'quadric':
a, b, c, d, e, f, g, h, j, k = coeffs
@ -324,13 +337,13 @@ class Plane(Surface):
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
A : float, optional
a : float, optional
The 'A' parameter for the plane. Defaults to 1.
B : float, optional
b : float, optional
The 'B' parameter for the plane. Defaults to 0.
C : float, optional
c : float, optional
The 'C' parameter for the plane. Defaults to 0.
D : float, optional
d : float, optional
The 'D' parameter for the plane. Defaults to 0.
name : str, optional
Name of the plane. If not specified, the name will be the empty string.
@ -363,56 +376,61 @@ class Plane(Surface):
"""
_type = 'plane'
_coeff_keys = ('A', 'B', 'C', 'D')
_coeff_keys = ('a', 'b', 'c', 'd')
def __init__(self, surface_id=None, boundary_type='transmission',
A=1., B=0., C=0., D=0., name=''):
a=1., b=0., c=0., d=0., name='', **kwargs):
super().__init__(surface_id, boundary_type, name=name)
self._periodic_surface = None
self.a = A
self.b = B
self.c = C
self.d = D
self.a = a
self.b = b
self.c = c
self.d = d
for k, v in kwargs.items():
if k in 'ABCD':
warn(_WARNING_UPPER.format(type(self).__name__, k.lower(), k),
FutureWarning)
setattr(self, k.lower(), v)
@property
def a(self):
return self.coefficients['A']
return self.coefficients['a']
@property
def b(self):
return self.coefficients['B']
return self.coefficients['b']
@property
def c(self):
return self.coefficients['C']
return self.coefficients['c']
@property
def d(self):
return self.coefficients['D']
return self.coefficients['d']
@property
def periodic_surface(self):
return self._periodic_surface
@a.setter
def a(self, A):
check_type('A coefficient', A, Real)
self._coefficients['A'] = A
def a(self, a):
check_type('A coefficient', a, Real)
self._coefficients['a'] = a
@b.setter
def b(self, B):
check_type('B coefficient', B, Real)
self._coefficients['B'] = B
def b(self, b):
check_type('B coefficient', b, Real)
self._coefficients['b'] = b
@c.setter
def c(self, C):
check_type('C coefficient', C, Real)
self._coefficients['C'] = C
def c(self, c):
check_type('C coefficient', c, Real)
self._coefficients['c'] = c
@d.setter
def d(self, D):
check_type('D coefficient', D, Real)
self._coefficients['D'] = D
def d(self, d):
check_type('D coefficient', d, Real)
self._coefficients['d'] = d
@periodic_surface.setter
def periodic_surface(self, periodic_surface):
@ -432,13 +450,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
@ -456,6 +495,38 @@ class Plane(Surface):
element.set("periodic_surface_id", str(self.periodic_surface.id))
return element
@classmethod
def from_points(cls, p1, p2, p3, **kwargs):
"""Return a plane given three points that pass through it.
Parameters
----------
p1, p2, p3 : 3-tuples
Points that pass through the plane
kwargs : dict
Keyword arguments passed to the :class:`Plane` constructor
Returns
-------
Plane
Plane that passes through the three points
"""
# Convert to numpy arrays
p1 = np.asarray(p1)
p2 = np.asarray(p2)
p3 = np.asarray(p3)
# Find normal vector to plane by taking cross product of two vectors
# connecting p1->p2 and p1->p3
n = np.cross(p2 - p1, p3 - p1)
# The equation of the plane will by n·(<x,y,z> - p1) = 0. Determine
# coefficients a, b, c, and d based on that
a, b, c = n
d = np.dot(n, p1)
return cls(a=a, b=b, c=c, d=d, **kwargs)
class XPlane(Plane):
"""A plane perpendicular to the x axis of the form :math:`x - x_0 = 0`
@ -561,6 +632,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`
@ -667,6 +758,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`
@ -773,8 +884,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
@ -786,7 +917,7 @@ class Cylinder(Surface, metaclass=ABCMeta):
Boundary condition that defines the behavior for particles hitting the
surface. Defaults to transmissive boundary condition where particles
freely pass through the surface.
R : float, optional
r : float, optional
Radius of the cylinder. Defaults to 1.
name : str, optional
Name of the cylinder. If not specified, the name will be the empty
@ -810,23 +941,23 @@ class Cylinder(Surface, metaclass=ABCMeta):
"""
def __init__(self, surface_id=None, boundary_type='transmission',
R=1., name=''):
r=1., name=''):
super().__init__(surface_id, boundary_type, name=name)
self.r = R
self.r = r
@property
def r(self):
return self.coefficients['R']
return self.coefficients['r']
@r.setter
def r(self, R):
check_type('R coefficient', R, Real)
self._coefficients['R'] = R
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
class XCylinder(Cylinder):
"""An infinite cylinder whose length is parallel to the x-axis of the form
:math:`(y - y_0)^2 + (z - z_0)^2 = R^2`.
:math:`(y - y_0)^2 + (z - z_0)^2 = r^2`.
Parameters
----------
@ -841,7 +972,7 @@ class XCylinder(Cylinder):
y-coordinate of the center of the cylinder. Defaults to 0.
z0 : float, optional
z-coordinate of the center of the cylinder. Defaults to 0.
R : float, optional
r : float, optional
Radius of the cylinder. Defaults to 0.
name : str, optional
Name of the cylinder. If not specified, the name will be the empty
@ -868,11 +999,14 @@ class XCylinder(Cylinder):
"""
_type = 'x-cylinder'
_coeff_keys = ('y0', 'z0', 'R')
_coeff_keys = ('y0', 'z0', 'r')
def __init__(self, surface_id=None, boundary_type='transmission',
y0=0., z0=0., R=1., name=''):
super().__init__(surface_id, boundary_type, R, name=name)
y0=0., z0=0., r=1., name='', *, R=None):
if R is not None:
warn(_WARNING_UPPER.format(type(self).__name__, 'r', 'R'), FutureWarning)
r = R
super().__init__(surface_id, boundary_type, r, name=name)
self.y0 = y0
self.z0 = z0
@ -938,17 +1072,39 @@ class XCylinder(Cylinder):
Returns
-------
float
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2`
: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
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
:math:`(x - x_0)^2 + (z - z_0)^2 = R^2`.
:math:`(x - x_0)^2 + (z - z_0)^2 = r^2`.
Parameters
----------
@ -963,7 +1119,7 @@ class YCylinder(Cylinder):
x-coordinate of the center of the cylinder. Defaults to 0.
z0 : float, optional
z-coordinate of the center of the cylinder. Defaults to 0.
R : float, optional
r : float, optional
Radius of the cylinder. Defaults to 1.
name : str, optional
Name of the cylinder. If not specified, the name will be the empty
@ -990,11 +1146,14 @@ class YCylinder(Cylinder):
"""
_type = 'y-cylinder'
_coeff_keys = ('x0', 'z0', 'R')
_coeff_keys = ('x0', 'z0', 'r')
def __init__(self, surface_id=None, boundary_type='transmission',
x0=0., z0=0., R=1., name=''):
super().__init__(surface_id, boundary_type, R, name=name)
x0=0., z0=0., r=1., name='', *, R=None):
if R is not None:
warn(_WARNING_UPPER.format(type(self).__name__, 'r', 'R'), FutureWarning)
r = R
super().__init__(surface_id, boundary_type, r, name=name)
self.x0 = x0
self.z0 = z0
@ -1060,17 +1219,39 @@ class YCylinder(Cylinder):
Returns
-------
float
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2`
: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
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
:math:`(x - x_0)^2 + (y - y_0)^2 = R^2`.
:math:`(x - x_0)^2 + (y - y_0)^2 = r^2`.
Parameters
----------
@ -1085,7 +1266,7 @@ class ZCylinder(Cylinder):
x-coordinate of the center of the cylinder. Defaults to 0.
y0 : float, optional
y-coordinate of the center of the cylinder. Defaults to 0.
R : float, optional
r : float, optional
Radius of the cylinder. Defaults to 1.
name : str, optional
Name of the cylinder. If not specified, the name will be the empty
@ -1112,11 +1293,14 @@ class ZCylinder(Cylinder):
"""
_type = 'z-cylinder'
_coeff_keys = ('x0', 'y0', 'R')
_coeff_keys = ('x0', 'y0', 'r')
def __init__(self, surface_id=None, boundary_type='transmission',
x0=0., y0=0., R=1., name=''):
super().__init__(surface_id, boundary_type, R, name=name)
x0=0., y0=0., r=1., name='', *, R=None):
if R is not None:
warn(_WARNING_UPPER.format(type(self).__name__, 'r', 'R'), FutureWarning)
r = R
super().__init__(surface_id, boundary_type, r, name=name)
self.x0 = x0
self.y0 = y0
@ -1182,16 +1366,38 @@ class ZCylinder(Cylinder):
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2`
: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
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`.
"""A sphere of the form :math:`(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 = r^2`.
Parameters
----------
@ -1208,7 +1414,7 @@ class Sphere(Surface):
y-coordinate of the center of the sphere. Defaults to 0.
z0 : float, optional
z-coordinate of the center of the sphere. Defaults to 0.
R : float, optional
r : float, optional
Radius of the sphere. Defaults to 1.
name : str, optional
Name of the sphere. If not specified, the name will be the empty string.
@ -1238,15 +1444,18 @@ class Sphere(Surface):
"""
_type = 'sphere'
_coeff_keys = ('x0', 'y0', 'z0', 'R')
_coeff_keys = ('x0', 'y0', 'z0', 'r')
def __init__(self, surface_id=None, boundary_type='transmission',
x0=0., y0=0., z0=0., R=1., name=''):
x0=0., y0=0., z0=0., r=1., name='', *, R=None):
if R is not None:
warn(_WARNING_UPPER.format(type(self).__name__, 'r', 'R'), FutureWarning)
r = R
super().__init__(surface_id, boundary_type, name=name)
self.x0 = x0
self.y0 = y0
self.z0 = z0
self.r = R
self.r = r
@property
def x0(self):
@ -1262,7 +1471,7 @@ class Sphere(Surface):
@property
def r(self):
return self.coefficients['R']
return self.coefficients['r']
@x0.setter
def x0(self, x0):
@ -1280,9 +1489,9 @@ class Sphere(Surface):
self._coefficients['z0'] = z0
@r.setter
def r(self, R):
check_type('R coefficient', R, Real)
self._coefficients['R'] = R
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
def bounding_box(self, side):
"""Determine an axis-aligned bounding box.
@ -1329,7 +1538,7 @@ class Sphere(Surface):
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 + (z' - z_0)^2 - R^2`
:math:`(x' - x_0)^2 + (y' - y_0)^2 + (z' - z_0)^2 - r^2`
"""
x = point[0] - self.x0
@ -1337,8 +1546,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
@ -1356,7 +1588,7 @@ class Cone(Surface, metaclass=ABCMeta):
y-coordinate of the apex. Defaults to 0.
z0 : float
z-coordinate of the apex. Defaults to 0.
R2 : float
r2 : float
Parameter related to the aperature. Defaults to 1.
name : str
Name of the cone. If not specified, the name will be the empty string.
@ -1385,15 +1617,18 @@ class Cone(Surface, metaclass=ABCMeta):
"""
_coeff_keys = ('x0', 'y0', 'z0', 'R2')
_coeff_keys = ('x0', 'y0', 'z0', 'r2')
def __init__(self, surface_id=None, boundary_type='transmission',
x0=0., y0=0., z0=0., R2=1., name=''):
x0=0., y0=0., z0=0., r2=1., name='', *, R2=None):
if R2 is not None:
warn(_WARNING_UPPER.format(type(self).__name__, 'r2', 'R2'), FutureWarning)
r2 = R2
super().__init__(surface_id, boundary_type, name=name)
self.x0 = x0
self.y0 = y0
self.z0 = z0
self.r2 = R2
self.r2 = r2
@property
def x0(self):
@ -1409,7 +1644,7 @@ class Cone(Surface, metaclass=ABCMeta):
@property
def r2(self):
return self.coefficients['R2']
return self.coefficients['r2']
@x0.setter
def x0(self, x0):
@ -1427,14 +1662,37 @@ class Cone(Surface, metaclass=ABCMeta):
self._coefficients['z0'] = z0
@r2.setter
def r2(self, R2):
check_type('R^2 coefficient', R2, Real)
self._coefficients['R2'] = R2
def r2(self, r2):
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 =
R^2 (x - x_0)^2`.
r^2 (x - x_0)^2`.
Parameters
----------
@ -1451,7 +1709,7 @@ class XCone(Cone):
y-coordinate of the apex. Defaults to 0.
z0 : float, optional
z-coordinate of the apex. Defaults to 0.
R2 : float, optional
r2 : float, optional
Parameter related to the aperature. Defaults to 1.
name : str, optional
Name of the cone. If not specified, the name will be the empty string.
@ -1464,7 +1722,7 @@ class XCone(Cone):
y-coordinate of the apex
z0 : float
z-coordinate of the apex
R2 : float
r2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
@ -1494,7 +1752,7 @@ class XCone(Cone):
Returns
-------
float
:math:`(y' - y_0)^2 + (z' - z_0)^2 - R^2(x' - x_0)^2`
:math:`(y' - y_0)^2 + (z' - z_0)^2 - r^2(x' - x_0)^2`
"""
x = point[0] - self.x0
@ -1505,7 +1763,7 @@ class XCone(Cone):
class YCone(Cone):
"""A cone parallel to the y-axis of the form :math:`(x - x_0)^2 + (z - z_0)^2 =
R^2 (y - y_0)^2`.
r^2 (y - y_0)^2`.
Parameters
----------
@ -1522,7 +1780,7 @@ class YCone(Cone):
y-coordinate of the apex. Defaults to 0.
z0 : float, optional
z-coordinate of the apex. Defaults to 0.
R2 : float, optional
r2 : float, optional
Parameter related to the aperature. Defaults to 1.
name : str, optional
Name of the cone. If not specified, the name will be the empty string.
@ -1535,7 +1793,7 @@ class YCone(Cone):
y-coordinate of the apex
z0 : float
z-coordinate of the apex
R2 : float
r2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
@ -1565,7 +1823,7 @@ class YCone(Cone):
Returns
-------
float
:math:`(x' - x_0)^2 + (z' - z_0)^2 - R^2(y' - y_0)^2`
:math:`(x' - x_0)^2 + (z' - z_0)^2 - r^2(y' - y_0)^2`
"""
x = point[0] - self.x0
@ -1576,7 +1834,7 @@ class YCone(Cone):
class ZCone(Cone):
"""A cone parallel to the x-axis of the form :math:`(x - x_0)^2 + (y - y_0)^2 =
R^2 (z - z_0)^2`.
r^2 (z - z_0)^2`.
Parameters
----------
@ -1593,7 +1851,7 @@ class ZCone(Cone):
y-coordinate of the apex. Defaults to 0.
z0 : float, optional
z-coordinate of the apex. Defaults to 0.
R2 : float, optional
r2 : float, optional
Parameter related to the aperature. Defaults to 1.
name : str, optional
Name of the cone. If not specified, the name will be the empty string.
@ -1606,7 +1864,7 @@ class ZCone(Cone):
y-coordinate of the apex
z0 : float
z-coordinate of the apex
R2 : float
r2 : float
Parameter related to the aperature
boundary_type : {'transmission, 'vacuum', 'reflective'}
Boundary condition that defines the behavior for particles hitting the
@ -1636,7 +1894,7 @@ class ZCone(Cone):
Returns
-------
float
:math:`(x' - x_0)^2 + (y' - y_0)^2 - R^2(z' - z_0)^2`
:math:`(x' - x_0)^2 + (y' - y_0)^2 - r^2(z' - z_0)^2`
"""
x = point[0] - self.x0
@ -1810,6 +2068,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.
@ -1824,7 +2106,7 @@ class Halfspace(Region):
can be created from an existing Surface through the __neg__ and __pos__
operators, as the following example demonstrates:
>>> sphere = openmc.Sphere(surface_id=1, R=10.0)
>>> sphere = openmc.Sphere(surface_id=1, r=10.0)
>>> inside_sphere = -sphere
>>> outside_sphere = +sphere
>>> type(inside_sphere)
@ -1955,3 +2237,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)

View file

@ -182,14 +182,14 @@ def segment_pin(n_rings, n_wedges, r_fuel, r_gap, r_clad):
theta = np.linspace(0, 2*math.pi, n_wedges + 1)
# Compute surfaces
fuel_rings = [openmc.ZCylinder(x0=0, y0=0, R=r_rings[i])
fuel_rings = [openmc.ZCylinder(x0=0, y0=0, r=r_rings[i])
for i in range(n_rings)]
fuel_wedges = [openmc.Plane(A=math.cos(theta[i]), B=math.sin(theta[i]))
fuel_wedges = [openmc.Plane(a=math.cos(theta[i]), b=math.sin(theta[i]))
for i in range(n_wedges)]
gap_ring = openmc.ZCylinder(x0=0, y0=0, R=r_gap)
clad_ring = openmc.ZCylinder(x0=0, y0=0, R=r_clad)
gap_ring = openmc.ZCylinder(x0=0, y0=0, r=r_gap)
clad_ring = openmc.ZCylinder(x0=0, y0=0, r=r_clad)
# Create cells
fuel_cells = []

View file

@ -32,7 +32,7 @@ class DistribmatTestHarness(PyAPITestHarness):
c1 = openmc.Cell(cell_id=1, fill=moderator)
mod_univ = openmc.Universe(universe_id=1, cells=[c1])
r0 = openmc.ZCylinder(R=0.3)
r0 = openmc.ZCylinder(r=0.3)
c11 = openmc.Cell(cell_id=11, region=-r0)
c11.fill = [dense_fuel, None, light_fuel, dense_fuel]
c12 = openmc.Cell(cell_id=12, region=+r0, fill=moderator)

View file

@ -37,7 +37,7 @@ def test_fixed_source():
mat.add_nuclide('U238', 0.0001)
mat.set_density('g/cc', 7.5)
surf = openmc.Sphere(R=10.0, boundary_type='vacuum')
surf = openmc.Sphere(r=10.0, boundary_type='vacuum')
cell = openmc.Cell(fill=mat, region=-surf)
model = openmc.model.Model()

View file

@ -27,7 +27,7 @@ def make_model():
c1 = openmc.Cell(cell_id=1, fill=moderator)
mod_univ = openmc.Universe(universe_id=1, cells=(c1,))
r0 = openmc.ZCylinder(R=0.3)
r0 = openmc.ZCylinder(r=0.3)
c11 = openmc.Cell(cell_id=11, fill=dense_fuel, region=-r0)
c11.temperature = [500, 700, 0, 800]
c12 = openmc.Cell(cell_id=12, fill=moderator, region=+r0)

View file

@ -30,7 +30,7 @@ class PeriodicTest(PyAPITestHarness):
z_min = openmc.ZPlane(5, z0=-5., boundary_type='reflective')
z_max = openmc.ZPlane(6, z0=5., boundary_type='reflective')
z_cyl = openmc.ZCylinder(7, x0=-2.5, y0=2.5, R=2.0)
z_cyl = openmc.ZCylinder(7, x0=-2.5, y0=2.5, r=2.0)
outside_cyl = openmc.Cell(1, fill=water, region=(
+x_min & -x_max & +y_min & -y_max & +z_min & -z_max & +z_cyl))

View file

@ -14,7 +14,7 @@ class SourceTestHarness(PyAPITestHarness):
materials = openmc.Materials([mat])
materials.export_to_xml()
cyl = openmc.XCylinder(boundary_type='vacuum', R=1.0e-6)
cyl = openmc.XCylinder(boundary_type='vacuum', r=1.0e-6)
x_plane_left = openmc.XPlane(boundary_type='vacuum', x0=-1.0)
x_plane_center = openmc.XPlane(boundary_type='transmission', x0=1.0)
x_plane_right = openmc.XPlane(boundary_type='vacuum', x0=1.0e9)

View file

@ -15,19 +15,19 @@ class SourceTestHarness(PyAPITestHarness):
materials = openmc.Materials([mat])
materials.export_to_xml()
sphere = openmc.Sphere(R=1.0e9, boundary_type='reflective')
sphere = openmc.Sphere(r=1.0e9, boundary_type='reflective')
inside_sphere = openmc.Cell()
inside_sphere.region = -sphere
inside_sphere.fill = mat
geometry = openmc.Geometry([inside_sphere])
geometry.export_to_xml()
source = openmc.Source()
source.space = openmc.stats.Point((0, 0, 0))
source.angle = openmc.stats.Isotropic()
source.energy = openmc.stats.Discrete([10.0e6], [1.0])
source.particle = 'photon'
settings = openmc.Settings()
settings.particles = 10000
settings.batches = 1
@ -37,7 +37,7 @@ class SourceTestHarness(PyAPITestHarness):
settings.run_mode = 'fixed source'
settings.source = source
settings.export_to_xml()
particle_filter = openmc.ParticleFilter('photon')
tally = openmc.Tally()
tally.filters = [particle_filter]

View file

@ -14,7 +14,7 @@ class SourceTestHarness(PyAPITestHarness):
materials = openmc.Materials([mat1])
materials.export_to_xml()
sphere = openmc.Sphere(surface_id=1, R=10.0, boundary_type='vacuum')
sphere = openmc.Sphere(surface_id=1, r=10.0, boundary_type='vacuum')
inside_sphere = openmc.Cell(cell_id=1)
inside_sphere.region = -sphere
inside_sphere.fill = mat1

View file

@ -25,8 +25,8 @@ class SurfaceTallyTestHarness(PyAPITestHarness):
materials_file.export_to_xml()
# Instantiate ZCylinder surfaces
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, R=1, \
name='Fuel OR')
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=1,
name='Fuel OR')
left = openmc.XPlane(surface_id=2, x0=-2, name='left')
right = openmc.XPlane(surface_id=3, x0=2, name='right')
bottom = openmc.YPlane(y0=-2, name='bottom')

View file

@ -44,7 +44,7 @@ class TRISOTestHarness(PyAPITestHarness):
graphite.add_s_alpha_beta('c_Graphite')
# Create TRISO particles
spheres = [openmc.Sphere(R=r*1e-4)
spheres = [openmc.Sphere(r=r*1e-4)
for r in [212.5, 312.5, 347.5, 382.5]]
c1 = openmc.Cell(fill=fuel, region=-spheres[0])
c2 = openmc.Cell(fill=porous_carbon, region=+spheres[0] & -spheres[1])

View file

@ -25,10 +25,10 @@ class VolumeTest(PyAPITestHarness):
materials = openmc.Materials((water, fuel))
materials.export_to_xml()
cyl = openmc.ZCylinder(1, R=1.0, boundary_type='vacuum')
top_sphere = openmc.Sphere(2, z0=5., R=1., boundary_type='vacuum')
cyl = openmc.ZCylinder(1, r=1.0, boundary_type='vacuum')
top_sphere = openmc.Sphere(2, z0=5., r=1., boundary_type='vacuum')
top_plane = openmc.ZPlane(3, z0=5.)
bottom_sphere = openmc.Sphere(4, z0=-5., R=1., boundary_type='vacuum')
bottom_sphere = openmc.Sphere(4, z0=-5., r=1., boundary_type='vacuum')
bottom_plane = openmc.ZPlane(5, z0=-5.)
# Define geometry

View file

@ -46,7 +46,7 @@ def cell_with_lattice():
m_inside = [openmc.Material(), openmc.Material(), None, openmc.Material()]
m_outside = openmc.Material()
cyl = openmc.ZCylinder(R=1.0)
cyl = openmc.ZCylinder(r=1.0)
inside_cyl = openmc.Cell(fill=m_inside, region=-cyl)
outside_cyl = openmc.Cell(fill=m_outside, region=+cyl)
univ = openmc.Universe(cells=[inside_cyl, outside_cyl])
@ -63,7 +63,7 @@ def cell_with_lattice():
@pytest.fixture
def mixed_lattice_model(uo2, water):
cyl = openmc.ZCylinder(R=0.4)
cyl = openmc.ZCylinder(r=0.4)
c1 = openmc.Cell(fill=uo2, region=-cyl)
c1.temperature = 600.0
c2 = openmc.Cell(fill=water, region=+cyl)

View file

@ -10,8 +10,8 @@ def test_volume(run_in_tmpdir, uo2):
# Create model with nested spheres
model = openmc.model.Model()
model.materials.append(uo2)
inner = openmc.Sphere(R=1.)
outer = openmc.Sphere(R=2., boundary_type='vacuum')
inner = openmc.Sphere(r=1.)
outer = openmc.Sphere(r=2., boundary_type='vacuum')
c1 = openmc.Cell(fill=uo2, region=-inner)
c2 = openmc.Cell(region=+inner & -outer)
u = openmc.Universe(cells=[c1, c2])
@ -44,8 +44,8 @@ def test_volume(run_in_tmpdir, uo2):
def test_export_xml(run_in_tmpdir, uo2):
s1 = openmc.Sphere(R=1.)
s2 = openmc.Sphere(R=2., boundary_type='reflective')
s1 = openmc.Sphere(r=1.)
s2 = openmc.Sphere(r=2., boundary_type='reflective')
c1 = openmc.Cell(fill=uo2, region=-s1)
c2 = openmc.Cell(fill=uo2, region=+s1 & -s2)
geom = openmc.Geometry([c1, c2])

View file

@ -7,7 +7,7 @@ import pytest
@pytest.fixture(scope='module')
def pincell1(uo2, water):
cyl = openmc.ZCylinder(R=0.35)
cyl = openmc.ZCylinder(r=0.35)
fuel = openmc.Cell(fill=uo2, region=-cyl)
moderator = openmc.Cell(fill=water, region=+cyl)
@ -19,7 +19,7 @@ def pincell1(uo2, water):
@pytest.fixture(scope='module')
def pincell2(uo2, water):
cyl = openmc.ZCylinder(R=0.4)
cyl = openmc.ZCylinder(r=0.4)
fuel = openmc.Cell(fill=uo2, region=-cyl)
moderator = openmc.Cell(fill=water, region=+cyl)

View file

@ -5,7 +5,7 @@ import pytest
@pytest.fixture(scope='module')
def pincell1(uo2, water):
cyl = openmc.ZCylinder(R=0.35)
cyl = openmc.ZCylinder(r=0.35)
fuel = openmc.Cell(fill=uo2, region=-cyl)
moderator = openmc.Cell(fill=water, region=+cyl)
@ -17,7 +17,7 @@ def pincell1(uo2, water):
@pytest.fixture(scope='module')
def pincell2(uo2, water):
cyl = openmc.ZCylinder(R=0.4)
cyl = openmc.ZCylinder(r=0.4)
fuel = openmc.Cell(fill=uo2, region=-cyl)
moderator = openmc.Cell(fill=water, region=+cyl)
@ -89,12 +89,12 @@ def rlat3(pincell1, pincell2, uo2, water, zr):
def test_mesh2d(rlat2):
shape = np.array(rlat2.shape)
width = shape*rlat2.pitch
mesh1 = openmc.Mesh.from_rect_lattice(rlat2)
assert np.array_equal(mesh1.dimension, (3, 3))
assert np.array_equal(mesh1.lower_left, rlat2.lower_left)
assert np.array_equal(mesh1.upper_right, rlat2.lower_left + width)
mesh2 = openmc.Mesh.from_rect_lattice(rlat2, division=3)
assert np.array_equal(mesh2.dimension, (9, 9))
assert np.array_equal(mesh2.lower_left, rlat2.lower_left)

View file

@ -47,7 +47,7 @@ def centers_rectangular_prism():
@pytest.fixture(scope='module')
def centers_x_cylinder():
cylinder = openmc.XCylinder(R=1, y0=1, z0=2)
cylinder = openmc.XCylinder(r=1, y0=1, z0=2)
min_x = openmc.XPlane(x0=0)
max_x = openmc.XPlane(x0=1)
region = +min_x & -max_x & -cylinder
@ -57,7 +57,7 @@ def centers_x_cylinder():
@pytest.fixture(scope='module')
def centers_y_cylinder():
cylinder = openmc.YCylinder(R=1, x0=1, z0=2)
cylinder = openmc.YCylinder(r=1, x0=1, z0=2)
min_y = openmc.YPlane(y0=0)
max_y = openmc.YPlane(y0=1)
region = +min_y & -max_y & -cylinder
@ -67,7 +67,7 @@ def centers_y_cylinder():
@pytest.fixture(scope='module')
def centers_z_cylinder():
cylinder = openmc.ZCylinder(R=1, x0=1, y0=2)
cylinder = openmc.ZCylinder(r=1, x0=1, y0=2)
min_z = openmc.ZPlane(z0=0)
max_z = openmc.ZPlane(z0=1)
region = +min_z & -max_z & -cylinder
@ -77,7 +77,7 @@ def centers_z_cylinder():
@pytest.fixture(scope='module')
def centers_sphere():
sphere = openmc.Sphere(R=1, x0=1, y0=2, z0=3)
sphere = openmc.Sphere(r=1, x0=1, y0=2, z0=3)
region = -sphere
return openmc.model.pack_spheres(radius=_RADIUS, region=region,
pf=_PACKING_FRACTION, initial_pf=0.2)
@ -85,8 +85,8 @@ def centers_sphere():
@pytest.fixture(scope='module')
def centers_spherical_shell():
sphere = openmc.Sphere(R=1, x0=1, y0=2, z0=3)
inner_sphere = openmc.Sphere(R=0.5, x0=1, y0=2, z0=3)
sphere = openmc.Sphere(r=1, x0=1, y0=2, z0=3)
inner_sphere = openmc.Sphere(r=0.5, x0=1, y0=2, z0=3)
region = -sphere & +inner_sphere
return openmc.model.pack_spheres(radius=_RADIUS, region=region,
pf=_PACKING_FRACTION, initial_pf=0.2)
@ -94,7 +94,7 @@ def centers_spherical_shell():
@pytest.fixture(scope='module')
def triso_universe():
sphere = openmc.Sphere(R=_RADIUS)
sphere = openmc.Sphere(r=_RADIUS)
cell = openmc.Cell(region=-sphere)
univ = openmc.Universe(cells=[cell])
return univ
@ -179,10 +179,10 @@ def test_packing_fraction(container, centers):
def test_num_spheres():
"""Check that the function returns the correct number of spheres"""
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=-openmc.Sphere(R=1), num_spheres=50
radius=_RADIUS, region=-openmc.Sphere(r=1), num_spheres=50
)
assert len(centers) == 50
def test_triso_lattice(triso_universe, centers_rectangular_prism):
trisos = [openmc.model.TRISO(_RADIUS, triso_universe, c)
@ -203,7 +203,7 @@ def test_container_input(triso_universe):
# Invalid container shape
with pytest.raises(ValueError):
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=+openmc.Sphere(R=1), num_spheres=100
radius=_RADIUS, region=+openmc.Sphere(r=1), num_spheres=100
)
@ -211,17 +211,17 @@ def test_packing_fraction_input():
# Provide neither packing fraction nor number of spheres
with pytest.raises(ValueError):
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=-openmc.Sphere(R=1)
radius=_RADIUS, region=-openmc.Sphere(r=1)
)
# Specify a packing fraction that is too high for CRP
with pytest.raises(ValueError):
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=-openmc.Sphere(R=1), pf=1
radius=_RADIUS, region=-openmc.Sphere(r=1), pf=1
)
# Specify a packing fraction that is too high for RSP
with pytest.raises(ValueError):
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=-openmc.Sphere(R=1), pf=0.5, initial_pf=0.4
radius=_RADIUS, region=-openmc.Sphere(r=1), pf=0.5, initial_pf=0.4
)

View file

@ -45,7 +45,7 @@ def test_repr(myplot):
def test_from_geometry():
width = 25.
s = openmc.Sphere(R=width/2, boundary_type='vacuum')
s = openmc.Sphere(r=width/2, boundary_type='vacuum')
c = openmc.Cell(region=-s)
univ = openmc.Universe(cells=[c])
geom = openmc.Geometry(univ)

View file

@ -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,9 +67,15 @@ 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.)
zcyl = openmc.ZCylinder(surface_id=1, r=1.)
z0 = openmc.ZPlane(surface_id=2, z0=-5.)
z1 = openmc.ZPlane(surface_id=3, z0=5.)
outside = +zcyl | -z0 | +z1
@ -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()

View file

@ -13,7 +13,7 @@ def assert_infinite_bb(s):
def test_plane():
s = openmc.Plane(A=1, B=2, C=-1, D=3, name='my plane')
s = openmc.Plane(a=1, b=2, c=-1, d=3, name='my plane')
assert s.a == 1
assert s.b == 2
assert s.c == -1
@ -29,10 +29,28 @@ 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)
def test_plane_from_points():
# Generate the plane x - y = 1 given three points
p1 = (0, -1, 0)
p2 = (1, 0, 0)
p3 = (1, 0, 1)
s = openmc.Plane.from_points(p1, p2, p3)
# Confirm correct coefficients
assert s.a == 1.0
assert s.b == -1.0
assert s.c == 0.0
assert s.d == 1.0
def test_xplane():
s = openmc.XPlane(x0=3., boundary_type='reflective')
assert s.x0 == 3.
@ -55,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)
@ -80,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.)
@ -102,13 +128,17 @@ 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)
def test_xcylinder():
y, z, r = 3, 5, 2
s = openmc.XCylinder(y0=y, z0=z, R=r)
s = openmc.XCylinder(y0=y, z0=z, r=r)
assert s.y0 == y
assert s.z0 == z
assert s.r == r
@ -124,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)
@ -139,7 +175,7 @@ def test_periodic():
def test_ycylinder():
x, z, r = 3, 5, 2
s = openmc.YCylinder(x0=x, z0=z, R=r)
s = openmc.YCylinder(x0=x, z0=z, r=r)
assert s.x0 == x
assert s.z0 == z
assert s.r == r
@ -155,10 +191,16 @@ 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
s = openmc.ZCylinder(x0=x, y0=y, R=r)
s = openmc.ZCylinder(x0=x, y0=y, r=r)
assert s.x0 == x
assert s.y0 == y
assert s.r == r
@ -174,13 +216,19 @@ 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)
def test_sphere():
x, y, z, r = -3, 5, 6, 2
s = openmc.Sphere(x0=x, y0=y, z0=z, R=r)
s = openmc.Sphere(x0=x, y0=y, z0=z, r=r)
assert s.x0 == x
assert s.y0 == y
assert s.z0 == z
@ -197,13 +245,20 @@ 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)
def cone_common(apex, r2, cls):
x, y, z = apex
s = cls(x0=x, y0=y, z0=z, R2=r2)
s = cls(x0=x, y0=y, z0=z, r2=r2)
assert s.x0 == x
assert s.y0 == y
assert s.z0 == z
@ -215,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)
@ -257,3 +319,24 @@ def test_quadric():
# evaluate method
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
p1 = (0, 0, 0)
p2 = (1, 1, 0)
s = openmc.model.cylinder_from_points(p1, p2, 1)
# Points p1 and p2 need to be inside cylinder
assert p1 in -s
assert p2 in -s
assert (-1, 1, 0) in +s
assert (1, -1, 0) in +s
assert (0, 0, 1.5) in +s

View file

@ -34,8 +34,8 @@ def test_basic():
def test_bounding_box():
cyl1 = openmc.ZCylinder(R=1.0)
cyl2 = openmc.ZCylinder(R=2.0)
cyl1 = openmc.ZCylinder(r=1.0)
cyl2 = openmc.ZCylinder(r=2.0)
c1 = openmc.Cell(region=-cyl1)
c2 = openmc.Cell(region=+cyl1 & -cyl2)