Ability to link volume calculation results in order to compute microscopic cross

sections from mgxs
This commit is contained in:
Paul Romano 2016-08-01 07:02:16 -05:00
parent 66d772da79
commit 08c9036bb6
9 changed files with 178 additions and 64 deletions

View file

@ -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

View file

@ -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

View file

@ -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

View file

@ -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

View file

@ -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:

View file

@ -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)

View file

@ -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):

View file

@ -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

View file

@ -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