mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-26 21:25:36 -04:00
Merge pull request #1216 from dryuri92/newhex
Support different orientations for hexagonal lattices
This commit is contained in:
Sterling Harper
GitHub
commit
2c0b16e73d
12 changed files with 1034 additions and 139 deletions
1
.gitignore
vendored
1
.gitignore
vendored
|
|
@ -9,6 +9,7 @@
|
|||
*.pyc
|
||||
|
||||
# Python distribution
|
||||
.settings/
|
||||
dist/
|
||||
openmc.egg-info/
|
||||
|
||||
|
|
|
|||
|
|
@ -318,6 +318,13 @@ the following attributes or sub-elements:
|
|||
|
||||
*Default*: None
|
||||
|
||||
:orientation:
|
||||
The orientation of the hexagonal lattice. The string "x" indicates that two
|
||||
sides of the lattice are parallel to the x-axis, whereas the string "y"
|
||||
indicates that two sides are parallel to the y-axis.
|
||||
|
||||
*Default*: "y"
|
||||
|
||||
:center:
|
||||
The coordinates of the center of the lattice. If the lattice does not have
|
||||
axial sections then only the x- and y-coordinates are specified.
|
||||
|
|
|
|||
|
|
@ -27,6 +27,7 @@ enum class LatticeType {
|
|||
rect, hex
|
||||
};
|
||||
|
||||
|
||||
//==============================================================================
|
||||
// Global variables
|
||||
//==============================================================================
|
||||
|
|
@ -265,8 +266,20 @@ public:
|
|||
void to_hdf5_inner(hid_t group_id) const;
|
||||
|
||||
private:
|
||||
enum class Orientation {
|
||||
y, //!< Flat side of lattice parallel to y-axis
|
||||
x //!< Flat side of lattice parallel to x-axis
|
||||
};
|
||||
|
||||
//! Fill universes_ vector for 'y' orientation
|
||||
void fill_lattice_y(const std::vector<std::string>& univ_words);
|
||||
|
||||
//! Fill universes_ vector for 'x' orientation
|
||||
void fill_lattice_x(const std::vector<std::string>& univ_words);
|
||||
|
||||
int n_rings_; //!< Number of radial tile positions
|
||||
int n_axial_; //!< Number of axial tile positions
|
||||
Orientation orientation_; //!< Orientation of lattice
|
||||
Position center_; //!< Global center of lattice
|
||||
std::array<double, 2> pitch_; //!< Lattice tile width and height
|
||||
};
|
||||
|
|
|
|||
|
|
@ -36,7 +36,7 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
|
|||
outer : openmc.Universe
|
||||
A universe to fill all space outside the lattice
|
||||
universes : Iterable of Iterable of openmc.Universe
|
||||
A two- or three-dimensional list/array of universes filling each element
|
||||
A two-or three-dimensional list/array of universes filling each element
|
||||
of the lattice
|
||||
|
||||
"""
|
||||
|
|
@ -141,6 +141,10 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
|
|||
center = group['center'][()]
|
||||
pitch = group['pitch'][()]
|
||||
outer = group['outer'][()]
|
||||
if 'orientation' in group:
|
||||
orientation = group['orientation'][()].decode()
|
||||
else:
|
||||
orientation = "y"
|
||||
|
||||
universe_ids = group['universes'][()]
|
||||
|
||||
|
|
@ -148,67 +152,124 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
|
|||
lattice = openmc.HexLattice(lattice_id, name)
|
||||
lattice.center = center
|
||||
lattice.pitch = pitch
|
||||
|
||||
lattice.orientation = orientation
|
||||
# If the Universe specified outer the Lattice is not void
|
||||
if outer >= 0:
|
||||
lattice.outer = universes[outer]
|
||||
if orientation == "y":
|
||||
# Build array of Universe pointers for the Lattice. Note that
|
||||
# we need to convert between the HDF5's square array of
|
||||
# (x, alpha, z) to the Python API's format of a ragged nested
|
||||
# list of (z, ring, theta).
|
||||
uarray = []
|
||||
for z in range(n_axial):
|
||||
# Add a list for this axial level.
|
||||
uarray.append([])
|
||||
x = n_rings - 1
|
||||
a = 2*n_rings - 2
|
||||
for r in range(n_rings - 1, 0, -1):
|
||||
# Add a list for this ring.
|
||||
uarray[-1].append([])
|
||||
|
||||
# Build array of Universe pointers for the Lattice. Note that
|
||||
# we need to convert between the HDF5's square array of
|
||||
# (x, alpha, z) to the Python API's format of a ragged nested
|
||||
# list of (z, ring, theta).
|
||||
uarray = []
|
||||
for z in range(n_axial):
|
||||
# Add a list for this axial level.
|
||||
uarray.append([])
|
||||
x = n_rings - 1
|
||||
a = 2*n_rings - 2
|
||||
for r in range(n_rings - 1, 0, -1):
|
||||
# Add a list for this ring.
|
||||
uarray[-1].append([])
|
||||
# Climb down the top-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x += 1
|
||||
a -= 1
|
||||
|
||||
# Climb down the top-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x += 1
|
||||
# Climb down the right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
a -= 1
|
||||
|
||||
# Climb down the bottom-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x -= 1
|
||||
|
||||
# Climb up the bottom-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x -= 1
|
||||
a += 1
|
||||
|
||||
# Climb up the left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
a += 1
|
||||
|
||||
# Climb up the top-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x += 1
|
||||
|
||||
# Move down to the next ring.
|
||||
a -= 1
|
||||
|
||||
# Climb down the right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
# Convert the ids into Universe objects.
|
||||
uarray[-1][-1] = [universes[u_id]
|
||||
for u_id in uarray[-1][-1]]
|
||||
|
||||
# Handle the degenerate center ring separately.
|
||||
u_id = universe_ids[z, a, x]
|
||||
uarray[-1].append([universes[u_id]])
|
||||
else:
|
||||
# Build array of Universe pointers for the Lattice. Note that
|
||||
# we need to convert between the HDF5's square array of
|
||||
# (alpha, y, z) to the Python API's format of a ragged nested
|
||||
# list of (z, ring, theta).
|
||||
uarray = []
|
||||
for z in range(n_axial):
|
||||
# Add a list for this axial level.
|
||||
uarray.append([])
|
||||
a = 2*n_rings - 2
|
||||
y = n_rings - 1
|
||||
for r in range(n_rings - 1, 0, -1):
|
||||
# Add a list for this ring.
|
||||
uarray[-1].append([])
|
||||
|
||||
# Climb down the bottom-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
y -= 1
|
||||
|
||||
# Climb across the bottom.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
a -= 1
|
||||
|
||||
# Climb up the bottom-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
a -= 1
|
||||
y +=1
|
||||
|
||||
# Climb up the top-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
y += 1
|
||||
|
||||
# Climb across the top.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
a += 1
|
||||
|
||||
# Climb down the top-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, y, a])
|
||||
a += 1
|
||||
y -= 1
|
||||
|
||||
# Move down to the next ring.
|
||||
a -= 1
|
||||
|
||||
# Climb down the bottom-right.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x -= 1
|
||||
# Convert the ids into Universe objects.
|
||||
uarray[-1][-1] = [universes[u_id]
|
||||
for u_id in uarray[-1][-1]]
|
||||
|
||||
# Climb up the bottom-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x -= 1
|
||||
a += 1
|
||||
|
||||
# Climb up the left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
a += 1
|
||||
|
||||
# Climb up the top-left.
|
||||
for i in range(r):
|
||||
uarray[-1][-1].append(universe_ids[z, a, x])
|
||||
x += 1
|
||||
|
||||
# Move down to the next ring.
|
||||
a -= 1
|
||||
|
||||
# Convert the ids into Universe objects.
|
||||
uarray[-1][-1] = [universes[u_id]
|
||||
for u_id in uarray[-1][-1]]
|
||||
|
||||
# Handle the degenerate center ring separately.
|
||||
u_id = universe_ids[z, a, x]
|
||||
uarray[-1].append([universes[u_id]])
|
||||
# Handle the degenerate center ring separately.
|
||||
u_id = universe_ids[z, y, a]
|
||||
uarray[-1].append([universes[u_id]])
|
||||
|
||||
# Add the universes to the lattice.
|
||||
if len(pitch) == 2:
|
||||
|
|
@ -346,7 +407,9 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
|
|||
Lattice element indices. For a rectangular lattice, the indices are
|
||||
given in the :math:`(x,y)` or :math:`(x,y,z)` coordinate system. For
|
||||
hexagonal lattices, they are given in the :math:`x,\alpha` or
|
||||
:math:`x,\alpha,z` coordinate systems.
|
||||
:math:`x,\alpha,z` coordinate systems for "y" orientations and
|
||||
:math:`\alpha,y` or :math:`\alpha,y,z` coordinate systems for "x"
|
||||
orientations.
|
||||
|
||||
Returns
|
||||
-------
|
||||
|
|
@ -646,7 +709,8 @@ class RectLattice(Lattice):
|
|||
return (x, y, z)
|
||||
|
||||
def get_universe_index(self, idx):
|
||||
"""Return index in the universes array corresponding to a lattice element index
|
||||
"""Return index in the universes array corresponding
|
||||
to a lattice element index
|
||||
|
||||
Parameters
|
||||
----------
|
||||
|
|
@ -790,7 +854,8 @@ class RectLattice(Lattice):
|
|||
lat_id = int(get_text(elem, 'id'))
|
||||
name = get_text(elem, 'name')
|
||||
lat = cls(lat_id, name)
|
||||
lat.lower_left = [float(i) for i in get_text(elem, 'lower_left').split()]
|
||||
lat.lower_left = [float(i)
|
||||
for i in get_text(elem, 'lower_left').split()]
|
||||
lat.pitch = [float(i) for i in get_text(elem, 'pitch').split()]
|
||||
outer = get_text(elem, 'outer')
|
||||
if outer is not None:
|
||||
|
|
@ -800,7 +865,7 @@ class RectLattice(Lattice):
|
|||
dimension = get_text(elem, 'dimension').split()
|
||||
shape = np.array(dimension, dtype=int)[::-1]
|
||||
uarray = np.array([get_universe(int(i)) for i in
|
||||
get_text(elem, 'universes').split()])
|
||||
get_text(elem, 'universes').split()])
|
||||
uarray.shape = shape
|
||||
lat.universes = uarray
|
||||
return lat
|
||||
|
|
@ -815,10 +880,15 @@ class HexLattice(Lattice):
|
|||
|
||||
Most methods for this class use a natural indexing scheme wherein elements
|
||||
are assigned an index corresponding to their position relative to skewed
|
||||
:math:`(x,\alpha,z)` axes as described fully in
|
||||
:ref:`hexagonal_indexing`. However, note that when universes are assigned to
|
||||
lattice elements using the :attr:`HexLattice.universes` property, the array
|
||||
indices do not correspond to natural indices.
|
||||
:math:`(x,\alpha,z)` or :math:`(\alpha,y,z)` bases, depending on the lattice
|
||||
orientation, as described fully in :ref:`hexagonal_indexing`. However, note
|
||||
that when universes are assigned to lattice elements using the
|
||||
:attr:`HexLattice.universes` property, the array indices do not correspond
|
||||
to natural indices.
|
||||
|
||||
.. versionchanged:: 0.11
|
||||
The orientation of the lattice can now be changed with the
|
||||
:attr:`orientation` attribute.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
|
|
@ -855,6 +925,9 @@ class HexLattice(Lattice):
|
|||
possible, where z is the axial index, r is in the ring index (starting
|
||||
from the outermost ring), and i is the index with a ring starting from
|
||||
the top and proceeding clockwise.
|
||||
orientation : {'x', 'y'}
|
||||
str by default 'y' orientation of main lattice diagonal another option
|
||||
- 'x'
|
||||
num_rings : int
|
||||
Number of radial ring positions in the xy-plane
|
||||
num_axial : int
|
||||
|
|
@ -869,11 +942,14 @@ class HexLattice(Lattice):
|
|||
self._num_rings = None
|
||||
self._num_axial = None
|
||||
self._center = None
|
||||
self._orientation = 'y'
|
||||
|
||||
def __repr__(self):
|
||||
string = 'HexLattice\n'
|
||||
string += '{0: <16}{1}{2}\n'.format('\tID', '=\t', self._id)
|
||||
string += '{0: <16}{1}{2}\n'.format('\tName', '=\t', self._name)
|
||||
string += '{0: <16}{1}{2}\n'.format('\tOrientation', '=\t',
|
||||
self._orientation)
|
||||
string += '{0: <16}{1}{2}\n'.format('\t# Rings', '=\t', self._num_rings)
|
||||
string += '{0: <16}{1}{2}\n'.format('\t# Axial', '=\t', self._num_axial)
|
||||
string += '{0: <16}{1}{2}\n'.format('\tCenter', '=\t',
|
||||
|
|
@ -902,6 +978,10 @@ class HexLattice(Lattice):
|
|||
def num_rings(self):
|
||||
return self._num_rings
|
||||
|
||||
@property
|
||||
def orientation(self):
|
||||
return self._orientation
|
||||
|
||||
@property
|
||||
def num_axial(self):
|
||||
return self._num_axial
|
||||
|
|
@ -955,6 +1035,11 @@ class HexLattice(Lattice):
|
|||
cv.check_length('lattice center', center, 2, 3)
|
||||
self._center = center
|
||||
|
||||
@orientation.setter
|
||||
def orientation(self, orientation):
|
||||
cv.check_value('orientation', orientation.lower(), ('x', 'y'))
|
||||
self._orientation = orientation.lower()
|
||||
|
||||
@Lattice.pitch.setter
|
||||
def pitch(self, pitch):
|
||||
cv.check_type('lattice pitch', pitch, Iterable, Real)
|
||||
|
|
@ -1000,7 +1085,7 @@ class HexLattice(Lattice):
|
|||
# Check the center ring.
|
||||
if len(axial_slice[-1]) != 1:
|
||||
msg = 'HexLattice ID={0:d} has the wrong number of ' \
|
||||
'elements in the innermost ring. Only 1 element is ' \
|
||||
'elements in the innermost ring. Only 1 element is ' \
|
||||
'allowed in the innermost ring.'.format(self._id)
|
||||
raise ValueError(msg)
|
||||
|
||||
|
|
@ -1009,7 +1094,7 @@ class HexLattice(Lattice):
|
|||
if len(axial_slice[r]) != 6*(self._num_rings - 1 - r):
|
||||
msg = 'HexLattice ID={0:d} has the wrong number of ' \
|
||||
'elements in ring number {1:d} (counting from the '\
|
||||
'outermost ring). This ring should have {2:d} ' \
|
||||
'outermost ring). This ring should have {2:d} ' \
|
||||
'elements.'.format(self._id, r,
|
||||
6*(self._num_rings - 1 - r))
|
||||
raise ValueError(msg)
|
||||
|
|
@ -1045,7 +1130,7 @@ class HexLattice(Lattice):
|
|||
-------
|
||||
3-tuple of int
|
||||
Indices of corresponding lattice element in :math:`(x,\alpha,z)`
|
||||
bases
|
||||
or :math:`(\alpha,y,z)` bases
|
||||
numpy.ndarray
|
||||
Carestian coordinates of the point in the corresponding lattice
|
||||
element coordinate system
|
||||
|
|
@ -1059,15 +1144,21 @@ class HexLattice(Lattice):
|
|||
else:
|
||||
z = point[2] - self.center[2]
|
||||
iz = floor(z/self.pitch[1] + 0.5*self.num_axial)
|
||||
alpha = y - x/sqrt(3.)
|
||||
ix = floor(x/(sqrt(0.75) * self.pitch[0]))
|
||||
ia = floor(alpha/self.pitch[0])
|
||||
|
||||
if self._orientation == 'x':
|
||||
alpha = y - x*sqrt(3.)
|
||||
i1 = floor(-alpha/(sqrt(3.0) * self.pitch[0]))
|
||||
i2 = floor(y/(sqrt(0.75) * self.pitch[0]))
|
||||
else:
|
||||
alpha = y - x/sqrt(3.)
|
||||
i1 = floor(x/(sqrt(0.75) * self.pitch[0]))
|
||||
i2 = floor(alpha/self.pitch[0])
|
||||
# Check four lattice elements to see which one is closest based on local
|
||||
# coordinates
|
||||
indices = [(i1, i2, iz), (i1 + 1, i2, iz), (i1, i2 + 1, iz),
|
||||
(i1 + 1, i2 + 1, iz)]
|
||||
d_min = np.inf
|
||||
for idx in [(ix, ia, iz), (ix + 1, ia, iz), (ix, ia + 1, iz),
|
||||
(ix + 1, ia + 1, iz)]:
|
||||
|
||||
for idx in indices:
|
||||
p = self.get_local_coordinates(point, idx)
|
||||
d = p[0]**2 + p[1]**2
|
||||
if d < d_min:
|
||||
|
|
@ -1085,7 +1176,8 @@ class HexLattice(Lattice):
|
|||
point : Iterable of float
|
||||
Cartesian coordinates of point
|
||||
idx : Iterable of int
|
||||
Indices of lattice element in :math:`(x,\alpha,z)` bases
|
||||
Indices of lattice element in :math:`(x,\alpha,z)`
|
||||
or :math:`(\alpha,y,z)` bases
|
||||
|
||||
Returns
|
||||
-------
|
||||
|
|
@ -1094,8 +1186,17 @@ class HexLattice(Lattice):
|
|||
system
|
||||
|
||||
"""
|
||||
x = point[0] - (self.center[0] + sqrt(0.75)*self.pitch[0]*idx[0])
|
||||
y = point[1] - (self.center[1] + (0.5*idx[0] + idx[1])*self.pitch[0])
|
||||
if self._orientation == 'x':
|
||||
x = point[0] - (self.center[0] + self.pitch[0]*idx[0] +
|
||||
0.5*self.pitch[0]*idx[1])
|
||||
y = point[1] - (self.center[1] +
|
||||
sqrt(0.75)*self.pitch[0]*idx[1])
|
||||
else:
|
||||
x = point[0] - (self.center[0]
|
||||
+ sqrt(0.75)*self.pitch[0]*idx[0])
|
||||
y = point[1] - (self.center[1]
|
||||
+ (0.5*idx[0] + idx[1])*self.pitch[0])
|
||||
|
||||
if self._num_axial is None:
|
||||
z = point[2]
|
||||
else:
|
||||
|
|
@ -1104,18 +1205,21 @@ class HexLattice(Lattice):
|
|||
return (x, y, z)
|
||||
|
||||
def get_universe_index(self, idx):
|
||||
r"""Return index in the universes array corresponding to a lattice element index
|
||||
r"""Return index in the universes array corresponding
|
||||
to a lattice element index
|
||||
|
||||
Parameters
|
||||
----------
|
||||
idx : Iterable of int
|
||||
Lattice element indices in the :math:`(x,\alpha,z)` coordinate
|
||||
system
|
||||
system in 'y' orientation case, or indices in the
|
||||
:math:`(\alpha,y,z)` coordinate system in 'x' one
|
||||
|
||||
Returns
|
||||
-------
|
||||
2- or 3-tuple of int
|
||||
Indices used when setting the :attr:`HexLattice.universes` property
|
||||
2- or 3-tuple of int
|
||||
Indices used when setting the :attr:`HexLattice.universes`
|
||||
property
|
||||
|
||||
"""
|
||||
|
||||
|
|
@ -1138,6 +1242,9 @@ class HexLattice(Lattice):
|
|||
else:
|
||||
i_within = 5*g - z
|
||||
|
||||
if self._orientation == 'x' and g > 0:
|
||||
i_within = (i_within + 5*g) % (6*g)
|
||||
|
||||
if self.num_axial is None:
|
||||
return (i_ring, i_within)
|
||||
else:
|
||||
|
|
@ -1149,8 +1256,8 @@ class HexLattice(Lattice):
|
|||
Parameters
|
||||
----------
|
||||
idx : Iterable of int
|
||||
Lattice element indices in the :math:`(x,\alpha,z)` coordinate
|
||||
system
|
||||
Lattice element indices in the both :math:`(x,\alpha,z)`
|
||||
and :math:`(\alpha,y,z)` coordinate system
|
||||
|
||||
Returns
|
||||
-------
|
||||
|
|
@ -1193,6 +1300,9 @@ class HexLattice(Lattice):
|
|||
self._outer.create_xml_subelement(xml_element)
|
||||
|
||||
lattice_subelement.set("n_rings", str(self._num_rings))
|
||||
# If orientation is "x" export it to XML
|
||||
if self._orientation == 'x':
|
||||
lattice_subelement.set("orientation", "x")
|
||||
|
||||
if self._num_axial is not None:
|
||||
lattice_subelement.set("n_axial", str(self._num_axial))
|
||||
|
|
@ -1267,6 +1377,7 @@ class HexLattice(Lattice):
|
|||
lat = cls(lat_id, name)
|
||||
lat.center = [float(i) for i in get_text(elem, 'center').split()]
|
||||
lat.pitch = [float(i) for i in get_text(elem, 'pitch').split()]
|
||||
lat.orientation = get_text(elem, 'orientation', 'y')
|
||||
outer = get_text(elem, 'outer')
|
||||
if outer is not None:
|
||||
lat.outer = get_universe(int(outer))
|
||||
|
|
@ -1277,13 +1388,13 @@ class HexLattice(Lattice):
|
|||
|
||||
# Create empty nested lists for one axial level
|
||||
univs = [[None for _ in range(max(6*(n_rings - 1 - r), 1))]
|
||||
for r in range(n_rings)]
|
||||
for r in range(n_rings)]
|
||||
if n_axial > 1:
|
||||
univs = [deepcopy(univs) for i in range(n_axial)]
|
||||
|
||||
# Get flat array of universes numbers
|
||||
# Get flat array of universes
|
||||
uarray = np.array([get_universe(int(i)) for i in
|
||||
get_text(elem, 'universes').split()])
|
||||
get_text(elem, 'universes').split()])
|
||||
|
||||
# Fill nested lists
|
||||
j = 0
|
||||
|
|
@ -1291,34 +1402,174 @@ class HexLattice(Lattice):
|
|||
# Get list for a single axial level
|
||||
axial_level = univs[z] if n_axial > 1 else univs
|
||||
|
||||
# Start iterating from top
|
||||
x, alpha = 0, n_rings - 1
|
||||
while True:
|
||||
# Set entry in list based on (x,alpha,z) coordinates
|
||||
_, i_ring, i_within = lat.get_universe_index((x, alpha, z))
|
||||
axial_level[i_ring][i_within] = uarray[j]
|
||||
if lat.orientation == 'y':
|
||||
# Start iterating from top
|
||||
x, alpha = 0, n_rings - 1
|
||||
while True:
|
||||
# Set entry in list based on (x,alpha,z) coordinates
|
||||
_, i_ring, i_within = lat.get_universe_index((x, alpha, z))
|
||||
axial_level[i_ring][i_within] = uarray[j]
|
||||
|
||||
# Move to the right
|
||||
x += 2
|
||||
alpha -= 1
|
||||
if not lat.is_valid_index((x, alpha, z)):
|
||||
# Move down in y direction
|
||||
alpha += x - 1
|
||||
x = 1 - x
|
||||
# Move to the right
|
||||
x += 2
|
||||
alpha -= 1
|
||||
if not lat.is_valid_index((x, alpha, z)):
|
||||
# Move to the right
|
||||
x += 2
|
||||
alpha -= 1
|
||||
# Move down in y direction
|
||||
alpha += x - 1
|
||||
x = 1 - x
|
||||
if not lat.is_valid_index((x, alpha, z)):
|
||||
# Reached the bottom
|
||||
# Move to the right
|
||||
x += 2
|
||||
alpha -= 1
|
||||
if not lat.is_valid_index((x, alpha, z)):
|
||||
# Reached the bottom
|
||||
break
|
||||
j += 1
|
||||
else:
|
||||
# Start iterating from top
|
||||
alpha, y = 1 - n_rings, n_rings - 1
|
||||
while True:
|
||||
# Set entry in list based on (alpha,y,z) coordinates
|
||||
_, i_ring, i_within = lat.get_universe_index((alpha, y, z))
|
||||
axial_level[i_ring][i_within] = uarray[j]
|
||||
|
||||
# Move to the right
|
||||
alpha += 1
|
||||
if not lat.is_valid_index((alpha, y, z)):
|
||||
# Move down to next row
|
||||
alpha = 1 - n_rings
|
||||
y -= 1
|
||||
|
||||
# Check if we've reached the bottom
|
||||
if y == -n_rings:
|
||||
break
|
||||
j += 1
|
||||
|
||||
while not lat.is_valid_index((alpha, y, z)):
|
||||
# Move to the right
|
||||
alpha += 1
|
||||
j += 1
|
||||
|
||||
lat.universes = univs
|
||||
return lat
|
||||
|
||||
def _repr_axial_slice(self, universes):
|
||||
"""Return string representation for the given 2D group of universes.
|
||||
|
||||
The 'universes' argument should be a list of lists of universes where
|
||||
each sub-list represents a single ring. The first list should be the
|
||||
outer ring.
|
||||
"""
|
||||
if self._orientation == 'x':
|
||||
return self._repr_axial_slice_x(universes)
|
||||
else:
|
||||
return self._repr_axial_slice_y(universes)
|
||||
|
||||
def _repr_axial_slice_x(self, universes):
|
||||
"""Return string representation for the given 2D group of universes
|
||||
in 'x' orientation case.
|
||||
|
||||
The 'universes' argument should be a list of lists of universes where
|
||||
each sub-list represents a single ring. The first list should be the
|
||||
outer ring.
|
||||
"""
|
||||
|
||||
# Find the largest universe ID and count the number of digits so we can
|
||||
# properly pad the output string later.
|
||||
largest_id = max([max([univ._id for univ in ring])
|
||||
for ring in universes])
|
||||
n_digits = len(str(largest_id))
|
||||
pad = ' '*n_digits
|
||||
id_form = '{: ^' + str(n_digits) + 'd}'
|
||||
|
||||
# Initialize the list for each row.
|
||||
rows = [[] for i in range(2*self._num_rings - 1)]
|
||||
middle = self._num_rings - 1
|
||||
|
||||
# Start with the degenerate first ring.
|
||||
universe = universes[-1][0]
|
||||
rows[middle] = [id_form.format(universe._id)]
|
||||
|
||||
# Add universes one ring at a time.
|
||||
for r in range(1, self._num_rings):
|
||||
# r_prime increments down while r increments up.
|
||||
r_prime = self._num_rings - 1 - r
|
||||
theta = 0
|
||||
y = middle
|
||||
|
||||
# Climb down the bottom-right
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].append(id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
y += 1
|
||||
theta += 1
|
||||
|
||||
# Climb left across the bottom
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].insert(0, id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
theta += 1
|
||||
|
||||
# Climb up the bottom-left
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].insert(0, id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
y -= 1
|
||||
theta += 1
|
||||
|
||||
# Climb up the top-left
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].insert(0, id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
y -= 1
|
||||
theta += 1
|
||||
|
||||
# Climb right across the top
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].append(id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
theta += 1
|
||||
|
||||
# Climb down the top-right
|
||||
for i in range(r):
|
||||
# Add the universe.
|
||||
universe = universes[r_prime][theta]
|
||||
rows[y].append(id_form.format(universe._id))
|
||||
|
||||
# Translate the indices.
|
||||
y += 1
|
||||
theta += 1
|
||||
|
||||
# Flip the rows and join each row into a single string.
|
||||
rows = [pad.join(x) for x in rows]
|
||||
|
||||
# Pad the beginning of the rows so they line up properly.
|
||||
for y in range(self._num_rings - 1):
|
||||
rows[y] = (self._num_rings - 1 - y)*pad + rows[y]
|
||||
rows[-1 - y] = (self._num_rings - 1 - y)*pad + rows[-1 - y]
|
||||
|
||||
# Join the rows together and return the string.
|
||||
universe_ids = '\n'.join(rows)
|
||||
return universe_ids
|
||||
|
||||
def _repr_axial_slice_y(self, universes):
|
||||
"""Return string representation for the given 2D group of universes in
|
||||
'y' orientation case..
|
||||
|
||||
The 'universes' argument should be a list of lists of universes where
|
||||
each sub-list represents a single ring. The first list should be the
|
||||
outer ring.
|
||||
|
|
@ -1425,7 +1676,7 @@ class HexLattice(Lattice):
|
|||
return universe_ids
|
||||
|
||||
@staticmethod
|
||||
def show_indices(num_rings):
|
||||
def _show_indices_y(num_rings):
|
||||
"""Return a diagram of the hexagonal lattice layout with indices.
|
||||
|
||||
This method can be used to show the proper indices to be used when
|
||||
|
|
@ -1527,3 +1778,124 @@ class HexLattice(Lattice):
|
|||
|
||||
# Join the rows together and return the string.
|
||||
return '\n'.join(rows)
|
||||
|
||||
@staticmethod
|
||||
def _show_indices_x(num_rings):
|
||||
"""Return a diagram of the hexagonal lattice with x orientation
|
||||
layout with indices.
|
||||
|
||||
This method can be used to show the proper indices to be used when
|
||||
setting the :attr:`HexLattice.universes` property. For example,running
|
||||
this method with num_rings=3 will return the similar diagram::
|
||||
|
||||
(0, 8) (0, 9) (0,10)
|
||||
|
||||
(0, 7) (1, 4) (1, 5) (0,11)
|
||||
|
||||
(0, 6) (1, 3) (2, 0) (1, 0) (0, 0)
|
||||
|
||||
(0, 5) (1, 2) (1, 1) (0, 1)
|
||||
|
||||
(0, 4) (0, 3) (0, 2)
|
||||
|
||||
Parameters
|
||||
----------
|
||||
num_rings : int
|
||||
Number of rings in the hexagonal lattice
|
||||
|
||||
Returns
|
||||
-------
|
||||
str
|
||||
Diagram of the hexagonal lattice showing indices in OX orientation
|
||||
|
||||
"""
|
||||
|
||||
# Find the largest string and count the number of digits so we can
|
||||
# properly pad the output string later
|
||||
largest_index = 6*(num_rings - 1)
|
||||
n_digits_index = len(str(largest_index))
|
||||
n_digits_ring = len(str(num_rings - 1))
|
||||
str_form = '({{:{}}},{{:{}}})'.format(n_digits_ring, n_digits_index)
|
||||
pad = ' '*(n_digits_index + n_digits_ring + 3)
|
||||
|
||||
# Initialize the list for each row.
|
||||
rows = [[] for i in range(2*num_rings - 1)]
|
||||
middle = num_rings - 1
|
||||
|
||||
# Start with the degenerate first ring.
|
||||
rows[middle] = [str_form.format(num_rings - 1, 0)]
|
||||
|
||||
# Add universes one ring at a time.
|
||||
for r in range(1, num_rings):
|
||||
# r_prime increments down while r increments up.
|
||||
r_prime = num_rings - 1 - r
|
||||
theta = 0
|
||||
y = middle
|
||||
|
||||
for i in range(r):
|
||||
# Climb down the bottom-right
|
||||
rows[y].append(str_form.format(r_prime, theta))
|
||||
y += 1
|
||||
theta += 1
|
||||
|
||||
for i in range(r):
|
||||
# Climb left across the bottom
|
||||
rows[y].insert(0, str_form.format(r_prime, theta))
|
||||
theta += 1
|
||||
|
||||
for i in range(r):
|
||||
# Climb up the bottom-left
|
||||
rows[y].insert(0, str_form.format(r_prime, theta))
|
||||
y -= 1
|
||||
theta += 1
|
||||
|
||||
for i in range(r):
|
||||
# Climb up the top-left
|
||||
rows[y].insert(0, str_form.format(r_prime, theta))
|
||||
y -= 1
|
||||
theta += 1
|
||||
|
||||
for i in range(r):
|
||||
# Climb right across the top
|
||||
rows[y].append(str_form.format(r_prime, theta))
|
||||
theta += 1
|
||||
|
||||
for i in range(r):
|
||||
# Climb down the top-right
|
||||
rows[y].append(str_form.format(r_prime, theta))
|
||||
y += 1
|
||||
theta += 1
|
||||
|
||||
# Flip the rows and join each row into a single string.
|
||||
rows = [pad.join(x) for x in rows]
|
||||
|
||||
# Pad the beginning of the rows so they line up properly.
|
||||
for y in range(num_rings - 1):
|
||||
rows[y] = (num_rings - 1 - y)*pad + rows[y]
|
||||
rows[-1 - y] = (num_rings - 1 - y)*pad + rows[-1 - y]
|
||||
|
||||
# Join the rows together and return the string.
|
||||
return '\n\n'.join(rows)
|
||||
|
||||
@staticmethod
|
||||
def show_indices(num_rings, orientation="y"):
|
||||
"""Return a diagram of the hexagonal lattice layout with indices.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
num_rings : int
|
||||
Number of rings in the hexagonal lattice
|
||||
orientation : {"x", "y"}
|
||||
Orientation of the hexagonal lattice
|
||||
|
||||
Returns
|
||||
-------
|
||||
str
|
||||
Diagram of the hexagonal lattice showing indices
|
||||
|
||||
"""
|
||||
|
||||
if orientation == 'x':
|
||||
return HexLattice._show_indices_x(num_rings)
|
||||
else:
|
||||
return HexLattice._show_indices_y(num_rings)
|
||||
|
|
|
|||
223
src/lattice.cpp
223
src/lattice.cpp
|
|
@ -440,6 +440,21 @@ HexLattice::HexLattice(pugi::xml_node lat_node)
|
|||
is_3d_ = false;
|
||||
}
|
||||
|
||||
// Read the lattice orientation. Default to 'y'.
|
||||
if (check_for_node(lat_node, "orientation")) {
|
||||
std::string orientation = get_node_value(lat_node, "orientation");
|
||||
if (orientation == "y") {
|
||||
orientation_ = Orientation::y;
|
||||
} else if (orientation == "x") {
|
||||
orientation_ = Orientation::x;
|
||||
} else {
|
||||
fatal_error("Unrecognized orientation '" + orientation
|
||||
+ "' for lattice " + std::to_string(id_));
|
||||
}
|
||||
} else {
|
||||
orientation_ = Orientation::y;
|
||||
}
|
||||
|
||||
// Read the lattice center.
|
||||
std::string center_str {get_node_value(lat_node, "center")};
|
||||
std::vector<std::string> center_words {split(center_str)};
|
||||
|
|
@ -482,14 +497,79 @@ HexLattice::HexLattice(pugi::xml_node lat_node)
|
|||
|
||||
// Parse the universes.
|
||||
// Universes in hexagonal lattices are stored in a manner that represents
|
||||
// a skewed coordinate system: (x, alpha) rather than (x, y). There is
|
||||
// a skewed coordinate system: (x, alpha) in case of 'y' orientation
|
||||
// and (alpha,y) in 'x' one rather than (x, y). There is
|
||||
// no obvious, direct relationship between the order of universes in the
|
||||
// input and the order that they will be stored in the skewed array so
|
||||
// the following code walks a set of index values across the skewed array
|
||||
// in a manner that matches the input order. Note that i_x = 0, i_a = 0
|
||||
// corresponds to the center of the hexagonal lattice.
|
||||
|
||||
// or i_a = 0, i_y = 0 corresponds to the center of the hexagonal lattice.
|
||||
universes_.resize((2*n_rings_-1) * (2*n_rings_-1) * n_axial_, C_NONE);
|
||||
if (orientation_ == Orientation::y) {
|
||||
fill_lattice_y(univ_words);
|
||||
} else {
|
||||
fill_lattice_x(univ_words);
|
||||
}
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
void
|
||||
HexLattice::fill_lattice_x(const std::vector<std::string>& univ_words)
|
||||
{
|
||||
int input_index = 0;
|
||||
for (int m = 0; m < n_axial_; m++) {
|
||||
// Initialize lattice indecies.
|
||||
int i_a = -(n_rings_ - 1);
|
||||
int i_y = n_rings_ - 1;
|
||||
|
||||
// Map upper region of hexagonal lattice which is found in the
|
||||
// first n_rings-1 rows of the input.
|
||||
for (int k = 0; k < n_rings_-1; k++) {
|
||||
|
||||
// Iterate over the input columns.
|
||||
for (int j = 0; j < k+n_rings_; j++) {
|
||||
int indx = (2*n_rings_-1)*(2*n_rings_-1) * m
|
||||
+ (2*n_rings_-1) * (i_y+n_rings_-1)
|
||||
+ (i_a+n_rings_-1);
|
||||
universes_[indx] = std::stoi(univ_words[input_index]);
|
||||
input_index++;
|
||||
// Move to the next right neighbour cell
|
||||
i_a += 1;
|
||||
}
|
||||
|
||||
// Return the lattice index to the start of the current row.
|
||||
i_a = -(n_rings_ - 1);
|
||||
i_y -= 1;
|
||||
}
|
||||
|
||||
// Map the lower region from the centerline of cart to down side
|
||||
for (int k = 0; k < n_rings_; k++) {
|
||||
// Walk the index to the lower-right neighbor of the last row start.
|
||||
i_a = -(n_rings_ - 1) + k;
|
||||
|
||||
// Iterate over the input columns.
|
||||
for (int j = 0; j < 2*n_rings_-k-1; j++) {
|
||||
int indx = (2*n_rings_-1)*(2*n_rings_-1) * m
|
||||
+ (2*n_rings_-1) * (i_y+n_rings_-1)
|
||||
+ (i_a+n_rings_-1);
|
||||
universes_[indx] = std::stoi(univ_words[input_index]);
|
||||
input_index++;
|
||||
// Move to the next right neighbour cell
|
||||
i_a += 1;
|
||||
}
|
||||
|
||||
// Return lattice index to start of current row.
|
||||
i_y -= 1;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
||||
void
|
||||
HexLattice::fill_lattice_y(const std::vector<std::string>& univ_words)
|
||||
{
|
||||
int input_index = 0;
|
||||
for (int m = 0; m < n_axial_; m++) {
|
||||
// Initialize lattice indecies.
|
||||
|
|
@ -610,9 +690,33 @@ std::pair<double, std::array<int, 3>>
|
|||
HexLattice::distance(Position r, Direction u, const std::array<int, 3>& i_xyz)
|
||||
const
|
||||
{
|
||||
// Compute the direction on the hexagonal basis.
|
||||
double beta_dir = u.x * std::sqrt(3.0) / 2.0 + u.y / 2.0;
|
||||
double gamma_dir = u.x * std::sqrt(3.0) / 2.0 - u.y / 2.0;
|
||||
// Short description of the direction vectors used here. The beta, gamma, and
|
||||
// delta vectors point towards the flat sides of each hexagonal tile.
|
||||
// Y - orientation:
|
||||
// basis0 = (1, 0)
|
||||
// basis1 = (-1/sqrt(3), 1) = +120 degrees from basis0
|
||||
// beta = (sqrt(3)/2, 1/2) = +30 degrees from basis0
|
||||
// gamma = (sqrt(3)/2, -1/2) = -60 degrees from beta
|
||||
// delta = (0, 1) = +60 degrees from beta
|
||||
// X - orientation:
|
||||
// basis0 = (1/sqrt(3), -1)
|
||||
// basis1 = (0, 1) = +120 degrees from basis0
|
||||
// beta = (1, 0) = +30 degrees from basis0
|
||||
// gamma = (1/2, -sqrt(3)/2) = -60 degrees from beta
|
||||
// delta = (1/2, sqrt(3)/2) = +60 degrees from beta
|
||||
// The z-axis is considered separately.
|
||||
double beta_dir;
|
||||
double gamma_dir;
|
||||
double delta_dir;
|
||||
if (orientation_ == Orientation::y) {
|
||||
beta_dir = u.x * std::sqrt(3.0) / 2.0 + u.y / 2.0;
|
||||
gamma_dir = u.x * std::sqrt(3.0) / 2.0 - u.y / 2.0;
|
||||
delta_dir = u.y;
|
||||
} else {
|
||||
beta_dir = u.x;
|
||||
gamma_dir = u.x / 2.0 - u.y * std::sqrt(3.0) / 2.0;
|
||||
delta_dir = u.x / 2.0 + u.y * std::sqrt(3.0) / 2.0;
|
||||
}
|
||||
|
||||
// Note that hexagonal lattice distance calculations are performed
|
||||
// using the particle's coordinates relative to the neighbor lattice
|
||||
|
|
@ -620,7 +724,7 @@ const
|
|||
// because there is significant disagreement between neighboring cells
|
||||
// on where the lattice boundary is due to finite precision issues.
|
||||
|
||||
// Upper-right and lower-left sides.
|
||||
// beta direction
|
||||
double d {INFTY};
|
||||
std::array<int, 3> lattice_trans;
|
||||
double edge = -copysign(0.5*pitch_[0], beta_dir); // Oncoming edge
|
||||
|
|
@ -632,7 +736,12 @@ const
|
|||
const std::array<int, 3> i_xyz_t {i_xyz[0]-1, i_xyz[1], i_xyz[2]};
|
||||
r_t = get_local_position(r, i_xyz_t);
|
||||
}
|
||||
double beta = r_t.x * std::sqrt(3.0) / 2.0 + r_t.y / 2.0;
|
||||
double beta;
|
||||
if (orientation_ == Orientation::y) {
|
||||
beta = r_t.x * std::sqrt(3.0) / 2.0 + r_t.y / 2.0;
|
||||
} else {
|
||||
beta = r_t.x;
|
||||
}
|
||||
if ((std::abs(beta - edge) > FP_PRECISION) && beta_dir != 0) {
|
||||
d = (edge - beta) / beta_dir;
|
||||
if (beta_dir > 0) {
|
||||
|
|
@ -642,7 +751,7 @@ const
|
|||
}
|
||||
}
|
||||
|
||||
// Lower-right and upper-left sides.
|
||||
// gamma direction
|
||||
edge = -copysign(0.5*pitch_[0], gamma_dir);
|
||||
if (gamma_dir > 0) {
|
||||
const std::array<int, 3> i_xyz_t {i_xyz[0]+1, i_xyz[1]-1, i_xyz[2]};
|
||||
|
|
@ -651,7 +760,12 @@ const
|
|||
const std::array<int, 3> i_xyz_t {i_xyz[0]-1, i_xyz[1]+1, i_xyz[2]};
|
||||
r_t = get_local_position(r, i_xyz_t);
|
||||
}
|
||||
double gamma = r_t.x * std::sqrt(3.0) / 2.0 - r_t.y / 2.0;
|
||||
double gamma;
|
||||
if (orientation_ == Orientation::y) {
|
||||
gamma = r_t.x * std::sqrt(3.0) / 2.0 - r_t.y / 2.0;
|
||||
} else {
|
||||
gamma = r_t.x / 2.0 - r_t.y * std::sqrt(3.0) / 2.0;
|
||||
}
|
||||
if ((std::abs(gamma - edge) > FP_PRECISION) && gamma_dir != 0) {
|
||||
double this_d = (edge - gamma) / gamma_dir;
|
||||
if (this_d < d) {
|
||||
|
|
@ -664,19 +778,25 @@ const
|
|||
}
|
||||
}
|
||||
|
||||
// Upper and lower sides.
|
||||
edge = -copysign(0.5*pitch_[0], u.y);
|
||||
if (u.y > 0) {
|
||||
// delta direction
|
||||
edge = -copysign(0.5*pitch_[0], delta_dir);
|
||||
if (delta_dir > 0) {
|
||||
const std::array<int, 3> i_xyz_t {i_xyz[0], i_xyz[1]+1, i_xyz[2]};
|
||||
r_t = get_local_position(r, i_xyz_t);
|
||||
} else {
|
||||
const std::array<int, 3> i_xyz_t {i_xyz[0], i_xyz[1]-1, i_xyz[2]};
|
||||
r_t = get_local_position(r, i_xyz_t);
|
||||
}
|
||||
if ((std::abs(r_t.y - edge) > FP_PRECISION) && u.y != 0) {
|
||||
double this_d = (edge - r_t.y) / u.y;
|
||||
double delta;
|
||||
if (orientation_ == Orientation::y) {
|
||||
delta = r_t.y;
|
||||
} else {
|
||||
delta = r_t.x / 2.0 + r_t.y * std::sqrt(3.0) / 2.0;
|
||||
}
|
||||
if ((std::abs(delta - edge) > FP_PRECISION) && delta_dir != 0) {
|
||||
double this_d = (edge - delta) / delta_dir;
|
||||
if (this_d < d) {
|
||||
if (u.y > 0) {
|
||||
if (delta_dir > 0) {
|
||||
lattice_trans = {0, 1, 0};
|
||||
} else {
|
||||
lattice_trans = {0, -1, 0};
|
||||
|
|
@ -727,16 +847,25 @@ HexLattice::get_indices(Position r, Direction u) const
|
|||
}
|
||||
}
|
||||
|
||||
// Convert coordinates into skewed bases. The (x, alpha) basis is used to
|
||||
// find the index of the global coordinates to within 4 cells.
|
||||
double alpha = r_o.y - r_o.x / std::sqrt(3.0);
|
||||
int ix = std::floor(r_o.x / (0.5*std::sqrt(3.0) * pitch_[0]));
|
||||
int ia = std::floor(alpha / pitch_[0]);
|
||||
int i1, i2;
|
||||
if (orientation_ == Orientation::y) {
|
||||
// Convert coordinates into skewed bases. The (x, alpha) basis is used to
|
||||
// find the index of the global coordinates to within 4 cells.
|
||||
double alpha = r_o.y - r_o.x / std::sqrt(3.0);
|
||||
i1 = std::floor(r_o.x / (0.5*std::sqrt(3.0) * pitch_[0]));
|
||||
i2 = std::floor(alpha / pitch_[0]);
|
||||
} else {
|
||||
// Convert coordinates into skewed bases. The (alpha, y) basis is used to
|
||||
// find the index of the global coordinates to within 4 cells.
|
||||
double alpha = r_o.y - r_o.x * std::sqrt(3.0);
|
||||
i1 = std::floor(-alpha / (std::sqrt(3.0) * pitch_[0]));
|
||||
i2 = std::floor(r_o.y / (0.5*std::sqrt(3.0) * pitch_[0]));
|
||||
}
|
||||
|
||||
// Add offset to indices (the center cell is (i_x, i_alpha) = (0, 0) but
|
||||
// Add offset to indices (the center cell is (i1, i2) = (0, 0) but
|
||||
// the array is offset so that the indices never go below 0).
|
||||
ix += n_rings_-1;
|
||||
ia += n_rings_-1;
|
||||
i1 += n_rings_-1;
|
||||
i2 += n_rings_-1;
|
||||
|
||||
// Calculate the (squared) distance between the particle and the centers of
|
||||
// the four possible cells. Regular hexagonal tiles form a Voronoi
|
||||
|
|
@ -755,14 +884,14 @@ HexLattice::get_indices(Position r, Direction u) const
|
|||
// is kept (i.e. the cell with the lowest dot product as the vectors will be
|
||||
// completely opposed if the particle is moving directly toward the center of
|
||||
// the cell).
|
||||
int ix_chg {};
|
||||
int ia_chg {};
|
||||
int i1_chg {};
|
||||
int i2_chg {};
|
||||
double d_min {INFTY};
|
||||
double dp_min {INFTY};
|
||||
for (int i = 0; i < 2; i++) {
|
||||
for (int j = 0; j < 2; j++) {
|
||||
// get local coordinates
|
||||
const std::array<int, 3> i_xyz {ix + j, ia + i, 0};
|
||||
const std::array<int, 3> i_xyz {i1 + j, i2 + i, 0};
|
||||
Position r_t = get_local_position(r, i_xyz);
|
||||
// calculate distance
|
||||
double d = r_t.x*r_t.x + r_t.y*r_t.y;
|
||||
|
|
@ -777,18 +906,18 @@ HexLattice::get_indices(Position r, Direction u) const
|
|||
if (on_boundary && dp > dp_min) continue;
|
||||
// update values
|
||||
d_min = d;
|
||||
ix_chg = j;
|
||||
ia_chg = i;
|
||||
i1_chg = j;
|
||||
i2_chg = i;
|
||||
dp_min = dp;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
// update outgoing indices
|
||||
ix += ix_chg;
|
||||
ia += ia_chg;
|
||||
i1 += i1_chg;
|
||||
i2 += i2_chg;
|
||||
|
||||
return {ix, ia, iz};
|
||||
return {i1, i2, iz};
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
|
@ -797,15 +926,26 @@ Position
|
|||
HexLattice::get_local_position(Position r, const std::array<int, 3> i_xyz)
|
||||
const
|
||||
{
|
||||
// x_l = x_g - (center + pitch_x*cos(30)*index_x)
|
||||
r.x -= center_.x + std::sqrt(3.0)/2.0 * (i_xyz[0] - n_rings_ + 1) * pitch_[0];
|
||||
// y_l = y_g - (center + pitch_x*index_x + pitch_y*sin(30)*index_y)
|
||||
r.y -= (center_.y + (i_xyz[1] - n_rings_ + 1) * pitch_[0]
|
||||
+ (i_xyz[0] - n_rings_ + 1) * pitch_[0] / 2.0);
|
||||
if (is_3d_) {
|
||||
r.z -= center_.z - (0.5 * n_axial_ - i_xyz[2] - 0.5) * pitch_[1];
|
||||
if (orientation_ == Orientation::y) {
|
||||
// x_l = x_g - (center + pitch_x*cos(30)*index_x)
|
||||
r.x -= center_.x
|
||||
+ std::sqrt(3.0)/2.0 * (i_xyz[0] - n_rings_ + 1) * pitch_[0];
|
||||
// y_l = y_g - (center + pitch_x*index_x + pitch_y*sin(30)*index_y)
|
||||
r.y -= (center_.y + (i_xyz[1] - n_rings_ + 1) * pitch_[0]
|
||||
+ (i_xyz[0] - n_rings_ + 1) * pitch_[0] / 2.0);
|
||||
} else {
|
||||
// x_l = x_g - (center + pitch_x*index_a + pitch_y*sin(30)*index_y)
|
||||
r.x -= (center_.x + (i_xyz[0] - n_rings_ + 1) * pitch_[0]
|
||||
+ (i_xyz[1] - n_rings_ + 1) * pitch_[0] / 2.0);
|
||||
// y_l = y_g - (center + pitch_y*cos(30)*index_y)
|
||||
r.y -= center_.y
|
||||
+ std::sqrt(3.0)/2.0 * (i_xyz[1] - n_rings_ + 1) * pitch_[0];
|
||||
}
|
||||
|
||||
if (is_3d_) {
|
||||
r.z -= center_.z - (0.5 * n_axial_ - i_xyz[2] - 0.5) * pitch_[1];
|
||||
}
|
||||
|
||||
return r;
|
||||
}
|
||||
|
||||
|
|
@ -863,6 +1003,11 @@ HexLattice::to_hdf5_inner(hid_t lat_group) const
|
|||
write_string(lat_group, "type", "hexagonal", false);
|
||||
write_dataset(lat_group, "n_rings", n_rings_);
|
||||
write_dataset(lat_group, "n_axial", n_axial_);
|
||||
if (orientation_ == Orientation::y) {
|
||||
write_string(lat_group, "orientation", "y", false);
|
||||
} else {
|
||||
write_string(lat_group, "orientation", "x", false);
|
||||
}
|
||||
if (is_3d_) {
|
||||
write_dataset(lat_group, "pitch", pitch_);
|
||||
write_dataset(lat_group, "center", center_);
|
||||
|
|
|
|||
|
|
@ -48,6 +48,7 @@ element geometry {
|
|||
(element n_axial { xsd:int } | attribute n_axial { xsd:int })? &
|
||||
(element center { list { xsd:double+ } } | attribute center { list { xsd:double+ } }) &
|
||||
(element pitch { list { xsd:double+ } } | attribute pitch { list { xsd:double+ } }) &
|
||||
(element orientation { ( "x" | "y" ) } | attribute orientation { ( "x" | "y" ) })? &
|
||||
(element universes { list { xsd:int+ } } | attribute universes { list { xsd:int+ } }) &
|
||||
(element outer { xsd:int } | attribute outer { xsd:int })?
|
||||
}*
|
||||
|
|
|
|||
|
|
@ -398,6 +398,22 @@
|
|||
</list>
|
||||
</attribute>
|
||||
</choice>
|
||||
<optional>
|
||||
<choice>
|
||||
<element name="orientation">
|
||||
<choice>
|
||||
<value>x</value>
|
||||
<value>y</value>
|
||||
</choice>
|
||||
</element>
|
||||
<attribute name="orientation">
|
||||
<choice>
|
||||
<value>x</value>
|
||||
<value>y</value>
|
||||
</choice>
|
||||
</attribute>
|
||||
</choice>
|
||||
</optional>
|
||||
<choice>
|
||||
<element name="universes">
|
||||
<list>
|
||||
|
|
|
|||
0
tests/regression_tests/lattice_hex_x/__init__.py
Normal file
0
tests/regression_tests/lattice_hex_x/__init__.py
Normal file
123
tests/regression_tests/lattice_hex_x/inputs_true.dat
Normal file
123
tests/regression_tests/lattice_hex_x/inputs_true.dat
Normal file
|
|
@ -0,0 +1,123 @@
|
|||
<?xml version='1.0' encoding='utf-8'?>
|
||||
<geometry>
|
||||
<cell id="1" material="1" region="-1" universe="1" />
|
||||
<cell id="2" material="4" region="-2 1" universe="1" />
|
||||
<cell id="3" material="2" region="2" universe="1" />
|
||||
<cell id="4" material="2" region="-3" universe="2" />
|
||||
<cell id="5" material="4" region="-4 3" universe="2" />
|
||||
<cell id="6" material="2" region="4" universe="2" />
|
||||
<cell id="7" material="3" region="-5" universe="3" />
|
||||
<cell id="8" material="4" region="-6 5" universe="3" />
|
||||
<cell id="9" material="2" region="-7 6" universe="3" />
|
||||
<cell id="10" material="4" region="-8 7" universe="3" />
|
||||
<cell id="11" material="2" region="8" universe="3" />
|
||||
<cell id="12" material="2" universe="4" />
|
||||
<cell fill="9" id="13" name="container assembly cell" region="-11 12 -13 14 15 -16 9 -10" universe="5" />
|
||||
<hex_lattice id="9" n_axial="2" n_rings="11" name="regular fuel assembly" orientation="x">
|
||||
<pitch>1.235 5.0</pitch>
|
||||
<outer>4</outer>
|
||||
<center>0.0 0.0 5.0</center>
|
||||
<universes>
|
||||
1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 3 1 1 1 1 3 1 1 1 3 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 2 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 3 1 1 1 3 1 1 1 1 3 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 3 1 1 1 1 3 1 1 1 3 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 2 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 3 1 1 1 3 1 1 1 1 3 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 3 1 1 1 1 3 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 3 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1 1
|
||||
1 1 1 1 1 1 1 1 1 1 1</universes>
|
||||
</hex_lattice>
|
||||
<surface coeffs="0.0 0.0 0.386" id="1" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.4582" id="2" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.45" id="3" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.5177" id="4" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.35" id="5" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.41" id="6" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.545" id="7" type="z-cylinder" />
|
||||
<surface coeffs="0.0 0.0 0.6323" id="8" type="z-cylinder" />
|
||||
<surface boundary="reflective" coeffs="0.0" id="9" type="z-plane" />
|
||||
<surface boundary="reflective" coeffs="10.0" id="10" type="z-plane" />
|
||||
<surface boundary="reflective" coeffs="11.8" id="11" type="y-plane" />
|
||||
<surface boundary="reflective" coeffs="-11.8" id="12" type="y-plane" />
|
||||
<surface boundary="reflective" coeffs="1.7320508075688772 1.0 0.0 23.6" id="13" type="plane" />
|
||||
<surface boundary="reflective" coeffs="-1.7320508075688772 1.0 0.0 -23.6" id="14" type="plane" />
|
||||
<surface boundary="reflective" coeffs="1.7320508075688772 1.0 0.0 -23.6" id="15" type="plane" />
|
||||
<surface boundary="reflective" coeffs="-1.7320508075688772 1.0 0.0 23.6" id="16" type="plane" />
|
||||
</geometry>
|
||||
<?xml version='1.0' encoding='utf-8'?>
|
||||
<materials>
|
||||
<material depletable="true" id="1" name="UO2">
|
||||
<density units="sum" />
|
||||
<nuclide ao="0.0008737" name="U235" />
|
||||
<nuclide ao="0.018744" name="U238" />
|
||||
<nuclide ao="0.039235" name="O16" />
|
||||
</material>
|
||||
<material id="2" name="borated H2O">
|
||||
<density units="sum" />
|
||||
<nuclide ao="0.06694" name="H1" />
|
||||
<nuclide ao="0.03347" name="O16" />
|
||||
<nuclide ao="6.6262e-06" name="B10" />
|
||||
<nuclide ao="2.6839e-05" name="B11" />
|
||||
</material>
|
||||
<material id="3" name="pellet B4C">
|
||||
<density units="sum" />
|
||||
<nuclide ao="0.01966" name="C0" />
|
||||
<nuclide ao="4.7344e-06" name="B11" />
|
||||
<nuclide ao="1.9177e-05" name="B10" />
|
||||
</material>
|
||||
<material id="4" name="Zirc4">
|
||||
<density units="sum" />
|
||||
<nuclide ao="0.021763349999999997" name="Zr90" />
|
||||
<nuclide ao="0.00474606" name="Zr91" />
|
||||
<nuclide ao="0.00725445" name="Zr92" />
|
||||
<nuclide ao="0.00735174" name="Zr94" />
|
||||
<nuclide ao="0.0011844" name="Zr96" />
|
||||
</material>
|
||||
</materials>
|
||||
<?xml version='1.0' encoding='utf-8'?>
|
||||
<settings>
|
||||
<run_mode>eigenvalue</run_mode>
|
||||
<particles>1000</particles>
|
||||
<batches>10</batches>
|
||||
<inactive>5</inactive>
|
||||
<source strength="1.0">
|
||||
<space type="box">
|
||||
<parameters>-13.62546635287517 -13.62546635287517 0.0 13.62546635287517 13.62546635287517 10.0</parameters>
|
||||
</space>
|
||||
</source>
|
||||
<seed>22</seed>
|
||||
</settings>
|
||||
2
tests/regression_tests/lattice_hex_x/results_true.dat
Normal file
2
tests/regression_tests/lattice_hex_x/results_true.dat
Normal file
|
|
@ -0,0 +1,2 @@
|
|||
k-combined:
|
||||
1.355663E+00 2.896562E-02
|
||||
206
tests/regression_tests/lattice_hex_x/test.py
Normal file
206
tests/regression_tests/lattice_hex_x/test.py
Normal file
|
|
@ -0,0 +1,206 @@
|
|||
from tests.testing_harness import PyAPITestHarness
|
||||
import openmc
|
||||
import numpy as np
|
||||
|
||||
|
||||
class HexLatticeOXTestHarness(PyAPITestHarness):
|
||||
|
||||
def _build_inputs(self):
|
||||
materials = openmc.Materials()
|
||||
|
||||
fuel_mat = openmc.Material(material_id=1, name="UO2")
|
||||
fuel_mat.set_density('sum')
|
||||
fuel_mat.add_nuclide('U235', 0.87370e-03)
|
||||
fuel_mat.add_nuclide('U238', 1.87440e-02)
|
||||
fuel_mat.add_nuclide('O16', 3.92350e-02)
|
||||
materials.append(fuel_mat)
|
||||
|
||||
coolant = openmc.Material(material_id=2, name="borated H2O")
|
||||
coolant.set_density('sum')
|
||||
coolant.add_nuclide('H1', 0.06694)
|
||||
coolant.add_nuclide('O16', 0.03347)
|
||||
coolant.add_nuclide('B10', 6.6262e-6)
|
||||
coolant.add_nuclide('B11', 2.6839e-5)
|
||||
materials.append(coolant)
|
||||
|
||||
absorber = openmc.Material(material_id=3, name="pellet B4C")
|
||||
absorber.set_density('sum')
|
||||
absorber.add_nuclide('C0', 0.01966)
|
||||
absorber.add_nuclide('B11', 4.7344e-6)
|
||||
absorber.add_nuclide('B10', 1.9177e-5)
|
||||
materials.append(absorber)
|
||||
|
||||
zirc = openmc.Material(material_id=4, name="Zirc4")
|
||||
zirc.set_density('sum')
|
||||
zirc.add_element('Zr', 4.23e-2)
|
||||
materials.append(zirc)
|
||||
|
||||
materials.export_to_xml()
|
||||
|
||||
# Geometry #
|
||||
|
||||
pin_rad = 0.7 # cm
|
||||
assembly_pitch = 1.235 # cm
|
||||
hexagonal_pitch = 23.6 # cm
|
||||
length = 10.0 # cm
|
||||
|
||||
# Fuel pin surfaces
|
||||
|
||||
cylfuelin = openmc.ZCylinder(surface_id=1, r=0.386)
|
||||
cylfuelout = openmc.ZCylinder(surface_id=2, r=0.4582)
|
||||
|
||||
# Fuel cells
|
||||
|
||||
infcell = openmc.Cell(cell_id=1)
|
||||
infcell.region = -cylfuelin
|
||||
infcell.fill = fuel_mat
|
||||
|
||||
clfcell = openmc.Cell(cell_id=2)
|
||||
clfcell.region = -cylfuelout & +cylfuelin
|
||||
clfcell.fill = zirc
|
||||
|
||||
outfcell = openmc.Cell(cell_id=3)
|
||||
outfcell.region = +cylfuelout
|
||||
outfcell.fill = coolant
|
||||
|
||||
# Fuel universe
|
||||
|
||||
fuel_ch_univ = openmc.Universe(universe_id=1, name="Fuel channel",
|
||||
cells=[infcell, clfcell, outfcell])
|
||||
|
||||
# Central tube surfaces
|
||||
|
||||
cyltubein = openmc.ZCylinder(surface_id=3, r=0.45)
|
||||
cyltubeout = openmc.ZCylinder(surface_id=4, r=0.5177)
|
||||
|
||||
# Central tube cells
|
||||
|
||||
inctcell = openmc.Cell(cell_id=4)
|
||||
inctcell.region = -cyltubein
|
||||
inctcell.fill = coolant
|
||||
|
||||
clctcell = openmc.Cell(cell_id=5)
|
||||
clctcell.region = -cyltubeout & +cyltubein
|
||||
clctcell.fill = zirc
|
||||
|
||||
outctcell = openmc.Cell(cell_id=6)
|
||||
outctcell.region = +cyltubeout
|
||||
outctcell.fill = coolant
|
||||
|
||||
# Central tubel universe
|
||||
|
||||
tube_ch_univ = openmc.Universe(universe_id=2,
|
||||
name="Central tube channel",
|
||||
cells=[inctcell, clctcell, outctcell])
|
||||
|
||||
# Absorber tube surfaces
|
||||
|
||||
cylabsin = openmc.ZCylinder(surface_id=5, r=0.35)
|
||||
cylabsout = openmc.ZCylinder(surface_id=6, r=0.41)
|
||||
cylabsclin = openmc.ZCylinder(surface_id=7, r=0.545)
|
||||
cylabsclout = openmc.ZCylinder(surface_id=8, r=0.6323)
|
||||
|
||||
# Absorber tube cells
|
||||
|
||||
inabscell = openmc.Cell(cell_id=7)
|
||||
inabscell.region = -cylabsin
|
||||
inabscell.fill = absorber
|
||||
|
||||
clabscell = openmc.Cell(cell_id=8)
|
||||
clabscell.region = -cylabsout & +cylabsin
|
||||
clabscell.fill = zirc
|
||||
|
||||
interabscell = openmc.Cell(cell_id=9)
|
||||
interabscell.region = -cylabsclin & +cylabsout
|
||||
interabscell.fill = coolant
|
||||
|
||||
clatcell = openmc.Cell(cell_id=10)
|
||||
clatcell.region = -cylabsclout & +cylabsclin
|
||||
clatcell.fill = zirc
|
||||
|
||||
outabscell = openmc.Cell(cell_id=11)
|
||||
outabscell.region = +cylabsclout
|
||||
outabscell.fill = coolant
|
||||
|
||||
# Absorber tube universe
|
||||
|
||||
abs_ch_univ = openmc.Universe(universe_id=3,
|
||||
name="Central tube channel",
|
||||
cells=[inabscell, clabscell,
|
||||
interabscell,
|
||||
clatcell, outabscell])
|
||||
# Assembly surfaces
|
||||
|
||||
edge_length = (1./np.sqrt(3.0)) * hexagonal_pitch
|
||||
fuel_bottom = openmc.ZPlane(surface_id=9, z0=0.0,
|
||||
boundary_type='reflective')
|
||||
fuel_top = openmc.ZPlane(surface_id=10, z0=length,
|
||||
boundary_type='reflective')
|
||||
|
||||
# a hex surface for the core to go inside of
|
||||
|
||||
hexprism = openmc.model.get_hexagonal_prism(edge_length=edge_length,
|
||||
origin=(0.0, 0.0),
|
||||
boundary_type='reflective',
|
||||
orientation='x')
|
||||
region = hexprism & +fuel_bottom & -fuel_top
|
||||
|
||||
inf_mat = openmc.Cell(cell_id=12)
|
||||
inf_mat.fill = coolant
|
||||
inf_mat_univ = openmc.Universe(universe_id=4, cells=[inf_mat])
|
||||
|
||||
# Fill lattice by channels
|
||||
|
||||
nring = 11
|
||||
universes = []
|
||||
for ring in range(nring - 1, -1, -1):
|
||||
arr = []
|
||||
arr.append(fuel_ch_univ)
|
||||
for cell in range(ring * 6 - 1):
|
||||
arr.append(fuel_ch_univ)
|
||||
universes.append(arr)
|
||||
universes[-1] = [tube_ch_univ]
|
||||
channels = [(7, 2), (7, 5), (7, 8), (7, 11), (7, 14), (7, 17), (5, 0),
|
||||
(4, 3), (5, 5), (4, 9), (5, 10), (4, 15), (5, 15),
|
||||
(4, 21), (5, 20), (4, 27), (5, 25), (4, 33)]
|
||||
for i, j in channels:
|
||||
universes[i][j] = abs_ch_univ
|
||||
lattice = openmc.HexLattice(name="regular fuel assembly")
|
||||
lattice.orientation = "x"
|
||||
lattice.center = (0., 0., length/2.0)
|
||||
lattice.pitch = (assembly_pitch, length/2.0)
|
||||
lattice.universes = 2*[universes]
|
||||
lattice.outer = inf_mat_univ
|
||||
|
||||
assembly_cell = openmc.Cell(cell_id=13,
|
||||
name="container assembly cell")
|
||||
assembly_cell.region = region
|
||||
assembly_cell.fill = lattice
|
||||
|
||||
root_univ = openmc.Universe(universe_id=5, name="root universe",
|
||||
cells=[assembly_cell])
|
||||
|
||||
geom = openmc.Geometry(root_univ)
|
||||
geom.export_to_xml()
|
||||
|
||||
# Settings #
|
||||
|
||||
settings = openmc.Settings()
|
||||
settings.run_mode = 'eigenvalue'
|
||||
|
||||
source = openmc.Source()
|
||||
ll = [-edge_length, -edge_length, 0.0]
|
||||
ur = [edge_length, edge_length, 10.0]
|
||||
source.space = openmc.stats.Box(ll, ur)
|
||||
source.strength = 1.0
|
||||
settings.source = source
|
||||
settings.batches = 10
|
||||
settings.inactive = 5
|
||||
settings.particles = 1000
|
||||
settings.seed = 22
|
||||
settings.export_to_xml()
|
||||
|
||||
|
||||
def test_lattice_hex_ox_surf():
|
||||
harness = HexLatticeOXTestHarness('statepoint.10.h5')
|
||||
harness.main()
|
||||
|
|
@ -202,6 +202,13 @@ def test_get_universe(rlat2, rlat3, hlat2, hlat3):
|
|||
assert hlat2.get_universe((1, 0)) == u1
|
||||
assert hlat2.get_universe((-2, 2)) == u1
|
||||
|
||||
hlat2.orientation = 'x'
|
||||
assert hlat2.get_universe((2, 0)) == u2
|
||||
assert hlat2.get_universe((1, 0)) == u2
|
||||
assert hlat2.get_universe((1, 1)) == u1
|
||||
assert hlat2.get_universe((-1, 1)) == u1
|
||||
hlat2.orientation = 'y'
|
||||
|
||||
u1, u2, u3, outer = hlat3.univs
|
||||
assert hlat3.get_universe((0, 0, 0)) == u2
|
||||
assert hlat3.get_universe((0, 0, 1)) == u3
|
||||
|
|
@ -342,3 +349,5 @@ def test_show_indices():
|
|||
for i in range(1, 11):
|
||||
lines = openmc.HexLattice.show_indices(i).split('\n')
|
||||
assert len(lines) == 4*i - 3
|
||||
lines_x = openmc.HexLattice.show_indices(i, 'x').split('\n')
|
||||
assert len(lines) == 4*i - 3
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue