diff --git a/openmc/cell.py b/openmc/cell.py index 9055f10e59..23246a0f68 100644 --- a/openmc/cell.py +++ b/openmc/cell.py @@ -85,6 +85,10 @@ class Cell(object): Array of offsets used for distributed cell searches distribcell_index : int Index of this cell in distribcell arrays + volume_information : dict + Estimate of the volume and total number of atoms of each nuclide from a + stochastic volume calculation. This information is set with the + :meth:`Cell.add_volume_information` method. """ @@ -100,6 +104,7 @@ class Cell(object): self._translation = None self._offsets = None self._distribcell_index = None + self._volume_information = None def __contains__(self, point): if self.region is None: @@ -212,6 +217,10 @@ class Cell(object): def distribcell_index(self): return self._distribcell_index + @property + def volume_information(self): + return self._volume_information + @id.setter def id(self, cell_id): if cell_id is None: @@ -351,6 +360,22 @@ class Cell(object): else: self.region = Intersection(self.region, region) + def add_volume_information(self, volume_calc): + """Add volume information to a cell. + + Parameters + ---------- + volume_calc : openmc.VolumeCalculation + Results from a stochastic volume calculation + + """ + for cell_id in volume_calc.results: + if cell_id == self.id: + self._volume_information = volume_calc.results[cell_id] + break + else: + raise ValueError('No volume information found for this cell.') + def get_cell_instance(self, path, distribcell_index): # If the Cell is filled by a Material @@ -368,8 +393,19 @@ class Cell(object): return offset - def get_all_nuclides(self): - """Return all nuclides contained in the cell + def get_nuclides(self): + """Returns all nuclides in the cell + + Returns + ------- + nuclides : list + List of nuclide names + + """ + return self.fill.get_nuclides() if self.fill_type != 'void' else [] + + def get_nuclide_densities(self): + """Return all nuclides contained in the cell and their densities Returns ------- @@ -381,8 +417,24 @@ class Cell(object): nuclides = OrderedDict() - if self.fill_type != 'void': - nuclides.update(self.fill.get_all_nuclides()) + if self.fill_type == 'material': + nuclides.update(self.fill.get_nuclide_densities()) + elif self.fill_type == 'void': + pass + else: + if self.volume_information is not None: + volume = self.volume_information['volume'][0] + for full_name, atoms in self.volume_information['atoms']: + name, xs = full_name.split('.') + nuclide = openmc.Nuclide(name, xs) + density = 1.0e-24 * atoms[0]/volume # density in atoms/b-cm + nuclides[name] = (nuclide, density) + else: + raise RuntimeError( + 'Volume information is needed to calculate microscopic cross ' + 'sections for cell {}. This can be done by running a ' + 'stochastic volume calculation via the ' + 'openmc.VolumeCalculation object'.format(self.id)) return nuclides diff --git a/openmc/geometry.py b/openmc/geometry.py index 14fd48fb1e..8314a9049c 100644 --- a/openmc/geometry.py +++ b/openmc/geometry.py @@ -50,6 +50,19 @@ class Geometry(object): self._root_universe = root_universe + def add_volume_information(self, volume_calc): + """Add volume information to from a stochastic volume calculation. + + Parameters + ---------- + volume_calc : openmc.VolumeCalculation + Results from a stochastic volume calculation + + """ + for cell in self.get_all_cells(): + if cell.id in volume_calc.results: + cell.add_volume_information(volume_calc) + def export_to_xml(self): """Create a geometry.xml file that can be used for a simulation. @@ -166,24 +179,6 @@ class Geometry(object): universes.sort(key=lambda x: x.id) return universes - def get_all_nuclides(self): - """Return all nuclides assigned to a material in the geometry - - Returns - ------- - list of openmc.Nuclide - Nuclides in the geometry - - """ - - nuclides = OrderedDict() - materials = self.get_all_materials() - - for material in materials: - nuclides.update(material.get_all_nuclides()) - - return nuclides - def get_all_materials(self): """Return all materials assigned to a cell diff --git a/openmc/lattice.py b/openmc/lattice.py index 81144e4d81..b26aab4c67 100644 --- a/openmc/lattice.py +++ b/openmc/lattice.py @@ -144,25 +144,26 @@ class Lattice(object): return univs - def get_all_nuclides(self): - """Return all nuclides contained in the lattice + def get_nuclides(self): + """Returns all nuclides in the lattice Returns ------- - nuclides : collections.OrderedDict - Dictionary whose keys are nuclide names and values are 2-tuples of - (nuclide, density) + nuclides : list + List of nuclide names """ - nuclides = OrderedDict() + nuclides = [] # Get all unique Universes contained in each of the lattice cells unique_universes = self.get_unique_universes() # Append all Universes containing each cell to the dictionary - for universe_id, universe in unique_universes.items(): - nuclides.update(universe.get_all_nuclides()) + for universe in unique_universes.values(): + for nuclide in universe.get_nuclides(): + if nuclide not in nuclides: + nuclides.append(nuclide) return nuclides diff --git a/openmc/material.py b/openmc/material.py index d7ffd04e94..5603d9d25c 100644 --- a/openmc/material.py +++ b/openmc/material.py @@ -492,9 +492,31 @@ class Material(object): for element, percent, percent_type in self._elements: element.scattering = 'iso-in-lab' - def get_all_nuclides(self): + def get_nuclides(self): """Returns all nuclides in the material + Returns + ------- + nuclides : list + List of nuclide names + + """ + + nuclides = [] + + for nuclide, density, density_type in self._nuclides: + nuclides.append(nuclide.name) + + for element, density, density_type in self._elements: + # Expand natural element into isotopes + for isotope, abundance in element.expand(): + nuclides.append(isotope.name) + + return nuclides + + def get_nuclide_densities(self): + """Returns all nuclides in the material and their densities + Returns ------- nuclides : dict diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index 84571ce236..ee11d0ef64 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -1008,7 +1008,7 @@ class Library(object): # Create the xsdata object and add it to the mgxs_file for i, domain in enumerate(self.domains): if self.by_nuclide: - nuclides = list(domain.get_all_nuclides().keys()) + nuclides = domain.get_nuclides() else: nuclides = ['total'] for nuclide in nuclides: diff --git a/openmc/mgxs/mgxs.py b/openmc/mgxs/mgxs.py index e1c2e62217..7aa12d0418 100644 --- a/openmc/mgxs/mgxs.py +++ b/openmc/mgxs/mgxs.py @@ -288,9 +288,7 @@ class MGXS(object): # If this is a by-nuclide cross-section, add nuclides to Tally if self.by_nuclide and score != 'flux': - all_nuclides = self.get_all_nuclides() - for nuclide in all_nuclides: - self._tallies[key].nuclides.append(nuclide) + self._tallies[key].nuclides += self.get_nuclides() else: self._tallies[key].nuclides.append('total') @@ -329,14 +327,14 @@ class MGXS(object): @property def num_nuclides(self): if self.by_nuclide: - return len(self.get_all_nuclides()) + return len(self.get_nuclides()) else: return 1 @property def nuclides(self): if self.by_nuclide: - return self.get_all_nuclides() + return self.get_nuclides() else: return 'sum' @@ -503,7 +501,7 @@ class MGXS(object): mgxs.name = name return mgxs - def get_all_nuclides(self): + def get_nuclides(self): """Get all nuclides in the cross section's spatial domain. Returns @@ -528,8 +526,7 @@ class MGXS(object): # Otherwise, return all nuclides in the spatial domain else: - nuclides = self.domain.get_all_nuclides() - return list(nuclides.keys()) + return self.domain.get_nuclides() def get_nuclide_density(self, nuclide): """Get the atomic number density in units of atoms/b-cm for a nuclide @@ -556,7 +553,7 @@ class MGXS(object): cv.check_type('nuclide', nuclide, basestring) # Get list of all nuclides in the spatial domain - nuclides = self.domain.get_all_nuclides() + nuclides = self.domain.get_nuclide_densities() if nuclide not in nuclides: msg = 'Unable to get density for nuclide "{0}" which is not in ' \ @@ -597,14 +594,14 @@ class MGXS(object): # Sum the atomic number densities for all nuclides if nuclides == 'sum': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() densities = np.zeros(1, dtype=np.float) for nuclide in nuclides: densities[0] += self.get_nuclide_density(nuclide) # Tabulate the atomic number densities for all nuclides elif nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() densities = np.zeros(self.num_nuclides, dtype=np.float) for i, nuclide in enumerate(nuclides): densities[i] += self.get_nuclide_density(nuclide) @@ -635,7 +632,7 @@ class MGXS(object): # If computing xs for each nuclide, replace CrossNuclides with originals if self.by_nuclide: self.xs_tally._nuclides = [] - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() for nuclide in nuclides: self.xs_tally.nuclides.append(openmc.Nuclide(nuclide)) @@ -796,7 +793,7 @@ class MGXS(object): # Construct a collection of the nuclides to retrieve from the xs tally if self.by_nuclide: if nuclides == 'all' or nuclides == 'sum' or nuclides == ['sum']: - query_nuclides = self.get_all_nuclides() + query_nuclides = self.get_nuclides() else: query_nuclides = nuclides else: @@ -1164,7 +1161,7 @@ class MGXS(object): # Construct a collection of the nuclides to report if self.by_nuclide: if nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() elif nuclides == 'sum': nuclides = ['sum'] else: @@ -1302,7 +1299,7 @@ class MGXS(object): # Construct a collection of the nuclides to report if self.by_nuclide: if nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() densities = np.zeros(len(nuclides), dtype=np.float) elif nuclides == 'sum': nuclides = ['sum'] @@ -1484,7 +1481,7 @@ class MGXS(object): if self.by_nuclide and nuclides == 'sum': # Use tally summation to sum across all nuclides - query_nuclides = self.get_all_nuclides() + query_nuclides = self.get_nuclides() xs_tally = self.xs_tally.summation(nuclides=query_nuclides) df = xs_tally.get_pandas_dataframe( distribcell_paths=distribcell_paths) @@ -1790,7 +1787,7 @@ class MatrixMGXS(MGXS): # Construct a collection of the nuclides to retrieve from the xs tally if self.by_nuclide: if nuclides == 'all' or nuclides == 'sum' or nuclides == ['sum']: - query_nuclides = self.get_all_nuclides() + query_nuclides = self.get_nuclides() else: query_nuclides = nuclides else: @@ -1937,7 +1934,7 @@ class MatrixMGXS(MGXS): # Construct a collection of the nuclides to report if self.by_nuclide: if nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() if nuclides == 'sum': nuclides = ['sum'] else: @@ -3622,7 +3619,7 @@ class ScatterMatrixXS(MatrixMGXS): # Construct a collection of the nuclides to retrieve from the xs tally if self.by_nuclide: if nuclides == 'all' or nuclides == 'sum' or nuclides == ['sum']: - query_nuclides = self.get_all_nuclides() + query_nuclides = self.get_nuclides() else: query_nuclides = nuclides else: @@ -3789,7 +3786,7 @@ class ScatterMatrixXS(MatrixMGXS): # Construct a collection of the nuclides to report if self.by_nuclide: if nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() if nuclides == 'sum': nuclides = ['sum'] else: @@ -4595,7 +4592,7 @@ class Chi(MGXS): nu_fission_out = self.tallies['nu-fission-out'] # Sum out all nuclides - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() nu_fission_in = nu_fission_in.summation(nuclides=nuclides) nu_fission_out = nu_fission_out.summation(nuclides=nuclides) @@ -4615,7 +4612,7 @@ class Chi(MGXS): # Get chi for all nuclides in the domain elif nuclides == 'all': - nuclides = self.get_all_nuclides() + nuclides = self.get_nuclides() xs = self.xs_tally.get_values(filters=filters, filter_bins=filter_bins, nuclides=nuclides, value=value) diff --git a/openmc/statepoint.py b/openmc/statepoint.py index 6337746650..58a79fa219 100644 --- a/openmc/statepoint.py +++ b/openmc/statepoint.py @@ -2,6 +2,7 @@ import sys import re import os import warnings +import glob import numpy as np @@ -22,8 +23,9 @@ class StatePoint(object): filename : str Path to file to load autolink : bool, optional - Whether to automatically link in metadata from a summary.h5 - file. Defaults to True. + Whether to automatically link in metadata from a summary.h5 file and + stochastic volume calculation results from volume_*.h5 files. Defaults + to True. Attributes ---------- @@ -143,6 +145,12 @@ class StatePoint(object): su = openmc.Summary(path_summary) self.link_with_summary(su) + path_volume = os.path.join(os.path.dirname(filename), 'volume_*.h5') + for path_i in glob.glob(path_volume): + if re.search(r'volume_\d+\.h5', path_i): + vol = openmc.VolumeCalculation.from_hdf5(path_i) + self.add_volume_information(vol) + def close(self): self._f.close() @@ -501,6 +509,19 @@ class StatePoint(object): for tally_id in self.tallies: self.tallies[tally_id].sparse = self.sparse + def add_volume_information(self, volume_calc): + """Add volume information to the geometry within the file + + Parameters + ---------- + volume_calc : openmc.VolumeCalculation + Results from a stochastic volume calculation + + """ + if self.summary is not None: + self.summary.add_volume_information(volume_calc) + + def get_tally(self, scores=[], filters=[], nuclides=[], name=None, id=None, estimator=None, exact_filters=False, exact_nuclides=False, exact_scores=False): diff --git a/openmc/summary.py b/openmc/summary.py index 8c626940e5..6fbace72ea 100644 --- a/openmc/summary.py +++ b/openmc/summary.py @@ -594,6 +594,17 @@ class Summary(object): # Add Tally to the global dictionary of all Tallies self.tallies[tally_id] = tally + def add_volume_information(self, volume_calc): + """Add volume information to the geometry within the summary file + + Parameters + ---------- + volume_calc : openmc.VolumeCalculation + Results from a stochastic volume calculation + + """ + self.openmc_geometry.add_volume_information(volume_calc) + def get_material_by_id(self, material_id): """Return a Material object given the material id diff --git a/openmc/universe.py b/openmc/universe.py index c8e7fcab1e..f2b8641699 100644 --- a/openmc/universe.py +++ b/openmc/universe.py @@ -349,7 +349,27 @@ class Universe(object): # Return the offset computed at all nested Universe levels return offset - def get_all_nuclides(self): + def get_nuclides(self): + """Returns all nuclides in the universe + + Returns + ------- + nuclides : list + List of nuclide names + + """ + + nuclides = [] + + # Append all Nuclides in each Cell in the Universe to the dictionary + for cell in self.cells.values(): + for nuclide in cell.get_nuclides(): + if nuclide not in nuclides: + nuclides.append(nuclide) + + return nuclides + + def get_nuclide_densities(self): """Return all nuclides contained in the universe Returns @@ -360,13 +380,8 @@ class Universe(object): """ - nuclides = OrderedDict() - - # Append all Nuclides in each Cell in the Universe to the dictionary - for cell in self._cells.values(): - nuclides.update(cell.get_all_nuclides()) - - return nuclides + raise NotImplementedError('Determining average nuclide densities over ' + 'an entire universe not yet supported.') def get_all_cells(self): """Return all cells that are contained within the universe