From 4ebeb737ab7310d8f36bda7d4e2f371b999d0b31 Mon Sep 17 00:00:00 2001 From: amandalund Date: Thu, 5 Jul 2018 23:25:42 -0500 Subject: [PATCH] Address @paulromano comments in #1023 --- openmc/data/endf.py | 5 ++- openmc/data/photon.py | 35 ++++++++------- scripts/openmc-get-nndc-data | 5 +-- scripts/openmc-get-photon-data | 43 ++++++++++++++++--- scripts/openmc-make-compton | 5 +-- src/material_header.F90 | 11 +---- src/math_functions.cpp | 13 +++--- src/math_functions.h | 22 +++++----- .../photon_source/inputs_true.dat | 4 ++ tests/regression_tests/photon_source/test.py | 20 ++++----- tools/ci/travis-before-script.sh | 2 +- 11 files changed, 95 insertions(+), 70 deletions(-) diff --git a/openmc/data/endf.py b/openmc/data/endf.py index aac40b3c9..bab3190ce 100644 --- a/openmc/data/endf.py +++ b/openmc/data/endf.py @@ -24,8 +24,9 @@ from openmc.stats.univariate import Uniform, Tabular, Legendre _LIBRARY = {0: 'ENDF/B', 1: 'ENDF/A', 2: 'JEFF', 3: 'EFF', 4: 'ENDF/B High Energy', 5: 'CENDL', 6: 'JENDL', - 31: 'INDL/V', 32: 'INDL/A', 33: 'FENDL', 34: 'IRDF', - 35: 'BROND', 36: 'INGDB-90', 37: 'FENDL/A', 41: 'BROND'} + 17: 'TENDL', 18: 'ROSFOND', 21: 'SG-21', 31: 'INDL/V', + 32: 'INDL/A', 33: 'FENDL', 34: 'IRDF', 35: 'BROND', + 36: 'INGDB-90', 37: 'FENDL/A', 41: 'BROND'} _SUBLIBRARY = { 0: 'Photo-nuclear data', diff --git a/openmc/data/photon.py b/openmc/data/photon.py index 4fdc60579..4ed6ef8c8 100644 --- a/openmc/data/photon.py +++ b/openmc/data/photon.py @@ -105,6 +105,7 @@ _STOPPING_POWERS = {} # for each element are in a 2D array with shape (n, k) stored on the key 'Z'. _BREMSSTRAHLUNG = {} + class AtomicRelaxation(EqualityMixin): """Atomic relaxation data. @@ -328,7 +329,7 @@ class AtomicRelaxation(EqualityMixin): class IncidentPhoton(EqualityMixin): - """Photon interaction data. + r"""Photon interaction data. This class stores photo-atomic, photo-nuclear, atomic relaxation, Compton profile, stopping power, and bremsstrahlung data assembled from @@ -350,9 +351,9 @@ class IncidentPhoton(EqualityMixin): Atomic relaxation data bremsstrahlung : dict Dictionary of bremsstrahlung DCS data with keys 'electron_energy' - (incident electron kinetic energy values in eV), 'photon_energy' + (incident electron kinetic energy values in [eV]), 'photon_energy' (ratio of the energy of the emitted photon to the incident electron - kinetic energy), and 'dcs' (cross sectin values in mb). The cross + kinetic energy), and 'dcs' (cross sectin values in [b]). The cross sections are in scaled form: :math:`(\beta^2/Z^2) E_k (d\sigma/dE_k)`, where :math:`E_k` is the energy of the emitted photon. compton_profiles : dict @@ -366,10 +367,10 @@ class IncidentPhoton(EqualityMixin): Contains the cross sections for each photon reaction. The keys are MT values and the values are instances of :class:`PhotonReaction`. stopping_powers : dict - Dictionary of stopping power data with keys 'energy' (in eV), 'I' (mean + Dictionary of stopping power data with keys 'energy' (in [eV]), 'I' (mean excitation energy), 's_collision' (collision stopping power in - eV cm:sup:`2`/g), and 's_radiative' (radiative stopping power in - eV cm:sup:`2`/g) + [eV cm\ :sup:`2`/g]), and 's_radiative' (radiative stopping power in + [eV cm\ :sup:`2`/g]) summed_reactions : collections.OrderedDict Contains summed cross sections. The keys are MT values and the values are instances of :class:`PhotonReaction`. @@ -587,13 +588,13 @@ class IncidentPhoton(EqualityMixin): _STOPPING_POWERS['energy'] = f['energy'].value*EV_PER_MEV for i in range(1, 99): group = f['{:03}'.format(i)] - _STOPPING_POWERS[i] = {'I': group.attrs['I'], - 's_collision': group['s_collision'].value, - 's_radiative': group['s_radiative'].value} # Units are in MeV cm^2/g; convert to eV cm^2/g - _STOPPING_POWERS[i]['s_collision'] *= EV_PER_MEV - _STOPPING_POWERS[i]['s_radiative'] *= EV_PER_MEV + _STOPPING_POWERS[i] = { + 'I': group.attrs['I'], + 's_collision': group['s_collision'].value*EV_PER_MEV, + 's_radiative': group['s_radiative'].value*EV_PER_MEV + } # Add stopping power data if Z < 99: @@ -616,8 +617,9 @@ class IncidentPhoton(EqualityMixin): # Index in data p = 39 - # Get log of incident electron kinetic energy values, used for cubic - # spline interpolation in log energy. Units are in MeV, so convert to eV. + # Get log of incident electron kinetic energy values, used for + # cubic spline interpolation in log energy. Units are in MeV, so + # convert to eV. logx = np.log(np.fromiter(brem[p:p+n], float, n)*EV_PER_MEV) p += n @@ -628,9 +630,10 @@ class IncidentPhoton(EqualityMixin): for i in range(1, 101): dcs = np.empty([len(log_energy), k]) - # Get the scaled cross section values for each electron energy and - # reduced photon energy for this Z - y = np.reshape(np.fromiter(brem[p:p+n*k], float, n*k), (n, k)) + # Get the scaled cross section values for each electron energy + # and reduced photon energy for this Z. Units are in mb, so + # convert to b. + y = np.reshape(np.fromiter(brem[p:p+n*k], float, n*k), (n, k))*1.0e-3 p += k*n for j in range(k): diff --git a/scripts/openmc-get-nndc-data b/scripts/openmc-get-nndc-data index 7a378b7fb..a84a9b282 100755 --- a/scripts/openmc-get-nndc-data +++ b/scripts/openmc-get-nndc-data @@ -54,10 +54,7 @@ for f in files: req = urlopen(url) # Get file size from header - if sys.version_info[0] < 3: - file_size = int(req.info().getheaders('Content-Length')[0]) - else: - file_size = req.length + file_size = req.length downloaded = 0 # Check if file already downloaded diff --git a/scripts/openmc-get-photon-data b/scripts/openmc-get-photon-data index 182fdf589..ae268a4cc 100755 --- a/scripts/openmc-get-photon-data +++ b/scripts/openmc-get-photon-data @@ -12,8 +12,7 @@ import shutil import zipfile import argparse from io import BytesIO - -import requests +from urllib.request import urlopen import openmc.data @@ -32,6 +31,7 @@ args = parser.parse_args() base_url = 'http://www.nndc.bnl.gov/endf/b7.1/zips/' files = ['ENDF-B-VII.1-photoat.zip', 'ENDF-B-VII.1-atomic_relax.zip'] +block_size = 16384 # ============================================================================== # DOWNLOAD FILES FROM NNDC SITE @@ -39,13 +39,44 @@ files = ['ENDF-B-VII.1-photoat.zip', 'ENDF-B-VII.1-atomic_relax.zip'] if not os.path.exists('photon_hdf5'): os.mkdir('photon_hdf5') - for f in files: # Establish connection to URL - print('Downloading {}...'.format(f)) url = base_url + f - r = requests.get(url, stream=True) - zipfile.ZipFile(BytesIO(r.content)).extractall() + req = urlopen(url) + + # Get file size from header + file_size = req.length + downloaded = 0 + + # Check if file already downloaded + if os.path.exists(f): + if os.path.getsize(f) == file_size: + print('Skipping ' + f) + continue + else: + overwrite = input('Overwrite {}? ([y]/n) '.format(f)) + if overwrite.lower().startswith('n'): + continue + + # Copy file to disk + print('Downloading {}... '.format(f), end='') + with open(f, 'wb') as fh: + while True: + chunk = req.read(block_size) + if not chunk: break + fh.write(chunk) + downloaded += len(chunk) + status = '{0:10} [{1:3.2f}%]'.format( + downloaded, downloaded * 100. / file_size) + print(status + chr(8)*len(status), end='') + print('') + +# ============================================================================== +# EXTRACT FILES + +for f in files: + print('Extracting {0}...'.format(f)) + zipfile.ZipFile(f).extractall() # ============================================================================== # GENERATE HDF5 DATA LIBRARY diff --git a/scripts/openmc-make-compton b/scripts/openmc-make-compton index 2591dce94..c4bd5b06b 100755 --- a/scripts/openmc-make-compton +++ b/scripts/openmc-make-compton @@ -20,10 +20,7 @@ block_size = 16384 req = urlopen(base_url + filename) # Get file size from header -if sys.version_info[0] < 3: - file_size = int(req.info().getheaders('Content-Length')[0]) -else: - file_size = req.length +file_size = req.length downloaded = 0 # Check if file already downloaded diff --git a/src/material_header.F90 b/src/material_header.F90 index 1664f75e7..9efe932d4 100644 --- a/src/material_header.F90 +++ b/src/material_header.F90 @@ -735,11 +735,7 @@ contains mat => materials(i_material) ! Determine whether we are generating electron or positron data - if (particle == POSITRON) then - positron_ = .true. - else - positron_ = .false. - end if + positron_ = (particle == POSITRON) ! Get the size of the energy grids n_k = size(ttb_k_grid) @@ -853,7 +849,7 @@ contains beta = sqrt(e*(e + TWO*MASS_ELECTRON)) / (e + MASS_ELECTRON) ! Compute the integrand of the PDF - f(j) = (1.0e-3_8 * x) / (beta**2 * stopping_power(j) * w) + f(j) = x / (beta**2 * stopping_power(j) * w) end do ! Number of points to integrate @@ -913,9 +909,6 @@ contains this % yield = -500.0_8 end where - deallocate(stopping_power_collision, stopping_power_radiative, & - stopping_power, dcs, f, z) - end subroutine bremsstrahlung_init end module material_header diff --git a/src/math_functions.cpp b/src/math_functions.cpp index d89cd1152..3342e3a91 100644 --- a/src/math_functions.cpp +++ b/src/math_functions.cpp @@ -691,7 +691,7 @@ void broaden_wmp_polynomials_c(double E, double dopp, int n, double factors[]) { } -void spline_c(int n, double x[], double y[], double z[]) +void spline_c(int n, const double x[], const double y[], double z[]) { double c_new[n-1]; @@ -718,7 +718,8 @@ void spline_c(int n, double x[], double y[], double z[]) } -double spline_interpolate_c(int n, double x[], double y[], double z[], double xint) +double spline_interpolate_c(int n, const double x[], const double y[], + const double z[], double xint) { // Find the lower bounding index in x of xint int i = n - 1; @@ -734,14 +735,12 @@ double spline_interpolate_c(int n, double x[], double y[], double z[], double xi double c = z[i]/2.0; double d = (z[i+1] - z[i])/(h*6.0); - double yint = y[i] + b*r + c*r*r + d*r*r*r; - - return yint; + return y[i] + b*r + c*r*r + d*r*r*r; } -double spline_integrate_c(int n, double x[], double y[], double z[], double xa, - double xb) +double spline_integrate_c(int n, const double x[], const double y[], + const double z[], double xa, double xb) { // Find the lower bounding index in x of the lower limit of integration. int ia = n - 1; diff --git a/src/math_functions.h b/src/math_functions.h index 2e92039c0..566e68ea6 100644 --- a/src/math_functions.h +++ b/src/math_functions.h @@ -157,15 +157,15 @@ extern "C" void broaden_wmp_polynomials_c(double E, double dopp, int n, //! used in any subsequent calls to spline_interpolate or spline_integrate for //! the same set of x and y values. //! -//! @param n Number of points -//! @param x Values of the independent variable, which must be strictly +//! @param n Number of points +//! @param x Values of the independent variable, which must be strictly //! increasing. -//! @param y Values of the dependent variable. -//! @param z The second derivative of the interpolating function at each value -//! of x. +//! @param y Values of the dependent variable. +//! @param[out] z The second derivative of the interpolating function at each +//! value of x. //============================================================================== -extern "C" void spline_c(int n, double x[], double y[], double z[]); +extern "C" void spline_c(int n, const double x[], const double y[], double z[]); //============================================================================== //! Determine the cubic spline interpolated y-value for a given x-value. @@ -180,8 +180,8 @@ extern "C" void spline_c(int n, double x[], double y[], double z[]); //! @result Interpolated value //============================================================================== -extern "C" double spline_interpolate_c(int n, double x[], double y[], double z[], - double xint); +extern "C" double spline_interpolate_c(int n, const double x[], const double y[], + const double z[], double xint); //============================================================================== //! Evaluate the definite integral of the interpolating cubic spline between @@ -198,8 +198,8 @@ extern "C" double spline_interpolate_c(int n, double x[], double y[], double z[] //! @result Integral //============================================================================== -extern "C" double spline_integrate_c(int n, double x[], double y[], double z[], - double xa, double xb); +extern "C" double spline_integrate_c(int n, const double x[], const double y[], + const double z[], double xa, double xb); } // namespace openmc -#endif // MATH_FUNCTIONS_H \ No newline at end of file +#endif // MATH_FUNCTIONS_H diff --git a/tests/regression_tests/photon_source/inputs_true.dat b/tests/regression_tests/photon_source/inputs_true.dat index 56ef025e7..425738a47 100644 --- a/tests/regression_tests/photon_source/inputs_true.dat +++ b/tests/regression_tests/photon_source/inputs_true.dat @@ -35,7 +35,11 @@ + + 2 + + 1 flux diff --git a/tests/regression_tests/photon_source/test.py b/tests/regression_tests/photon_source/test.py index f278e0c27..30961ab5a 100644 --- a/tests/regression_tests/photon_source/test.py +++ b/tests/regression_tests/photon_source/test.py @@ -19,9 +19,7 @@ class SourceTestHarness(PyAPITestHarness): inside_sphere = openmc.Cell() inside_sphere.region = -sphere inside_sphere.fill = mat - root = openmc.Universe() - root.add_cell(inside_sphere) - geometry = openmc.Geometry(root) + geometry = openmc.Geometry([inside_sphere]) geometry.export_to_xml() source = openmc.Source() @@ -40,20 +38,22 @@ class SourceTestHarness(PyAPITestHarness): settings.source = source settings.export_to_xml() + particle_filter = openmc.ParticleFilter('photon') tally = openmc.Tally() + tally.filters = [particle_filter] tally.scores = ['flux'] tallies = openmc.Tallies([tally]) tallies.export_to_xml() def _get_results(self): - sp = openmc.StatePoint(self._sp_name) - outstr = '' - t = sp.get_tally() - outstr += 'tally {}:\n'.format(t.id) - outstr += 'sum = {:12.6E}\n'.format(t.sum[0, 0, 0]) - outstr += 'sum_sq = {:12.6E}\n'.format(t.sum_sq[0, 0, 0]) + with openmc.StatePoint(self._sp_name) as sp: + outstr = '' + t = sp.get_tally() + outstr += 'tally {}:\n'.format(t.id) + outstr += 'sum = {:12.6E}\n'.format(t.sum[0, 0, 0]) + outstr += 'sum_sq = {:12.6E}\n'.format(t.sum_sq[0, 0, 0]) - return outstr + return outstr def test_source(): diff --git a/tools/ci/travis-before-script.sh b/tools/ci/travis-before-script.sh index bbb34358b..64611f794 100755 --- a/tools/ci/travis-before-script.sh +++ b/tools/ci/travis-before-script.sh @@ -13,7 +13,7 @@ fi # Download ENDF/B-VII.1 distribution ENDF=$HOME/endf-b-vii.1/ if [[ ! -d $ENDF/neutrons || ! -d $ENDF/photoat || ! -d $ENDF/atomic_relax ]]; then - wget https://anl.box.com/shared/static/4kd2gxnf4gtk4w1c8eua5fsua22kvgjb.xz -O - | tar -C $HOME -xvJ + wget https://anl.box.com/shared/static/yw7xe3k9gbps0e6muyf0cg134tq0punw.xz -O - | tar -C $HOME -xvJ fi # Download multipole library