diff --git a/openmc/data/endf.c b/openmc/data/endf.c index fee6af179e..6e9c1f5b4f 100644 --- a/openmc/data/endf.c +++ b/openmc/data/endf.c @@ -1,14 +1,31 @@ #include +//! Convert string representation of a floating point number into a double +// +//! This function handles converting floating point numbers from an ENDF 11 +//! character field into a double, covering all the corner cases. Floating point +//! numbers are allowed to contain whitespace (which is ignored). Also, in +//! exponential notation, it allows the 'e' to be omitted. A field containing +//! only whitespace is to be interpreted as a zero. +// +//! \param buffer character input from an ENDF file +//! \param n Length of character input +//! \return Floating point number + double cfloat_endf(const char* buffer, int n) { char arr[12]; // 11 characters plus a null terminator int j = 0; // current position in arr int found_significand = 0; int found_exponent = 0; + + // limit n to 11 characters + n = n > 11 ? 11 : n; + for (int i = 0; i < n; ++i) { - // Skip whitespace characters char c = buffer[i]; + + // Skip whitespace characters if (c == ' ') continue; if (found_significand) { diff --git a/openmc/data/grid.py b/openmc/data/grid.py index 6aec569dc8..33a5b7b00e 100644 --- a/openmc/data/grid.py +++ b/openmc/data/grid.py @@ -37,12 +37,18 @@ def linearize(x, f, tolerance=0.001): y_stack.insert(0, f(x[i + 1])) while True: + # Get the bounding points currently on the stack x_high, x_low = x_stack[-2:] y_high, y_low = y_stack[-2:] + + # Evaluate the function at the midpoint x_mid = 0.5*(x_low + x_high) y_mid = f(x_mid) + # Linearly interpolate between the bounding points y_interp = y_low + (y_high - y_low)/(x_high - x_low)*(x_mid - x_low) + + # Check the error on the interpolated point and compare to tolerance error = abs((y_interp - y_mid)/y_mid) if error > tolerance: x_stack.insert(-1, x_mid) diff --git a/openmc/data/kalbach_mann.py b/openmc/data/kalbach_mann.py index c7df015a75..f4c914d909 100644 --- a/openmc/data/kalbach_mann.py +++ b/openmc/data/kalbach_mann.py @@ -245,6 +245,7 @@ class KalbachMann(AngleEnergy): eout_i = Mixture([p_discrete, 1. - p_discrete], [eout_discrete, eout_continuous]) + # Precompound factor and slope are on rows 3 and 4, respectively km_r = Tabulated1D(data[0, j:j+n], data[3, j:j+n]) km_a = Tabulated1D(data[0, j:j+n], data[4, j:j+n]) diff --git a/openmc/data/multipole.py b/openmc/data/multipole.py index 4f096b6653..c7390d14db 100644 --- a/openmc/data/multipole.py +++ b/openmc/data/multipole.py @@ -5,7 +5,8 @@ import h5py import numpy as np import openmc.checkvalue as cv -from openmc.mixin import EqualityMixin +from ..exceptions import DataError +from ..mixin import EqualityMixin from . import WMP_VERSION, WMP_VERSION_MAJOR from .data import K_BOLTZMANN @@ -332,19 +333,21 @@ class WindowedMultipole(EqualityMixin): if isinstance(group_or_filename, h5py.Group): group = group_or_filename + need_to_close = False else: h5file = h5py.File(str(group_or_filename), 'r') + need_to_close = True # Make sure version matches if 'version' in h5file.attrs: major, minor = h5file.attrs['version'] if major != WMP_VERSION_MAJOR: - raise IOError( + raise DataError( 'WMP data format uses version {}. {} whereas your ' 'installation of the OpenMC Python API expects version ' '{}.x.'.format(major, minor, WMP_VERSION_MAJOR)) else: - raise IOError( + raise DataError( 'WMP data does not indicate a version. Your installation of ' 'the OpenMC Python API expects version {}.x data.' .format(WMP_VERSION_MAJOR)) @@ -382,6 +385,10 @@ class WindowedMultipole(EqualityMixin): raise ValueError("Windowed multipole is only supported for " "curvefits with 3 or more terms.") + # If HDF5 file was opened here, make sure it gets closed + if need_to_close: + h5file.close() + return out def _evaluate(self, E, T): @@ -505,7 +512,7 @@ class WindowedMultipole(EqualityMixin): ---------- path : str Path to write HDF5 file to - mode : {'r', r+', 'w', 'x', 'a'} + mode : {'r+', 'w', 'x', 'a'} Mode that is used to open the HDF5 file. This is the second argument to the :class:`h5py.File` constructor. libver : {'earliest', 'latest'} diff --git a/openmc/data/neutron.py b/openmc/data/neutron.py index 06ba86f4d8..47ebfabf83 100644 --- a/openmc/data/neutron.py +++ b/openmc/data/neutron.py @@ -20,7 +20,7 @@ from .function import Tabulated1D, Sum, ResonancesWithBackground from .grid import linearize, thin from .njoy import make_ace from .product import Product -from .reaction import Reaction, _get_photon_products_ace +from .reaction import Reaction, _get_photon_products_ace, FISSION_MTS from . import resonance as res from . import resonance_covariance as res_cov from .urr import ProbabilityTables @@ -345,7 +345,8 @@ class IncidentNeutron(EqualityMixin): # Add grid around each resonance that includes the peak +/- the # width times each value in _RESONANCE_ENERGY_GRID. Values are # constrained so that points around one resonance don't overlap - # with points around another. This algorithm is from Fudge. + # with points around another. This algorithm is from Fudge + # (https://doi.org/10.1063/1.1945057). energies = [] for e, g, e_lower, e_upper in zip(e_peak, gamma, e_mid[:-1], e_mid[1:]): @@ -410,7 +411,7 @@ class IncidentNeutron(EqualityMixin): ---------- path : str Path to write HDF5 file to - mode : {'r', 'r+', 'w', 'x', 'a'} + mode : {'r+', 'w', 'x', 'a'} Mode that is used to open the HDF5 file. This is the second argument to the :class:`h5py.File` constructor. libver : {'earliest', 'latest'} @@ -424,64 +425,62 @@ class IncidentNeutron(EqualityMixin): 'originated from an ENDF file.') # Open file and write version - f = h5py.File(str(path), mode, libver=libver) - f.attrs['filetype'] = np.string_('data_neutron') - f.attrs['version'] = np.array(HDF5_VERSION) + with h5py.File(str(path), mode, libver=libver) as f: + f.attrs['filetype'] = np.string_('data_neutron') + f.attrs['version'] = np.array(HDF5_VERSION) - # Write basic data - g = f.create_group(self.name) - g.attrs['Z'] = self.atomic_number - g.attrs['A'] = self.mass_number - g.attrs['metastable'] = self.metastable - g.attrs['atomic_weight_ratio'] = self.atomic_weight_ratio - ktg = g.create_group('kTs') - for i, temperature in enumerate(self.temperatures): - ktg.create_dataset(temperature, data=self.kTs[i]) + # Write basic data + g = f.create_group(self.name) + g.attrs['Z'] = self.atomic_number + g.attrs['A'] = self.mass_number + g.attrs['metastable'] = self.metastable + g.attrs['atomic_weight_ratio'] = self.atomic_weight_ratio + ktg = g.create_group('kTs') + for i, temperature in enumerate(self.temperatures): + ktg.create_dataset(temperature, data=self.kTs[i]) - # Write energy grid - eg = g.create_group('energy') - for temperature in self.temperatures: - eg.create_dataset(temperature, data=self.energy[temperature]) + # Write energy grid + eg = g.create_group('energy') + for temperature in self.temperatures: + eg.create_dataset(temperature, data=self.energy[temperature]) - # Write 0K energy grid if needed - if '0K' in self.energy and '0K' not in eg: - eg.create_dataset('0K', data=self.energy['0K']) + # Write 0K energy grid if needed + if '0K' in self.energy and '0K' not in eg: + eg.create_dataset('0K', data=self.energy['0K']) - # Write reaction data - rxs_group = g.create_group('reactions') - for rx in self.reactions.values(): - # Skip writing redundant reaction if it doesn't have photon - # production or is a summed transmutation reaction. MT=4 is also - # sometimes needed for probability tables. Also write gas - # production, heating, and damage energy production. - if rx.redundant: - photon_rx = any(p.particle == 'photon' for p in rx.products) - keep_mts = (4, 16, 103, 104, 105, 106, 107, - 203, 204, 205, 206, 207, 301, 444, 901) - if not (photon_rx or rx.mt in keep_mts): - continue + # Write reaction data + rxs_group = g.create_group('reactions') + for rx in self.reactions.values(): + # Skip writing redundant reaction if it doesn't have photon + # production or is a summed transmutation reaction. MT=4 is also + # sometimes needed for probability tables. Also write gas + # production, heating, and damage energy production. + if rx.redundant: + photon_rx = any(p.particle == 'photon' for p in rx.products) + keep_mts = (4, 16, 103, 104, 105, 106, 107, + 203, 204, 205, 206, 207, 301, 444, 901) + if not (photon_rx or rx.mt in keep_mts): + continue - rx_group = rxs_group.create_group('reaction_{:03}'.format(rx.mt)) - rx.to_hdf5(rx_group) + rx_group = rxs_group.create_group('reaction_{:03}'.format(rx.mt)) + rx.to_hdf5(rx_group) - # Write total nu data if available - if len(rx.derived_products) > 0 and 'total_nu' not in g: - tgroup = g.create_group('total_nu') - rx.derived_products[0].to_hdf5(tgroup) + # Write total nu data if available + if len(rx.derived_products) > 0 and 'total_nu' not in g: + tgroup = g.create_group('total_nu') + rx.derived_products[0].to_hdf5(tgroup) - # Write unresolved resonance probability tables - if self.urr: - urr_group = g.create_group('urr') - for temperature, urr in self.urr.items(): - tgroup = urr_group.create_group(temperature) - urr.to_hdf5(tgroup) + # Write unresolved resonance probability tables + if self.urr: + urr_group = g.create_group('urr') + for temperature, urr in self.urr.items(): + tgroup = urr_group.create_group(temperature) + urr.to_hdf5(tgroup) - # Write fission energy release data - if self.fission_energy is not None: - fer_group = g.create_group('fission_energy_release') - self.fission_energy.to_hdf5(fer_group) - - f.close() + # Write fission energy release data + if self.fission_energy is not None: + fer_group = g.create_group('fission_energy_release') + self.fission_energy.to_hdf5(fer_group) @classmethod def from_hdf5(cls, group_or_filename): @@ -543,7 +542,7 @@ class IncidentNeutron(EqualityMixin): data.reactions[rx.mt] = rx # Read total nu data if available - if rx.mt in (18, 19, 20, 21, 38) and 'total_nu' in group: + if rx.mt in FISSION_MTS and 'total_nu' in group: tgroup = group['total_nu'] rx.derived_products.append(Product.from_hdf5(tgroup)) @@ -617,18 +616,20 @@ class IncidentNeutron(EqualityMixin): absorption_xs = ace.xss[i + 2*n_energy : i + 3*n_energy] heating_number = ace.xss[i + 4*n_energy : i + 5*n_energy]*EV_PER_MEV - # Create redundant reactions (total, absorption, and heating) + # Create redundant reaction for total (MT=1) total = Reaction(1) total.xs[strT] = Tabulated1D(energy, total_xs) total.redundant = True data.reactions[1] = total + # Create redundant reaction for absorption (MT=101) if np.count_nonzero(absorption_xs) > 0: absorption = Reaction(101) absorption.xs[strT] = Tabulated1D(energy, absorption_xs) absorption.redundant = True data.reactions[101] = absorption + # Create redundant reaction for heating (MT=301) heating = Reaction(301) heating.xs[strT] = Tabulated1D(energy, heating_number*total_xs) heating.redundant = True diff --git a/openmc/data/photon.py b/openmc/data/photon.py index 4be2417b16..a85bff5448 100644 --- a/openmc/data/photon.py +++ b/openmc/data/photon.py @@ -676,8 +676,10 @@ class IncidentPhoton(EqualityMixin): """ if isinstance(group_or_filename, h5py.Group): group = group_or_filename + need_to_close = False else: h5file = h5py.File(str(group_or_filename), 'r') + need_to_close = True # Make sure version matches if 'version' in h5file.attrs: @@ -738,6 +740,10 @@ class IncidentPhoton(EqualityMixin): 'num_electrons', 'photon_energy'): data.bremsstrahlung[key] = rgroup[key][()] + # If HDF5 file was opened here, make sure it gets closed + if need_to_close: + h5file.close() + return data def export_to_hdf5(self, path, mode='a', libver='earliest'): @@ -747,7 +753,7 @@ class IncidentPhoton(EqualityMixin): ---------- path : str Path to write HDF5 file to - mode : {'r', 'r+', 'w', 'x', 'a'} + mode : {'r+', 'w', 'x', 'a'} Mode that is used to open the HDF5 file. This is the second argument to the :class:`h5py.File` constructor. libver : {'earliest', 'latest'} @@ -755,69 +761,69 @@ class IncidentPhoton(EqualityMixin): that are less backwards compatible but have performance benefits. """ - # Open file and write version - f = h5py.File(str(path), mode, libver=libver) - f.attrs['filetype'] = np.string_('data_photon') - if 'version' not in f.attrs: - f.attrs['version'] = np.array(HDF5_VERSION) + with h5py.File(str(path), mode, libver=libver) as f: + # Write filetype and version + f.attrs['filetype'] = np.string_('data_photon') + if 'version' not in f.attrs: + f.attrs['version'] = np.array(HDF5_VERSION) - group = f.create_group(self.name) - group.attrs['Z'] = Z = self.atomic_number + group = f.create_group(self.name) + group.attrs['Z'] = Z = self.atomic_number - # Determine union energy grid - union_grid = np.array([]) - for rx in self: - union_grid = np.union1d(union_grid, rx.xs.x) - group.create_dataset('energy', data=union_grid) + # Determine union energy grid + union_grid = np.array([]) + for rx in self: + union_grid = np.union1d(union_grid, rx.xs.x) + group.create_dataset('energy', data=union_grid) - # Write cross sections - shell_group = group.create_group('subshells') - designators = [] - for mt, rx in self.reactions.items(): - name, key = _REACTION_NAME[mt] - if mt in (502, 504, 515, 517, 522, 525): - sub_group = group.create_group(key) - elif mt >= 534 and mt <= 572: - # Subshell - designators.append(key) - sub_group = shell_group.create_group(key) + # Write cross sections + shell_group = group.create_group('subshells') + designators = [] + for mt, rx in self.reactions.items(): + name, key = _REACTION_NAME[mt] + if mt in (502, 504, 515, 517, 522, 525): + sub_group = group.create_group(key) + elif mt >= 534 and mt <= 572: + # Subshell + designators.append(key) + sub_group = shell_group.create_group(key) - # Write atomic relaxation - if self.atomic_relaxation is not None: - if key in self.atomic_relaxation.subshells: - self.atomic_relaxation.to_hdf5(sub_group, key) - else: - continue - - rx.to_hdf5(sub_group, union_grid, Z) - - shell_group.attrs['designators'] = np.array(designators, dtype='S') - - # Write Compton profiles - if self.compton_profiles: - compton_group = group.create_group('compton_profiles') - - profile = self.compton_profiles - compton_group.create_dataset('num_electrons', - data=profile['num_electrons']) - compton_group.create_dataset('binding_energy', - data=profile['binding_energy']) - - # Get electron momentum values - compton_group.create_dataset('pz', data=profile['J'][0].x) - - # Create/write 2D array of profiles - J = np.array([Jk.y for Jk in profile['J']]) - compton_group.create_dataset('J', data=J) - - # Write bremsstrahlung - if self.bremsstrahlung: - brem_group = group.create_group('bremsstrahlung') - for key, value in self.bremsstrahlung.items(): - if key == 'I': - brem_group.attrs[key] = value + # Write atomic relaxation + if self.atomic_relaxation is not None: + if key in self.atomic_relaxation.subshells: + self.atomic_relaxation.to_hdf5(sub_group, key) else: - brem_group.create_dataset(key, data=value) + continue + + rx.to_hdf5(sub_group, union_grid, Z) + + shell_group.attrs['designators'] = np.array(designators, dtype='S') + + # Write Compton profiles + if self.compton_profiles: + compton_group = group.create_group('compton_profiles') + + profile = self.compton_profiles + compton_group.create_dataset('num_electrons', + data=profile['num_electrons']) + compton_group.create_dataset('binding_energy', + data=profile['binding_energy']) + + # Get electron momentum values + compton_group.create_dataset('pz', data=profile['J'][0].x) + + # Create/write 2D array of profiles + J = np.array([Jk.y for Jk in profile['J']]) + compton_group.create_dataset('J', data=J) + + # Write bremsstrahlung + if self.bremsstrahlung: + brem_group = group.create_group('bremsstrahlung') + for key, value in self.bremsstrahlung.items(): + if key == 'I': + brem_group.attrs[key] = value + else: + brem_group.create_dataset(key, data=value) def _add_bremsstrahlung(self): """Add the data used in the thick-target bremsstrahlung approximation @@ -837,7 +843,8 @@ class IncidentPhoton(EqualityMixin): } filename = os.path.join(os.path.dirname(__file__), 'BREMX.DAT') - brem = open(filename, 'r').read().split() + with open(filename, 'r') as fh: + brem = fh.read().split() # Incident electron kinetic energy grid in eV _BREMSSTRAHLUNG['electron_energy'] = np.logspace(3, 9, 200) diff --git a/openmc/data/product.py b/openmc/data/product.py index f7778cfa96..a6b2fd89e5 100644 --- a/openmc/data/product.py +++ b/openmc/data/product.py @@ -31,25 +31,22 @@ class Product(EqualityMixin): delayed neutron precursor). A special value of 'total' is used when the yield represents particles from prompt and delayed sources. particle : str - What particle the reaction product is. + The particle type of the reaction product yield_ : openmc.data.Function1D Yield of secondary particle in the reaction. """ def __init__(self, particle='neutron'): - self.particle = particle - self.decay_rate = 0.0 - self.emission_mode = 'prompt' - self.distribution = [] self.applicability = [] - self.yield_ = Polynomial((1,)) # 0-order polynomial i.e. a constant + self.decay_rate = 0.0 + self.distribution = [] + self.emission_mode = 'prompt' + self.particle = particle + self.yield_ = Polynomial((1,)) # 0-order polynomial, i.e., a constant def __repr__(self): - if isinstance(self.yield_, Real): - return "".format( - self.particle, self.emission_mode, self.yield_) - elif isinstance(self.yield_, Tabulated1D): + if isinstance(self.yield_, Tabulated1D): if np.all(self.yield_.y == self.yield_.y[0]): return "".format( self.particle, self.emission_mode, self.yield_.y[0]) diff --git a/openmc/data/reaction.py b/openmc/data/reaction.py index 2283b084af..2e88021ecb 100644 --- a/openmc/data/reaction.py +++ b/openmc/data/reaction.py @@ -64,6 +64,8 @@ REACTION_NAME.update({i: '(n,3He{})'.format(i - 750) for i in range(750, 799)}) REACTION_NAME.update({i: '(n,a{})'.format(i - 800) for i in range(800, 849)}) REACTION_NAME.update({i: '(n,2n{})'.format(i - 875) for i in range(875, 891)}) +FISSION_MTS = (18, 19, 20, 21, 38) + def _get_products(ev, mt): """Generate products from MF=6 in an ENDF evaluation @@ -492,7 +494,7 @@ def _get_activation_products(ev, rx): # Determine if file 9/10 are present present = {9: False, 10: False} - for i in range(n_states): + for _ in range(n_states): if decay_sublib: items = get_cont_record(file_obj) else: @@ -1021,7 +1023,7 @@ class Reaction(EqualityMixin): neutron.yield_ = yield_ rx.products.append(neutron) else: - assert mt in (18, 19, 20, 21, 38) + assert mt in FISSION_MTS rx.products, rx.derived_products = _get_fission_products_ace(ace) for p in rx.products: @@ -1126,13 +1128,13 @@ class Reaction(EqualityMixin): # Get fission product yields (nu) as well as delayed neutron energy # distributions - if mt in (18, 19, 20, 21, 38): + if mt in FISSION_MTS: rx.products, rx.derived_products = _get_fission_products_endf(ev) if (6, mt) in ev.section: # Product angle-energy distribution for product in _get_products(ev, mt): - if mt in (18, 19, 20, 21, 38) and product.particle == 'neutron': + if mt in FISSION_MTS and product.particle == 'neutron': rx.products[0].applicability = product.applicability rx.products[0].distribution = product.distribution else: @@ -1177,7 +1179,7 @@ class Reaction(EqualityMixin): for dist in neutron.distribution: dist.angle = AngleDistribution.from_endf(ev, mt) - if mt in (18, 19, 20, 21, 38) and (5, mt) in ev.section: + if mt in FISSION_MTS and (5, mt) in ev.section: # For fission reactions, rx.products[0].applicability = neutron.applicability rx.products[0].distribution = neutron.distribution diff --git a/openmc/data/thermal.py b/openmc/data/thermal.py index 5f939e59e8..b940994e0e 100644 --- a/openmc/data/thermal.py +++ b/openmc/data/thermal.py @@ -444,7 +444,7 @@ class ThermalScattering(EqualityMixin): ---------- path : str Path to write HDF5 file to - mode : {'r', r+', 'w', 'x', 'a'} + mode : {'r+', 'w', 'x', 'a'} Mode that is used to open the HDF5 file. This is the second argument to the :class:`h5py.File` constructor. libver : {'earliest', 'latest'} diff --git a/openmc/deplete/chain.py b/openmc/deplete/chain.py index 960cf06456..b78590bc25 100644 --- a/openmc/deplete/chain.py +++ b/openmc/deplete/chain.py @@ -95,9 +95,11 @@ def replace_missing(product, decay_data): # Iterate until we find an existing nuclide while product not in decay_data: if Z > 98: + # Assume alpha decay occurs for Z=99 and above Z -= 2 A -= 4 else: + # Otherwise assume a beta- or beta+ if beta_minus: Z += 1 else: @@ -314,7 +316,7 @@ class Chain: nuclide.reactions.append(ReactionTuple( name, daughter, q_value, 1.0)) - if any(mt in reactions_available for mt in [18, 19, 20, 21, 38]): + if any(mt in reactions_available for mt in openmc.data.FISSION_MTS): if parent in fpy_data: q_value = reactions[parent][18] nuclide.reactions.append( @@ -334,10 +336,10 @@ class Chain: yield_energies = [0.0] yield_data = {} - for E, table in zip(yield_energies, fpy.independent): + for E, yield_table in zip(yield_energies, fpy.independent): yield_replace = 0.0 yields = defaultdict(float) - for product, y in table.items(): + for product, y in yield_table.items(): # Handle fission products that have no decay data if product not in decay_data: daughter = replace_missing(product, decay_data) @@ -735,26 +737,27 @@ class Chain: rxn_Q = parent.reactions[rxn_index[0]].Q # Remove existing reactions - for ix in reversed(rxn_index): parent.reactions.pop(ix) + # Add new reactions all_meta = True - - for tgt, br in new_ratios.items(): - all_meta = all_meta and ("_m" in tgt) + for target, br in new_ratios.items(): + all_meta = all_meta and ("_m" in target) parent.reactions.append(ReactionTuple( - reaction, tgt, rxn_Q, br)) + reaction, target, rxn_Q, br)) + # If branching ratios don't add to unity, add reaction to ground + # with remainder of branching ratio if all_meta and sums[parent_name] != 1.0: ground_br = 1.0 - sums[parent_name] - ground_tgt = grounds.get(parent_name) - if ground_tgt is None: + ground_target = grounds.get(parent_name) + if ground_target is None: pz, pa, pm = zam(parent_name) - ground_tgt = gnd_name(pz, pa + 1, 0) - new_ratios[ground_tgt] = ground_br + ground_target = gnd_name(pz, pa + 1, 0) + new_ratios[ground_target] = ground_br parent.reactions.append(ReactionTuple( - reaction, ground_tgt, rxn_Q, ground_br)) + reaction, ground_target, rxn_Q, ground_br)) @property def fission_yields(self): diff --git a/openmc/deplete/nuclide.py b/openmc/deplete/nuclide.py index 06f217f222..c7e05b235b 100644 --- a/openmc/deplete/nuclide.py +++ b/openmc/deplete/nuclide.py @@ -232,9 +232,9 @@ class Nuclide: elem.set('half_life', str(self.half_life)) elem.set('decay_modes', str(len(self.decay_modes))) elem.set('decay_energy', str(self.decay_energy)) - for mode, daughter, br in self.decay_modes: + for mode_type, daughter, br in self.decay_modes: mode_elem = ET.SubElement(elem, 'decay') - mode_elem.set('type', mode) + mode_elem.set('type', mode_type) mode_elem.set('target', daughter or "Nothing") mode_elem.set('branching_ratio', str(br)) diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py index 5f86dd5e8a..19895b4d9c 100644 --- a/openmc/deplete/results.py +++ b/openmc/deplete/results.py @@ -12,7 +12,7 @@ import numpy as np from . import comm, have_mpi, MPI from .reaction_rates import ReactionRates -_VERSION_RESULTS = (1, 0) +VERSION_RESULTS = (1, 0) __all__ = ["Results"] @@ -174,8 +174,8 @@ class Results: """ new = Results() new.volume = {lm: self.volume[lm] for lm in local_materials} - new.mat_to_ind = dict(zip( - local_materials, range(len(local_materials)))) + new.mat_to_ind = {mat: idx for (idx, mat) in enumerate(local_materials)} + # Direct transfer direct_attrs = ("time", "k", "power", "nuc_to_ind", "mat_to_hdf5_ind", "proc_time") @@ -228,7 +228,7 @@ class Results: # Store concentration mat and nuclide dictionaries (along with volumes) - handle.attrs['version'] = np.array(_VERSION_RESULTS) + handle.attrs['version'] = np.array(VERSION_RESULTS) handle.attrs['filetype'] = np.string_('depletion results') mat_list = sorted(self.mat_to_hdf5_ind, key=int) diff --git a/openmc/deplete/results_list.py b/openmc/deplete/results_list.py index 6267a0fccc..f9f50b8220 100644 --- a/openmc/deplete/results_list.py +++ b/openmc/deplete/results_list.py @@ -1,7 +1,7 @@ import h5py import numpy as np -from .results import Results, _VERSION_RESULTS +from .results import Results, VERSION_RESULTS from openmc.checkvalue import check_filetype_version, check_value @@ -30,7 +30,7 @@ class ResultsList(list): New instance of depletion results """ with h5py.File(str(filename), "r") as fh: - check_filetype_version(fh, 'depletion results', _VERSION_RESULTS[0]) + check_filetype_version(fh, 'depletion results', VERSION_RESULTS[0]) new = cls() # Get number of results stored @@ -68,9 +68,9 @@ class ResultsList(list): Returns ------- - time : numpy.ndarray + times : numpy.ndarray Array of times in units of ``time_units`` - concentration : numpy.ndarray + concentrations : numpy.ndarray Concentration of specified nuclide in units of ``nuc_units`` """ @@ -78,30 +78,30 @@ class ResultsList(list): check_value("nuc_units", nuc_units, {"atoms", "atom/b-cm", "atom/cm3"}) - time = np.empty_like(self, dtype=float) - concentration = np.empty_like(self, dtype=float) + times = np.empty_like(self, dtype=float) + concentrations = np.empty_like(self, dtype=float) # Evaluate value in each region for i, result in enumerate(self): - time[i] = result.time[0] - concentration[i] = result[0, mat, nuc] + times[i] = result.time[0] + concentrations[i] = result[0, mat, nuc] # Unit conversions if time_units == "d": - time /= (60 * 60 * 24) + times /= (60 * 60 * 24) elif time_units == "h": - time /= (60 * 60) + times /= (60 * 60) elif time_units == "min": - time /= 60 + times /= 60 if nuc_units != "atoms": # Divide by volume to get density - concentration /= self[0].volume[mat] + concentrations /= self[0].volume[mat] if nuc_units == "atom/b-cm": # 1 barn = 1e-24 cm^2 - concentration *= 1e-24 + concentrations *= 1e-24 - return time, concentration + return times, concentrations def get_reaction_rate(self, mat, nuc, rx): """Get reaction rate in a single material/nuclide over time @@ -125,44 +125,44 @@ class ResultsList(list): Returns ------- - time : numpy.ndarray + times : numpy.ndarray Array of times in [s] - rate : numpy.ndarray + rates : numpy.ndarray Array of reaction rates """ - time = np.empty_like(self, dtype=float) - rate = np.empty_like(self, dtype=float) + times = np.empty_like(self, dtype=float) + rates = np.empty_like(self, dtype=float) # Evaluate value in each region for i, result in enumerate(self): - time[i] = result.time[0] - rate[i] = result.rates[0].get(mat, nuc, rx) * result[0, mat, nuc] + times[i] = result.time[0] + rates[i] = result.rates[0].get(mat, nuc, rx) * result[0, mat, nuc] - return time, rate + return times, rates def get_eigenvalue(self): """Evaluates the eigenvalue from a results list. Returns ------- - time : numpy.ndarray + times : numpy.ndarray Array of times in [s] - eigenvalue : numpy.ndarray + eigenvalues : numpy.ndarray k-eigenvalue at each time. Column 0 contains the eigenvalue, while column 1 contains the associated uncertainty """ - time = np.empty_like(self, dtype=float) - eigenvalue = np.empty((len(self), 2), dtype=float) + times = np.empty_like(self, dtype=float) + eigenvalues = np.empty((len(self), 2), dtype=float) # Get time/eigenvalue at each point for i, result in enumerate(self): - time[i] = result.time[0] - eigenvalue[i] = result.k[0] + times[i] = result.time[0] + eigenvalues[i] = result.k[0] - return time, eigenvalue + return times, eigenvalues def get_depletion_time(self): """Return an array of the average time to deplete a material @@ -175,7 +175,7 @@ class ResultsList(list): Returns ------- - times : :class:`numpy.ndarray` + times : numpy.ndarray Vector of average time to deplete a single material across all processes and materials. diff --git a/tests/unit_tests/test_endf.py b/tests/unit_tests/test_endf.py index 11a8db4432..9e69708673 100644 --- a/tests/unit_tests/test_endf.py +++ b/tests/unit_tests/test_endf.py @@ -22,6 +22,7 @@ def test_float_endf(): assert endf.float_endf('1.+2') == approx(100.0) assert endf.float_endf('-1.+2') == approx(-100.0) assert endf.float_endf(' ') == 0.0 + assert endf.float_endf('9.876540000000000') == approx(9.87654) def test_int_endf():