crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 05:35:49 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Merge branch 'develop' into temp_interp_tol

Browse source
This commit is contained in:
josh 2022-10-30 03:10:57 +00:00
commit 0792d987d0
41 changed files with 835 additions and 479 deletions
Download patch file Download diff file Expand all files Collapse all files

20
.github/workflows/ci.yml vendored
View file

@ -25,7 +25,7 @@ jobs:
runs-on: ubuntu-20.04
strategy:
matrix:
python-version: [3.8]
python-version: ['3.10']
mpi: [n, y]
omp: [n, y]
dagmc: [n]
@ -37,24 +37,30 @@ jobs:
- python-version: 3.7
omp: n
mpi: n
- python-version: 3.8
omp: n
mpi: n
- python-version: 3.9
omp: n
mpi: n
- dagmc: y
python-version: 3.8
python-version: '3.10'
mpi: y
omp: y
- libmesh: y
python-version: 3.8
python-version: '3.10'
mpi: y
omp: y
- libmesh: y
python-version: 3.8
python-version: '3.10'
mpi: n
omp: y
- event: y
python-version: 3.8
python-version: '3.10'
omp: y
mpi: n
- vectfit: y
python-version: 3.8
python-version: '3.10'
omp: n
mpi: y
name: "Python ${{ matrix.python-version }} (omp=${{ matrix.omp }},
@ -75,7 +81,7 @@ jobs:
- uses: actions/checkout@v2
- name: Set up Python ${{ matrix.python-version }}
uses: actions/setup-python@v2
uses: actions/setup-python@v4
with:
python-version: ${{ matrix.python-version }}

17
CMakeLists.txt
View file

@ -3,8 +3,8 @@ project(openmc C CXX)
# Set version numbers
set(OPENMC_VERSION_MAJOR 0)
set(OPENMC_VERSION_MINOR 14)
set(OPENMC_VERSION_RELEASE 0)
set(OPENMC_VERSION_MINOR 13)
set(OPENMC_VERSION_RELEASE 3)
set(OPENMC_VERSION ${OPENMC_VERSION_MAJOR}.${OPENMC_VERSION_MINOR}.${OPENMC_VERSION_RELEASE})
configure_file(include/openmc/version.h.in "${CMAKE_BINARY_DIR}/include/openmc/version.h" @ONLY)
@ -482,6 +482,19 @@ if (OPENMC_USE_MPI)
target_link_libraries(libopenmc MPI::MPI_CXX)
endif()
#===============================================================================
# Log build info that this executable can report later
#===============================================================================
target_compile_definitions(libopenmc PRIVATE BUILD_TYPE=${CMAKE_BUILD_TYPE})
target_compile_definitions(libopenmc PRIVATE COMPILER_ID=${CMAKE_CXX_COMPILER_ID})
target_compile_definitions(libopenmc PRIVATE COMPILER_VERSION=${CMAKE_CXX_COMPILER_VERSION})
if (OPENMC_ENABLE_PROFILE)
target_compile_definitions(libopenmc PRIVATE PROFILINGBUILD)
endif()
if (OPENMC_ENABLE_COVERAGE)
target_compile_definitions(libopenmc PRIVATE COVERAGEBUILD)
endif()
#===============================================================================
# openmc executable
#===============================================================================

4
docs/source/conf.py
View file

@ -69,9 +69,9 @@ copyright = '2011-2022, Massachusetts Institute of Technology, UChicago Argonne
# built documents.
#
# The short X.Y version.
version = "0.14"
version = "0.13"
# The full version, including alpha/beta/rc tags.
release = "0.14.0"
release = "0.13.3"
# The language for content autogenerated by Sphinx. Refer to documentation
# for a list of supported languages.

3
docs/source/io_formats/statepoint.rst
View file

@ -109,7 +109,8 @@ The current version of the statepoint file format is 17.0.
- **y** (*double[]*) -- Interpolant values for energyfunction
interpolation. Only used for 'energyfunction' filters.
:Attributes: - **interpolation** (*int*) -- Interpolation type. Only used for
:Attributes:
- **interpolation** (*int*) -- Interpolation type. Only used for
'energyfunction' filters.
**/tallies/derivatives/derivative <id>/**

1
docs/source/pythonapi/model.rst
View file

@ -11,7 +11,6 @@ Convenience Functions
:template: myfunction.rst
openmc.model.borated_water
openmc.model.cylinder_from_points
openmc.model.hexagonal_prism
openmc.model.rectangular_prism
openmc.model.subdivide

6
docs/source/pythonapi/stats.rst
View file

@ -22,6 +22,12 @@ Univariate Probability Distributions
openmc.stats.Legendre
openmc.stats.Mixture
openmc.stats.Normal
.. autosummary::
:toctree: generated
:nosignatures:
:template: myfunction.rst
openmc.stats.muir
Angular Distributions

99
docs/source/releasenotes/0.13.2.rst Normal file
View file

@ -0,0 +1,99 @@
====================
What's New in 0.13.2
====================
.. currentmodule:: openmc
-------
Summary
-------
This release of OpenMC includes several bug fixes, performance improvements for
complex geometries and depletion simulations, and other general enhancements.
Notably, a capability has been added to compute the photon spectra from decay of
unstable nuclides. Alongside that, a new :data:`openmc.config` configuration
variable has been introduced that allows easier configuration of data sources.
Additionally, users can now perform cell or material rejection when sampling
external source distributions. Finally,
------------------------------------
Compatibility Notes and Deprecations
------------------------------------
- If you are building against libMesh for unstructured mesh tally support,
version 1.6 or higher is now required.
- The ``openmc.stats.Muir`` class has been replaced by a
:func:`openmc.stats.muir` function that returns an instance of
:class:`openmc.stats.Normal`.
------------
New Features
------------
- The :meth:`openmc.Material.get_nuclide_atom_densities` method now takes an
optional ``nuclide`` argument.
- Functions/methods in the :mod:`openmc.deplete` module now accept times in
Julian years (``'a'``).
- The :meth:`openmc.Universe.plot` method now allows a pre-existing axes object
to be passed in.
- Performance optimization for geometries with many complex regions.
- Performance optimization for depletion by avoiding deepcopies and caching
reaction rates.
- The :class:`openmc.RegularMesh` class now has a
:meth:`~openmc.RegularMesh.from_domain` classmethod.
- The :class:`openmc.CylindricalMesh` class now has a
:meth:`~openmc.CylindricalMesh.from_domain` classmethod.
- Improved method to condense diffusion coefficients from the :mod:`openmc.mgxs`
module.
- A new :data:`openmc.config` configuration variable has been introduced that
allows data sources to be specified at runtime or via environment variables.
- The :class:`openmc.EnergyFunctionFilter` class now supports multiple
interpolation schemes, not just linear-linear interpolation.
- The :class:`openmc.DAGMCUniverse` class now has ``material_names``,
``n_cells``, and ``n_surfaces`` attributes.
- A new :func:`openmc.data.decay_photon_energy` function has been added that
returns the energy spectrum of photons emitted from the decay of an unstable
nuclide.
- The :class:`openmc.Material` class also has a new
:attr:`~openmc.Material.decay_photon_energy` attribute that gives the decay
photon energy spectrum from the material based on its constituent nuclides.
- The :class:`openmc.deplete.StepResult` now has a
:meth:`~openmc.deplete.StepResult.get_material` method.
- The :class:`openmc.Source` class now takes a ``domains`` argument that
specifies a list of cells, materials, or universes that is used to reject
source sites (i.e., if the sampled sites are not within the specified domain,
they are rejected).
---------
Bug Fixes
---------
- `Delay call to Tally::set_strides <https://github.com/openmc-dev/openmc/pull/2183>`_
- `Fix reading reference direction from XML for angular distributions <https://github.com/openmc-dev/openmc/pull/2204>`_
- `Fix erroneous behavior in Material.add_components <https://github.com/openmc-dev/openmc/pull/2205>`_
- `Fix reading thermal elastic data from ACE <https://github.com/openmc-dev/openmc/pull/2226>`_
- `Fix reading source file with time attribute <https://github.com/openmc-dev/openmc/pull/2228>`_
- `Fix conversion of multiple thermal scattering data files from ACE <https://github.com/openmc-dev/openmc/pull/2232>`_
- `Fix reading values from wwinp file <https://github.com/openmc-dev/openmc/pull/2242>`_
- `Handle possibility of .ppm file in Universe.plot <https://github.com/openmc-dev/openmc/pull/2251>`_
- `Update volume calc types to mitigate overflow issues <https://github.com/openmc-dev/openmc/pull/2270>`_
------------
Contributors
------------
- `Lewis Gross <https://github.com/lewisgross1296>`_
- `Andrew Johnson <https://github.com/drewejohnson>`_
- `Miriam Kreher <https://github.com/mkreher13>`_
- `James Logan <https://github.com/jlogan03>`_
- `Jose Ignacio Marquez Damien <https://github.com/marquezj>`_
- `Josh May <https://github.com/joshmay1>`_
- `Patrick Myers <https://github.com/myerspat>`_
- `Adam Nelson <https://github.com/nelsonag>`_
- `April Novak <https://github.com/aprilnovak>`_
- `Ethan Peterson <https://github.com/eepeterson>`_
- `Gavin Ridley <https://github.com/gridley>`_
- `Paul Romano <https://github.com/paulromano>`_
- `Patrick Shriwise <https://github.com/pshriwise>`_
- `Jonathan Shimwell <https://github.com/Shimwell>`_
- `Olek Yardas <https://github.com/yardasol>`_

1
docs/source/releasenotes/index.rst
View file

@ -7,6 +7,7 @@ Release Notes
.. toctree::
:maxdepth: 1
0.13.2
0.13.1
0.13.0
0.12.2

3
include/openmc/capi.h
View file

@ -38,6 +38,9 @@ int openmc_energyfunc_filter_get_energy(
int openmc_energyfunc_filter_get_y(int32_t index, size_t* n, const double** y);
int openmc_energyfunc_filter_set_data(
int32_t index, size_t n, const double* energies, const double* y);
int openmc_energyfunc_filter_set_interpolation(
int32_t index, const char* interp);
int openmc_energyfunc_filter_get_interpolation(int32_t index, int* interp);
int openmc_extend_cells(int32_t n, int32_t* index_start, int32_t* index_end);
int openmc_extend_filters(int32_t n, int32_t* index_start, int32_t* index_end);
int openmc_extend_materials(

1
include/openmc/constants.h
View file

@ -5,6 +5,7 @@
#define OPENMC_CONSTANTS_H
#include <cmath>
#include <cstdint>
#include <limits>
#include "openmc/array.h"

3
include/openmc/output.h
View file

@ -40,6 +40,9 @@ void print_usage();
//! Display current version and copright/license information
void print_version();
//! Display compile flags employed, etc
void print_build_info();
//! Display header listing what physical values will displayed
void print_columns();

6
include/openmc/tallies/filter_energyfunc.h
View file

@ -40,8 +40,9 @@ public:
const vector<double>& energy() const { return energy_; }
const vector<double>& y() const { return y_; }
Interpolation interpolation_;
Interpolation interpolation() const { return interpolation_; }
void set_data(gsl::span<const double> energy, gsl::span<const double> y);
void set_interpolation(const std::string& interpolation);
private:
//----------------------------------------------------------------------------
@ -52,6 +53,9 @@ private:
//! Interpolant values.
vector<double> y_;
//! Interpolation scheme
Interpolation interpolation_ {Interpolation::lin_lin};
};
} // namespace openmc

7
include/openmc/volume_calc.h
View file

@ -1,6 +1,9 @@
#ifndef OPENMC_VOLUME_CALC_H
#define OPENMC_VOLUME_CALC_H
#include <cstdint>
#include <string>
#include "openmc/array.h"
#include "openmc/position.h"
#include "openmc/tallies/trigger.h"
@ -10,7 +13,6 @@
#include "xtensor/xtensor.hpp"
#include <gsl/gsl-lite.hpp>
#include <string>
namespace openmc {
@ -69,7 +71,8 @@ private:
//! \param[in] i_material Index in global materials vector
//! \param[in,out] indices Vector of material indices
//! \param[in,out] hits Number of hits corresponding to each material
void check_hit(int i_material, vector<int>& indices, vector<int>& hits) const;
void check_hit(
int i_material, vector<uint64_t>& indices, vector<uint64_t>& hits) const;
};
//==============================================================================

2
openmc/__init__.py
View file

@ -38,4 +38,4 @@ from .config import *
from openmc.model import rectangular_prism, hexagonal_prism, Model
__version__ = '0.14.0-dev'
__version__ = '0.13.3-dev'

6
openmc/cell.py
View file

@ -634,6 +634,9 @@ class Cell(IDManagerMixin):
if self.rotation is not None:
element.set("rotation", ' '.join(map(str, self.rotation.ravel())))
if self.volume is not None:
element.set("volume", str(self.volume))
return element
@classmethod
@ -687,6 +690,9 @@ class Cell(IDManagerMixin):
c.temperature = [float(t_i) for t_i in t.split()]
else:
c.temperature = float(t)
v = get_text(elem, 'volume')
if v is not None:
c.volume = float(v)
for key in ('temperature', 'rotation', 'translation'):
value = get_text(elem, key)
if value is not None:

6
openmc/data/decay.py
View file

@ -558,7 +558,9 @@ class Decay(EqualityMixin):
raise NotImplementedError("Multiple interpolation regions: {name}, {particle}")
interpolation = INTERPOLATION_SCHEME[f.interpolation[0]]
if interpolation not in ('histogram', 'linear-linear'):
raise NotImplementedError("Continuous spectra with {interpolation} interpolation ({name}, {particle}) not supported")
raise NotImplementedError(
f"Continuous spectra with {interpolation} interpolation "
f"({name}, {particle}) not supported")
intensity = spectra['continuous_normalization'].n
rates = decay_constant * intensity * f.y
@ -585,7 +587,7 @@ def decay_photon_energy(nuclide: str) -> Optional[Univariate]:
for the first time, you need to ensure that a depletion chain has been
specified in openmc.config['chain_file'].
.. versionadded:: 0.14.0
.. versionadded:: 0.13.2
Parameters
----------

1
openmc/deplete/chain.py
View file

@ -1025,6 +1025,7 @@ class Chain:
new_nuclide = Nuclide(previous.name)
new_nuclide.half_life = previous.half_life
new_nuclide.decay_energy = previous.decay_energy
new_nuclide.sources = previous.sources.copy()
if hasattr(previous, '_fpy'):
new_nuclide._fpy = previous._fpy

2
openmc/deplete/stepresult.py
View file

@ -202,7 +202,7 @@ class StepResult:
def get_material(self, mat_id):
"""Return material object for given depleted composition
.. versionadded:: 0.14.0
.. versionadded:: 0.13.2
Parameters
----------

2
openmc/filter.py
View file

@ -2095,7 +2095,7 @@ class EnergyFunctionFilter(Filter):
y = [float(x) for x in get_text(elem, 'y').split()]
out = cls(energy, y, filter_id=filter_id)
if elem.find('interpolation') is not None:
out.interpolation = elem.find('interpolation')
out.interpolation = elem.find('interpolation').text
return out
def can_merge(self, other):

26
openmc/lib/filter.py
View file

@ -7,6 +7,7 @@ import numpy as np
from numpy.ctypeslib import as_array
from openmc.exceptions import AllocationError, InvalidIDError
from openmc.data.function import INTERPOLATION_SCHEME
from . import _dll
from .core import _FortranObjectWithID
from .error import _error_handler
@ -16,9 +17,9 @@ from .mesh import _get_mesh
__all__ = [
'Filter', 'AzimuthalFilter', 'CellFilter', 'CellbornFilter', 'CellfromFilter',
'CellInstanceFilter', 'CollisionFilter', 'DistribcellFilter', 'DelayedGroupFilter',
'EnergyFilter', 'EnergyoutFilter', 'EnergyFunctionFilter', 'LegendreFilter',
'MaterialFilter', 'MeshFilter', 'MeshSurfaceFilter', 'MuFilter', 'ParticleFilter',
'CellInstanceFilter', 'CollisionFilter', 'DistribcellFilter', 'DelayedGroupFilter',
'EnergyFilter', 'EnergyoutFilter', 'EnergyFunctionFilter', 'LegendreFilter',
'MaterialFilter', 'MeshFilter', 'MeshSurfaceFilter', 'MuFilter', 'ParticleFilter',
'PolarFilter', 'SphericalHarmonicsFilter', 'SpatialLegendreFilter', 'SurfaceFilter',
'UniverseFilter', 'ZernikeFilter', 'ZernikeRadialFilter', 'filters'
]
@ -47,6 +48,12 @@ _dll.openmc_energyfunc_filter_get_y.resttpe = c_int
_dll.openmc_energyfunc_filter_get_y.errcheck = _error_handler
_dll.openmc_energyfunc_filter_get_y.argtypes = [
c_int32, POINTER(c_size_t), POINTER(POINTER(c_double))]
_dll.openmc_energyfunc_filter_get_interpolation.resttpe = c_int
_dll.openmc_energyfunc_filter_get_interpolation.errcheck = _error_handler
_dll.openmc_energyfunc_filter_get_interpolation.argtypes = [c_int32, POINTER(c_int)]
_dll.openmc_energyfunc_filter_set_interpolation.resttpe = c_int
_dll.openmc_energyfunc_filter_set_interpolation.errcheck = _error_handler
_dll.openmc_energyfunc_filter_set_interpolation.argtypes = [c_int32, c_char_p]
_dll.openmc_filter_get_id.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_filter_get_id.restype = c_int
_dll.openmc_filter_get_id.errcheck = _error_handler
@ -247,6 +254,8 @@ class EnergyFunctionFilter(Filter):
Independent variable for the interpolation
y : numpy.ndarray
Dependent variable for the interpolation
interpolation : {'histogram', 'linear-linear', 'linear-log', 'log-linear', 'log-log', 'quadratic', 'cubic'}
Interpolation scheme
"""
energy_array = np.asarray(energy)
y_array = np.asarray(y)
@ -264,6 +273,17 @@ class EnergyFunctionFilter(Filter):
def y(self):
return self._get_attr(_dll.openmc_energyfunc_filter_get_y)
@property
def interpolation(self) -> str:
interp = c_int()
_dll.openmc_energyfunc_filter_get_interpolation(self._index, interp)
return INTERPOLATION_SCHEME[interp.value]
@interpolation.setter
def interpolation(self, interp: str):
interp_ptr = c_char_p(interp.encode())
_dll.openmc_energyfunc_filter_set_interpolation(self._index, interp_ptr)
def _get_attr(self, cfunc):
array_p = POINTER(c_double)()
n = c_size_t()

6
openmc/material.py
View file

@ -98,7 +98,7 @@ class Material(IDManagerMixin):
this distribution is the total intensity of the photon source in
[decay/sec].
.. versionadded:: 0.14.0
.. versionadded:: 0.13.2
"""
@ -825,6 +825,8 @@ class Material(IDManagerMixin):
element : str
Specifies the element to match when searching through the nuclides
.. versionadded:: 0.13.2
Returns
-------
nuclides : list of str
@ -877,6 +879,8 @@ class Material(IDManagerMixin):
Nuclide for which atom density is desired. If not specified, the
atom density for each nuclide in the material is given.
.. versionadded:: 0.13.2
Returns
-------
nuclides : dict

4
openmc/plots.py
View file

@ -171,7 +171,9 @@ def _get_plot_image(plot, cwd):
stem = plot.filename if plot.filename is not None else f'plot_{plot.id}'
png_file = Path(cwd) / f'{stem}.png'
if not png_file.exists():
raise FileNotFoundError(f"Could not find .png image for plot {plot.id}")
raise FileNotFoundError(
f"Could not find .png image for plot {plot.id}. Your version of "
"OpenMC may not be built against libpng.")
return Image(str(png_file))

705
openmc/settings.py
View file

File diff suppressed because it is too large Load diff

2
openmc/stats/univariate.py
View file

@ -729,7 +729,7 @@ def muir(e0, m_rat, kt):
distribution: the mean energy of particles ``e0``, the mass of reactants
``m_rat``, and the ion temperature ``kt``.
.. versionadded:: 0.14.0
.. versionadded:: 0.13.2
Parameters
----------

14
openmc/universe.py
View file

@ -358,7 +358,10 @@ class Universe(UniverseBase):
model.plot_geometry(False, cwd=tmpdir, openmc_exec=openmc_exec)
# Read image from file
img = mpimg.imread(Path(tmpdir) / f'plot_{plot.id}.png')
img_path = Path(tmpdir) / f'plot_{plot.id}.png'
if not img_path.is_file():
img_path = img_path.with_suffix('.ppm')
img = mpimg.imread(img_path)
# Create a figure sized such that the size of the axes within
# exactly matches the number of pixels specified
@ -647,19 +650,26 @@ class DAGMCUniverse(UniverseBase):
bounding_box : 2-tuple of numpy.array
Lower-left and upper-right coordinates of an axis-aligned bounding box
of the universe.
.. versionadded:: 0.13.1
material_names : list of str
Return a sorted list of materials names that are contained within the
DAGMC h5m file. This is useful when naming openmc.Material() objects
as each material name present in the DAGMC h5m file must have a
matching openmc.Material() with the same name.
.. versionadded:: 0.13.2
n_cells : int
The number of cells in the DAGMC model. This is the number of cells at
runtime and accounts for the implicit complement whether or not is it
present in the DAGMC file.
.. versionadded:: 0.13.2
n_surfaces : int
The number of surfaces in the model.
.. versionadded:: 0.13.1
.. versionadded:: 0.13.2
"""
def __init__(self,

12
openmc/weight_windows.py
View file

@ -206,9 +206,9 @@ class WeightWindows(IDManagerMixin):
# reshape data according to mesh and energy bins
bounds = np.asarray(bounds)
if isinstance(self.mesh, UnstructuredMesh):
bounds.reshape(-1, self.num_energy_bins)
bounds = bounds.reshape(-1, self.num_energy_bins)
else:
bounds.reshape(*self.mesh.dimension, self.num_energy_bins)
bounds = bounds.reshape(*self.mesh.dimension, self.num_energy_bins)
self._lower_ww_bounds = bounds
@property
@ -225,9 +225,9 @@ class WeightWindows(IDManagerMixin):
# reshape data according to mesh and energy bins
bounds = np.asarray(bounds)
if isinstance(self.mesh, UnstructuredMesh):
bounds.reshape(-1, self.num_energy_bins)
bounds = bounds.reshape(-1, self.num_energy_bins)
else:
bounds.reshape(*self.mesh.dimension, self.num_energy_bins)
bounds = bounds.reshape(*self.mesh.dimension, self.num_energy_bins)
self._upper_ww_bounds = bounds
@property
@ -341,6 +341,10 @@ class WeightWindows(IDManagerMixin):
particle_type = get_text(elem, 'particle_type')
survival_ratio = float(get_text(elem, 'survival_ratio'))
ww_shape = (len(e_bounds) - 1,) + mesh.dimension[::-1]
lower_ww_bounds = np.array(lower_ww_bounds).reshape(ww_shape).T
upper_ww_bounds = np.array(upper_ww_bounds).reshape(ww_shape).T
max_lower_bound_ratio = None
if get_text(elem, 'max_lower_bound_ratio'):
max_lower_bound_ratio = float(get_text(elem, 'max_lower_bound_ratio'))

1
setup.py
View file

@ -57,6 +57,7 @@ kwargs = {
'Programming Language :: Python :: 3.7',
'Programming Language :: Python :: 3.8',
'Programming Language :: Python :: 3.9',
'Programming Language :: Python :: 3.10',
],
# Dependencies

2
src/distribution.cpp
View file

@ -179,7 +179,7 @@ Tabular::Tabular(pugi::xml_node node)
interp_ = Interpolation::lin_lin;
} else {
openmc::fatal_error(
"Unknown interpolation type for distribution: " + temp);
"Unsupported interpolation type for distribution: " + temp);
}
} else {
interp_ = Interpolation::histogram;

1
src/initialize.cpp
View file

@ -260,6 +260,7 @@ int parse_command_line(int argc, char* argv[])
} else if (arg == "-v" || arg == "--version") {
print_version();
print_build_info();
return OPENMC_E_UNASSIGNED;
} else if (arg == "-t" || arg == "--track") {

60
src/output.cpp
View file

@ -1,8 +1,8 @@
#include "openmc/output.h"
#include <algorithm> // for transform, max
#include <cstring> // for strlen
#include <cstdio> // for stdout
#include <cstring> // for strlen
#include <ctime> // for time, localtime
#include <fstream>
#include <iomanip> // for setw, setprecision, put_time
@ -134,9 +134,10 @@ void header(const char* msg, int level)
auto out = header(msg);
// Print header based on verbosity level.
if (settings::verbosity >= level)
if (settings::verbosity >= level) {
fmt::print("\n{}\n\n", out);
std::fflush(stdout);
}
}
//==============================================================================
@ -347,6 +348,61 @@ void print_version()
//==============================================================================
void print_build_info()
{
const std::string n("no");
const std::string y("yes");
std::string mpi(n);
std::string phdf5(n);
std::string dagmc(n);
std::string libmesh(n);
std::string png(n);
std::string profiling(n);
std::string coverage(n);
#ifdef PHDF5
phdf5 = y;
#endif
#ifdef OPENMC_MPI
mpi = y;
#endif
#ifdef DAGMC
dagmc = y;
#endif
#ifdef LIBMESH
libmesh = y;
#endif
#ifdef USE_LIBPNG
png = y;
#endif
#ifdef PROFILINGBUILD
profiling = y;
#endif
#ifdef COVERAGEBUILD
coverage = y;
#endif
// Wraps macro variables in quotes
#define STRINGIFY(x) STRINGIFY2(x)
#define STRINGIFY2(x) #x
if (mpi::master) {
fmt::print("Build type: {}\n", STRINGIFY(BUILD_TYPE));
fmt::print("Compiler ID: {} {}\n", STRINGIFY(COMPILER_ID),
STRINGIFY(COMPILER_VERSION));
fmt::print("MPI enabled: {}\n", mpi);
fmt::print("Parallel HDF5 enabled: {}\n", phdf5);
fmt::print("PNG support: {}\n", png);
fmt::print("DAGMC support: {}\n", dagmc);
fmt::print("libMesh support: {}\n", libmesh);
fmt::print("Coverage testing: {}\n", coverage);
fmt::print("Profiling flags: {}\n", profiling);
}
}
//==============================================================================
void print_columns()
{
if (settings::entropy_on) {

3
src/source.cpp
View file

@ -224,7 +224,8 @@ SourceSite IndependentSource::sample(uint64_t* seed) const
auto id = (domain_type_ == DomainType::CELL)
? model::cells[coord.cell]->id_
: model::universes[coord.universe]->id_;
if (found = contains(domain_ids_, id)) break;
if ((found = contains(domain_ids_, id)))
break;
}
}
}

114
src/tallies/filter_energyfunc.cpp
View file

@ -25,42 +25,14 @@ void EnergyFunctionFilter::from_xml(pugi::xml_node node)
fatal_error("y values not specified for EnergyFunction filter.");
auto y = get_node_array<double>(node, "y");
this->set_data(energy, y);
// default to linear-linear interpolation
interpolation_ = Interpolation::lin_lin;
if (check_for_node(node, "interpolation")) {
std::string interpolation = get_node_value(node, "interpolation");
if (interpolation == "histogram") {
interpolation_ = Interpolation::histogram;
} else if (interpolation == "linear-linear") {
interpolation_ = Interpolation::lin_lin;
} else if (interpolation == "linear-log") {
interpolation_ = Interpolation::lin_log;
} else if (interpolation == "log-linear") {
interpolation_ = Interpolation::log_lin;
} else if (interpolation == "log-log") {
interpolation_ = Interpolation::log_log;
} else if (interpolation == "quadratic") {
if (energy.size() < 3)
fatal_error(
fmt::format("Quadratic interpolation on EnergyFunctionFilter {} "
"requires at least 3 data points.",
id()));
interpolation_ = Interpolation::quadratic;
} else if (interpolation == "cubic") {
if (energy.size() < 4)
fatal_error(fmt::format("Cubic interpolation on EnergyFunctionFilter "
"{} requires at least 4 data points.",
id()));
interpolation_ = Interpolation::cubic;
} else {
fatal_error(fmt::format(
"Found invalid interpolation type '{}' on EnergyFunctionFilter {}.",
interpolation, id()));
}
this->set_interpolation(interpolation);
}
this->set_data(energy, y);
}
void EnergyFunctionFilter::set_data(
@ -86,6 +58,38 @@ void EnergyFunctionFilter::set_data(
}
}
void EnergyFunctionFilter::set_interpolation(const std::string& interpolation)
{
if (interpolation == "histogram") {
interpolation_ = Interpolation::histogram;
} else if (interpolation == "linear-linear") {
interpolation_ = Interpolation::lin_lin;
} else if (interpolation == "linear-log") {
interpolation_ = Interpolation::lin_log;
} else if (interpolation == "log-linear") {
interpolation_ = Interpolation::log_lin;
} else if (interpolation == "log-log") {
interpolation_ = Interpolation::log_log;
} else if (interpolation == "quadratic") {
if (energy_.size() < 3)
fatal_error(
fmt::format("Quadratic interpolation on EnergyFunctionFilter {} "
"requires at least 3 data points.",
this->id()));
interpolation_ = Interpolation::quadratic;
} else if (interpolation == "cubic") {
if (energy_.size() < 4)
fatal_error(fmt::format("Cubic interpolation on EnergyFunctionFilter "
"{} requires at least 4 data points.",
this->id()));
interpolation_ = Interpolation::cubic;
} else {
fatal_error(fmt::format(
"Found invalid interpolation type '{}' on EnergyFunctionFilter {}.",
interpolation, this->id()));
}
}
void EnergyFunctionFilter::get_all_bins(
const Particle& p, TallyEstimator estimator, FilterMatch& match) const
{
@ -105,7 +109,8 @@ void EnergyFunctionFilter::to_statepoint(hid_t filter_group) const
write_dataset(filter_group, "energy", energy_);
write_dataset(filter_group, "y", y_);
hid_t y_dataset = open_dataset(filter_group, "y");
write_attribute<int>(y_dataset, "interpolation", static_cast<int>(interpolation_));
write_attribute<int>(
y_dataset, "interpolation", static_cast<int>(interpolation_));
close_dataset(y_dataset);
}
@ -189,4 +194,51 @@ extern "C" int openmc_energyfunc_filter_get_y(
return 0;
}
extern "C" int openmc_energyfunc_filter_set_interpolation(
int32_t index, const char* interp)
{
// ensure this is a valid index to allocated filter
if (int err = verify_filter(index))
return err;
// get a pointer to the filter
const auto& filt_base = model::tally_filters[index].get();
// downcast to EnergyFunctionFilter
auto* filt = dynamic_cast<EnergyFunctionFilter*>(filt_base);
// check if a valid filter was produced
if (!filt) {
set_errmsg(
"Tried to set interpolation data for non-energy function filter.");
return OPENMC_E_INVALID_TYPE;
}
// Set interpolation
filt->set_interpolation(interp);
return 0;
}
extern "C" int openmc_energyfunc_filter_get_interpolation(
int32_t index, int* interp)
{
// ensure this is a valid index to allocated filter
if (int err = verify_filter(index))
return err;
// get a pointer to the filter
const auto& filt_base = model::tally_filters[index].get();
// downcast to EnergyFunctionFilter
auto* filt = dynamic_cast<EnergyFunctionFilter*>(filt_base);
// check if a valid filter was produced
if (!filt) {
set_errmsg(
"Tried to set interpolation data for non-energy function filter.");
return OPENMC_E_INVALID_TYPE;
}
*interp = static_cast<int>(filt->interpolation());
return 0;
}
} // namespace openmc

46
src/volume_calc.cpp
View file

@ -99,16 +99,16 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
{
// Shared data that is collected from all threads
int n = domain_ids_.size();
vector<vector<int>> master_indices(
vector<vector<uint64_t>> master_indices(
n); // List of material indices for each domain
vector<vector<int>> master_hits(
vector<vector<uint64_t>> master_hits(
n); // Number of hits for each material in each domain
int iterations = 0;
// Divide work over MPI processes
size_t min_samples = n_samples_ / mpi::n_procs;
size_t remainder = n_samples_ % mpi::n_procs;
size_t i_start, i_end;
uint64_t min_samples = n_samples_ / mpi::n_procs;
uint64_t remainder = n_samples_ % mpi::n_procs;
uint64_t i_start, i_end;
if (mpi::rank < remainder) {
i_start = (min_samples + 1) * mpi::rank;
i_end = i_start + min_samples + 1;
@ -123,14 +123,14 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
#pragma omp parallel
{
// Variables that are private to each thread
vector<vector<int>> indices(n);
vector<vector<int>> hits(n);
vector<vector<uint64_t>> indices(n);
vector<vector<uint64_t>> hits(n);
Particle p;
// Sample locations and count hits
#pragma omp for
for (size_t i = i_start; i < i_end; i++) {
int64_t id = iterations * n_samples_ + i;
uint64_t id = iterations * n_samples_ + i;
uint64_t seed = init_seed(id, STREAM_VOLUME);
p.n_coord() = 1;
@ -223,7 +223,14 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
// bump iteration counter and get total number
// of samples at this point
iterations++;
size_t total_samples = iterations * n_samples_;
uint64_t total_samples = iterations * n_samples_;
// warn user if total sample size is greater than what the uin64_t type can
// represent
if (total_samples == std::numeric_limits<uint64_t>::max()) {
warning("The number of samples has exceeded the type used to track hits. "
"Volume results may be inaccurate.");
}
// reset
double trigger_val = -INFTY;
@ -238,18 +245,19 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
// Create 2D array to store atoms/uncertainty for each nuclide. Later this
// is compressed into vectors storing only those nuclides that are
// non-zero
auto n_nuc = settings::run_CE ? data::nuclides.size()
: data::mg.nuclides_.size();
auto n_nuc =
settings::run_CE ? data::nuclides.size() : data::mg.nuclides_.size();
xt::xtensor<double, 2> atoms({n_nuc, 2}, 0.0);
#ifdef OPENMC_MPI
if (mpi::master) {
for (int j = 1; j < mpi::n_procs; j++) {
int q;
// retrieve results
MPI_Recv(
&q, 1, MPI_INTEGER, j, 2 * j, mpi::intracomm, MPI_STATUS_IGNORE);
vector<int> buffer(2 * q);
MPI_Recv(buffer.data(), 2 * q, MPI_INTEGER, j, 2 * j + 1,
&q, 1, MPI_UINT64_T, j, 2 * j, mpi::intracomm, MPI_STATUS_IGNORE);
vector<uint64_t> buffer(2 * q);
MPI_Recv(buffer.data(), 2 * q, MPI_UINT64_T, j, 2 * j + 1,
mpi::intracomm, MPI_STATUS_IGNORE);
for (int k = 0; k < q; ++k) {
bool already_added = false;
@ -268,20 +276,20 @@ vector<VolumeCalculation::Result> VolumeCalculation::execute() const
}
} else {
int q = master_indices[i_domain].size();
vector<int> buffer(2 * q);
vector<uint64_t> buffer(2 * q);
for (int k = 0; k < q; ++k) {
buffer[2 * k] = master_indices[i_domain][k];
buffer[2 * k + 1] = master_hits[i_domain][k];
}
MPI_Send(&q, 1, MPI_INTEGER, 0, 2 * mpi::rank, mpi::intracomm);
MPI_Send(buffer.data(), 2 * q, MPI_INTEGER, 0, 2 * mpi::rank + 1,
MPI_Send(&q, 1, MPI_UINT64_T, 0, 2 * mpi::rank, mpi::intracomm);
MPI_Send(buffer.data(), 2 * q, MPI_UINT64_T, 0, 2 * mpi::rank + 1,
mpi::intracomm);
}
#endif
if (mpi::master) {
int total_hits = 0;
size_t total_hits = 0;
for (int j = 0; j < master_indices[i_domain].size(); ++j) {
total_hits += master_hits[i_domain][j];
double f =
@ -464,7 +472,7 @@ void VolumeCalculation::to_hdf5(
}
void VolumeCalculation::check_hit(
int i_material, vector<int>& indices, vector<int>& hits) const
int i_material, vector<uint64_t>& indices, vector<uint64_t>& hits) const
{
// Check if this material was previously hit and if so, increment count

8
tests/regression_tests/weightwindows/inputs_true.dat
View file

File diff suppressed because one or more lines are too long

2
tests/regression_tests/weightwindows/results_true.dat
View file

@ -1 +1 @@
f10c722e27d1f0a69f700bc72c4b7751f375dc0a74e049c9abb4b261d2265155c0e516e7e38f9751b83a24eac07e944e51740d398b720524cf8cd6bf0b8c51fc
ebc761815175b25fc95a226174928c226a3ab5dbf3b2a2abf09e079a0b87dee1dee74aa9b9eaec35acd58b8c481d264be7e9b3f052905bcc14ccd76f36b01549

2
tests/unit_tests/test_cell.py
View file

@ -294,9 +294,11 @@ def test_to_xml_element(cell_with_lattice):
c = cells[0]
c.temperature = 900.0
c.volume = 1.0
elem = c.create_xml_subelement(root)
assert elem.get('region') == str(c.region)
assert elem.get('temperature') == str(c.temperature)
assert elem.get('volume') == str(c.volume)
@pytest.mark.parametrize("rotation", [

2
tests/unit_tests/test_deplete_chain.py
View file

@ -500,6 +500,7 @@ def test_reduce(gnd_simple_chain, endf_chain):
assert u5_round0.n_decay_modes == ref_U5.n_decay_modes
assert u5_round0.half_life == ref_U5.half_life
assert u5_round0.decay_energy == ref_U5.decay_energy
assert u5_round0.sources == ref_U5.sources
for newmode, refmode in zip(u5_round0.decay_modes, ref_U5.decay_modes):
assert newmode.target is None
assert newmode.type == refmode.type
@ -522,6 +523,7 @@ def test_reduce(gnd_simple_chain, endf_chain):
assert bareI5.n_decay_modes == ref_iodine.n_decay_modes
assert bareI5.half_life == ref_iodine.half_life
assert bareI5.decay_energy == ref_iodine.decay_energy
assert bareI5.sources == ref_iodine.sources
for newmode, refmode in zip(bareI5.decay_modes, ref_iodine.decay_modes):
assert newmode.target is None
assert newmode.type == refmode.type

15
tests/unit_tests/test_filters.py
View file

@ -254,3 +254,18 @@ def test_lethargy_bin_width():
energy_bins = openmc.mgxs.GROUP_STRUCTURES['VITAMIN-J-175']
assert f.lethargy_bin_width[0] == np.log10(energy_bins[1]/energy_bins[0])
assert f.lethargy_bin_width[-1] == np.log10(energy_bins[-1]/energy_bins[-2])
def test_energyfunc():
f = openmc.EnergyFunctionFilter(
[0.0, 10.0, 2.0e3, 1.0e6, 20.0e6],
[1.0, 0.9, 0.8, 0.7, 0.6],
'histogram'
)
# Make sure XML roundtrip works
elem = f.to_xml_element()
new_f = openmc.EnergyFunctionFilter.from_xml_element(elem)
np.testing.assert_allclose(f.energy, new_f.energy)
np.testing.assert_allclose(f.y, new_f.y)
assert f.interpolation == new_f.interpolation

5
tests/unit_tests/test_lib.py
View file

@ -283,6 +283,11 @@ def test_energy_function_filter(lib_init):
assert len(efunc.y) == 2
assert (efunc.y == [0.0, 2.0]).all()
# Default should be lin-lin
assert efunc.interpolation == 'linear-linear'
efunc.interpolation = 'histogram'
assert efunc.interpolation == 'histogram'
def test_tally(lib_init):
t = openmc.lib.tallies[1]

78
tests/unit_tests/test_stats.py
View file

@ -7,9 +7,12 @@ import openmc.stats
def assert_sample_mean(samples, expected_mean):
std_dev = samples.std() / np.sqrt(samples.size)
assert np.abs(expected_mean - samples.mean()) < 3*std_dev
# Calculate sample standard deviation
std_dev = samples.std() / np.sqrt(samples.size - 1)
# Means should agree within 4 sigma 99.993% of the time. Note that this is
# expected to fail about 1 out of 16,000 times
assert np.abs(expected_mean - samples.mean()) < 4*std_dev
def test_discrete():
@ -40,9 +43,9 @@ def test_discrete():
d3 = openmc.stats.Discrete(vals, probs)
# sample discrete distribution and check that the mean of the samples is
# within 3 std. dev. of the expected mean
# within 4 std. dev. of the expected mean
n_samples = 1_000_000
samples = d3.sample(n_samples, seed=100)
samples = d3.sample(n_samples)
assert_sample_mean(samples, exp_mean)
@ -86,11 +89,11 @@ def test_uniform():
np.testing.assert_array_equal(t.p, [1/(b-a), 1/(b-a)])
assert t.interpolation == 'histogram'
# Sample distribution and check that the mean of the samples is within 3
# Sample distribution and check that the mean of the samples is within 4
# std. dev. of the expected mean
exp_mean = 0.5 * (a + b)
n_samples = 1_000_000
samples = d.sample(n_samples, seed=100)
samples = d.sample(n_samples)
assert_sample_mean(samples, exp_mean)
@ -105,12 +108,13 @@ def test_powerlaw():
assert d.n == n
assert len(d) == 3
exp_mean = 100.0 * (n+1) / (n+2)
# Determine mean of distribution
exp_mean = (n+1)*(b**(n+2) - a**(n+2))/((n+2)*(b**(n+1) - a**(n+1)))
# sample power law distribution and check that the mean of the samples is
# within 3 std. dev. of the expected mean
# within 4 std. dev. of the expected mean
n_samples = 1_000_000
samples = d.sample(n_samples, seed=100)
samples = d.sample(n_samples)
assert_sample_mean(samples, exp_mean)
@ -126,13 +130,13 @@ def test_maxwell():
exp_mean = 3/2 * theta
# sample maxwell distribution and check that the mean of the samples is
# within 3 std. dev. of the expected mean
# within 4 std. dev. of the expected mean
n_samples = 1_000_000
samples = d.sample(n_samples, seed=100)
samples = d.sample(n_samples)
assert_sample_mean(samples, exp_mean)
# A second sample with a different seed
samples_2 = d.sample(n_samples, seed=200)
# A second sample starting from a different seed
samples_2 = d.sample(n_samples)
assert_sample_mean(samples_2, exp_mean)
assert samples_2.mean() != samples.mean()
@ -154,9 +158,9 @@ def test_watt():
exp_mean = 3/2 * a + a**2 * b / 4
# sample Watt distribution and check that the mean of the samples is within
# 3 std. dev. of the expected mean
# 4 std. dev. of the expected mean
n_samples = 1_000_000
samples = d.sample(n_samples, seed=100)
samples = d.sample(n_samples)
assert_sample_mean(samples, exp_mean)
@ -176,28 +180,16 @@ def test_tabular():
d = openmc.stats.Tabular(x, p)
n_samples = 100_000
samples = d.sample(n_samples, seed=100)
diff = np.abs(samples - d.mean())
# within_1_sigma = np.count_nonzero(diff < samples.std())
# assert within_1_sigma / n_samples >= 0.68
within_2_sigma = np.count_nonzero(diff < 2*samples.std())
assert within_2_sigma / n_samples >= 0.95
within_3_sigma = np.count_nonzero(diff < 3*samples.std())
assert within_3_sigma / n_samples >= 0.99
samples = d.sample(n_samples)
assert_sample_mean(samples, d.mean())
# test histogram sampling
d = openmc.stats.Tabular(x, p, interpolation='histogram')
d.normalize()
assert d.integral() == pytest.approx(1.0)
samples = d.sample(n_samples, seed=100)
diff = np.abs(samples - d.mean())
# within_1_sigma = np.count_nonzero(diff < samples.std())
# assert within_1_sigma / n_samples >= 0.68
within_2_sigma = np.count_nonzero(diff < 2*samples.std())
assert within_2_sigma / n_samples >= 0.95
within_3_sigma = np.count_nonzero(diff < 3*samples.std())
assert within_3_sigma / n_samples >= 0.99
samples = d.sample(n_samples)
assert_sample_mean(samples, d.mean())
def test_legendre():
@ -361,15 +353,9 @@ def test_normal():
assert len(d) == 2
# sample normal distribution
n_samples = 10000
samples = d.sample(n_samples, seed=100)
samples = np.abs(samples - mean)
within_1_sigma = np.count_nonzero(samples < std_dev)
assert within_1_sigma / n_samples >= 0.68
within_2_sigma = np.count_nonzero(samples < 2*std_dev)
assert within_2_sigma / n_samples >= 0.95
within_3_sigma = np.count_nonzero(samples < 3*std_dev)
assert within_3_sigma / n_samples >= 0.99
n_samples = 100_000
samples = d.sample(n_samples)
assert_sample_mean(samples, mean)
def test_muir():
@ -386,15 +372,9 @@ def test_muir():
assert isinstance(d, openmc.stats.Normal)
# sample muir distribution
n_samples = 10000
samples = d.sample(n_samples, seed=100)
samples = np.abs(samples - mean)
within_1_sigma = np.count_nonzero(samples < d.std_dev)
assert within_1_sigma / n_samples >= 0.68
within_2_sigma = np.count_nonzero(samples < 2*d.std_dev)
assert within_2_sigma / n_samples >= 0.95
within_3_sigma = np.count_nonzero(samples < 3*d.std_dev)
assert within_3_sigma / n_samples >= 0.99
n_samples = 100_000
samples = d.sample(n_samples)
assert_sample_mean(samples, mean)
def test_combine_distributions():

16
tests/unit_tests/weightwindows/test.py
View file

@ -195,3 +195,19 @@ def test_weightwindows(model):
compare_results('neutron', analog_tally, ww_tally)
compare_results('photon', analog_tally, ww_tally)
def test_lower_ww_bounds_shape():
"""checks that lower_ww_bounds is reshaped to the mesh dimension when set"""
ww_mesh = openmc.RegularMesh()
ww_mesh.lower_left = (-10, -10, -10)
ww_mesh.upper_right = (10, 10, 10)
ww_mesh.dimension = (2, 3, 4)
ww = openmc.WeightWindows(
mesh=ww_mesh,
lower_ww_bounds=[1]*24,
upper_bound_ratio=5,
energy_bounds=(1, 1e40)
)
assert ww.lower_ww_bounds.shape == (2, 3, 4, 1)