From 022d8d3a44c0151dfad75c2baf72aa6f199dd08f Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Mon, 28 Dec 2015 23:14:05 -0500 Subject: [PATCH 1/5] Replace filter % offset w/ cell % distribcell_ind --- openmc/filter.py | 14 ---------- openmc/statepoint.py | 4 --- src/geometry_header.F90 | 9 ++++--- src/initialize.F90 | 58 ++++++++++++++--------------------------- src/output.F90 | 23 ++++++++-------- src/state_point.F90 | 1 - src/summary.F90 | 3 ++- src/tally.F90 | 6 +++-- src/tally_header.F90 | 1 - 9 files changed, 44 insertions(+), 75 deletions(-) diff --git a/openmc/filter.py b/openmc/filter.py index 04935b8ed..dce3d4924 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -42,8 +42,6 @@ class Filter(object): The number of filter bins mesh : Mesh or None A Mesh object for 'mesh' type filters. - offset : Integral - A value used to index tally bins for 'distribcell' tallies. stride : Integral The number of filter, nuclide and score bins within each of this filter's bins. @@ -57,7 +55,6 @@ class Filter(object): self._num_bins = 0 self._bins = None self._mesh = None - self._offset = -1 self._stride = None if type is not None: @@ -93,7 +90,6 @@ class Filter(object): clone._bins = copy.deepcopy(self.bins, memo) clone._num_bins = self.num_bins clone._mesh = copy.deepcopy(self.mesh, memo) - clone._offset = self.offset clone._stride = self.stride memo[id(self)] = clone @@ -108,7 +104,6 @@ class Filter(object): string = 'Filter\n' string += '{0: <16}{1}{2}\n'.format('\tType', '=\t', self.type) string += '{0: <16}{1}{2}\n'.format('\tBins', '=\t', self.bins) - string += '{0: <16}{1}{2}\n'.format('\tOffset', '=\t', self.offset) return string @property @@ -134,10 +129,6 @@ class Filter(object): def mesh(self): return self._mesh - @property - def offset(self): - return self._offset - @property def stride(self): return self._stride @@ -226,11 +217,6 @@ class Filter(object): self.type = 'mesh' self.bins = self.mesh.id - @offset.setter - def offset(self, offset): - cv.check_type('filter offset', offset, Integral) - self._offset = offset - @stride.setter def stride(self, stride): cv.check_type('filter stride', stride, Integral) diff --git a/openmc/statepoint.py b/openmc/statepoint.py index 2edf9badd..bc020cfbc 100644 --- a/openmc/statepoint.py +++ b/openmc/statepoint.py @@ -360,9 +360,6 @@ class StatePoint(object): # Read the Filter type filter_type = self._f['{0}{1}/type'.format(subbase, j)].value.decode() - # Read the Filter offset - offset = self._f['{0}{1}/offset'.format(subbase, j)].value - n_bins = self._f['{0}{1}/n_bins'.format(subbase, j)].value # Read the bin values @@ -370,7 +367,6 @@ class StatePoint(object): # Create Filter object filter = openmc.Filter(filter_type, bins) - filter.offset = offset filter.num_bins = n_bins if filter_type == 'mesh': diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 index 6ca13740b..684de883c 100644 --- a/src/geometry_header.F90 +++ b/src/geometry_header.F90 @@ -131,10 +131,13 @@ module geometry_header integer, allocatable :: offset (:) ! Distribcell offset for tally ! counter integer, allocatable :: region(:) ! Definition of spatial region as - ! Boolean expression of half-spaces + ! Boolean expression of half-spaces integer, allocatable :: rpn(:) ! Reverse Polish notation for region - ! expression - logical :: simple ! Is the region simple (intersections only) + ! expression + logical :: simple ! Is the region simple (intersections + ! only) + integer :: distribcell_ind ! Index corresponding to this cell in + ! distribcell arrays ! Rotation matrix and translation vector real(8), allocatable :: translation(:) diff --git a/src/initialize.F90 b/src/initialize.F90 index 9d63eb4e9..acfe92b92 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -1050,12 +1050,12 @@ contains do i = 1, n_tallies t => tallies(i) - do j = 1, t%n_filters - filter => t%filters(j) + do j = 1, t % n_filters + filter => t % filters(j) - if (filter%type == FILTER_DISTRIBCELL) then - if (.not. cell_list%contains(filter%int_bins(1))) then - call cell_list%add(filter%int_bins(1)) + if (filter % type == FILTER_DISTRIBCELL) then + if (.not. cell_list % contains(filter % int_bins(1))) then + call cell_list % add(filter % int_bins(1)) end if end if @@ -1066,8 +1066,8 @@ contains ! to determine the number of offset tables to allocate do i = 1, n_universes univ => universes(i) - do j = 1, univ%n_cells - if (cell_list%contains(univ%cells(j))) then + do j = 1, univ % n_cells + if (cell_list % contains(univ % cells(j))) then n_maps = n_maps + 1 end if end do @@ -1089,56 +1089,38 @@ contains do i = 1, n_universes univ => universes(i) - do j = 1, univ%n_cells + do j = 1, univ % n_cells + if (.not. cell_list % contains(univ % cells(j))) cycle - if (cell_list%contains(univ%cells(j))) then + cells(univ % cells(j)) % distribcell_ind = k - ! Loop over all tallies - do l = 1, n_tallies - t => tallies(l) - - do m = 1, t%n_filters - filter => t%filters(m) - - ! Loop over only distribcell filters - ! If filter points to cell we just found, set offset index - if (filter%type == FILTER_DISTRIBCELL) then - if (filter%int_bins(1) == univ%cells(j)) then - filter%offset = k - end if - end if - - end do - end do - - univ_list(k) = univ%id - k = k + 1 - end if + univ_list(k) = univ % id + k = k + 1 end do end do ! Allocate the offset tables for lattices do i = 1, n_lattices - lat => lattices(i)%obj + lat => lattices(i) % obj select type(lat) type is (RectLattice) - allocate(lat%offset(n_maps, lat%n_cells(1), lat%n_cells(2), & - lat%n_cells(3))) + allocate(lat % offset(n_maps, lat % n_cells(1), lat % n_cells(2), & + lat % n_cells(3))) type is (HexLattice) - allocate(lat%offset(n_maps, 2 * lat%n_rings - 1, & - 2 * lat%n_rings - 1, lat%n_axial)) + allocate(lat % offset(n_maps, 2 * lat % n_rings - 1, & + 2 * lat % n_rings - 1, lat % n_axial)) end select - lat%offset(:, :, :, :) = 0 + lat % offset(:, :, :, :) = 0 end do ! Allocate offset table for fill cells do i = 1, n_cells - if (cells(i)%material == NONE) then - allocate(cells(i)%offset(n_maps)) + if (cells(i) % material == NONE) then + allocate(cells(i) % offset(n_maps)) end if end do diff --git a/src/output.F90 b/src/output.F90 index c31b20b70..8be6c09da 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1378,8 +1378,7 @@ contains label = '' univ => universes(BASE_UNIVERSE) offset = 0 - call find_offset(t % filters(i_filter) % offset, & - t % filters(i_filter) % int_bins(1), & + call find_offset(t % filters(i_filter) % int_bins(1), & univ, bin-1, offset, label) case (FILTER_SURFACE) i = t % filters(i_filter) % int_bins(bin) @@ -1413,15 +1412,15 @@ contains ! with the given offset !=============================================================================== - recursive subroutine find_offset(map, goal, univ, final, offset, path) + recursive subroutine find_offset(goal, univ, final, offset, path) - integer, intent(in) :: map ! Index in maps vector - integer, intent(in) :: goal ! The target cell ID + integer, intent(in) :: goal ! The target cell index type(Universe), intent(in) :: univ ! Universe to begin search integer, intent(in) :: final ! Target offset integer, intent(inout) :: offset ! Current offset character(*), intent(inout) :: path ! Path to offset + integer :: map ! Index in maps vector integer :: i, j ! Index over cells integer :: k, l, m ! Indices in lattice integer :: old_k, old_l, old_m ! Previous indices in lattice @@ -1436,6 +1435,9 @@ contains type(Universe), pointer :: next_univ ! Next universe to loop through class(Lattice), pointer :: lat ! Pointer to current lattice + ! Get the distribcell index for this cell + map = cells(goal) % distribcell_ind + n = univ % n_cells ! Write to the geometry stack @@ -1537,7 +1539,7 @@ contains offset = c % offset(map) + offset next_univ => universes(c % fill) - call find_offset(map, goal, next_univ, final, offset, path) + call find_offset(goal, next_univ, final, offset, path) return ! ==================================================================== @@ -1577,7 +1579,7 @@ contains path = trim(path) // "(" // trim(to_str(k)) // & "," // trim(to_str(l)) // "," // & trim(to_str(m)) // ")" - call find_offset(map, goal, next_univ, final, offset, path) + call find_offset(goal, next_univ, final, offset, path) return else old_m = m @@ -1593,7 +1595,7 @@ contains path = trim(path) // "(" // trim(to_str(old_k)) // & "," // trim(to_str(old_l)) // "," // & trim(to_str(old_m)) // ")" - call find_offset(map, goal, next_univ, final, offset, path) + call find_offset(goal, next_univ, final, offset, path) return end if @@ -1638,8 +1640,7 @@ contains trim(to_str(k - lat % n_rings)) // "," // & trim(to_str(l - lat % n_rings)) // "," // & trim(to_str(m)) // ")" - call find_offset(map, goal, next_univ, final, offset, & - path) + call find_offset(goal, next_univ, final, offset, path) return else old_m = m @@ -1656,7 +1657,7 @@ contains trim(to_str(old_k - lat % n_rings)) // "," // & trim(to_str(old_l - lat % n_rings)) // "," // & trim(to_str(old_m)) // ")" - call find_offset(map, goal, next_univ, final, offset, path) + call find_offset(goal, next_univ, final, offset, path) return end if diff --git a/src/state_point.F90 b/src/state_point.F90 index f9f4b5b75..046c857fd 100644 --- a/src/state_point.F90 +++ b/src/state_point.F90 @@ -268,7 +268,6 @@ contains call write_dataset(filter_group, "type", "delayedgroup") end select - call write_dataset(filter_group, "offset", tally%filters(j)%offset) call write_dataset(filter_group, "n_bins", tally%filters(j)%n_bins) if (tally % filters(j) % type == FILTER_ENERGYIN .or. & tally % filters(j) % type == FILTER_ENERGYOUT .or. & diff --git a/src/summary.F90 b/src/summary.F90 index cb005d3be..11c4cd434 100644 --- a/src/summary.F90 +++ b/src/summary.F90 @@ -196,6 +196,8 @@ contains end do call write_dataset(cell_group, "region", adjustl(region_spec)) + call write_dataset(cell_group, "distribcell_ind", c % distribcell_ind) + call close_group(cell_group) end do CELL_LOOP @@ -540,7 +542,6 @@ contains filter_group = create_group(tally_group, "filter " // trim(to_str(j))) ! Write number of bins for this filter - call write_dataset(filter_group, "offset", t%filters(j)%offset) call write_dataset(filter_group, "n_bins", t%filters(j)%n_bins) ! Write filter bins diff --git a/src/tally.F90 b/src/tally.F90 index 09e30a242..2ccd59a58 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -1712,6 +1712,7 @@ contains integer :: j integer :: n ! number of bins for single filter integer :: offset ! offset for distribcell + integer :: distribcell_ind ! index in distribcell arrays real(8) :: E ! particle energy real(8) :: theta, phi ! Polar and Azimuthal Angles, respectively type(TallyObject), pointer :: t @@ -1756,12 +1757,13 @@ contains case (FILTER_DISTRIBCELL) ! determine next distribcell bin + distribcell_ind = cells(t % filters(i) % int_bins(1)) % distribcell_ind matching_bins(i) = NO_BIN_FOUND offset = 0 do j = 1, p % n_coord if (cells(p % coord(j) % cell) % type == CELL_FILL) then offset = offset + cells(p % coord(j) % cell) % & - offset(t % filters(i) % offset) + offset(distribcell_ind) elseif(cells(p % coord(j) % cell) % type == CELL_LATTICE) then if (lattices(p % coord(j + 1) % lattice) % obj & % are_valid_indices([& @@ -1769,7 +1771,7 @@ contains p % coord(j + 1) % lattice_y, & p % coord(j + 1) % lattice_z])) then offset = offset + lattices(p % coord(j + 1) % lattice) % obj % & - offset(t % filters(i) % offset, & + offset(distribcell_ind, & p % coord(j + 1) % lattice_x, & p % coord(j + 1) % lattice_y, & p % coord(j + 1) % lattice_z) diff --git a/src/tally_header.F90 b/src/tally_header.F90 index 18b921952..dcbba3d89 100644 --- a/src/tally_header.F90 +++ b/src/tally_header.F90 @@ -55,7 +55,6 @@ module tally_header type TallyFilter integer :: type = NONE integer :: n_bins = 0 - integer :: offset = 0 ! Only used for distribcell filters integer, allocatable :: int_bins(:) real(8), allocatable :: real_bins(:) ! Only used for energy filters end type TallyFilter From 56e8e7d6818937cdcb57fc208532025e52768a68 Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Tue, 29 Dec 2015 16:25:29 -0500 Subject: [PATCH 2/5] Update PyAPI with distribcell changes --- openmc/filter.py | 2 +- openmc/geometry.py | 20 +++++++++++++------- openmc/summary.py | 5 +++++ openmc/universe.py | 37 ++++++++++++++++++++++++++----------- src/input_xml.F90 | 3 ++- src/output.F90 | 14 +++++--------- 6 files changed, 52 insertions(+), 29 deletions(-) diff --git a/openmc/filter.py b/openmc/filter.py index dce3d4924..ea92f5257 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -609,7 +609,7 @@ class Filter(object): # If this region is in Cell corresponding to the # distribcell filter bin, store it in dictionary if cell_id == self.bins[0]: - offset = openmc_geometry.get_offset(path, self.offset) + offset = openmc_geometry.get_offset(path) offsets_to_coords[offset] = coords # Each distribcell offset is a DataFrame bin diff --git a/openmc/geometry.py b/openmc/geometry.py index e848e0cdd..16f6a45e5 100644 --- a/openmc/geometry.py +++ b/openmc/geometry.py @@ -42,7 +42,7 @@ class Geometry(object): self._root_universe = root_universe - def get_offset(self, path, filter_offset): + def get_offset(self, path): """Returns the corresponding location in the results array for a given path and filter number. This is primarily intended to post-processing result when a distribcell filter is used. @@ -55,8 +55,6 @@ class Geometry(object): lattice passed through. For the case of the lattice, a tuple should be provided to indicate which coordinates in the lattice should be entered. This should be in the form: (lat_id, i_x, i_y, i_z) - filter_offset : int - An integer that specifies which offset map the filter is using Returns ------- @@ -65,14 +63,22 @@ class Geometry(object): """ + # Find the distribcell index of the cell. + cells = self.get_all_cells() + if path[-1] in cells: + distribcell_ind = cells[path[-1]].distribcell_ind + else: + raise RuntimeError('Could not find cell {} specified in a \ + distribcell filter'.format(path[-1])) + # Return memoize'd offset if possible - if (path, filter_offset) in self._offsets: - offset = self._offsets[(path, filter_offset)] + if (path, distribcell_ind) in self._offsets: + offset = self._offsets[(path, distribcell_ind)] # Begin recursive call to compute offset starting with the base Universe else: - offset = self._root_universe.get_offset(path, filter_offset) - self._offsets[(path, filter_offset)] = offset + offset = self._root_universe.get_offset(path, distribcell_ind) + self._offsets[(path, distribcell_ind)] = offset # Return the final offset return offset diff --git a/openmc/summary.py b/openmc/summary.py index bc6551e7c..9e2145b34 100644 --- a/openmc/summary.py +++ b/openmc/summary.py @@ -271,6 +271,11 @@ class Summary(object): cell.region = Region.from_expression( region, {s.id: s for s in self.surfaces.values()}) + # Get the distribcell index + ind = self._f['geometry/cells'][key]['distribcell_ind'].value + if ind != 0: + cell.distribcell_ind = ind + # Add the Cell to the global dictionary of all Cells self.cells[index] = cell diff --git a/openmc/universe.py b/openmc/universe.py index 98367c381..b346ffd00 100644 --- a/openmc/universe.py +++ b/openmc/universe.py @@ -63,6 +63,8 @@ class Cell(object): that is used to translate (shift) the universe. offsets : ndarray Array of offsets used for distributed cell searches + distribcell_ind : int + Index of this cell in distribcell arrays """ @@ -76,6 +78,7 @@ class Cell(object): self._rotation = None self._translation = None self._offsets = None + self._distribcell_ind = None def __eq__(self, other): if not isinstance(other, Cell): @@ -122,6 +125,8 @@ class Cell(object): string += '{0: <16}{1}{2}\n'.format('\tTranslation', '=\t', self._translation) string += '{0: <16}{1}{2}\n'.format('\tOffset', '=\t', self._offsets) + string += '{0: <16}{1}{2}\n'.format('\tDistribcell index', '=\t', + self._distribcell_ind) return string @@ -164,6 +169,10 @@ class Cell(object): def offsets(self): return self._offsets + @property + def distribcell_ind(self): + return self._distribcell_ind + @id.setter def id(self, cell_id): if cell_id is None: @@ -231,6 +240,11 @@ class Cell(object): cv.check_type('cell region', region, Region) self._region = region + @distribcell_ind.setter + def distribcell_ind(self, ind): + cv.check_type('distribcell index', ind, Integral) + self._distribcell_ind = ind + def add_surface(self, surface, halfspace): """Add a half-space to the list of half-spaces whose intersection defines the cell. @@ -271,7 +285,7 @@ class Cell(object): else: self.region = Intersection(self.region, region) - def get_offset(self, path, filter_offset): + def get_offset(self, path, distribcell_ind): # Get the current element and remove it from the list cell_id = path[0] path = path[1:] @@ -282,12 +296,12 @@ class Cell(object): # If the Cell is filled by a Universe elif self._type == 'fill': - offset = self._offsets[filter_offset-1] - offset += self._fill.get_offset(path, filter_offset) + offset = self._offsets[distribcell_ind-1] + offset += self._fill.get_offset(path, distribcell_ind) # If the Cell is filled by a Lattice else: - offset = self._fill.get_offset(path, filter_offset) + offset = self._fill.get_offset(path, distribcell_ind) return offset @@ -591,7 +605,7 @@ class Universe(object): self._cells.clear() - def get_offset(self, path, filter_offset): + def get_offset(self, path, distribcell_ind): # Get the current element and remove it from the list path = path[1:] @@ -599,7 +613,7 @@ class Universe(object): cell_id = path[0] # Make a recursive call to the Cell within this Universe - offset = self._cells[cell_id].get_offset(path, filter_offset) + offset = self._cells[cell_id].get_offset(path, distribcell_ind) # Return the offset computed at all nested Universe levels return offset @@ -1059,21 +1073,22 @@ class RectLattice(Lattice): cv.check_greater_than('lattice pitch', dim, 0.0) self._pitch = pitch - def get_offset(self, path, filter_offset): + def get_offset(self, path, distribcell_ind): # Get the current element and remove it from the list i = path[0] path = path[1:] # For 2D Lattices if len(self._dimension) == 2: - offset = self._offsets[i[1]-1, i[2]-1, 0, filter_offset-1] - offset += self._universes[i[1]][i[2]].get_offset(path, filter_offset) + offset = self._offsets[i[1]-1, i[2]-1, 0, distribcell_ind-1] + offset += self._universes[i[1]][i[2]].get_offset(path, + distribcell_ind) # For 3D Lattices else: - offset = self._offsets[i[1]-1, i[2]-1, i[3]-1, filter_offset-1] + offset = self._offsets[i[1]-1, i[2]-1, i[3]-1, distribcell_ind-1] offset += self._universes[i[1]-1][i[2]-1][i[3]-1].get_offset(path, - filter_offset) + distribcell_ind) return offset diff --git a/src/input_xml.F90 b/src/input_xml.F90 index abd6e2909..a4317ed8e 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -1049,8 +1049,9 @@ contains do i = 1, n_cells c => cells(i) - ! Initialize the number of cell instances - this is a base case for distribcells + ! Initialize distribcell instances and distribcell index c % instances = 0 + c % distribcell_ind = NONE ! Get pointer to i-th cell node call get_list_item(node_cell_list, i, node_cell) diff --git a/src/output.F90 b/src/output.F90 index 8be6c09da..5b938e98d 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1449,17 +1449,13 @@ contains ! Look through all cells in this universe do i = 1, n - - cell_index = univ % cells(i) - c => cells(cell_index) - - ! If the cell ID matches the goal and the offset matches final, - ! write to the geometry stack - if (cell_dict % get_key(c % id) == goal .AND. offset == final) then - path = trim(path) // "->" // to_str(c%id) + ! If the cell matches the goal and the offset matches final, write to the + ! geometry stack + if (univ % cells(i) == goal .AND. offset == final) then + c => cells(univ % cells(i)) + path = trim(path) // "->" // to_str(c % id) return end if - end do ! Find the fill cell or lattice cell that we need to enter From 33375dfa30f80de5fe69dcfb47425d21642e479a Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Wed, 30 Dec 2015 15:57:44 -0500 Subject: [PATCH 3/5] Small fixes for #538 --- openmc/filter.py | 2 +- openmc/geometry.py | 5 +++-- openmc/universe.py | 22 +++++++++++----------- src/initialize.F90 | 14 +++++++------- 4 files changed, 22 insertions(+), 21 deletions(-) diff --git a/openmc/filter.py b/openmc/filter.py index ea92f5257..d61fb648f 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -609,7 +609,7 @@ class Filter(object): # If this region is in Cell corresponding to the # distribcell filter bin, store it in dictionary if cell_id == self.bins[0]: - offset = openmc_geometry.get_offset(path) + offset = openmc_geometry.get_cell_instance(path) offsets_to_coords[offset] = coords # Each distribcell offset is a DataFrame bin diff --git a/openmc/geometry.py b/openmc/geometry.py index 16f6a45e5..baadfce3c 100644 --- a/openmc/geometry.py +++ b/openmc/geometry.py @@ -42,7 +42,7 @@ class Geometry(object): self._root_universe = root_universe - def get_offset(self, path): + def get_cell_instance(self, path): """Returns the corresponding location in the results array for a given path and filter number. This is primarily intended to post-processing result when a distribcell filter is used. @@ -77,7 +77,8 @@ class Geometry(object): # Begin recursive call to compute offset starting with the base Universe else: - offset = self._root_universe.get_offset(path, distribcell_ind) + offset = self._root_universe.get_cell_instance(path, + distribcell_ind) self._offsets[(path, distribcell_ind)] = offset # Return the final offset diff --git a/openmc/universe.py b/openmc/universe.py index b346ffd00..e7b2151a3 100644 --- a/openmc/universe.py +++ b/openmc/universe.py @@ -285,7 +285,7 @@ class Cell(object): else: self.region = Intersection(self.region, region) - def get_offset(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_ind): # Get the current element and remove it from the list cell_id = path[0] path = path[1:] @@ -296,12 +296,12 @@ class Cell(object): # If the Cell is filled by a Universe elif self._type == 'fill': - offset = self._offsets[distribcell_ind-1] - offset += self._fill.get_offset(path, distribcell_ind) + offset = self.offsets[distribcell_ind-1] + offset += self.fill.get_cell_instance(path, distribcell_ind) # If the Cell is filled by a Lattice else: - offset = self._fill.get_offset(path, distribcell_ind) + offset = self.fill.get_cell_instance(path, distribcell_ind) return offset @@ -605,7 +605,7 @@ class Universe(object): self._cells.clear() - def get_offset(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_ind): # Get the current element and remove it from the list path = path[1:] @@ -613,7 +613,7 @@ class Universe(object): cell_id = path[0] # Make a recursive call to the Cell within this Universe - offset = self._cells[cell_id].get_offset(path, distribcell_ind) + offset = self.cells[cell_id].get_cell_instance(path, distribcell_ind) # Return the offset computed at all nested Universe levels return offset @@ -1073,7 +1073,7 @@ class RectLattice(Lattice): cv.check_greater_than('lattice pitch', dim, 0.0) self._pitch = pitch - def get_offset(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_ind): # Get the current element and remove it from the list i = path[0] path = path[1:] @@ -1081,14 +1081,14 @@ class RectLattice(Lattice): # For 2D Lattices if len(self._dimension) == 2: offset = self._offsets[i[1]-1, i[2]-1, 0, distribcell_ind-1] - offset += self._universes[i[1]][i[2]].get_offset(path, - distribcell_ind) + offset += self._universes[i[1]][i[2]].get_cell_instance(path, + distribcell_ind) # For 3D Lattices else: offset = self._offsets[i[1]-1, i[2]-1, i[3]-1, distribcell_ind-1] - offset += self._universes[i[1]-1][i[2]-1][i[3]-1].get_offset(path, - distribcell_ind) + offset += self._universes[i[1]-1][i[2]-1][i[3]-1].get_cell_instance( + path, distribcell_ind) return offset diff --git a/src/initialize.F90 b/src/initialize.F90 index acfe92b92..1134040e9 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -1086,16 +1086,16 @@ contains found(:,:) = .false. k = 1 + ! Search through universes for distributed cells and assign each one a + ! unique distribcell array index. do i = 1, n_universes univ => universes(i) - do j = 1, univ % n_cells - if (.not. cell_list % contains(univ % cells(j))) cycle - - cells(univ % cells(j)) % distribcell_ind = k - - univ_list(k) = univ % id - k = k + 1 + if (cell_list % contains(univ % cells(j))) then + cells(univ % cells(j)) % distribcell_ind = k + univ_list(k) = univ % id + k = k + 1 + end if end do end do From 7fea70650f576674f7663fbbbdee2f543956293b Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Sun, 3 Jan 2016 13:44:42 -0500 Subject: [PATCH 4/5] Improve geometry.get_cell_instance docstring --- openmc/geometry.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/openmc/geometry.py b/openmc/geometry.py index baadfce3c..eb2b9db55 100644 --- a/openmc/geometry.py +++ b/openmc/geometry.py @@ -43,9 +43,9 @@ class Geometry(object): self._root_universe = root_universe def get_cell_instance(self, path): - """Returns the corresponding location in the results array for a given path and - filter number. This is primarily intended to post-processing result when - a distribcell filter is used. + """Return the instance number for the final cell in a geometry path. + + The instance is an index into tally distribcell filter arrays. Parameters ---------- @@ -58,8 +58,8 @@ class Geometry(object): Returns ------- - offset : int - Location in the results array for the path and filter + instance : int + Index in tally results array for distribcell filters """ From bae0ad713683bff2e6bc2bc234c0ec91f66709eb Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Mon, 4 Jan 2016 18:25:59 -0500 Subject: [PATCH 5/5] Small fixes for #538 --- docs/source/usersguide/output/statepoint.rst | 4 --- docs/source/usersguide/output/summary.rst | 4 +++ openmc/geometry.py | 10 +++--- openmc/summary.py | 4 +-- openmc/universe.py | 38 ++++++++++---------- src/geometry_header.F90 | 2 +- src/initialize.F90 | 2 +- src/input_xml.F90 | 2 +- src/output.F90 | 2 +- src/summary.F90 | 2 +- src/tally.F90 | 9 ++--- 11 files changed, 40 insertions(+), 39 deletions(-) diff --git a/docs/source/usersguide/output/statepoint.rst b/docs/source/usersguide/output/statepoint.rst index d3c1729af..15bc79f73 100644 --- a/docs/source/usersguide/output/statepoint.rst +++ b/docs/source/usersguide/output/statepoint.rst @@ -185,10 +185,6 @@ if run_mode == 'k-eigenvalue': Type of the j-th filter. Can be 'universe', 'material', 'cell', 'cellborn', 'surface', 'mesh', 'energy', 'energyout', or 'distribcell'. -**/tallies/tally /filter /offset** (*int*) - - Filter offset (used for distribcell filter). - **/tallies/tally /filter /n_bins** (*int*) Number of bins for the j-th filter. diff --git a/docs/source/usersguide/output/summary.rst b/docs/source/usersguide/output/summary.rst index f87f60c4a..302f34079 100644 --- a/docs/source/usersguide/output/summary.rst +++ b/docs/source/usersguide/output/summary.rst @@ -121,6 +121,10 @@ The current revision of the summary file format is 1. Region specification for the cell. +**/geometry/cells/cell /distribcell_index** (*int*) + + Index of this cell in distribcell filter arrays. + **/geometry/surfaces/surface /index** (*int*) Index in surfaces array used internally in OpenMC. diff --git a/openmc/geometry.py b/openmc/geometry.py index eb2b9db55..fc8e19f07 100644 --- a/openmc/geometry.py +++ b/openmc/geometry.py @@ -66,20 +66,20 @@ class Geometry(object): # Find the distribcell index of the cell. cells = self.get_all_cells() if path[-1] in cells: - distribcell_ind = cells[path[-1]].distribcell_ind + distribcell_index = cells[path[-1]].distribcell_index else: raise RuntimeError('Could not find cell {} specified in a \ distribcell filter'.format(path[-1])) # Return memoize'd offset if possible - if (path, distribcell_ind) in self._offsets: - offset = self._offsets[(path, distribcell_ind)] + if (path, distribcell_index) in self._offsets: + offset = self._offsets[(path, distribcell_index)] # Begin recursive call to compute offset starting with the base Universe else: offset = self._root_universe.get_cell_instance(path, - distribcell_ind) - self._offsets[(path, distribcell_ind)] = offset + distribcell_index) + self._offsets[(path, distribcell_index)] = offset # Return the final offset return offset diff --git a/openmc/summary.py b/openmc/summary.py index 9e2145b34..1c4c69aab 100644 --- a/openmc/summary.py +++ b/openmc/summary.py @@ -272,9 +272,9 @@ class Summary(object): region, {s.id: s for s in self.surfaces.values()}) # Get the distribcell index - ind = self._f['geometry/cells'][key]['distribcell_ind'].value + ind = self._f['geometry/cells'][key]['distribcell_index'].value if ind != 0: - cell.distribcell_ind = ind + cell.distribcell_index = ind # Add the Cell to the global dictionary of all Cells self.cells[index] = cell diff --git a/openmc/universe.py b/openmc/universe.py index e7b2151a3..8c50ddbc7 100644 --- a/openmc/universe.py +++ b/openmc/universe.py @@ -63,7 +63,7 @@ class Cell(object): that is used to translate (shift) the universe. offsets : ndarray Array of offsets used for distributed cell searches - distribcell_ind : int + distribcell_index : int Index of this cell in distribcell arrays """ @@ -78,7 +78,7 @@ class Cell(object): self._rotation = None self._translation = None self._offsets = None - self._distribcell_ind = None + self._distribcell_index = None def __eq__(self, other): if not isinstance(other, Cell): @@ -126,7 +126,7 @@ class Cell(object): self._translation) string += '{0: <16}{1}{2}\n'.format('\tOffset', '=\t', self._offsets) string += '{0: <16}{1}{2}\n'.format('\tDistribcell index', '=\t', - self._distribcell_ind) + self._distribcell_index) return string @@ -170,8 +170,8 @@ class Cell(object): return self._offsets @property - def distribcell_ind(self): - return self._distribcell_ind + def distribcell_index(self): + return self._distribcell_index @id.setter def id(self, cell_id): @@ -240,10 +240,10 @@ class Cell(object): cv.check_type('cell region', region, Region) self._region = region - @distribcell_ind.setter - def distribcell_ind(self, ind): + @distribcell_index.setter + def distribcell_index(self, ind): cv.check_type('distribcell index', ind, Integral) - self._distribcell_ind = ind + self._distribcell_index = ind def add_surface(self, surface, halfspace): """Add a half-space to the list of half-spaces whose intersection defines the @@ -285,7 +285,7 @@ class Cell(object): else: self.region = Intersection(self.region, region) - def get_cell_instance(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_index): # Get the current element and remove it from the list cell_id = path[0] path = path[1:] @@ -296,12 +296,12 @@ class Cell(object): # If the Cell is filled by a Universe elif self._type == 'fill': - offset = self.offsets[distribcell_ind-1] - offset += self.fill.get_cell_instance(path, distribcell_ind) + offset = self.offsets[distribcell_index-1] + offset += self.fill.get_cell_instance(path, distribcell_index) # If the Cell is filled by a Lattice else: - offset = self.fill.get_cell_instance(path, distribcell_ind) + offset = self.fill.get_cell_instance(path, distribcell_index) return offset @@ -605,7 +605,7 @@ class Universe(object): self._cells.clear() - def get_cell_instance(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_index): # Get the current element and remove it from the list path = path[1:] @@ -613,7 +613,7 @@ class Universe(object): cell_id = path[0] # Make a recursive call to the Cell within this Universe - offset = self.cells[cell_id].get_cell_instance(path, distribcell_ind) + offset = self.cells[cell_id].get_cell_instance(path, distribcell_index) # Return the offset computed at all nested Universe levels return offset @@ -1073,22 +1073,22 @@ class RectLattice(Lattice): cv.check_greater_than('lattice pitch', dim, 0.0) self._pitch = pitch - def get_cell_instance(self, path, distribcell_ind): + def get_cell_instance(self, path, distribcell_index): # Get the current element and remove it from the list i = path[0] path = path[1:] # For 2D Lattices if len(self._dimension) == 2: - offset = self._offsets[i[1]-1, i[2]-1, 0, distribcell_ind-1] + offset = self._offsets[i[1]-1, i[2]-1, 0, distribcell_index-1] offset += self._universes[i[1]][i[2]].get_cell_instance(path, - distribcell_ind) + distribcell_index) # For 3D Lattices else: - offset = self._offsets[i[1]-1, i[2]-1, i[3]-1, distribcell_ind-1] + offset = self._offsets[i[1]-1, i[2]-1, i[3]-1, distribcell_index-1] offset += self._universes[i[1]-1][i[2]-1][i[3]-1].get_cell_instance( - path, distribcell_ind) + path, distribcell_index) return offset diff --git a/src/geometry_header.F90 b/src/geometry_header.F90 index 684de883c..3dad3ba39 100644 --- a/src/geometry_header.F90 +++ b/src/geometry_header.F90 @@ -136,7 +136,7 @@ module geometry_header ! expression logical :: simple ! Is the region simple (intersections ! only) - integer :: distribcell_ind ! Index corresponding to this cell in + integer :: distribcell_index ! Index corresponding to this cell in ! distribcell arrays ! Rotation matrix and translation vector diff --git a/src/initialize.F90 b/src/initialize.F90 index 1134040e9..30654be9d 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -1092,7 +1092,7 @@ contains univ => universes(i) do j = 1, univ % n_cells if (cell_list % contains(univ % cells(j))) then - cells(univ % cells(j)) % distribcell_ind = k + cells(univ % cells(j)) % distribcell_index = k univ_list(k) = univ % id k = k + 1 end if diff --git a/src/input_xml.F90 b/src/input_xml.F90 index a4317ed8e..07d739b88 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -1051,7 +1051,7 @@ contains ! Initialize distribcell instances and distribcell index c % instances = 0 - c % distribcell_ind = NONE + c % distribcell_index = NONE ! Get pointer to i-th cell node call get_list_item(node_cell_list, i, node_cell) diff --git a/src/output.F90 b/src/output.F90 index 5b938e98d..66e683137 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1436,7 +1436,7 @@ contains class(Lattice), pointer :: lat ! Pointer to current lattice ! Get the distribcell index for this cell - map = cells(goal) % distribcell_ind + map = cells(goal) % distribcell_index n = univ % n_cells diff --git a/src/summary.F90 b/src/summary.F90 index 11c4cd434..7a430cf74 100644 --- a/src/summary.F90 +++ b/src/summary.F90 @@ -196,7 +196,7 @@ contains end do call write_dataset(cell_group, "region", adjustl(region_spec)) - call write_dataset(cell_group, "distribcell_ind", c % distribcell_ind) + call write_dataset(cell_group, "distribcell_index", c % distribcell_index) call close_group(cell_group) end do CELL_LOOP diff --git a/src/tally.F90 b/src/tally.F90 index 2ccd59a58..7210efc5c 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -1712,7 +1712,7 @@ contains integer :: j integer :: n ! number of bins for single filter integer :: offset ! offset for distribcell - integer :: distribcell_ind ! index in distribcell arrays + integer :: distribcell_index ! index in distribcell arrays real(8) :: E ! particle energy real(8) :: theta, phi ! Polar and Azimuthal Angles, respectively type(TallyObject), pointer :: t @@ -1757,13 +1757,14 @@ contains case (FILTER_DISTRIBCELL) ! determine next distribcell bin - distribcell_ind = cells(t % filters(i) % int_bins(1)) % distribcell_ind + distribcell_index = cells(t % filters(i) % int_bins(1)) & + % distribcell_index matching_bins(i) = NO_BIN_FOUND offset = 0 do j = 1, p % n_coord if (cells(p % coord(j) % cell) % type == CELL_FILL) then offset = offset + cells(p % coord(j) % cell) % & - offset(distribcell_ind) + offset(distribcell_index) elseif(cells(p % coord(j) % cell) % type == CELL_LATTICE) then if (lattices(p % coord(j + 1) % lattice) % obj & % are_valid_indices([& @@ -1771,7 +1772,7 @@ contains p % coord(j + 1) % lattice_y, & p % coord(j + 1) % lattice_z])) then offset = offset + lattices(p % coord(j + 1) % lattice) % obj % & - offset(distribcell_ind, & + offset(distribcell_index, & p % coord(j + 1) % lattice_x, & p % coord(j + 1) % lattice_y, & p % coord(j + 1) % lattice_z)