diff --git a/openmc/data/library.py b/openmc/data/library.py index 9f6d93146..3beb25543 100644 --- a/openmc/data/library.py +++ b/openmc/data/library.py @@ -5,6 +5,7 @@ import h5py from openmc.mixin import EqualityMixin from openmc.clean_xml import clean_xml_indentation +from openmc.checkvalue import check_type class DataLibrary(EqualityMixin): @@ -95,13 +96,14 @@ class DataLibrary(EqualityMixin): method='xml') @classmethod - def from_xml(cls, path): + def from_xml(cls, path=None): """Read cross section data library from an XML file. Parameters ---------- - path : str - Path to XML file to read. + path : str, optional + Path to XML file to read. If not provided, the + `OPENMC_CROSS_SECTIONS` environment variable will be used. Returns ------- @@ -112,6 +114,18 @@ class DataLibrary(EqualityMixin): data = cls() + # If cross_sections is None, get the cross sections from the + # OPENMC_CROSS_SECTIONS environment variable + if path is None: + path = os.environ.get('OPENMC_CROSS_SECTIONS') + + # Check to make sure there was an environmental variable. + if path is None: + raise ValueError("Either path or OPENMC_CROSS_SECTIONS " + "environmental variable must be set") + + check_type('path', path, str) + tree = ET.parse(path) root = tree.getroot() if root.find('directory') is not None: diff --git a/openmc/element.py b/openmc/element.py index 6fe988195..905cdc6c9 100644 --- a/openmc/element.py +++ b/openmc/element.py @@ -262,7 +262,7 @@ class Element(object): return isotopes def plot_xs(self, types, divisor_types=None, temperature=294., axis=None, - Erange=(1.E-5, 20.E6), sab_name=None, cross_sections=None, + energy_range=(1.E-5, 20.E6), sab_name=None, cross_sections=None, enrichment=None, **kwargs): """Creates a figure of continuous-energy microscopic cross sections for this element @@ -283,7 +283,7 @@ class Element(object): axis : matplotlib.axes, optional A previously generated axis to use for plotting. If not specified, a new axis and figure will be generated. - Erange : tuple of floats + energy_range : tuple of floats Energy range (in eV) to plot the cross section within sab_name : str, optional Name of S(a,b) library to apply to MT=2 data when applicable. @@ -324,12 +324,12 @@ class Element(object): E = np.union1d(Enum, Ediv) data_new = np.zeros((len(types), len(E))) - for l in range(len(types)): - data_new[l, :] = \ - np.divide(np.interp(E, Enum, data[l, :]), - np.interp(E, Ediv, data_div[l, :])) - if divisor_types[l] != 'unity': - types[l] = types[l] + ' / ' + divisor_types[l] + for line in range(len(types)): + data_new[line, :] = \ + np.divide(np.interp(E, Enum, data[line, :]), + np.interp(E, Ediv, data_div[line, :])) + if divisor_types[line] != 'unity': + types[line] = types[line] + ' / ' + divisor_types[line] data = data_new # Generate the plot @@ -339,13 +339,14 @@ class Element(object): else: fig = None ax = axis + # Set to loglog or semilogx depending on if we are plotting a data + # type which we expect to vary linearly + if set(types).issubset(PLOT_TYPES_LINEAR): + plot_func = ax.semilogx + else: + plot_func = ax.loglog + # Plot the data for i in range(len(data)): - # Set to loglog or semilogx depending on if we are plotting a data - # type which we expect to vary linearly - if types[i] in PLOT_TYPES_LINEAR: - plot_func = ax.semilogx - else: - plot_func = ax.loglog if np.sum(data[i, :]) > 0.: print('max = {:.2E}'.format(np.max(data[i, :]))) plot_func(E, data[i, :], label=types[i]) @@ -356,7 +357,7 @@ class Element(object): else: ax.set_ylabel('Elemental Cross Section [b]') ax.legend(loc='best') - ax.set_xlim(Erange) + ax.set_xlim(energy_range) if self.name is not None: title = 'Cross Section for ' + self.name ax.set_title(title) @@ -388,11 +389,11 @@ class Element(object): Returns ------- - unionE : numpy.array + energy_grid : numpy.array Energies at which cross sections are calculated, in units of eV data : numpy.ndarray Macroscopic cross sections calculated at the energy grid described - by unionE + by energy_grid """ @@ -402,18 +403,8 @@ class Element(object): cv.check_iterable_type('types', types, str) cv.check_type('temperature', temperature, Real) - # If cross_sections is None, get the cross sections from the - # OPENMC_CROSS_SECTIONS environment variable - if cross_sections is None: - cross_sections = os.environ.get('OPENMC_CROSS_SECTIONS') - - # If a cross_sections library is present, check natural nuclides - # against the nuclides in the library - if cross_sections is not None: - library = openmc.data.DataLibrary.from_xml(cross_sections) - else: - raise ValueError("cross_sections or OPENMC_CROSS_SECTIONS " - "environmental variable must be set") + # Load the library + library = openmc.data.DataLibrary.from_xml(cross_sections) # Expand elements in to nuclides with atomic densities nuclides = self.expand(100., 'ao', enrichment=enrichment, @@ -443,17 +434,17 @@ class Element(object): # Condense the data for every nuclide # First create a union energy grid - unionE = E[0] + energy_grid = E[0] for n in range(1, len(E)): - unionE = np.union1d(unionE, E[n]) + energy_grid = np.union1d(energy_grid, E[n]) # Now we can combine all the nuclidic data - data = np.zeros((len(types), len(unionE))) - for l in range(len(types)): - if types[l] == 'unity': - data[l, :] = 1. + data = np.zeros((len(types), len(energy_grid))) + for line in range(len(types)): + if types[line] == 'unity': + data[line, :] = 1. else: for n in range(len(nuclides)): - data[l, :] += nuc_fractions[n] * xs[n][l](unionE) + data[line, :] += nuc_fractions[n] * xs[n][line](energy_grid) - return unionE, data + return energy_grid, data diff --git a/openmc/material.py b/openmc/material.py index 1889a4765..cd80b2089 100644 --- a/openmc/material.py +++ b/openmc/material.py @@ -627,20 +627,8 @@ class Material(object): import scipy.constants as sc - cv.check_type('cross_sections', cross_sections, str) - - # If cross_sections is None, get the cross sections from the - # OPENMC_CROSS_SECTIONS environment variable - if cross_sections is None: - cross_sections = os.environ.get('OPENMC_CROSS_SECTIONS') - - # If a cross_sections library is present, check natural nuclides - # against the nuclides in the library - if cross_sections is not None: - library = openmc.data.DataLibrary.from_xml(cross_sections) - else: - raise ValueError("cross_sections or OPENMC_CROSS_SECTIONS " - "environmental variable must be set") + # Load the library + library = openmc.data.DataLibrary.from_xml(cross_sections) # Expand elements in to nuclides nuclides = self.get_nuclide_densities() @@ -704,15 +692,13 @@ class Material(object): nuc_densities = density * nuc_densities nuclides = OrderedDict() - n = -1 - for nuc in nucs: - n += 1 + for n, nuc in enumerate(nucs): nuclides[nuc] = (nuc, nuc_densities[n]) return nuclides def plot_xs(self, types, divisor_types=None, temperature=294., - axis=None, Erange=(1.E-5, 20.E6), cross_sections=None, + axis=None, energy_range=(1.E-5, 20.E6), cross_sections=None, **kwargs): """Creates a figure of continuous-energy macroscopic cross sections for this material @@ -733,7 +719,7 @@ class Material(object): axis : matplotlib.axes, optional A previously generated axis to use for plotting. If not specified, a new axis and figure will be generated. - Erange : tuple of floats, optional + energy_range : tuple of floats, optional Energy range (in eV) to plot the cross section within cross_sections : str, optional Location of cross_sections.xml file. Default is None. @@ -767,12 +753,12 @@ class Material(object): E = np.union1d(Enum, Ediv) data_new = np.zeros((len(types), len(E))) - for l in range(len(types)): - data_new[l, :] = \ - np.divide(np.interp(E, Enum, data[l, :]), - np.interp(E, Ediv, data_div[l, :])) - if divisor_types[l] != 'unity': - types[l] = types[l] + ' / ' + divisor_types[l] + for line in range(len(types)): + data_new[line, :] = \ + np.divide(np.interp(E, Enum, data[line, :]), + np.interp(E, Ediv, data_div[line, :])) + if divisor_types[line] != 'unity': + types[line] = types[line] + ' / ' + divisor_types[line] data = data_new # Generate the plot @@ -782,13 +768,14 @@ class Material(object): else: fig = None ax = axis + # Set to loglog or semilogx depending on if we are plotting a data + # type which we expect to vary linearly + if set(types).issubset(PLOT_TYPES_LINEAR): + plot_func = ax.semilogx + else: + plot_func = ax.loglog + # Plot the data for i in range(len(data)): - # Set to loglog or semilogx depending on if we are plotting a data - # type which we expect to vary linearly - if types[i] in PLOT_TYPES_LINEAR: - plot_func = ax.semilogx - else: - plot_func = ax.loglog if np.sum(data[i, :]) > 0.: plot_func(E, data[i, :], label=types[i]) @@ -798,7 +785,7 @@ class Material(object): else: ax.set_ylabel('Macroscopic Cross Section [1/cm]') ax.legend(loc='best') - ax.set_xlim(Erange) + ax.set_xlim(energy_range) if self.name is not None: title = 'Macroscopic Cross Section for ' + self.name ax.set_title(title) @@ -823,11 +810,11 @@ class Material(object): Returns ------- - unionE : numpy.array + energy_grid : numpy.array Energies at which cross sections are calculated, in units of eV data : numpy.ndarray Macroscopic cross sections calculated at the energy grid described - by unionE + by energy_grid """ @@ -840,18 +827,8 @@ class Material(object): cv.check_type('temperature', temperature, Real) T = temperature - # If cross_sections is None, get the cross sections from the - # OPENMC_CROSS_SECTIONS environment variable - if cross_sections is None: - cross_sections = os.environ.get('OPENMC_CROSS_SECTIONS') - - # If a cross_sections library is present, check natural nuclides - # against the nuclides in the library - if cross_sections is not None: - library = openmc.data.DataLibrary.from_xml(cross_sections) - else: - raise ValueError("cross_sections or OPENMC_CROSS_SECTIONS " - "environmental variable must be set") + # Load the library + library = openmc.data.DataLibrary.from_xml(cross_sections) # Expand elements in to nuclides with atomic densities nuclides = self.get_nuclide_atom_densities(cross_sections) @@ -881,20 +858,20 @@ class Material(object): # Condense the data for every nuclide # First create a union energy grid - unionE = E[0] + energy_grid = E[0] for n in range(1, len(E)): - unionE = np.union1d(unionE, E[n]) + energy_grid = np.union1d(energy_grid, E[n]) # Now we can combine all the nuclidic data - data = np.zeros((len(types), len(unionE))) - for l in range(len(types)): - if types[l] == 'unity': - data[l, :] = 1. + data = np.zeros((len(types), len(energy_grid))) + for line in range(len(types)): + if types[line] == 'unity': + data[line, :] = 1. else: for n in range(len(nuclides)): - data[l, :] += nuc_densities[n] * xs[n][l](unionE) + data[line, :] += nuc_densities[n] * xs[n][line](energy_grid) - return unionE, data + return energy_grid, data def _get_nuclide_xml(self, nuclide, distrib=False): xml_element = ET.Element("nuclide") diff --git a/openmc/nuclide.py b/openmc/nuclide.py index 439aedf48..c42d4f8d5 100644 --- a/openmc/nuclide.py +++ b/openmc/nuclide.py @@ -98,10 +98,9 @@ class Nuclide(object): self._scattering = scattering def plot_xs(self, types, divisor_types=None, temperature=294., axis=None, - Erange=(1.E-5, 20.E6), sab_name=None, cross_sections=None, + energy_range=(1.E-5, 20.E6), sab_name=None, cross_sections=None, **kwargs): - """Creates a figure of continuous-energy cross sections for this - nuclide + """Creates a figure of continuous-energy cross sections for this nuclide Parameters ---------- @@ -119,7 +118,7 @@ class Nuclide(object): axis : matplotlib.axes, optional A previously generated axis to use for plotting. If not specified, a new axis and figure will be generated. - Erange : tuple of floats + energy_range : tuple of floats Energy range (in eV) to plot the cross section within sab_name : str, optional Name of S(a,b) library to apply to MT=2 data when applicable. @@ -156,11 +155,12 @@ class Nuclide(object): E = np.union1d(Enum, Ediv) data_new = [] - for l in range(len(types)): - data_new.append(openmc.data.Combination([data[l], data_div[l]], + for line in range(len(types)): + data_new.append(openmc.data.Combination([data[line], + data_div[line]], [np.divide])) - if divisor_types[l] != 'unity': - types[l] = types[l] + ' / ' + divisor_types[l] + if divisor_types[line] != 'unity': + types[line] = types[line] + ' / ' + divisor_types[line] data = data_new # Generate the plot @@ -170,14 +170,15 @@ class Nuclide(object): else: fig = None ax = axis + # Set to loglog or semilogx depending on if we are plotting a data + # type which we expect to vary linearly + if set(types).issubset(PLOT_TYPES_LINEAR): + plot_func = ax.semilogx + else: + plot_func = ax.loglog + # Plot the data for i in range(len(data)): to_plot = data[i](E) - # Set to loglog or semilogx depending on if we are plotting a data - # type which we expect to vary linearly - if types[i] in PLOT_TYPES_LINEAR: - plot_func = ax.semilogx - else: - plot_func = ax.loglog if np.sum(to_plot) > 0.: plot_func(E, to_plot, label=types[i]) @@ -187,7 +188,7 @@ class Nuclide(object): else: ax.set_ylabel('Microscopic Cross Section [b]') ax.legend(loc='best') - ax.set_xlim(Erange) + ax.set_xlim(energy_range) if self.name is not None: title = 'Microscopic Cross Section for ' + self.name ax.set_title(title) @@ -216,10 +217,11 @@ class Nuclide(object): Returns ------- - E : numpy.array + energy_grid : numpy.array Energies at which cross sections are calculated, in units of eV data : numpy.ndarray - Cross sections calculated at the energy grid described by unionE + Cross sections calculated at the energy grid described by + energy_grid """ @@ -246,18 +248,8 @@ class Nuclide(object): yields.append((False,)) ops.append(()) - # If cross_sections is None, get the cross sections from the - # OPENMC_CROSS_SECTIONS environment variable - if cross_sections is None: - cross_sections = os.environ.get('OPENMC_CROSS_SECTIONS') - - # If a cross_sections library is present, check natural nuclides - # against the nuclides in the library - if cross_sections is not None: - library = openmc.data.DataLibrary.from_xml(cross_sections) - else: - raise ValueError("cross_sections or OPENMC_CROSS_SECTIONS " - "environmental variable must be set") + # Load the library + library = openmc.data.DataLibrary.from_xml(cross_sections) # Convert temperature to format needed for access in the library cv.check_type('temperature', temperature, Real) @@ -265,7 +257,7 @@ class Nuclide(object): T = temperature # Now we can create the data sets to be plotted - E = [] + energy_grid = [] xs = [] lib = library.get_by_material(self.name) if lib is not None: @@ -325,9 +317,9 @@ class Nuclide(object): if inelastic.x[-1] > sab_Emax: sab_Emax = inelastic.x[-1] sab_funcs.append(inelastic) - E = grid + energy_grid = grid else: - E = nuc.energy[nucT] + energy_grid = nuc.energy[nucT] for i, mt_set in enumerate(mts): # Get the reaction xs data from the nuclide @@ -385,4 +377,4 @@ class Nuclide(object): else: raise ValueError(nuclide[0] + " not in library") - return E, xs + return energy_grid, xs diff --git a/openmc/plot_data.py b/openmc/plot_data.py index ed87f2f00..1efeda00c 100644 --- a/openmc/plot_data.py +++ b/openmc/plot_data.py @@ -3,7 +3,7 @@ import numpy as np # Supported keywords for material xs plotting PLOT_TYPES = ['total', 'scatter', 'elastic', 'inelastic', 'fission', 'absorption', 'capture', 'nu-fission', 'nu-scatter', 'unity', - 'slowing-down power'] + 'slowing-down power', 'damage'] # Special MT values UNITY_MT = -1 @@ -44,7 +44,8 @@ PLOT_TYPES_MT = {'total': (2, 3,), 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 183, 184, 190, 194, 196, 198, 199, 200, 875, - 891, XI_MT)} + 891, XI_MT), + 'damage': (444,)} # Operations to use when combining MTs the first np.add is used in reference # to zero PLOT_TYPES_OP = {'total': (np.add,), @@ -56,7 +57,8 @@ PLOT_TYPES_OP = {'total': (np.add,), 'nu-scatter': (np.add,) * (len(PLOT_TYPES_MT['nu-scatter']) - 1), 'unity': (), 'slowing-down power': - (np.add,) * (len(PLOT_TYPES_MT['slowing-down power']) - 2) + (np.multiply,)} + (np.add,) * (len(PLOT_TYPES_MT['slowing-down power']) - 2) + (np.multiply,), + 'damage': ()} # Whether or not to multiply the reaction by the yield as well PLOT_TYPES_YIELD = {'total': (False, False), @@ -68,7 +70,9 @@ PLOT_TYPES_YIELD = {'total': (False, False), 'nu-scatter': (True,) * len(PLOT_TYPES_MT['nu-scatter']), 'unity': (False,), 'slowing-down power': - (True,) * len(PLOT_TYPES_MT['slowing-down power'])} + (True,) * len(PLOT_TYPES_MT['slowing-down power']), + 'damage': (False,)} # Types of plots to plot linearly in y -PLOT_TYPES_LINEAR = ['nu-fission / fission', 'nu-scatter / scatter'] +PLOT_TYPES_LINEAR = {'nu-fission / fission', 'nu-scatter / scatter', + 'nu-fission / absorption', 'fission / absorption'}