From 3b8d1ffb4b21a0b64a9577f9dcd81e86303695c0 Mon Sep 17 00:00:00 2001 From: jingang Date: Wed, 15 Aug 2018 22:08:12 -0400 Subject: [PATCH] update for new wmp library v1.0 - python api and docs --- docs/source/io_formats/data_wmp.rst | 60 +--- docs/source/methods/cross_sections.rst | 5 +- openmc/data/__init__.py | 2 +- openmc/data/multipole.py | 364 +++++++------------------ 4 files changed, 111 insertions(+), 320 deletions(-) diff --git a/docs/source/io_formats/data_wmp.rst b/docs/source/io_formats/data_wmp.rst index 6174e85ff7..25916c5229 100644 --- a/docs/source/io_formats/data_wmp.rst +++ b/docs/source/io_formats/data_wmp.rst @@ -5,7 +5,7 @@ Windowed Multipole Library Format ================================= **/version** (*char[]*) - The format version of the file. The current version is "v0.2" + The format version of the file. The current version is "v1.0" **/nuclide/** - **broaden_poly** (*int[]*) @@ -23,55 +23,25 @@ Windowed Multipole Library Format \text{data}[:,i] = [\text{pole},~\text{residue}_1,~\text{residue}_2, ~\ldots] - The residues are in the order: total, competitive if present, - absorption, fission. Complex numbers are stored by forming a type with - ":math:`r`" and ":math:`i`" identifiers, similar to how `h5py`_ does it. - - **end_E** (*double*) + The residues are in the order: scattering, absorption, fission. Complex + numbers are stored by forming a type with ":math:`r`" and ":math:`i`" + identifiers, similar to how `h5py`_ does it. + - **E_max** (*double*) Highest energy the windowed multipole part of the library is valid for. - - **formalism** (*int*) - The formalism of the underlying data. Uses the `ENDF-6`_ format - formalism numbers. - - .. table:: Table of supported formalisms. - - +-------------+------------------+ - | Formalism | Formalism number | - +=============+==================+ - | MLBW | 2 | - +-------------+------------------+ - | Reich-Moore | 3 | - +-------------+------------------+ - - - **l_value** (*int[]*) - The index for a corresponding pole. Equivalent to the :math:`l` quantum - number of the resonance the pole comes from :math:`+1`. - - **pseudo_K0RS** (*double[]*) - :math:`l` dependent value of - - .. math:: - \sqrt{\frac{2 m_n}{\hbar}}\frac{AWR}{AWR + 1} r_{s,l} - - Where :math:`m_n` is mass of neutron, :math:`AWR` is the atomic weight - ratio of the target to the neutron, and :math:`r_{s,l}` is the - scattering radius for a given :math:`l`. + - **E_min** (*double*) + Lowest energy the windowed multipole part of the library is valid for. - **spacing** (*double*) .. math:: - \frac{\sqrt{E_{max}}- \sqrt{E_{min}}}{n_w} + \frac{\sqrt{E_{max}} - \sqrt{E_{min}}}{n_w} - Where :math:`E_{max}` is the maximum energy the windows go up to. This - is not equivalent to the maximum energy for which the windowed multipole - data is valid for. It is slightly higher to ensure an integer number of - windows. :math:`E_{min}` is the minimum energy and equivalent to - ``start_E``, and :math:`n_w` is the number of windows, given by - ``windows``. + Where :math:`E_{max}` is the maximum energy the windows go up to. + :math:`E_{min}` is the minimum energy, and :math:`n_w` is the number of + windows, given by ``windows``. - **sqrtAWR** (*double*) Square root of the atomic weight ratio. - - **start_E** (*double*) - Lowest energy the windowed multipole part of the library is valid for. - - **w_start** (*int[]*) - The pole to start from for each window. - - **w_end** (*int[]*) - The pole to end at for each window. + - **windows** (*int[][]*) + The poles to start from and end at for each window. windows[i, 0] and + windows[i, 1] are, respectively, the indexes (1-based) of the first and + last pole in window i. .. _h5py: http://docs.h5py.org/en/latest/ -.. _ENDF-6: https://www.oecd-nea.org/dbdata/data/manual-endf/endf102.pdf diff --git a/docs/source/methods/cross_sections.rst b/docs/source/methods/cross_sections.rst index 9f677bac43..4940275d3a 100644 --- a/docs/source/methods/cross_sections.rst +++ b/docs/source/methods/cross_sections.rst @@ -90,7 +90,7 @@ Assuming free-gas thermal motion, cross sections in the multipole form can be analytically Doppler broadened to give the form: .. math:: - \sigma(E, T) = \frac{1}{2 E \sqrt{\xi}} \sum_j \text{Re} \left[i r_j + \sigma(E, T) = \frac{1}{2 E \sqrt{\xi}} \sum_j \text{Re} \left[r_j \sqrt{\pi} W_i(z) - \frac{r_j}{\sqrt{\pi}} C \left(\frac{p_j}{\sqrt{\xi}}, \frac{u}{2 \sqrt{\xi}}\right)\right] .. math:: @@ -141,7 +141,7 @@ scattering does not occur in the resolved resonance region. This is usually, but not always the case. Future library versions may eliminate this issue. The data format used by OpenMC to represent windowed multipole data is specified -in :ref:`io_data_wmp`. +in :ref:`io_data_wmp` with a publicly available `WMP library`_. .. _temperature_treatment: @@ -270,6 +270,7 @@ or even isotropic scattering. https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-14-24530.pdf .. _Hwang: http://www.ans.org/pubs/journals/nse/a_16381 .. _Josey: http://dx.doi.org/10.1016/j.jcp.2015.08.013 +.. _WMP Library: https://github.com/mit-crpg/WMP_Library .. _MCNP: http://mcnp.lanl.gov .. _Serpent: http://montecarlo.vtt.fi .. _NJOY: http://t2.lanl.gov/codes.shtml diff --git a/openmc/data/__init__.py b/openmc/data/__init__.py index 390a67c9b9..dec11cd3a9 100644 --- a/openmc/data/__init__.py +++ b/openmc/data/__init__.py @@ -4,7 +4,7 @@ HDF5_VERSION_MINOR = 0 HDF5_VERSION = (HDF5_VERSION_MAJOR, HDF5_VERSION_MINOR) # Version of WMP nuclear data format -WMP_VERSION = 'v0.2' +WMP_VERSION = 'v1.0' from .data import * diff --git a/openmc/data/multipole.py b/openmc/data/multipole.py index f7a78d9532..c6e4f123c9 100644 --- a/openmc/data/multipole.py +++ b/openmc/data/multipole.py @@ -10,26 +10,16 @@ import openmc.checkvalue as cv from openmc.mixin import EqualityMixin -# Formalisms -_FORM_MLBW = 2 -_FORM_RM = 3 - # Constants that determine which value to access _MP_EA = 0 # Pole -# Reich-Moore indices -_RM_RT = 1 # Residue total -_RM_RA = 2 # Residue absorption -_RM_RF = 3 # Residue fission - -# Multi-level Breit Wigner indices -_MLBW_RT = 1 # Residue total -_MLBW_RX = 2 # Residue competitive -_MLBW_RA = 3 # Residue absorption -_MLBW_RF = 4 # Residue fission +# Residue indices +_MP_RS = 1 # Residue scattering +_MP_RA = 2 # Residue absorption +_MP_RF = 3 # Residue fission # Polynomial fit indices -_FIT_T = 0 # Total +_FIT_S = 0 # Scattering _FIT_A = 1 # Absorption _FIT_F = 2 # Fission @@ -143,98 +133,61 @@ class WindowedMultipole(EqualityMixin): Parameters ---------- - formalism : {'MLBW', 'RM'} - The R-matrix formalism used to reconstruct resonances. Either 'MLBW' - for multi-level Breit Wigner or 'RM' for Reich-Moore. Attributes ---------- - num_l : Integral - Number of possible l quantum states for this nuclide. fit_order : Integral Order of the windowed curvefit. fissionable : bool Whether or not the target nuclide has fission data. - formalism : {'MLBW', 'RM'} - The R-matrix formalism used to reconstruct resonances. Either 'MLBW' - for multi-level Breit Wigner or 'RM' for Reich-Moore. spacing : Real The width of each window in sqrt(E)-space. For example, the frst window - will end at (sqrt(start_E) + spacing)**2 and the second window at - (sqrt(start_E) + 2*spacing)**2. + will end at (sqrt(E_min) + spacing)**2 and the second window at + (sqrt(E_min) + 2*spacing)**2. sqrtAWR : Real Square root of the atomic weight ratio of the target nuclide. - start_E : Real + E_min : Real Lowest energy in eV the library is valid for. - end_E : Real + E_max : Real Highest energy in eV the library is valid for. data : np.ndarray A 2D array of complex poles and residues. data[i, 0] gives the energy at which pole i is located. data[i, 1:] gives the residues associated - with the i-th pole. There are 3 residues for Reich-Moore data, one each - for the total, absorption, and fission channels. Multi-level - Breit Wigner data has an additional residue for the competitive channel. - pseudo_k0RS : np.ndarray - A 1D array of Real values. There is one value for each valid l - quantum number. The values are equal to - sqrt(2 m / hbar) * AWR / (AWR + 1) * r - where m is the neutron mass, AWR is the atomic weight ratio, and r - is the l-dependent scattering radius. - l_value : np.ndarray - A 1D array of Integral values equal to the l quantum number for each - pole + 1. - w_start : np.ndarray - A 1D array of Integral values. w_start[i] - 1 is the index of the first - pole in window i. - w_end : np.ndarray - A 1D array of Integral values. w_end[i] - 1 is the index of the last - pole in window i. + with the i-th pole. There are 3 residues, one each for the scattering, + absorption, and fission channels. + windows : np.ndarray + A 2D array of Integral values. windows[i, 0] - 1 is the index of the + first pole in window i. windows[i, 1] - 1 is the index of the last pole + in window i. broaden_poly : np.ndarray A 1D array of boolean values indicating whether or not the polynomial curvefit in that window should be Doppler broadened. curvefit : np.ndarray A 3D array of Real curvefit polynomial coefficients. curvefit[i, 0, :] - gives coefficients for the total cross section in window i. + gives coefficients for the scattering cross section in window i. curvefit[i, 1, :] gives absorption coefficients and curvefit[i, 2, :] gives fission coefficients. The polynomial terms are increasing powers of sqrt(E) starting with 1/E e.g: a/E + b/sqrt(E) + c + d sqrt(E) + ... """ - def __init__(self, formalism): - self._num_l = None - self.formalism = formalism + def __init__(self): self.spacing = None self.sqrtAWR = None - self.start_E = None - self.end_E = None + self.E_min = None + self.E_max = None self.data = None - self.pseudo_k0RS = None - self.l_value = None - self.w_start = None - self.w_end = None + self.windows = None self.broaden_poly = None self.curvefit = None - @property - def num_l(self): - return self._num_l - @property def fit_order(self): return self.curvefit.shape[1] - 1 @property def fissionable(self): - if self.formalism == 'RM': - return self.data.shape[1] == 4 - else: - # Assume self.formalism == 'MLBW' - return self.data.shape[1] == 5 - - @property - def formalism(self): - return self._formalism + return self.data.shape[1] == 4 @property def spacing(self): @@ -245,32 +198,24 @@ class WindowedMultipole(EqualityMixin): return self._sqrtAWR @property - def start_E(self): - return self._start_E + def E_min(self): + return self._E_min @property - def end_E(self): - return self._end_E + def E_max(self): + return self._E_max @property def data(self): return self._data - @property - def pseudo_k0RS(self): - return self._pseudo_k0RS - @property def l_value(self): return self._l_value @property - def w_start(self): - return self._w_start - - @property - def w_end(self): - return self._w_end + def windows(self): + return self._windows @property def broaden_poly(self): @@ -280,116 +225,64 @@ class WindowedMultipole(EqualityMixin): def curvefit(self): return self._curvefit - @formalism.setter - def formalism(self, formalism): - cv.check_type('formalism', formalism, str) - cv.check_value('formalism', formalism, ('MLBW', 'RM')) - self._formalism = formalism - @spacing.setter def spacing(self, spacing): if spacing is not None: - cv.check_type('spacing', spacing, Real) - cv.check_greater_than('spacing', spacing, 0.0, equality=False) + check_type('spacing', spacing, Real) + check_greater_than('spacing', spacing, 0.0, equality=False) self._spacing = spacing @sqrtAWR.setter def sqrtAWR(self, sqrtAWR): if sqrtAWR is not None: - cv.check_type('sqrtAWR', sqrtAWR, Real) - cv.check_greater_than('sqrtAWR', sqrtAWR, 0.0, equality=False) + check_type('sqrtAWR', sqrtAWR, Real) + check_greater_than('sqrtAWR', sqrtAWR, 0.0, equality=False) self._sqrtAWR = sqrtAWR - @start_E.setter - def start_E(self, start_E): - if start_E is not None: - cv.check_type('start_E', start_E, Real) - cv.check_greater_than('start_E', start_E, 0.0, equality=True) - self._start_E = start_E + @E_min.setter + def E_min(self, E_min): + if E_min is not None: + check_type('E_min', E_min, Real) + check_greater_than('E_min', E_min, 0.0, equality=True) + self._E_min = E_min - @end_E.setter - def end_E(self, end_E): - if end_E is not None: - cv.check_type('end_E', end_E, Real) - cv.check_greater_than('end_E', end_E, 0.0, equality=False) - self._end_E = end_E + @E_max.setter + def E_max(self, E_max): + if E_max is not None: + check_type('E_max', E_max, Real) + check_greater_than('E_max', E_max, 0.0, equality=False) + self._E_max = E_max @data.setter def data(self, data): if data is not None: - cv.check_type('data', data, np.ndarray) + check_type('data', data, np.ndarray) if len(data.shape) != 2: raise ValueError('Multipole data arrays must be 2D') - if self.formalism == 'RM': - if data.shape[1] not in (3, 4): - raise ValueError('For the Reich-Moore formalism, ' - 'data.shape[1] must be 3 or 4. One value for the pole.' - ' One each for the total and absorption residues. ' - 'Possibly one more for a fission residue.') - else: - # Assume self.formalism == 'MLBW' - if data.shape[1] not in (4, 5): - raise ValueError('For the Multi-level Breit-Wigner ' - 'formalism, data.shape[1] must be 4 or 5. One value ' - 'for the pole. One each for the total, competitive, ' - 'and absorption residues. Possibly one more for a ' - 'fission residue.') + if data.shape[1] not in (3, 4): + raise ValueError( + 'data.shape[1] must be 3 or 4. One value for the pole.' + ' One each for the scattering and absorption residues. ' + 'Possibly one more for a fission residue.') if not np.issubdtype(data.dtype, complex): raise TypeError('Multipole data arrays must be complex dtype') self._data = data - @pseudo_k0RS.setter - def pseudo_k0RS(self, pseudo_k0RS): - if pseudo_k0RS is not None: - cv.check_type('pseudo_k0RS', pseudo_k0RS, np.ndarray) - if len(pseudo_k0RS.shape) != 1: - raise ValueError('Multipole pseudo_k0RS arrays must be 1D') - if not np.issubdtype(pseudo_k0RS.dtype, float): - raise TypeError('Multipole data arrays must be float dtype') - self._pseudo_k0RS = pseudo_k0RS - - @l_value.setter - def l_value(self, l_value): - if l_value is not None: - cv.check_type('l_value', l_value, np.ndarray) - if len(l_value.shape) != 1: - raise ValueError('Multipole l_value arrays must be 1D') - if not np.issubdtype(l_value.dtype, int): - raise TypeError('Multipole l_value arrays must be integer' + @windows.setter + def windows(self, windows): + if windows is not None: + check_type('windows', windows, np.ndarray) + if len(windows.shape) != 2: + raise ValueError('Multipole windows arrays must be 2D') + if not np.issubdtype(windows.dtype, int): + raise TypeError('Multipole windows arrays must be integer' ' dtype') - - self._num_l = len(np.unique(l_value)) - - else: - self._num_l = None - - self._l_value = l_value - - @w_start.setter - def w_start(self, w_start): - if w_start is not None: - cv.check_type('w_start', w_start, np.ndarray) - if len(w_start.shape) != 1: - raise ValueError('Multipole w_start arrays must be 1D') - if not np.issubdtype(w_start.dtype, int): - raise TypeError('Multipole w_start arrays must be integer' - ' dtype') - self._w_start = w_start - - @w_end.setter - def w_end(self, w_end): - if w_end is not None: - cv.check_type('w_end', w_end, np.ndarray) - if len(w_end.shape) != 1: - raise ValueError('Multipole w_end arrays must be 1D') - if not np.issubdtype(w_end.dtype, int): - raise TypeError('Multipole w_end arrays must be integer dtype') - self._w_end = w_end + self._windows = windows @broaden_poly.setter def broaden_poly(self, broaden_poly): if broaden_poly is not None: - cv.check_type('broaden_poly', broaden_poly, np.ndarray) + check_type('broaden_poly', broaden_poly, np.ndarray) if len(broaden_poly.shape) != 1: raise ValueError('Multipole broaden_poly arrays must be 1D') if not np.issubdtype(broaden_poly.dtype, bool): @@ -400,10 +293,10 @@ class WindowedMultipole(EqualityMixin): @curvefit.setter def curvefit(self, curvefit): if curvefit is not None: - cv.check_type('curvefit', curvefit, np.ndarray) + check_type('curvefit', curvefit, np.ndarray) if len(curvefit.shape) != 3: raise ValueError('Multipole curvefit arrays must be 3D') - if curvefit.shape[2] not in (2, 3): # sig_t, sig_a (maybe sig_f) + if curvefit.shape[2] not in (2, 3): # sig_s, sig_a (maybe sig_f) raise ValueError('The third dimension of multipole curvefit' ' arrays must have a length of 2 or 3') if not np.issubdtype(curvefit.dtype, float): @@ -443,19 +336,14 @@ class WindowedMultipole(EqualityMixin): 'Python API expects version ' + WMP_VERSION) group = h5file['nuclide'] - # Read scalars. + out = cls() - if group['formalism'].value == _FORM_MLBW: - out = cls('MLBW') - elif group['formalism'].value == _FORM_RM: - out = cls('RM') - else: - raise ValueError('Unrecognized/Unsupported R-matrix formalism') + # Read scalars. out.spacing = group['spacing'].value out.sqrtAWR = group['sqrtAWR'].value - out.start_E = group['start_E'].value - out.end_E = group['end_E'].value + out.E_min = group['E_min'].value + out.E_max = group['E_max'].value # Read arrays. @@ -463,27 +351,15 @@ class WindowedMultipole(EqualityMixin): out.data = group['data'].value - out.l_value = group['l_value'].value - if out.l_value.shape[0] != out.data.shape[0]: - raise ValueError(err.format('l_value', 'data')) - - out.pseudo_k0RS = group['pseudo_K0RS'].value - if out.pseudo_k0RS.shape[0] != out.num_l: - raise ValueError(err.format('pseudo_k0RS', 'l_value')) - - out.w_start = group['w_start'].value - - out.w_end = group['w_end'].value - if out.w_end.shape[0] != out.w_start.shape[0]: - raise ValueError(err.format('w_end', 'w_start')) + out.windows = group['windows'].value out.broaden_poly = group['broaden_poly'].value.astype(np.bool) - if out.broaden_poly.shape[0] != out.w_start.shape[0]: - raise ValueError(err.format('broaden_poly', 'w_start')) + if out.broaden_poly.shape[0] != out.windows.shape[0]: + raise ValueError(err.format('broaden_poly', 'windows')) out.curvefit = group['curvefit'].value - if out.curvefit.shape[0] != out.w_start.shape[0]: - raise ValueError(err.format('curvefit', 'w_start')) + if out.curvefit.shape[0] != out.windows.shape[0]: + raise ValueError(err.format('curvefit', 'windows')) # _broaden_wmp_polynomials assumes the curve fit has at least 3 terms. if out.fit_order < 2: @@ -493,7 +369,7 @@ class WindowedMultipole(EqualityMixin): return out def _evaluate(self, E, T): - """Compute total, absorption, and fission cross sections. + """Compute scattering, absorption, and fission cross sections. Parameters ---------- @@ -510,8 +386,8 @@ class WindowedMultipole(EqualityMixin): """ - if E < self.start_E: return (0, 0, 0) - if E > self.end_E: return (0, 0, 0) + if E < self.E_min: return (0, 0, 0) + if E > self.E_max: return (0, 0, 0) # ====================================================================== # Bookkeeping @@ -525,34 +401,12 @@ class WindowedMultipole(EqualityMixin): # the 1-based vs. 0-based indexing. Similarly startw needs to be # decreased by 1. endw does not need to be decreased because # range(startw, endw) does not include endw. - i_window = int(np.floor((sqrtE - sqrt(self.start_E)) / self.spacing)) - startw = self.w_start[i_window] - 1 - endw = self.w_end[i_window] - - # Fill in factors. Because of the unique interference dips in scatering - # resonances, the total cross section has a special "factor" that does - # not appear in the absorption and fission equations. - if startw <= endw: - twophi = np.zeros(self.num_l, dtype=np.float) - sig_t_factor = np.zeros(self.num_l, dtype=np.cfloat) - - for iL in range(self.num_l): - twophi[iL] = self.pseudo_k0RS[iL] * sqrtE - if iL == 1: - twophi[iL] = twophi[iL] - np.arctan(twophi[iL]) - elif iL == 2: - arg = 3.0 * twophi[iL] / (3.0 - twophi[iL]**2) - twophi[iL] = twophi[iL] - np.arctan(arg) - elif iL == 3: - arg = (twophi[iL] * (15.0 - twophi[iL]**2) - / (15.0 - 6.0 * twophi[iL]**2)) - twophi[iL] = twophi[iL] - np.arctan(arg) - - twophi = 2.0 * twophi - sig_t_factor = np.cos(twophi) - 1j*np.sin(twophi) + i_window = int(np.floor((sqrtE - sqrt(self.E_min)) / self.spacing)) + startw = self.windows[i_window, 0] - 1 + endw = self.windows[i_window, 1] # Initialize the ouptut cross sections. - sig_t = 0.0 + sig_s = 0.0 sig_a = 0.0 sig_f = 0.0 @@ -565,7 +419,7 @@ class WindowedMultipole(EqualityMixin): broadened_polynomials = _broaden_wmp_polynomials(E, dopp, self.fit_order + 1) for i_poly in range(self.fit_order+1): - sig_t += (self.curvefit[i_window, i_poly, _FIT_T] + sig_s += (self.curvefit[i_window, i_poly, _FIT_S] * broadened_polynomials[i_poly]) sig_a += (self.curvefit[i_window, i_poly, _FIT_A] * broadened_polynomials[i_poly]) @@ -575,7 +429,7 @@ class WindowedMultipole(EqualityMixin): else: temp = invE for i_poly in range(self.fit_order+1): - sig_t += self.curvefit[i_window, i_poly, _FIT_T] * temp + sig_s += self.curvefit[i_window, i_poly, _FIT_S] * temp sig_a += self.curvefit[i_window, i_poly, _FIT_A] * temp if self.fissionable: sig_f += self.curvefit[i_window, i_poly, _FIT_F] * temp @@ -589,22 +443,10 @@ class WindowedMultipole(EqualityMixin): for i_pole in range(startw, endw): psi_chi = -1j / (self.data[i_pole, _MP_EA] - sqrtE) c_temp = psi_chi / E - if self.formalism == 'MLBW': - sig_t += ((self.data[i_pole, _MLBW_RT] * c_temp * - sig_t_factor[self.l_value[i_pole]-1]).real - + (self.data[i_pole, _MLBW_RX] * c_temp).real) - sig_a += (self.data[i_pole, _MLBW_RA] * c_temp).real - if self.fissionable: - sig_f += (self.data[i_pole, _MLBW_RF] * c_temp).real - elif self.formalism == 'RM': - sig_t += (self.data[i_pole, _RM_RT] * c_temp * - sig_t_factor[self.l_value[i_pole]-1]).real - sig_a += (self.data[i_pole, _RM_RA] * c_temp).real - if self.fissionable: - sig_f += (self.data[i_pole, _RM_RF] * c_temp).real - else: - raise ValueError('Unrecognized/Unsupported R-matrix' - ' formalism') + sig_s += (self.data[i_pole, _MP_RS] * c_temp).real + sig_a += (self.data[i_pole, _MP_RA] * c_temp).real + if self.fissionable: + sig_f += (self.data[i_pole, _MP_RF] * c_temp).real else: # At temperature, use Faddeeva function-based form. @@ -612,27 +454,15 @@ class WindowedMultipole(EqualityMixin): for i_pole in range(startw, endw): Z = (sqrtE - self.data[i_pole, _MP_EA]) * dopp w_val = _faddeeva(Z) * dopp * invE * sqrt(pi) - if self.formalism == 'MLBW': - sig_t += ((self.data[i_pole, _MLBW_RT] * - sig_t_factor[self.l_value[i_pole]-1] + - self.data[i_pole, _MLBW_RX]) * w_val).real - sig_a += (self.data[i_pole, _MLBW_RA] * w_val).real - if self.fissionable: - sig_f += (self.data[i_pole, _MLBW_RF] * w_val).real - elif self.formalism == 'RM': - sig_t += (self.data[i_pole, _RM_RT] * w_val * - sig_t_factor[self.l_value[i_pole]-1]).real - sig_a += (self.data[i_pole, _RM_RA] * w_val).real - if self.fissionable: - sig_f += (self.data[i_pole, _RM_RF] * w_val).real - else: - raise ValueError('Unrecognized/Unsupported R-matrix' - ' formalism') + sig_s += (self.data[i_pole, _MP_RS] * w_val).real + sig_a += (self.data[i_pole, _MP_RA] * w_val).real + if self.fissionable: + sig_f += (self.data[i_pole, _MP_RF] * w_val).real - return sig_t, sig_a, sig_f + return sig_s, sig_a, sig_f def __call__(self, E, T): - """Compute total, absorption, and fission cross sections. + """Compute scattering, absorption, and fission cross sections. Parameters ---------- @@ -674,24 +504,14 @@ class WindowedMultipole(EqualityMixin): g = f.create_group('nuclide') # Write scalars. - if self.formalism == 'MLBW': - g.create_dataset('formalism', - data=np.array(_FORM_MLBW, dtype=np.int32)) - else: - # Assume RM. - g.create_dataset('formalism', - data=np.array(_FORM_RM, dtype=np.int32)) g.create_dataset('spacing', data=np.array(self.spacing)) g.create_dataset('sqrtAWR', data=np.array(self.sqrtAWR)) - g.create_dataset('start_E', data=np.array(self.start_E)) - g.create_dataset('end_E', data=np.array(self.end_E)) + g.create_dataset('E_min', data=np.array(self.E_min)) + g.create_dataset('E_max', data=np.array(self.E_max)) # Write arrays. g.create_dataset('data', data=self.data) - g.create_dataset('l_value', data=self.l_value) - g.create_dataset('pseudo_K0RS', data=self.pseudo_k0RS) - g.create_dataset('w_start', data=self.w_start) - g.create_dataset('w_end', data=self.w_end) + g.create_dataset('windows', data=self.windows) g.create_dataset('broaden_poly', data=self.broaden_poly.astype(np.int8)) g.create_dataset('curvefit', data=self.curvefit)