Merge pull request #1381 from GiudGiud/PR_recursive_LNS

Lattice discretization with degenerate and recursive local neighbor symmetry algorithms
This commit is contained in:
Sterling Harper 2020-01-28 18:22:45 -05:00 committed by GitHub
commit f1bb96eaef
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
8 changed files with 425 additions and 13 deletions

View file

@ -422,12 +422,18 @@ class Cell(IDManagerMixin):
return universes
def clone(self, memo=None):
def clone(self, clone_materials=True, clone_regions=True, memo=None):
"""Create a copy of this cell with a new unique ID, and clones
the cell's region and fill.
Parameters
----------
clone_materials : bool
Whether to create separate copies of the materials filling cells
contained in this cell, or the material filling this cell.
clone_regions : bool
Whether to create separate copies of the regions bounding cells
contained in this cell, and the region bounding this cell.
memo : dict or None
A nested dictionary of previously cloned objects. This parameter
is used internally and should not be specified by the user.
@ -456,13 +462,25 @@ class Cell(IDManagerMixin):
self._paths = paths
if self.region is not None:
clone.region = self.region.clone(memo)
if clone_regions:
clone.region = self.region.clone(memo)
else:
clone.region = self.region
if self.fill is not None:
if self.fill_type == 'distribmat':
clone.fill = [fill.clone(memo) if fill is not None else None
for fill in self.fill]
if not clone_materials:
clone.fill = self.fill
else:
clone.fill = [fill.clone(memo) if fill is not None else
None for fill in self.fill]
elif self.fill_type == 'material':
if not clone_materials:
clone.fill = self.fill
else:
clone.fill = self.fill.clone(memo)
else:
clone.fill = self.fill.clone(memo)
clone.fill = self.fill.clone(clone_materials,
clone_regions, memo)
# Memoize the clone
memo[self] = clone

View file

@ -386,7 +386,8 @@ class Geometry(object):
surfaces = OrderedDict()
for cell in self.get_all_cells().values():
surfaces = cell.region.get_surfaces(surfaces)
if cell.region is not None:
surfaces = cell.region.get_surfaces(surfaces)
return surfaces
def get_redundant_surfaces(self):

View file

@ -7,6 +7,8 @@ from numbers import Real
from xml.etree import ElementTree as ET
import numpy as np
import warnings
import types
import openmc.checkvalue as cv
import openmc
@ -454,12 +456,18 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
return []
return [(self, idx)] + u.find(p)
def clone(self, memo=None):
def clone(self, clone_materials=True, clone_regions=True, memo=None):
"""Create a copy of this lattice with a new unique ID, and clones
all universes within this lattice.
Parameters
----------
clone_materials : bool
Whether to create separate copies of the materials filling cells
contained in this lattice and its outer universe.
clone_regions : bool
Whether to create separate copies of the regions bounding cells
contained in this lattice and its outer universe.
memo : dict or None
A nested dictionary of previously cloned objects. This parameter
is used internally and should not be specified by the user.
@ -480,19 +488,23 @@ class Lattice(IDManagerMixin, metaclass=ABCMeta):
clone.id = None
if self.outer is not None:
clone.outer = self.outer.clone(memo)
clone.outer = self.outer.clone(clone_materials, clone_regions,
memo)
# Assign universe clones to the lattice clone
for i in self.indices:
if isinstance(self, RectLattice):
clone.universes[i] = self.universes[i].clone(memo)
clone.universes[i] = self.universes[i].clone(
clone_materials, clone_regions, memo)
else:
if self.ndim == 2:
clone.universes[i[0]][i[1]] = \
self.universes[i[0]][i[1]].clone(memo)
self.universes[i[0]][i[1]].clone(clone_materials,
clone_regions, memo)
else:
clone.universes[i[0]][i[1]][i[2]] = \
self.universes[i[0]][i[1]][i[2]].clone(memo)
self.universes[i[0]][i[1]][i[2]].clone(
clone_materials, clone_regions, memo)
# Memoize the clone
memo[self] = clone
@ -759,6 +771,245 @@ class RectLattice(Lattice):
0 <= idx[1] < self.shape[1] and
0 <= idx[2] < self.shape[2])
def discretize(self, strategy="degenerate",
universes_to_ignore=[],
materials_to_clone=[],
lattice_neighbors=[], key=lambda univ: univ.id):
"""Discretize the lattice with either a degenerate or a local neighbor
symmetry strategy
'Degenerate' clones every universe in the lattice, thus making them all
uniquely defined. This is typically required if depletion or thermal
hydraulics will make every universe's environment unique.
'Local neighbor symmetry' groups universes with similar neighborhoods.
These clusters of cells and materials provide increased convergence
speed to multi-group cross sections tallies. The user can specify
the lattice's neighbors to discriminate between two sides of a
lattice for example.
Parameters
----------
strategy : {'degenerate', 'lns'}
Which strategy to adopt when discretizing the lattice
universes_to_ignore : Iterable of Universe
Lattice universes that need not be discretized
materials_to_clone : Iterable of Material
List of materials that should be cloned when discretizing
lattice_neighbors : Iterable of Universe
List of the lattice's neighbors. By default, if present, the
lattice outer universe will be used. The neighbors should be
ordered as follows [top left, top, top right, left, right,
bottom left, bottom, bottom right]
key : function
Function of argument a universe that is used to extract a
comparison key. This function will be called on each universe's
neighbors in the lattice to form a neighbor pattern. This pattern
is then used to identify unique neighbor symmetries.
"""
# Check routine inputs
if self.ndim != 2:
raise NotImplementedError("LNS discretization is not implemented "
"for 1D and 3D lattices")
cv.check_value('strategy', strategy, ('degenerate', 'lns'))
cv.check_type('universes_to_ignore', universes_to_ignore, Iterable,
openmc.Universe)
cv.check_type('materials_to_clone', materials_to_clone, Iterable,
openmc.Material)
cv.check_type('lattice_neighbors', lattice_neighbors, Iterable,
openmc.Universe)
cv.check_value('number of lattice_neighbors', len(lattice_neighbors),
(0, 8))
cv.check_type('key', key, types.FunctionType)
# Use outer universe if neighbors are missing and outer is defined
if self.outer is not None and len(lattice_neighbors) == 0:
lattice_neighbors = [key(self.outer) for i in range(8)]
elif len(lattice_neighbors) == 8:
lattice_neighbors = [key(universe) for universe in
lattice_neighbors]
# Dictionary that will keep track of where each pattern appears, how
# it was rotated and/or symmetrized
patterns = {}
# Initialize pattern array
pattern = np.empty(shape=(3, 3), dtype=type(key(self.universes[0][0])))
# Define an auxiliary function that returns a universe's neighbors
# that are outside the lattice
def find_edge_neighbors(pattern, i, j):
# If no neighbors have been specified, start with an empty array
if len(lattice_neighbors) == 0:
return
# Left edge
if i == 0:
pattern[:, 0] = lattice_neighbors[3]
if j == 0:
pattern[0, 0] = lattice_neighbors[0]
elif j == self.shape[1] - 1:
pattern[2, 0] = lattice_neighbors[5]
# Bottom edge
if j == 0:
pattern[0, 1] = lattice_neighbors[1]
if i != 0:
pattern[0, 0] = lattice_neighbors[1]
if i != self.shape[0] - 1:
pattern[0, 2] = lattice_neighbors[1]
# Right edge
if i == self.shape[0] - 1:
pattern[:, 2] = lattice_neighbors[4]
if j == 0:
pattern[0, 2] = lattice_neighbors[2]
elif j == self.shape[1] - 1:
pattern[2, 2] = lattice_neighbors[7]
# Top edge
if j == self.shape[1] - 1:
pattern[2, 1] = lattice_neighbors[6]
if i != 0:
pattern[2, 0] = lattice_neighbors[6]
if i != self.shape[0] - 1:
pattern[2, 2] = lattice_neighbors[6]
# Define an auxiliary function that returns a universe's neighbors
# among the universes inside the lattice
def find_lattice_neighbors(pattern, i, j):
# Away from left edge
if i != 0:
if j > 0:
pattern[0, 0] = key(self.universes[j-1][i-1])
pattern[1, 0] = key(self.universes[j][i-1])
if j < self.shape[1] - 1:
pattern[2, 0] = key(self.universes[j+1][i-1])
# Away from bottom edge
if j != 0:
if i > 0:
pattern[0, 0] = key(self.universes[j-1][i-1])
pattern[0, 1] = key(self.universes[j-1][i])
if i < self.shape[0] - 1:
pattern[0, 2] = key(self.universes[j-1][i+1])
# Away from right edge
if i != self.shape[0] - 1:
if j > 0:
pattern[0, 2] = key(self.universes[j-1][i+1])
pattern[1, 2] = key(self.universes[j][i+1])
if j < self.shape[1] - 1:
pattern[2, 2] = key(self.universes[j+1][i+1])
# Away from top edge
if j != self.shape[1] - 1:
if i > 0:
pattern[2, 0] = key(self.universes[j+1][i-1])
pattern[2, 1] = key(self.universes[j+1][i])
if i < self.shape[0] - 1:
pattern[2, 2] = key(self.universes[j+1][i+1])
# Analyze lattice, find unique patterns in groups of universes
for j in range(self.shape[1]):
for i in range(self.shape[0]):
# Skip universes to ignore
if self.universes[j][i] in universes_to_ignore:
continue
# Create a neighborhood pattern based on the universe's
# neighbors in the grid, and lattice's neighbors at the edges
# Degenerate discretization has all universes be different
if strategy == "degenerate":
patterns[(i, j)] = {'locations': [(i, j)]}
continue
# Find neighbors among lattice's neighbors at the edges
find_edge_neighbors(pattern, i, j)
# Find neighbors among the lattice's universes
find_lattice_neighbors(pattern, i, j)
pattern[1, 1] = key(self.universes[j][i])
# Look for pattern in dictionary of patterns found
found = False
for known_pattern, pattern_data in patterns.items():
# Look at all rotations of pattern
for rot in range(4):
if not found and tuple(map(tuple, pattern)) ==\
known_pattern:
found = True
# Save location of the pattern in the lattice
pattern_data['locations'].append((i, j))
# Rotate pattern
pattern = np.rot90(pattern)
# Look at transpose of pattern and its rotations
pattern = np.transpose(pattern)
for rot in range(4):
if not found and tuple(map(tuple, pattern)) ==\
known_pattern:
found = True
# Save location of the pattern in the lattice
pattern_data['locations'].append((i, j))
# Rotate pattern
pattern = np.rot90(pattern)
# Transpose pattern back for the next search
pattern = np.transpose(pattern)
# Create new pattern and add to the patterns dictionary
if not found:
patterns[tuple(map(tuple, pattern))] =\
{'locations': [(i, j)]}
# Discretize lattice
for pattern, pattern_data in patterns.items():
first_pos = pattern_data['locations'][0]
# Create a clone of the universe, without cloning materials
new_universe = self.universes[first_pos[1]][first_pos[0]].clone(
clone_materials=False, clone_regions=False)
# Replace only the materials in materials_to_clone
for material in materials_to_clone:
material_cloned = False
for cell in new_universe.get_all_cells().values():
if cell.fill_type == 'material':
if cell.fill.id == material.id:
# Only a single clone of each material is necessary
if not material_cloned:
material_clone = material.clone()
material_cloned = True
cell.fill = material_clone
elif cell.fill_type == 'distribmat':
raise(ValueError, "Lattice discretization should not "
"be used with distributed materials")
elif len(cell.temperature) > 1 or len(cell.fill) > 1:
raise(ValueError, "Lattice discretization should not "
"be used with distributed cells")
# Rebuild lattice from list of locations with this pattern
for index, location in enumerate(pattern_data['locations']):
self.universes[location[1]][location[0]] = new_universe
def create_xml_subelement(self, xml_element, memo=None):
"""Add the lattice xml representation to an incoming xml element

View file

@ -482,12 +482,18 @@ class Universe(IDManagerMixin):
return universes
def clone(self, memo=None):
def clone(self, clone_materials=True, clone_regions=True, memo=None):
"""Create a copy of this universe with a new unique ID, and clones
all cells within this universe.
Parameters
----------
clone_materials : bool
Whether to create separates copies of the materials filling cells
contained in this universe.
clone_regions : bool
Whether to create separates copies of the regions bounding cells
contained in this universe.
memo : dict or None
A nested dictionary of previously cloned objects. This parameter
is used internally and should not be specified by the user.
@ -510,7 +516,8 @@ class Universe(IDManagerMixin):
# Clone all cells for the universe clone
clone._cells = OrderedDict()
for cell in self._cells.values():
clone.add_cell(cell.clone(memo))
clone.add_cell(cell.clone(clone_materials, clone_regions,
memo))
# Memoize the clone
memo[self] = clone

View file

@ -54,6 +54,18 @@ def test_clone():
assert c2.region != c.region
assert c2.temperature == c.temperature
c3 = c.clone(clone_materials=False)
assert c3.id != c.id
assert c3.fill == c.fill
assert c3.region != c.region
assert c3.temperature == c.temperature
c4 = c.clone(clone_regions=False)
assert c4.id != c.id
assert c4.fill != c.fill
assert c4.region == c.region
assert c4.temperature == c.temperature
def test_temperature(cell_with_lattice):
# Make sure temperature propagates through universes

View file

@ -278,6 +278,16 @@ def test_clone(rlat2, hlat2, hlat3):
assert hlat_clone.center == hlat3.center
assert hlat_clone.pitch == hlat3.pitch
rlat_clone = rlat2.clone(clone_materials=False)
assert rlat_clone.get_all_materials() == rlat2.get_all_materials()
rlat_clone = rlat2.clone(clone_materials=False, clone_regions=False)
for c1 in rlat_clone.cells:
for c2 in rlat2.cells:
if c1.fill == c2.fill:
print(c1.fill)
assert c1.region == c2.region
def test_repr(rlat2, rlat3, hlat2, hlat3):
repr(rlat2)

View file

@ -0,0 +1,91 @@
from math import sqrt
import xml.etree.ElementTree as ET
import openmc
import pytest
from tests.unit_tests.test_lattice import zr, pincell1, pincell2, rlat2
def test_discretization_degenerate(rlat2):
rlat_clone = rlat2.clone()
rlat_clone.discretize()
assert rlat_clone.get_universe((0, 0)).id !=\
rlat_clone.get_universe((1, 0)).id
assert rlat_clone.get_universe((1, 0)).id !=\
rlat_clone.get_universe((0, 1)).id
assert rlat_clone.get_universe((0, 1)).id !=\
rlat_clone.get_universe((0, 0)).id
assert rlat_clone.get_universe((0, 2)).id !=\
rlat_clone.get_universe((2, 0)).id
assert rlat_clone.get_universe((2, 1)).id !=\
rlat_clone.get_universe((1, 2)).id
def test_discretization_skip_universe(rlat2):
rlat_clone = rlat2.clone()
rlat_clone.discretize(
universes_to_ignore=[rlat_clone.get_universe((0, 0))])
assert rlat_clone.get_universe((0, 0)) == rlat_clone.get_universe((0, 1))
assert rlat_clone.get_universe((0, 1)) == rlat_clone.get_universe((2, 1))
assert rlat_clone.get_universe((0, 1)) == rlat_clone.get_universe((1, 2))
def test_discretization_clone_only_some_materials(rlat2):
rlat_clone = rlat2.clone()
fuel1 = next(iter(rlat_clone.get_universe((0, 0)).cells.values())).fill
rlat_clone.discretize(materials_to_clone=[fuel1])
assert next(reversed(rlat_clone.get_universe((0, 0)).cells.values())).fill\
== next(reversed(rlat_clone.get_universe((1, 0)).cells.values())).fill
assert next(iter(rlat_clone.get_universe((0, 0)).cells.values())).fill\
!= next(iter(rlat_clone.get_universe((1, 0)).cells.values())).fill
def test_discretization_lns(rlat2):
rlat_clone = rlat2.clone()
rlat_clone.discretize(strategy="lns")
assert rlat_clone.get_universe((0, 2)) == rlat_clone.get_universe((2, 2))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((0, 0))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((0, 1))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((1, 0))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((1, 1))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((1, 2))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((2, 0))
assert rlat_clone.get_universe((0, 2)) != rlat_clone.get_universe((2, 1))
def test_discretization_lns_using_names(rlat2):
rlat_clone = rlat2.clone()
rlat_clone.discretize(strategy="lns", key=lambda univ: univ.name)
assert rlat_clone.get_universe((0, 2)) == rlat_clone.get_universe((2, 2))
rlat_clone = rlat2.clone()
rlat_clone.get_universe((0, 1)).name="u1"
rlat_clone.get_universe((1, 0)).name="u1"
rlat_clone.discretize(strategy="lns", key=lambda univ: univ.name)
assert rlat_clone.get_universe((0, 2)) == rlat_clone.get_universe((2, 0))
assert rlat_clone.get_universe((2, 2)) == rlat_clone.get_universe((0, 0))
assert rlat_clone.get_universe((0, 2)) == rlat_clone.get_universe((2, 2))
assert rlat_clone.get_universe((1, 2)) == rlat_clone.get_universe((1, 0))
assert rlat_clone.get_universe((2, 1)) == rlat_clone.get_universe((0, 1))
assert rlat_clone.get_universe((1, 2)) == rlat_clone.get_universe((0, 1))
assert rlat_clone.get_universe((0, 0)) != rlat_clone.get_universe((0, 1))
def test_discretization_lns_with_neighbor_list(rlat2):
rlat_clone = rlat2.clone()
u1 = rlat_clone.get_universe((1, 0))
u2 = rlat_clone.get_universe((0, 0))
rlat_clone.discretize(strategy="lns",
lattice_neighbors=[u1,u1,u2,u1,u2,u2,u1,u2])
assert rlat_clone.get_universe((0, 1)) == rlat_clone.get_universe((0, 0))
assert rlat_clone.get_universe((1, 1)) != rlat_clone.get_universe((0, 0))

View file

@ -99,6 +99,28 @@ def test_get_all_universes():
assert not (univs ^ {u1, u2})
def test_clone():
c1 = openmc.Cell(cell_id=1)
c1.region = -openmc.ZCylinder(r=1.0)
c2 = openmc.Cell(cell_id=2)
c2.fill = openmc.Material()
c3 = openmc.Cell()
u1 = openmc.Universe(name='cool', cells=(c1, c2, c3))
u2 = u1.clone()
assert u2.name == u1.name
assert u2.cells != u1.cells
assert u2.get_all_materials() != u1.get_all_materials()
u2 = u1.clone(clone_materials=False)
assert u2.get_all_materials() == u1.get_all_materials()
u3 = u1.clone(clone_regions=False)
assert next(iter(u3.cells.values())).region ==\
next(iter(u1.cells.values())).region
def test_create_xml(cell_with_lattice):
cells = [openmc.Cell() for i in range(5)]
u = openmc.Universe(cells=cells)