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

Merge branch 'develop' into centre_for_cylinder_spherical_meshes

Browse source
This commit is contained in:
RemDelaporteMathurin 2023-02-02 08:47:28 -05:00
commit 90995dc512
210 changed files with 9470 additions and 5064 deletions
Download patch file Download diff file Expand all files Collapse all files

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

@ -25,10 +25,11 @@ jobs:
runs-on: ubuntu-20.04
strategy:
matrix:
python-version: [3.8]
python-version: ['3.10']
mpi: [n, y]
omp: [n, y]
dagmc: [n]
ncrystal: [n]
libmesh: [n]
event: [n]
vectfit: [n]
@ -37,28 +38,38 @@ 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
- ncrystal: y
python-version: '3.10'
mpi: n
omp: n
- libmesh: y
python-version: '3.10'
mpi: y
omp: y
- libmesh: y
python-version: 3.8
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 }},
mpi=${{ matrix.mpi }}, dagmc=${{ matrix.dagmc }},
mpi=${{ matrix.mpi }}, dagmc=${{ matrix.dagmc }}, ncrystal=${{ matrix.ncrystal }},
libmesh=${{ matrix.libmesh }}, event=${{ matrix.event }}
vectfit=${{ matrix.vectfit }})"
@ -67,15 +78,16 @@ jobs:
PHDF5: ${{ matrix.mpi }}
OMP: ${{ matrix.omp }}
DAGMC: ${{ matrix.dagmc }}
NCRYSTAL: ${{ matrix.ncrystal }}
EVENT: ${{ matrix.event }}
VECTFIT: ${{ matrix.vectfit }}
LIBMESH: ${{ matrix.libmesh }}
steps:
- uses: actions/checkout@v2
- uses: actions/checkout@v3
- name: Set up Python ${{ matrix.python-version }}
uses: actions/setup-python@v2
uses: actions/setup-python@v4
with:
python-version: ${{ matrix.python-version }}
@ -112,7 +124,9 @@ jobs:
- name: test
shell: bash
run: $GITHUB_WORKSPACE/tools/ci/gha-script.sh
run: |
CTEST_OUTPUT_ON_FAILURE=1 make test -C $GITHUB_WORKSPACE/build/
$GITHUB_WORKSPACE/tools/ci/gha-script.sh
- name: after_success
shell: bash

2
.github/workflows/dockerhub-publish-release-dagmc-libmesh.yml vendored
View file

@ -8,7 +8,7 @@ jobs:
main:
runs-on: ubuntu-latest
steps:
- uses: actions/checkout@v2
- uses: actions/checkout@v3
- name: Set env
run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
-

2
.github/workflows/dockerhub-publish-release-dagmc.yml vendored
View file

@ -8,7 +8,7 @@ jobs:
main:
runs-on: ubuntu-latest
steps:
- uses: actions/checkout@v2
- uses: actions/checkout@v3
- name: Set env
run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
-

2
.github/workflows/dockerhub-publish-release-libmesh.yml vendored
View file

@ -8,7 +8,7 @@ jobs:
main:
runs-on: ubuntu-latest
steps:
- uses: actions/checkout@v2
- uses: actions/checkout@v3
- name: Set env
run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
-

2
.github/workflows/dockerhub-publish-release.yml vendored
View file

@ -8,7 +8,7 @@ jobs:
main:
runs-on: ubuntu-latest
steps:
- uses: actions/checkout@v2
- uses: actions/checkout@v3
- name: Set env
run: echo "RELEASE_VERSION=${GITHUB_REF#refs/*/}" >> $GITHUB_ENV
-

3
.gitmodules vendored
View file

@ -13,3 +13,6 @@
[submodule "vendor/fmt"]
path = vendor/fmt
url = https://github.com/fmtlib/fmt.git
[submodule "vendor/Catch2"]
path = vendor/Catch2
url = https://github.com/catchorg/Catch2.git

67
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)
@ -31,11 +31,14 @@ endif()
#===============================================================================
option(OPENMC_USE_OPENMP "Enable shared-memory parallelism with OpenMP" ON)
option(OPENMC_BUILD_TESTS "Build tests" ON)
option(OPENMC_ENABLE_PROFILE "Compile with profiling flags" OFF)
option(OPENMC_ENABLE_COVERAGE "Compile with coverage analysis flags" OFF)
option(OPENMC_USE_DAGMC "Enable support for DAGMC (CAD) geometry" OFF)
option(OPENMC_USE_LIBMESH "Enable support for libMesh unstructured mesh tallies" OFF)
option(OPENMC_USE_MPI "Enable MPI" OFF)
option(OPENMC_USE_MCPL "Enable MCPL" OFF)
option(OPENMC_USE_NCRYSTAL "Enable support for NCrystal scattering" OFF)
# Warnings for deprecated options
foreach(OLD_OPT IN ITEMS "openmp" "profile" "coverage" "dagmc" "libmesh")
@ -98,6 +101,15 @@ macro(find_package_write_status pkg)
endif()
endmacro()
#===============================================================================
# NCrystal Scattering Support
#===============================================================================
if(OPENMC_USE_NCRYSTAL)
find_package(NCrystal REQUIRED)
message(STATUS "Found NCrystal: ${NCrystal_DIR} (version ${NCrystal_VERSION})")
endif()
#===============================================================================
# DAGMC Geometry Support - need DAGMC/MOAB
#===============================================================================
@ -158,6 +170,15 @@ if(${HDF5_VERSION} VERSION_GREATER_EQUAL 1.12.0)
list(APPEND cxxflags -DH5Oget_info_by_idx_vers=1 -DH5O_info_t_vers=1)
endif()
#===============================================================================
# MCPL
#===============================================================================
if (OPENMC_USE_MCPL)
find_package(MCPL REQUIRED)
message(STATUS "Found MCPL: ${MCPL_DIR} (found version \"${MCPL_VERSION}\")")
endif()
#===============================================================================
# Set compile/link flags based on which compiler is being used
#===============================================================================
@ -265,6 +286,17 @@ if (NOT gsl-lite_FOUND)
target_compile_definitions(gsl-lite-v1 INTERFACE gsl_CONFIG_ALLOWS_NONSTRICT_SPAN_COMPARISON=1)
endif()
#===============================================================================
# Catch2 library
#===============================================================================
if(OPENMC_BUILD_TESTS)
find_package_write_status(Catch2)
if (NOT Catch2)
add_subdirectory(vendor/Catch2)
endif()
endif()
#===============================================================================
# RPATH information
#===============================================================================
@ -323,10 +355,12 @@ list(APPEND libopenmc_SOURCES
src/lattice.cpp
src/material.cpp
src/math_functions.cpp
src/mcpl_interface.cpp
src/mesh.cpp
src/message_passing.cpp
src/mgxs.cpp
src/mgxs_interface.cpp
src/ncrystal_interface.cpp
src/nuclide.cpp
src/output.cpp
src/particle.cpp
@ -482,6 +516,35 @@ if (OPENMC_USE_MPI)
target_link_libraries(libopenmc MPI::MPI_CXX)
endif()
if (OPENMC_BUILD_TESTS)
# Add cpp tests directory
include(CTest)
add_subdirectory(tests/cpp_unit_tests)
endif()
if (OPENMC_USE_MCPL)
target_compile_definitions(libopenmc PUBLIC OPENMC_MCPL)
target_link_libraries(libopenmc MCPL::mcpl)
endif()
if(OPENMC_USE_NCRYSTAL)
target_compile_definitions(libopenmc PRIVATE NCRYSTAL)
target_link_libraries(libopenmc NCrystal::NCrystal)
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
#===============================================================================

2
LICENSE
View file

@ -1,4 +1,4 @@
Copyright (c) 2011-2022 Massachusetts Institute of Technology, UChicago Argonne
Copyright (c) 2011-2023 Massachusetts Institute of Technology, UChicago Argonne
LLC, and OpenMC contributors
Permission is hereby granted, free of charge, to any person obtaining a copy of

9
cmake/OpenMCConfig.cmake.in
View file

@ -9,6 +9,11 @@ if(@OPENMC_USE_DAGMC@)
find_package(DAGMC REQUIRED HINTS @DAGMC_DIR@)
endif()
if(@OPENMC_USE_NCRYSTAL@)
find_package(NCrystal REQUIRED)
message(STATUS "Found NCrystal: ${NCrystal_DIR} (version ${NCrystal_VERSION})")
endif()
if(@OPENMC_USE_LIBMESH@)
include(FindPkgConfig)
list(APPEND CMAKE_PREFIX_PATH @LIBMESH_PREFIX@)
@ -25,3 +30,7 @@ endif()
if(@OPENMC_USE_MPI@)
find_package(MPI REQUIRED)
endif()
if(@OPENMC_USE_MCPL@)
find_package(MCPL REQUIRED)
endif()

10
docs/source/conf.py
View file

@ -62,16 +62,16 @@ master_doc = 'index'
# General information about the project.
project = 'OpenMC'
copyright = '2011-2022, Massachusetts Institute of Technology, UChicago Argonne LLC, and OpenMC contributors'
copyright = '2011-2023, Massachusetts Institute of Technology, UChicago Argonne LLC, and OpenMC contributors'
# The version info for the project you're documenting, acts as replacement for
# |version| and |release|, also used in various other places throughout the
# 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.
@ -247,7 +247,7 @@ napoleon_use_ivar = True
intersphinx_mapping = {
'python': ('https://docs.python.org/3', None),
'numpy': ('https://numpy.org/doc/stable/', None),
'scipy': ('https://docs.scipy.org/doc/scipy/reference', None),
'scipy': ('https://docs.scipy.org/doc/scipy/', None),
'pandas': ('https://pandas.pydata.org/pandas-docs/stable/', None),
'matplotlib': ('https://matplotlib.org/', None)
'matplotlib': ('https://matplotlib.org/stable/', None)
}

35
docs/source/devguide/tests.rst
View file

@ -31,9 +31,13 @@ Prerequisites
- The test suite requires a specific set of cross section data in order for
tests to pass. A download URL for the data that OpenMC expects can be found
within ``tools/ci/download-xs.sh``.
within ``tools/ci/download-xs.sh``. Once the tarball is downloaded and
unpacked, set the :envvar:`OPENMC_CROSS_SECTIONS` environment variable to the
path of the ``cross_sections.xml`` file within the unpacked data.
- In addition to the HDF5 data, some tests rely on ENDF files. A download URL
for those can also be found in ``tools/ci/download-xs.sh``.
for those can also be found in ``tools/ci/download-xs.sh``. Once the tarball
is downloaded and unpacked, set the :envvar:`OPENMC_ENDF_DATA` environment
variable to the top-level directory of the unpacked tarball.
- Some tests require `NJOY <https://www.njoy21.io/NJOY2016>`_ to preprocess
cross section data. The test suite assumes that you have an ``njoy``
executable available on your :envvar:`PATH`.
@ -41,7 +45,7 @@ Prerequisites
Running Tests
-------------
To execute the test suite, go to the ``tests/`` directory and run::
To execute the Python test suite, go to the ``tests/`` directory and run::
pytest
@ -51,6 +55,14 @@ installed and run::
pytest --cov=../openmc --cov-report=html
To execute the C++ test suite, go to your build directory and run::
ctest
If you want to view testing output on failure run::
ctest --output-on-failure
Generating XML Inputs
---------------------
@ -62,6 +74,23 @@ run::
pytest --build-inputs <name-of-test>
Adding C++ Unit Tests
---------------------
The C++ test suite uses Catch2 integrated with CTest. Each header file should
have a corresponding test file in ``tests/cpp_unit_tests/``. If the test file
does not exist run::
touch test_<name-of-header-file>.cpp
The file must be added to the CMake build system in
``tests/cpp_unit_tests/CMakeLists.txt``. ``test_<name-of-header-file>`` should
be added to ``TEST_NAMES``.
To add a test case to ``test_<name-of-header-file>.cpp`` ensure
``catch2/catch_test_macros.hpp`` is included. A unit test can then be added
using the ``TEST_CASE`` macro and the ``REQUIRE`` assertion from Catch2.
Adding Tests to the Regression Suite
------------------------------------

40
docs/source/io_formats/settings.rst
View file

@ -32,6 +32,17 @@ standard deviation.
*Default*: false
-------------------------------------
``<create_delayed_neutrons>`` Element
-------------------------------------
The ``<create_delayed_neutrons>`` element indicates whether delayed neutrons
are created in fission. If this element is set to "true", delayed neutrons
will be created in fission events; otherwise only prompt neutrons will be
created.
*Default*: true
-------------------------------------
``<create_fission_neutrons>`` Element
-------------------------------------
@ -732,6 +743,14 @@ attributes/sub-elements:
*Default*: false
:mcpl:
If this element is set to "true", the source point file containing the
source bank will be written as an MCPL_ file name ``source.mcpl`` instead of
an HDF5 file. This option is only applicable if the ``<separate>`` element
is set to true.
*Default*: false
------------------------------
``<surf_source_read>`` Element
------------------------------
@ -767,6 +786,16 @@ certain surfaces and write out the source bank in a separate file called
*Default*: None
:mcpl:
An optional boolean which indicates if the banked particles should be
written to a file in the MCPL_-format instead of the native HDF5-based
format. If activated the output file name is changed to
``surface_source.mcpl``.
*Default*: false
.. _MCPL: https://mctools.github.io/mcpl/mcpl.pdf
------------------------------
``<survival_biasing>`` Element
------------------------------
@ -832,7 +861,9 @@ cell, the nearest temperature at which cross sections are given is to be
applied, within a given tolerance (see :ref:`temperature_tolerance`). A value of
"interpolation" indicates that cross sections are to be linear-linear
interpolated between temperatures at which nuclear data are present (see
:ref:`temperature_treatment`).
:ref:`temperature_treatment`). With the "interpolation" method, temperatures
outside of the bounds of the nuclear data may be accepted, provided they still
fall within the tolerance (see :ref:`temperature_tolerance`).
*Default*: "nearest"
@ -871,7 +902,12 @@ The ``<temperature_tolerance>`` element specifies a tolerance in Kelvin that is
to be applied when the "nearest" temperature method is used. For example, if a
cell temperature is 340 K and the tolerance is 15 K, then the closest
temperature in the range of 325 K to 355 K will be used to evaluate cross
sections.
sections. If the ``<temperature_method>`` is "interpolation", the tolerance
specified applies to cell temperatures outside of the data bounds. For example,
if a cell is specified at 695K, a tolerance of 15K and data is only available
at 700K and 1000K, the cell's cross sections will be evaluated at 700K, since
the desired temperature of 695K is within the tolerance of the actual data
despite not being bounded on both sides.
*Default*: 10 K

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>/**

2
docs/source/license.rst
View file

@ -4,7 +4,7 @@
License Agreement
=================
Copyright © 2011-2022 Massachusetts Institute of Technology, UChicago Argonne
Copyright © 2011-2023 Massachusetts Institute of Technology, UChicago Argonne
LLC, and OpenMC contributors
Permission is hereby granted, free of charge, to any person obtaining a copy of

28
docs/source/methods/cross_sections.rst
View file

@ -178,6 +178,27 @@ been selected. There are three methods available:
section data is loaded for a single temperature and is used in the
unresolved resonance and fast energy ranges.
------------------
NCrystal materials
------------------
As an alternative of the standard thermal scattering treatment using
:math:`S(\alpha,\beta)` tables, OpenMC allows to create materials using
NCrystal_. In addition to the regular thermal elastic, and thermal inelastic
processes, NCrystal allows the generation of models for materials that cannot
currently included in ACE files such as oriented single crystals (see the
`NCrystal paper`_), and further extend the physics `using plugins`_. Thermal
scattering kernels are generated on the fly from dynamic and structural data, or
loaded from :math:`S(\alpha,\beta)` tables converted from ENDF6 evaluations.
These kernels are sampled in a direct way using a fast `rejection algorithm`_
that does not require previous processing. A `large library`_ of materials is
already included in the NCrystal distribution, and new materials can be easily
defined from scratch in the `NCMAT format`_ or `combining existing files`_.
The compositions of the materials defined in NCrystal are passed on to OpenMC
all other reactions except for thermal neutron scattering are handled by
continuous energy ACE libraries.
----------------
Multi-Group Data
----------------
@ -279,3 +300,10 @@ or even isotropic scattering.
.. _ENDF/B data: https://www.nndc.bnl.gov/endf-b8.0/
.. _Leppanen: https://doi.org/10.1016/j.anucene.2009.03.019
.. _algorithms: http://ab-initio.mit.edu/wiki/index.php/Faddeeva_Package
.. _NCrystal: https://github.com/mctools/ncrystal
.. _NCrystal paper: https://doi.org/10.1016/j.cpc.2019.07.015
.. _using plugins: https://doi.org/10.1016/j.cpc.2021.108082
.. _rejection algorithm: https://doi.org/10.1016/j.jcp.2018.11.043
.. _large library: https://github.com/mctools/ncrystal/wiki/Data-library
.. _NCMAT format: https://github.com/mctools/ncrystal/wiki/NCMAT-format
.. _combining existing files: https://github.com/mctools/ncrystal/wiki/Announcement-Release3.0.0#2-multiphase-materials

2
docs/source/pythonapi/base.rst
View file

@ -120,7 +120,7 @@ Constructing Tallies
openmc.MaterialFilter
openmc.CellFilter
openmc.CellFromFilter
openmc.CellbornFilter
openmc.CellBornFilter
openmc.CellInstanceFilter
openmc.CollisionFilter
openmc.SurfaceFilter

3
docs/source/pythonapi/data.rst
View file

@ -63,9 +63,10 @@ Core Functions
atomic_weight
combine_distributions
decay_constant
decay_energy
decay_photon_energy
dose_coefficients
gnd_name
gnds_name
half_life
isotopes
kalbach_slope

2
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
@ -32,6 +31,7 @@ Composite Surfaces
openmc.model.XConeOneSided
openmc.model.YConeOneSided
openmc.model.ZConeOneSided
openmc.model.Polygon
TRISO Fuel Modeling
-------------------

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

29
docs/source/usersguide/install.rst
View file

@ -271,6 +271,27 @@ Prerequisites
cmake -DOPENMC_USE_DAGMC=on -DCMAKE_PREFIX_PATH=/path/to/dagmc/installation ..
* MCPL_ library for reading and writing .mcpl files
This option allows OpenMC to read and write MCPL (Monte Carlo Particle
Lists) files instead of .h5 files for sources (external source
distribution, k-eigenvalue source distribution, and surface sources). To
turn this option on in the CMake configuration step, add the following
option:
cmake -DOPENMC_USE_MCPL=on ..
* NCrystal_ library for defining materials with enhanced thermal neutron transport
Adding this option allows the creation of materials from NCrystal, which
replaces the scattering kernel treatment of ACE files with a modular,
on-the-fly approach. To use it `install
<https://github.com/mctools/ncrystal/wiki/Get-NCrystal>`_ and `initialize
<https://github.com/mctools/ncrystal/wiki/Using-NCrystal#setting-up>`_
NCrystal and turn on the option in the CMake configuration step:
cmake -DOPENMC_USE_NCRYSTAL=on ..
* libMesh_ mesh library framework for numerical simulations of partial differential equations
This optional dependency enables support for unstructured mesh tally
@ -294,6 +315,8 @@ Prerequisites
.. _MOAB: https://bitbucket.org/fathomteam/moab
.. _libMesh: https://libmesh.github.io/
.. _libpng: http://www.libpng.org/pub/png/libpng.html
.. _MCPL: https://github.com/mctools/mcpl
.. _NCrystal: https://github.com/mctools/ncrystal
Obtaining the Source
--------------------
@ -362,6 +385,12 @@ OPENMC_USE_DAGMC
should also be defined as `DAGMC_ROOT` in the CMake configuration command.
(Default: off)
OPENMC_USE_NCRYSTAL
Turns on support for NCrystal materials. NCrystal must be
`installed <https://github.com/mctools/ncrystal/wiki/Get-NCrystal>`_ and
`initialized <https://github.com/mctools/ncrystal/wiki/Using-NCrystal#setting-up>`_.
(Default: off)
OPENMC_USE_LIBMESH
Enables the use of unstructured mesh tallies with libMesh_. (Default: off)

37
docs/source/usersguide/materials.rst
View file

@ -101,6 +101,42 @@ you would need to add hydrogen and oxygen to a material and then assign the
.. _usersguide_naming:
-------------------------
Adding NCrystal materials
-------------------------
Additional support for thermal scattering can be added by using NCrystal_. The
:meth:`Material.from_ncrystal` class method generates a :class:`openmc.Material`
object from an `NCrystal configuration string
<https://github.com/mctools/ncrystal/wiki/Using-NCrystal#uniform-material-configuration-syntax>`_.
Temperature, material composition, and density are passed from the configuration
string and the `NCMAT file
<https://github.com/mctools/ncrystal/wiki/NCMAT-format>`_ that define the
material, e.g.::
mat = openmc.Material.from_ncrystal('Al_sg225.ncmat;temp=300K')
defines a material containing polycrystalline alumnium,
::
mat = openmc.Material.from_ncrystal("""Ge_sg227.ncmat;dcutoff=0.5;mos=40arcsec;
dir1=@crys_hkl:5,1,1@lab:0,0,1;
dir2=@crys_hkl:0,-1,1@lab:0,1,0""")
defines an oriented germanium single crystal with 40 arcsec mosaicity.
NCrystal only handles low energy neutron interactions. Other interactions are
provided by standard ACE files. NCrystal_ comes with a `predefined library
<https://github.com/mctools/ncrystal/wiki/Data-library>`_ but more materials can
be added by creating NCMAT files or on-the-fly in the configuration string.
.. warning:: Currently, NCrystal_ materials cannot be modified after they are created.
Density, temperature and composition should be defined in the
configuration string or the NCMAT file.
.. _NCrystal: https://github.com/mctools/ncrystal
------------------
Naming Conventions
------------------
@ -222,3 +258,4 @@ been generated, you can tell OpenMC to use this file either by setting
materials.cross_sections = '/path/to/cross_sections.xml'
.. _MCNP: https://mcnp.lanl.gov/

4
docs/source/usersguide/scripts.rst
View file

@ -157,9 +157,9 @@ geometry.xml
added. Any 'surfaces' attributes/elements on a cell will be renamed 'region'.
materials.xml
Nuclide names will be changed from ACE aliases (e.g., Am-242m) to HDF5/GND
Nuclide names will be changed from ACE aliases (e.g., Am-242m) to HDF5/GNDS
names (e.g., Am242_m1). Thermal scattering table names will be changed from
ACE aliases (e.g., HH2O) to HDF5/GND names (e.g., c_H_in_H2O).
ACE aliases (e.g., HH2O) to HDF5/GNDS names (e.g., c_H_in_H2O).
----------------------
``openmc-update-mgxs``

1
docs/source/usersguide/settings.rst
View file

@ -474,7 +474,6 @@ selected::
Some features related to photon transport are not currently implemented,
including:
* Tallying photon energy deposition.
* Generating a photon source from a neutron calculation that can be used
for a later fixed source photon calculation.
* Photoneutron reactions.

5
docs/source/usersguide/tallies.rst
View file

@ -269,9 +269,8 @@ The following tables show all valid scores:
|heating |Total nuclear heating in units of eV per source |
| |particle. For neutrons, this corresponds to MT=301 |
| |produced by NJOY's HEATR module while for photons, |
| |this is tallied from either direct photon energy |
| |deposition (analog estimator) or pre-generated |
| |photon heating number. See :ref:`methods_heating` |
| |this is tallied from direct photon energy |
| |deposition. See :ref:`methods_heating`. |
+----------------------+---------------------------------------------------+
|heating-local |Total nuclear heating in units of eV per source |
| |particle assuming energy from secondary photons is |

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(

146
include/openmc/cell.h
View file

@ -52,6 +52,100 @@ extern vector<unique_ptr<Cell>> cells;
} // namespace model
//==============================================================================
class Region {
public:
//----------------------------------------------------------------------------
// Constructors
Region() {}
explicit Region(std::string region_spec, int32_t cell_id);
//----------------------------------------------------------------------------
// Methods
//! \brief Determine if a cell contains the particle at a given location.
//!
//! The bounds of the cell are determined by a logical expression involving
//! surface half-spaces. The expression used is given in infix notation
//!
//! The function is split into two cases, one for simple cells (those
//! involving only the intersection of half-spaces) and one for complex cells.
//! Both cases use short circuiting; however, in the case fo complex cells,
//! the complexity increases with the binary operators involved.
//! \param r The 3D Cartesian coordinate to check.
//! \param u A direction used to "break ties" the coordinates are very
//! close to a surface.
//! \param on_surface The signed index of a surface that the coordinate is
//! known to be on. This index takes precedence over surface sense
//! calculations.
bool contains(Position r, Direction u, int32_t on_surface) const;
//! Find the oncoming boundary of this cell.
std::pair<double, int32_t> distance(
Position r, Direction u, int32_t on_surface, Particle* p) const;
//! Get the BoundingBox for this cell.
BoundingBox bounding_box(int32_t cell_id) const;
//! Get the CSG expression as a string
std::string str() const;
//! Get a vector containing all the surfaces in the region expression
vector<int32_t> surfaces() const;
//----------------------------------------------------------------------------
// Accessors
//! Get Boolean of if the cell is simple or not
bool is_simple() const { return simple_; }
private:
//----------------------------------------------------------------------------
// Private Methods
//! Get a vector of the region expression in postfix notation
vector<int32_t> generate_postfix(int32_t cell_id) const;
//! Determine if a particle is inside the cell for a simple cell (only
//! intersection operators)
bool contains_simple(Position r, Direction u, int32_t on_surface) const;
//! Determine if a particle is inside the cell for a complex cell.
//!
//! Uses the comobination of half-spaces and binary operators to determine
//! if short circuiting can be used. Short cicuiting uses the relative and
//! absolute depth of parenthases in the expression.
bool contains_complex(Position r, Direction u, int32_t on_surface) const;
//! BoundingBox if the paritcle is in a simple cell.
BoundingBox bounding_box_simple() const;
//! BoundingBox if the particle is in a complex cell.
BoundingBox bounding_box_complex(vector<int32_t> postfix) const;
//! Enfource precedence: Parenthases, Complement, Intersection, Union
void add_precedence();
//! Add parenthesis to enforce precedence
std::vector<int32_t>::iterator add_parentheses(
std::vector<int32_t>::iterator start);
//! Remove complement operators from the expression
void remove_complement_ops();
//! Remove complement operators by using DeMorgan's laws
void apply_demorgan(
vector<int32_t>::iterator start, vector<int32_t>::iterator stop);
//----------------------------------------------------------------------------
// Private Data
//! Definition of spatial region as Boolean expression of half-spaces
// TODO: Should this be a vector of some other type
vector<int32_t> expression_;
bool simple_; //!< Does the region contain only intersections?
};
//==============================================================================
class Cell {
@ -108,6 +202,12 @@ public:
//! Get the BoundingBox for this cell.
virtual BoundingBox bounding_box() const = 0;
//! Get a vector of surfaces in the cell
virtual vector<int32_t> surfaces() const { return vector<int32_t>(); }
//! Check if the cell region expression is simple
virtual bool is_simple() const { return true; }
//----------------------------------------------------------------------------
// Accessors
@ -198,10 +298,6 @@ public:
//! T. The units are sqrt(eV).
vector<double> sqrtkT_;
//! Definition of spatial region as Boolean expression of half-spaces
vector<int32_t> region_;
bool simple_; //!< Does the region contain only intersections?
//! \brief Neighboring cells in the same universe.
NeighborList neighbors_;
@ -227,35 +323,36 @@ struct CellInstanceItem {
class CSGCell : public Cell {
public:
//----------------------------------------------------------------------------
// Constructors
CSGCell();
explicit CSGCell(pugi::xml_node cell_node);
bool contains(Position r, Direction u, int32_t on_surface) const override;
//----------------------------------------------------------------------------
// Methods
vector<int32_t> surfaces() const override { return region_.surfaces(); }
std::pair<double, int32_t> distance(
Position r, Direction u, int32_t on_surface, Particle* p) const override;
Position r, Direction u, int32_t on_surface, Particle* p) const override
{
return region_.distance(r, u, on_surface, p);
}
bool contains(Position r, Direction u, int32_t on_surface) const override
{
return region_.contains(r, u, on_surface);
}
BoundingBox bounding_box() const override
{
return region_.bounding_box(id_);
}
void to_hdf5_inner(hid_t group_id) const override;
BoundingBox bounding_box() const override;
bool is_simple() const override { return region_.is_simple(); }
protected:
bool contains_simple(Position r, Direction u, int32_t on_surface) const;
bool contains_complex(Position r, Direction u, int32_t on_surface) const;
BoundingBox bounding_box_simple() const;
static BoundingBox bounding_box_complex(vector<int32_t> postfix);
//! Applies DeMorgan's laws to a section of the RPN
//! \param start Starting point for token modification
//! \param stop Stopping point for token modification
static void apply_demorgan(
vector<int32_t>::iterator start, vector<int32_t>::iterator stop);
//! Removes complement operators from the RPN
//! \param rpn The rpn to remove complement operators from.
static void remove_complement_ops(vector<int32_t>& rpn);
//! Returns the beginning position of a parenthesis block (immediately before
//! two surface tokens) in the RPN given a starting position at the end of
//! that block (immediately after two surface tokens)
@ -263,6 +360,9 @@ protected:
//! \param rpn The rpn being searched
static vector<int32_t>::iterator find_left_parenthesis(
vector<int32_t>::iterator start, const vector<int32_t>& rpn);
private:
Region region_;
};
//==============================================================================

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"

5
include/openmc/cross_sections.h
View file

@ -62,6 +62,11 @@ extern vector<Library> libraries;
//! libraries
void read_cross_sections_xml();
//! Read cross sections file (either XML or multigroup H5) and populate data
//! libraries
//! \param[in] root node of the cross_sections.xml
void read_cross_sections_xml(pugi::xml_node root);
//! Load nuclide and thermal scattering data from HDF5 files
//
//! \param[in] nuc_temps Temperatures for each nuclide in [K]

12
include/openmc/distribution.h
View file

@ -47,14 +47,20 @@ public:
// Properties
const vector<double>& x() const { return x_; }
const vector<double>& p() const { return p_; }
const vector<double>& prob() const { return prob_; }
const vector<size_t>& alias() const { return alias_; }
private:
vector<double> x_; //!< Possible outcomes
vector<double> p_; //!< Probability of each outcome
vector<double> x_; //!< Possible outcomes
vector<double> prob_; //!< Probability of accepting the uniformly sampled bin,
//!< mapped to alias method table
vector<size_t> alias_; //!< Alias table
//! Normalize distribution so that probabilities sum to unity
void normalize();
//! Initialize alias tables for distribution
void init_alias(vector<double>& x, vector<double>& p);
};
//==============================================================================

18
include/openmc/file_utils.h
View file

@ -3,14 +3,32 @@
#include <fstream> // for ifstream
#include <string>
#include <sys/stat.h>
namespace openmc {
// TODO: replace with std::filesystem when switch to C++17 is made
//! Determine if a path is a directory
//! \param[in] path Path to check
//! \return Whether the path is a directory
inline bool dir_exists(const std::string& path)
{
struct stat s;
if (stat(path.c_str(), &s) != 0)
return false;
return s.st_mode & S_IFDIR;
}
//! Determine if a file exists
//! \param[in] filename Path to file
//! \return Whether file exists
inline bool file_exists(const std::string& filename)
{
// rule out file being a path to a directory
if (dir_exists(filename))
return false;
std::ifstream s {filename};
return s.good();
}

6
include/openmc/geometry_aux.h
View file

@ -10,6 +10,7 @@
#include <vector>
#include "openmc/vector.h"
#include "openmc/xml_interface.h"
namespace openmc {
@ -19,8 +20,13 @@ extern std::unordered_map<int32_t, std::unordered_map<int32_t, int32_t>>
extern std::unordered_map<int32_t, int32_t> universe_level_counts;
} // namespace model
//! Read geometry from XML file
void read_geometry_xml();
//! Read geometry from XML node
//! \param[in] root node of geometry XML element
void read_geometry_xml(pugi::xml_node root);
//==============================================================================
//! Replace Universe, Lattice, and Material IDs with indices.
//==============================================================================

10
include/openmc/initialize.h
View file

@ -1,6 +1,8 @@
#ifndef OPENMC_INITIALIZE_H
#define OPENMC_INITIALIZE_H
#include <string>
#ifdef OPENMC_MPI
#include "mpi.h"
#endif
@ -11,7 +13,13 @@ int parse_command_line(int argc, char* argv[]);
#ifdef OPENMC_MPI
void initialize_mpi(MPI_Comm intracomm);
#endif
void read_input_xml();
//! Read material, geometry, settings, and tallies from a single XML file
bool read_model_xml();
//! Read inputs from separate XML files
void read_separate_xml_files();
//! Write some output that occurs right after initialization
void initial_output();
} // namespace openmc

10
include/openmc/material.h
View file

@ -12,6 +12,7 @@
#include "openmc/bremsstrahlung.h"
#include "openmc/constants.h"
#include "openmc/memory.h" // for unique_ptr
#include "openmc/ncrystal_interface.h"
#include "openmc/particle.h"
#include "openmc/vector.h"
@ -151,12 +152,17 @@ public:
//! \return Temperature in [K]
double temperature() const;
//! Get pointer to NCrystal material object
//! \return Pointer to NCrystal material object
const NCrystalMat& ncrystal_mat() const { return ncrystal_mat_; };
//----------------------------------------------------------------------------
// Data
int32_t id_ {C_NONE}; //!< Unique ID
std::string name_; //!< Name of material
vector<int> nuclide_; //!< Indices in nuclides vector
vector<int> element_; //!< Indices in elements vector
NCrystalMat ncrystal_mat_; //!< NCrystal material object
xt::xtensor<double, 1> atom_density_; //!< Nuclide atom density in [atom/b-cm]
double density_; //!< Total atom density in [atom/b-cm]
double density_gpcc_; //!< Total atom density in [g/cm^3]
@ -221,6 +227,10 @@ double density_effect(const vector<double>& f, const vector<double>& e_b_sq,
//! Read material data from materials.xml
void read_materials_xml();
//! Read material data XML node
//! \param[in] root node of materials XML element
void read_materials_xml(pugi::xml_node root);
void free_memory_material();
} // namespace openmc

36
include/openmc/mcpl_interface.h Normal file
View file

@ -0,0 +1,36 @@
#ifndef OPENMC_MCPL_INTERFACE_H
#define OPENMC_MCPL_INTERFACE_H
#include "openmc/particle_data.h"
#include "openmc/vector.h"
#include <string>
namespace openmc {
//==============================================================================
// Constants
//==============================================================================
extern "C" const bool MCPL_ENABLED;
//==============================================================================
// Functions
//==============================================================================
//! Get a vector of source sites from an MCPL file
//
//! \param[in] path Path to MCPL file
//! \return Vector of source sites
vector<SourceSite> mcpl_source_sites(std::string path);
//! Write an MCPL source file
//
//! \param[in] filename Path to MCPL file
//! \param[in] surf_source_bank Whether to use the surface source bank
void write_mcpl_source_point(
const char* filename, bool surf_source_bank = false);
} // namespace openmc
#endif // OPENMC_MCPL_INTERFACE_H

94
include/openmc/ncrystal_interface.h Normal file
View file

@ -0,0 +1,94 @@
#ifndef OPENMC_NCRYSTAL_INTERFACE_H
#define OPENMC_NCRYSTAL_INTERFACE_H
#ifdef NCRYSTAL
#include "NCrystal/NCRNG.hh"
#include "NCrystal/NCrystal.hh"
#endif
#include "openmc/particle.h"
#include <cstdint> // for uint64_t
#include <limits> // for numeric_limits
#include <string>
namespace openmc {
//==============================================================================
// Constants
//==============================================================================
extern "C" const bool NCRYSTAL_ENABLED;
//! Energy in [eV] to switch between NCrystal and ENDF
constexpr double NCRYSTAL_MAX_ENERGY {5.0};
//==============================================================================
// Wrapper class an NCrystal material
//==============================================================================
class NCrystalMat {
public:
//----------------------------------------------------------------------------
// Constructors
NCrystalMat() = default;
explicit NCrystalMat(const std::string& cfg);
//----------------------------------------------------------------------------
// Methods
#ifdef NCRYSTAL
//! Return configuration string
std::string cfg() const;
//! Get cross section from NCrystal material
//
//! \param[in] p Particle object
//! \return Cross section in [b]
double xs(const Particle& p) const;
// Process scattering event
//
//! \param[in] p Particle object
void scatter(Particle& p) const;
//! Whether the object holds a valid NCrystal material
operator bool() const;
#else
//----------------------------------------------------------------------------
// Trivial methods when compiling without NCRYSTAL
std::string cfg() const
{
return "";
}
double xs(const Particle& p) const
{
return -1.0;
}
void scatter(Particle& p) const {}
operator bool() const
{
return false;
}
#endif
private:
//----------------------------------------------------------------------------
// Data members (only present when compiling with NCrystal support)
#ifdef NCRYSTAL
std::string cfg_; //!< NCrystal configuration string
std::shared_ptr<const NCrystal::ProcImpl::Process>
ptr_; //!< Pointer to NCrystal material object
#endif
};
//==============================================================================
// Functions
//==============================================================================
void ncrystal_update_micro(double xs, NuclideMicroXS& micro);
} // namespace openmc
#endif // OPENMC_NCRYSTAL_INTERFACE_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();

3
include/openmc/photon.h
View file

@ -100,6 +100,9 @@ public:
// Bremsstrahlung scaled DCS
xt::xtensor<double, 2> dcs_;
// Whether atomic relaxation data is present
bool has_atomic_relaxation_ {false};
// Constant data
static constexpr int MAX_STACK_SIZE =
7; //!< maximum possible size of atomic relaxation stack

5
include/openmc/physics.h
View file

@ -96,8 +96,9 @@ void inelastic_scatter(const Nuclide& nuc, const Reaction& rx, Particle& p);
void sample_secondary_photons(Particle& p, int i_nuclide);
//！Split or Roulette particles based their weight and the lower weight window
// bound.
//! Split or Roulette particles based their weight and the lower weight window
//! bound.
//
//! \param[in] p, particle to be split or rouletted with the weight window.
void split_particle(Particle& p);

4
include/openmc/plot.h
View file

@ -279,6 +279,10 @@ void voxel_finalize(hid_t dspace, hid_t dset, hid_t memspace);
//! Read plot specifications from a plots.xml file
void read_plots_xml();
//! Read plot specifications from an XML Node
//! \param[in] XML node containing plot info
void read_plots_xml(pugi::xml_node root);
//! Clear memory
void free_memory_plot();

8
include/openmc/settings.h
View file

@ -28,6 +28,7 @@ extern bool check_overlaps; //!< check overlaps in geometry?
extern bool confidence_intervals; //!< use confidence intervals for results?
extern bool
create_fission_neutrons; //!< create fission neutrons (fixed source)?
extern bool create_delayed_neutrons; //!< create delayed fission neutrons?
extern "C" bool cmfd_run; //!< is a CMFD run?
extern bool
delayed_photon_scaling; //!< Scale fission photon yield to include delayed
@ -47,7 +48,9 @@ extern "C" bool run_CE; //!< run with continuous-energy data?
extern bool source_latest; //!< write latest source at each batch?
extern bool source_separate; //!< write source to separate file?
extern bool source_write; //!< write source in HDF5 files?
extern bool source_mcpl_write; //!< write source in mcpl files?
extern bool surf_source_write; //!< write surface source file?
extern bool surf_mcpl_write; //!< write surface mcpl file?
extern bool surf_source_read; //!< read surface source file?
extern bool survival_biasing; //!< use survival biasing?
extern bool temperature_multipole; //!< use multipole data?
@ -127,9 +130,12 @@ extern double weight_survive; //!< Survival weight after Russian roulette
//==============================================================================
//! Read settings from XML file
//! \param[in] root XML node for <settings>
void read_settings_xml();
//! Read settings from XML node
//! \param[in] root XML node for <settings>
void read_settings_xml(pugi::xml_node root);
void free_memory_settings();
} // namespace openmc

1
include/openmc/source.h
View file

@ -102,6 +102,7 @@ class FileSource : public Source {
public:
// Constructors
explicit FileSource(std::string path);
explicit FileSource(const vector<SourceSite>& sites) : sites_ {sites} {}
// Methods
SourceSite sample(uint64_t* seed) const override;

33
include/openmc/tallies/filter.h
View file

@ -17,6 +17,34 @@
namespace openmc {
enum class FilterType {
AZIMUTHAL,
CELLBORN,
CELLFROM,
CELL,
CELL_INSTANCE,
COLLISION,
DELAYED_GROUP,
DISTRIBCELL,
ENERGY_FUNCTION,
ENERGY,
ENERGY_OUT,
LEGENDRE,
MATERIAL,
MESH,
MESH_SURFACE,
MU,
PARTICLE,
POLAR,
SPHERICAL_HARMONICS,
SPATIAL_LEGENDRE,
SURFACE,
TIME,
UNIVERSE,
ZERNIKE,
ZERNIKE_RADIAL
};
//==============================================================================
//! Modifies tally score events.
//==============================================================================
@ -58,7 +86,8 @@ public:
//----------------------------------------------------------------------------
// Methods
virtual std::string type() const = 0;
virtual std::string type_str() const = 0;
virtual FilterType type() const = 0;
//! Matches a tally event to a set of filter bins and weights.
//!
@ -72,7 +101,7 @@ public:
//! Writes data describing this filter to an HDF5 statepoint group.
virtual void to_statepoint(hid_t filter_group) const
{
write_dataset(filter_group, "type", type());
write_dataset(filter_group, "type", type_str());
write_dataset(filter_group, "n_bins", n_bins_);
}

3
include/openmc/tallies/filter_azimuthal.h
View file

@ -24,7 +24,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "azimuthal"; }
std::string type_str() const override { return "azimuthal"; }
FilterType type() const override { return FilterType::AZIMUTHAL; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_cell.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "cell"; }
std::string type_str() const override { return "cell"; }
FilterType type() const override { return FilterType::CELL; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_cell_instance.h
View file

@ -28,7 +28,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "cellinstance"; }
std::string type_str() const override { return "cellinstance"; }
FilterType type() const override { return FilterType::CELL_INSTANCE; }
void from_xml(pugi::xml_node node) override;

5
include/openmc/tallies/filter_cellborn.h
View file

@ -11,12 +11,13 @@ namespace openmc {
//! Specifies which cell the particle was born in.
//==============================================================================
class CellbornFilter : public CellFilter {
class CellBornFilter : public CellFilter {
public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "cellborn"; }
std::string type_str() const override { return "cellborn"; }
FilterType type() const override { return FilterType::CELLBORN; }
void get_all_bins(const Particle& p, TallyEstimator estimator,
FilterMatch& match) const override;

3
include/openmc/tallies/filter_cellfrom.h
View file

@ -16,7 +16,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "cellfrom"; }
std::string type_str() const override { return "cellfrom"; }
FilterType type() const override { return FilterType::CELLFROM; }
void get_all_bins(const Particle& p, TallyEstimator estimator,
FilterMatch& match) const override;

3
include/openmc/tallies/filter_collision.h
View file

@ -23,7 +23,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "collision"; }
std::string type_str() const override { return "collision"; }
FilterType type() const override { return FilterType::COLLISION; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_delayedgroup.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "delayedgroup"; }
std::string type_str() const override { return "delayedgroup"; }
FilterType type() const override { return FilterType::DELAYED_GROUP; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_distribcell.h
View file

@ -21,7 +21,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "distribcell"; }
std::string type_str() const override { return "distribcell"; }
FilterType type() const override { return FilterType::DISTRIBCELL; }
void from_xml(pugi::xml_node node) override;

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

@ -22,7 +22,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "energy"; }
std::string type_str() const override { return "energy"; }
FilterType type() const override { return FilterType::ENERGY; }
void from_xml(pugi::xml_node node) override;
@ -63,7 +64,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "energyout"; }
std::string type_str() const override { return "energyout"; }
FilterType type() const override { return FilterType::ENERGY_OUT; }
void get_all_bins(const Particle& p, TallyEstimator estimator,
FilterMatch& match) const override;

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

@ -24,7 +24,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "energyfunction"; }
std::string type_str() const override { return "energyfunction"; }
FilterType type() const override { return FilterType::ENERGY_FUNCTION; }
void from_xml(pugi::xml_node node) override;
@ -40,8 +41,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 +54,9 @@ private:
//! Interpolant values.
vector<double> y_;
//! Interpolation scheme
Interpolation interpolation_ {Interpolation::lin_lin};
};
} // namespace openmc

3
include/openmc/tallies/filter_legendre.h
View file

@ -21,7 +21,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "legendre"; }
std::string type_str() const override { return "legendre"; }
FilterType type() const override { return FilterType::LEGENDRE; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_material.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "material"; }
std::string type_str() const override { return "material"; }
FilterType type() const override { return FilterType::MATERIAL; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_mesh.h
View file

@ -24,7 +24,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "mesh"; }
std::string type_str() const override { return "mesh"; }
FilterType type() const override { return FilterType::MESH; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_meshsurface.h
View file

@ -10,7 +10,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "meshsurface"; }
std::string type_str() const override { return "meshsurface"; }
FilterType type() const override { return FilterType::MESH_SURFACE; }
void get_all_bins(const Particle& p, TallyEstimator estimator,
FilterMatch& match) const override;

3
include/openmc/tallies/filter_mu.h
View file

@ -23,7 +23,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "mu"; }
std::string type_str() const override { return "mu"; }
FilterType type() const override { return FilterType::MU; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_particle.h
View file

@ -21,7 +21,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "particle"; }
std::string type_str() const override { return "particle"; }
FilterType type() const override { return FilterType::PARTICLE; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_polar.h
View file

@ -24,7 +24,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "polar"; }
std::string type_str() const override { return "polar"; }
FilterType type() const override { return FilterType::POLAR; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_sph_harm.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "sphericalharmonics"; }
std::string type_str() const override { return "sphericalharmonics"; }
FilterType type() const override { return FilterType::SPHERICAL_HARMONICS; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_sptl_legendre.h
View file

@ -23,7 +23,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "spatiallegendre"; }
std::string type_str() const override { return "spatiallegendre"; }
FilterType type() const override { return FilterType::SPATIAL_LEGENDRE; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_surface.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "surface"; }
std::string type_str() const override { return "surface"; }
FilterType type() const override { return FilterType::SURFACE; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_time.h
View file

@ -22,7 +22,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "time"; }
std::string type_str() const override { return "time"; }
FilterType type() const override { return FilterType::TIME; }
void from_xml(pugi::xml_node node) override;

3
include/openmc/tallies/filter_universe.h
View file

@ -25,7 +25,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "universe"; }
std::string type_str() const override { return "universe"; }
FilterType type() const override { return FilterType::UNIVERSE; }
void from_xml(pugi::xml_node node) override;

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

@ -21,7 +21,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "zernike"; }
std::string type_str() const override { return "zernike"; }
FilterType type() const override { return FilterType::ZERNIKE; }
void from_xml(pugi::xml_node node) override;
@ -72,7 +73,8 @@ public:
//----------------------------------------------------------------------------
// Methods
std::string type() const override { return "zernikeradial"; }
std::string type_str() const override { return "zernikeradial"; }
FilterType type() const override { return FilterType::ZERNIKE_RADIAL; }
void get_all_bins(const Particle& p, TallyEstimator estimator,
FilterMatch& match) const override;

8
include/openmc/tallies/tally.h
View file

@ -48,10 +48,10 @@ public:
void set_nuclides(const vector<std::string>& nuclides);
//! returns vector of indices corresponding to the tally this is called on
const vector<int32_t>& filters() const { return filters_; }
const vector<int32_t>& filters() const { return filters_; }
//! \brief Returns the tally filter at index i
int32_t filters(int i) const { return filters_[i]; }
int32_t filters(int i) const { return filters_[i]; }
void set_filters(gsl::span<Filter*> filters);
@ -178,6 +178,10 @@ extern double global_tally_leakage;
//! Read tally specification from tallies.xml
void read_tallies_xml();
//! Read tally specification from an XML node
//! \param[in] root node of tallies XML element
void read_tallies_xml(pugi::xml_node root);
//! \brief Accumulate the sum of the contributions from each history within the
//! batch to a new random variable
void accumulate_tallies();

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
man/man1/openmc.1
View file

@ -57,7 +57,7 @@ Indicates the default path to an HDF5 file that contains multi-group cross
section libraries if the user has not specified the <cross_sections> tag in
.I materials.xml\fP.
.SH LICENSE
Copyright \(co 2011-2022 Massachusetts Institute of Technology, UChicago
Copyright \(co 2011-2023 Massachusetts Institute of Technology, UChicago
Argonne LLC, and OpenMC contributors.
.PP
Permission is hereby granted, free of charge, to any person obtaining a copy of

6
openmc/__init__.py
View file

@ -10,11 +10,11 @@ from openmc.material import *
from openmc.plots import *
from openmc.region import *
from openmc.volume import *
from openmc.source import *
from openmc.weight_windows import *
from openmc.settings import *
from openmc.surface import *
from openmc.universe import *
from openmc.source import *
from openmc.settings import *
from openmc.lattice import *
from openmc.filter import *
from openmc.filter_expansion import *
@ -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'

32
openmc/_xml.py
View file

@ -1,22 +1,40 @@
def clean_indentation(element, level=0, spaces_per_level=2):
"""
copy and paste from https://effbot.org/zone/element-lib.htm#prettyprint
it basically walks your tree and adds spaces and newlines so the tree is
printed in a nice way
def clean_indentation(element, level=0, spaces_per_level=2, trailing_indent=True):
"""Set indentation of XML element and its sub-elements.
Copied and pasted from https://effbot.org/zone/element-lib.htm#prettyprint.
It walks your tree and adds spaces and newlines so the tree is
printed in a nice way.
Parameters
----------
level : int
Indentation level for the element passed in (default 0)
spaces_per_level : int
Number of spaces per indentation level (default 2)
trailing_indent : bool
Whether or not to add indentation after closing the element
"""
i = "\n" + level*spaces_per_level*" "
# ensure there's always some tail for the element passed in
if not element.tail:
element.tail = ""
if len(element):
if not element.text or not element.text.strip():
element.text = i + spaces_per_level*" "
if not element.tail or not element.tail.strip():
if trailing_indent and (not element.tail or not element.tail.strip()):
element.tail = i
for sub_element in element:
# `trailing_indent` is intentionally not forwarded to the recursive
# call. Any child element of the topmost element should add
# indentation at the end to ensure its parent's indentation is
# correct.
clean_indentation(sub_element, level+1, spaces_per_level)
if not sub_element.tail or not sub_element.tail.strip():
sub_element.tail = i
else:
if level and (not element.tail or not element.tail.strip()):
if trailing_indent and level and (not element.tail or not element.tail.strip()):
element.tail = i

23
openmc/cell.py
View file

@ -1,6 +1,5 @@
from collections import OrderedDict
from collections.abc import Iterable
from copy import deepcopy
from math import cos, sin, pi
from numbers import Real
from xml.etree import ElementTree as ET
@ -318,12 +317,12 @@ class Cell(IDManagerMixin):
@temperature.setter
def temperature(self, temperature):
# Make sure temperatures are positive
cv.check_type('cell temperature', temperature, (Iterable, Real))
cv.check_type('cell temperature', temperature, (Iterable, Real), none_ok=True)
if isinstance(temperature, Iterable):
cv.check_type('cell temperature', temperature, Iterable, Real)
for T in temperature:
cv.check_greater_than('cell temperature', T, 0.0, True)
else:
elif isinstance(temperature, Real):
cv.check_greater_than('cell temperature', temperature, 0.0, True)
# If this cell is filled with a universe or lattice, propagate
@ -519,8 +518,14 @@ class Cell(IDManagerMixin):
paths = self._paths
self._paths = None
clone = deepcopy(self)
clone.id = None
clone = openmc.Cell(name=self.name)
clone.volume = self.volume
if self.temperature is not None:
clone.temperature = self.temperature
if self.translation is not None:
clone.translation = self.translation
if self.rotation is not None:
clone.rotation = self.rotation
clone._num_instances = None
# Restore paths on original instance
@ -634,6 +639,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
@ -647,7 +655,7 @@ class Cell(IDManagerMixin):
surfaces : dict
Dictionary mapping surface IDs to :class:`openmc.Surface` instances
materials : dict
Dictionary mapping material IDs to :class:`openmc.Material`
Dictionary mapping material ID strings to :class:`openmc.Material`
instances (defined in :math:`openmc.Geometry.from_xml`)
get_universe : function
Function returning universe (defined in
@ -687,6 +695,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:

4
openmc/checkvalue.py
View file

@ -1,12 +1,12 @@
import copy
import os
from typing import Union
import typing # required to prevent typing.Union namespace overwriting Union
from collections.abc import Iterable
import numpy as np
# Type for arguments that accept file paths
PathLike = Union[str, os.PathLike]
PathLike = typing.Union[str, os.PathLike]
def check_type(name, value, expected_type, expected_iter_type=None, *, none_ok=False):

4
openmc/config.py
View file

@ -42,7 +42,9 @@ class _Config(MutableMapping):
# Reset photon source data since it relies on chain file
_DECAY_PHOTON_ENERGY.clear()
else:
raise KeyError(f'Unrecognized config key: {key}')
raise KeyError(f'Unrecognized config key: {key}. Acceptable keys '
'are "cross_sections", "mg_cross_sections" and '
'"chain_file"')
def __iter__(self):
return iter(self._mapping)

4
openmc/data/ace.py
View file

@ -24,7 +24,7 @@ import numpy as np
import openmc.checkvalue as cv
from openmc.mixin import EqualityMixin
from .data import ATOMIC_SYMBOL, gnd_name, EV_PER_MEV, K_BOLTZMANN
from .data import ATOMIC_SYMBOL, gnds_name, EV_PER_MEV, K_BOLTZMANN
from .endf import ENDF_FLOAT_RE
@ -88,7 +88,7 @@ def get_metadata(zaid, metastable_scheme='nndc'):
# Determine name
element = ATOMIC_SYMBOL[Z]
name = gnd_name(Z, mass_number, metastable)
name = gnds_name(Z, mass_number, metastable)
return (name, element, Z, mass_number, metastable)

33
openmc/data/data.py
View file

@ -197,8 +197,8 @@ NEUTRON_MASS = 1.00866491595
# Used in atomic_mass function as a cache
_ATOMIC_MASS = {}
# Regex for GND nuclide names (used in zam function)
_GND_NAME_RE = re.compile(r'([A-Zn][a-z]*)(\d+)((?:_[em]\d+)?)')
# Regex for GNDS nuclide names (used in zam function)
_GNDS_NAME_RE = re.compile(r'([A-Zn][a-z]*)(\d+)((?:_[em]\d+)?)')
# Used in half_life function as a cache
_HALF_LIFE = {}
@ -223,19 +223,19 @@ def atomic_mass(isotope):
"""
if not _ATOMIC_MASS:
# Load data from AME2016 file
mass_file = os.path.join(os.path.dirname(__file__), 'mass16.txt')
# Load data from AME2020 file
mass_file = os.path.join(os.path.dirname(__file__), 'mass_1.mas20.txt')
with open(mass_file, 'r') as ame:
# Read lines in file starting at line 40
for line in itertools.islice(ame, 39, None):
# Read lines in file starting at line 37
for line in itertools.islice(ame, 36, None):
name = f'{line[20:22].strip()}{int(line[16:19])}'
mass = float(line[96:99]) + 1e-6*float(
line[100:106] + '.' + line[107:112])
mass = float(line[106:109]) + 1e-6*float(
line[110:116] + '.' + line[117:123])
_ATOMIC_MASS[name.lower()] = mass
# For isotopes found in some libraries that represent all natural
# isotopes of their element (e.g. C0), calculate the atomic mass as
# the sum of the atomic mass times the natural abudance of the isotopes
# the sum of the atomic mass times the natural abundance of the isotopes
# that make up the element.
for element in ['C', 'Zn', 'Pt', 'Os', 'Tl']:
isotope_zero = element.lower() + '0'
@ -436,8 +436,11 @@ def water_density(temperature, pressure=0.1013):
return coeff / pi / gamma1_pi
def gnd_name(Z, A, m=0):
"""Return nuclide name using GND convention
def gnds_name(Z, A, m=0):
"""Return nuclide name using GNDS convention
.. versionchanged:: 0.14.0
Function name changed from ``gnd_name`` to ``gnds_name``
Parameters
----------
@ -451,7 +454,7 @@ def gnd_name(Z, A, m=0):
Returns
-------
str
Nuclide name in GND convention, e.g., 'Am242_m1'
Nuclide name in GNDS convention, e.g., 'Am242_m1'
"""
if m > 0:
@ -502,7 +505,7 @@ def zam(name):
Parameters
----------
name : str
Name of nuclide using GND convention, e.g., 'Am242_m1'
Name of nuclide using GNDS convention, e.g., 'Am242_m1'
Returns
-------
@ -511,10 +514,10 @@ def zam(name):
"""
try:
symbol, A, state = _GND_NAME_RE.match(name).groups()
symbol, A, state = _GNDS_NAME_RE.match(name).groups()
except AttributeError:
raise ValueError(f"'{name}' does not appear to be a nuclide name in "
"GND format")
"GNDS format")
if symbol not in ATOMIC_NUMBER:
raise ValueError(f"'{symbol}' is not a recognized element symbol")

66
openmc/data/decay.py
View file

@ -1,10 +1,9 @@
from collections.abc import Iterable
from io import StringIO
from math import log
from pathlib import Path
import re
from typing import Optional
from warnings import warn
from xml.etree import ElementTree as ET
import numpy as np
from uncertainties import ufloat, UFloat
@ -145,7 +144,7 @@ class FissionProductYields(EqualityMixin):
# Assign basic nuclide properties
self.nuclide = {
'name': ev.gnd_name,
'name': ev.gnds_name,
'atomic_number': ev.target['atomic_number'],
'mass_number': ev.target['mass_number'],
'isomeric_state': ev.target['isomeric_state']
@ -559,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
@ -579,14 +580,14 @@ class Decay(EqualityMixin):
_DECAY_PHOTON_ENERGY = {}
def decay_photon_energy(nuclide: str) -> Univariate:
def decay_photon_energy(nuclide: str) -> Optional[Univariate]:
"""Get photon energy distribution resulting from the decay of a nuclide
This function relies on data stored in a depletion chain. Before calling it
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
----------
@ -595,9 +596,10 @@ def decay_photon_energy(nuclide: str) -> Univariate:
Returns
-------
openmc.stats.Univariate
Distribution of energies in [eV] of photons emitted from decay. Note
that the probabilities represent intensities, given as [decay/sec].
openmc.stats.Univariate or None
Distribution of energies in [eV] of photons emitted from decay, or None
if no photon source exists. Note that the probabilities represent
intensities, given as [decay/sec].
"""
if not _DECAY_PHOTON_ENERGY:
chain_file = openmc.config.get('chain_file')
@ -619,3 +621,49 @@ def decay_photon_energy(nuclide: str) -> Univariate:
"sources listed.")
return _DECAY_PHOTON_ENERGY.get(nuclide)
_DECAY_ENERGY = {}
def decay_energy(nuclide: str):
"""Get decay energy value resulting from the decay of a nuclide
This function relies on data stored in a depletion chain. Before calling it
for the first time, you need to ensure that a depletion chain has been
specified in openmc.config['chain_file'].
.. versionadded:: 0.13.3
Parameters
----------
nuclide : str
Name of nuclide, e.g., 'H3'
Returns
-------
float
Decay energy of nuclide in [eV]. If the nuclide is stable, a value of
0.0 is returned.
"""
if not _DECAY_ENERGY:
chain_file = openmc.config.get('chain_file')
if chain_file is None:
raise DataError(
"A depletion chain file must be specified with "
"openmc.config['chain_file'] in order to load decay data."
)
from openmc.deplete import Chain
chain = Chain.from_xml(chain_file)
for nuc in chain.nuclides:
if nuc.decay_energy:
_DECAY_ENERGY[nuc.name] = nuc.decay_energy
# If the chain file contained no decay energy, warn the user
if not _DECAY_ENERGY:
warn(f"Chain file '{chain_file}' does not have any decay energy.")
return _DECAY_ENERGY.get(nuclide, 0.0)

12
openmc/data/endf.py
View file

@ -12,7 +12,7 @@ import re
import numpy as np
from .data import gnd_name
from .data import gnds_name
from .function import Tabulated1D
try:
from ._endf import float_endf
@ -520,10 +520,10 @@ class Evaluation:
self.reaction_list.append((mf, mt, nc, mod))
@property
def gnd_name(self):
return gnd_name(self.target['atomic_number'],
self.target['mass_number'],
self.target['isomeric_state'])
def gnds_name(self):
return gnds_name(self.target['atomic_number'],
self.target['mass_number'],
self.target['isomeric_state'])
class Tabulated2D:
@ -531,7 +531,7 @@ class Tabulated2D:
This is a dummy class that is not really used other than to store the
interpolation information for a two-dimensional function. Once we refactor
to adopt GND-like data containers, this will probably be removed or
to adopt GNDS-like data containers, this will probably be removed or
extended.
Parameters

3475
openmc/data/mass16.txt
View file

File diff suppressed because it is too large Load diff

3594
openmc/data/mass_1.mas20.txt Normal file
View file

File diff suppressed because it is too large Load diff

12
openmc/data/multipole.py
View file

@ -748,12 +748,12 @@ class WindowedMultipole(EqualityMixin):
Parameters
----------
name : str
Name of the nuclide using the GND naming convention
Name of the nuclide using the GNDS naming convention
Attributes
----------
name : str
Name of the nuclide using the GND naming convention
Name of the nuclide using the GNDS naming convention
spacing : float
The width of each window in sqrt(E)-space. For example, the frst window
will end at (sqrt(E_min) + spacing)**2 and the second window at
@ -1154,8 +1154,8 @@ class WindowedMultipole(EqualityMixin):
Returns
-------
3-tuple of Real
Total, absorption, and fission microscopic cross sections at the
given energy and temperature.
Scattering, absorption, and fission microscopic cross sections
at the given energy and temperature.
"""
@ -1248,8 +1248,8 @@ class WindowedMultipole(EqualityMixin):
Returns
-------
3-tuple of Real or 3-tuple of numpy.ndarray
Total, absorption, and fission microscopic cross sections at the
given energy and temperature.
Scattering, absorption, and fission microscopic cross sections
at the given energy and temperature.
"""

7
openmc/data/neutron.py
View file

@ -44,7 +44,7 @@ class IncidentNeutron(EqualityMixin):
Parameters
----------
name : str
Name of the nuclide using the GND naming convention
Name of the nuclide using the GNDS naming convention
atomic_number : int
Number of protons in the target nucleus
mass_number : int
@ -75,7 +75,7 @@ class IncidentNeutron(EqualityMixin):
Metastable state of the target nucleus. A value of zero indicates ground
state.
name : str
Name of the nuclide using the GND naming convention
Name of the nuclide using the GNDS naming convention
reactions : collections.OrderedDict
Contains the cross sections, secondary angle and energy distributions,
and other associated data for each reaction. The keys are the MT values
@ -868,9 +868,8 @@ class IncidentNeutron(EqualityMixin):
heatr_evals = get_evaluations(kwargs["heatr"])
heatr_local_evals = get_evaluations(kwargs["heatr"] + "_local")
for ev, ev_local in zip(heatr_evals, heatr_local_evals):
temp = "{}K".format(round(ev.target["temperature"]))
for ev, ev_local, temp in zip(heatr_evals, heatr_local_evals, data.temperatures):
# Get total KERMA (originally from ACE file) and energy grid
kerma = data.reactions[301].xs[temp]
E = kerma.x

4
openmc/data/photon.py
View file

@ -13,7 +13,7 @@ from scipy.interpolate import CubicSpline
import openmc.checkvalue as cv
from openmc.mixin import EqualityMixin
from . import HDF5_VERSION
from . import HDF5_VERSION, HDF5_VERSION_MAJOR
from .ace import Table, get_metadata, get_table
from .data import ATOMIC_SYMBOL, EV_PER_MEV
from .endf import Evaluation, get_head_record, get_tab1_record, get_list_record
@ -143,6 +143,8 @@ class AtomicRelaxation(EqualityMixin):
Dictionary indicating the number of electrons in a subshell when neutral
(values) for given subshells (keys). The subshells should be given as
strings, e.g., 'K', 'L1', 'L2', etc.
subshells : list
List of subshells as strings, e.g. ``['K', 'L1', ...]``
transitions : pandas.DataFrame
Dictionary indicating allowed transitions and their probabilities
(values) for given subshells (keys). The subshells should be given as

2
openmc/data/reaction.py
View file

@ -555,7 +555,7 @@ def _get_activation_products(ev, rx):
Z, A = divmod(items[2], 1000)
excited_state = items[3]
# Get GND name for product
# Get GNDS name for product
symbol = ATOMIC_SYMBOL[Z]
if excited_state > 0:
name = '{}{}_e{}'.format(symbol, A, excited_state)

39
openmc/data/thermal.py
View file

@ -104,7 +104,7 @@ def get_thermal_name(name):
Returns
-------
str
GND-format thermal scattering name
GNDS-format thermal scattering name
"""
if name in _THERMAL_NAMES:
@ -396,7 +396,7 @@ class ThermalScattering(EqualityMixin):
Parameters
----------
name : str
Name of the material using GND convention, e.g. c_H_in_H2O
Name of the material using GNDS convention, e.g. c_H_in_H2O
atomic_weight_ratio : float
Atomic mass ratio of the target nuclide.
kTs : Iterable of float
@ -415,7 +415,7 @@ class ThermalScattering(EqualityMixin):
Inelastic scattering cross section derived in the incoherent
approximation
name : str
Name of the material using GND convention, e.g. c_H_in_H2O
Name of the material using GNDS convention, e.g. c_H_in_H2O
temperatures : Iterable of str
List of string representations the temperatures of the target nuclide
in the data set. The temperatures are strings of the temperature,
@ -491,7 +491,7 @@ class ThermalScattering(EqualityMixin):
ACE table to read from. If given as a string, it is assumed to be
the filename for the ACE file.
name : str
GND-conforming name of the material, e.g. c_H_in_H2O. If none is
GNDS-conforming name of the material, e.g. c_H_in_H2O. If none is
passed, the appropriate name is guessed based on the name of the ACE
table.
@ -596,7 +596,7 @@ class ThermalScattering(EqualityMixin):
ACE table to read from. If given as a string, it is assumed to be
the filename for the ACE file.
name : str
GND-conforming name of the material, e.g. c_H_in_H2O. If none is
GNDS-conforming name of the material, e.g. c_H_in_H2O. If none is
passed, the appropriate name is guessed based on the name of the ACE
table.
@ -670,12 +670,41 @@ class ThermalScattering(EqualityMixin):
mu_i = []
for j in range(n_energy_out[i]):
mu = ace.xss[idx + 4:idx + 4 + n_mu]
# The equiprobable angles produced by NJOY are not always
# sorted. This is problematic when the smearing algorithm
# is applied when sampling the angles. We sort the angles
# here, because they are equiprobable, so the order
# doesn't matter.
mu.sort()
p_mu = 1. / n_mu * np.ones(n_mu)
mu_ij = Discrete(mu, p_mu)
mu_ij.c = np.cumsum(p_mu)
mu_i.append(mu_ij)
idx += 3 + n_mu
# Check if the CDF for the outgoing energy distribution starts
# at 0. For NJOY and FRENDY evaluations, this is never the case,
# and can very rarely lead to negative energies when sampling
# the outgoing energy. From Eq. 7.6 of the ENDF manual, we can
# add an outgoing energy 0 eV that has a PDF of 0 (and of
# course, a CDF of 0 as well).
if eout_i.c[0] > 0.:
eout_i.x = np.insert(eout_i.x, 0, 0.)
eout_i.p = np.insert(eout_i.p, 0, 0.)
eout_i.c = np.insert(eout_i.c, 0, 0.)
# For this added outgoing energy (of 0 eV) we add a set of
# isotropic discrete angles.
dmu = 2. / n_mu
mu = np.linspace(-1. + 0.5*dmu, 1. - 0.5*dmu, n_mu)
p_mu = 1. / n_mu * np.ones(n_mu)
mu_0 = Discrete(mu, p_mu)
mu_0.c = np.cumsum(p_mu)
mu_i.insert(0, mu_0)
# We don't worry about renormalizing the outgoing energy PDF/CDF
# after this manipulation, because it never seems to be
# normalized to begin with (at least with NJOY).
energy_out.append(eout_i)
mu_out.append(mu_i)

10
openmc/deplete/abc.py
View file

@ -265,8 +265,8 @@ class ReactionRateHelper(ABC):
Ordering of reactions
"""
def divide_by_adens(self, number):
"""Normalize reaction rates by number of nuclides
def divide_by_atoms(self, number):
"""Normalize reaction rates by number of atoms
Acts on the current material examined by :meth:`get_material_rates`
@ -798,7 +798,7 @@ class Integrator(ABC):
t, self._i_res = self._get_start_data()
for i, (dt, source_rate) in enumerate(self):
if output:
if output and comm.rank == 0:
print(f"[openmc.deplete] t={t} s, dt={dt} s, source={source_rate}")
# Solve transport equation (or obtain result from restart)
@ -818,7 +818,7 @@ class Integrator(ABC):
conc = conc_list.pop()
StepResult.save(self.operator, conc_list, res_list, [t, t + dt],
source_rate, self._i_res + i, proc_time)
source_rate, self._i_res + i, proc_time)
t += dt
@ -826,7 +826,7 @@ class Integrator(ABC):
# source rate is passed to the transport operator (which knows to
# just return zero reaction rates without actually doing a transport
# solve)
if output and final_step:
if output and final_step and comm.rank == 0:
print(f"[openmc.deplete] t={t} (final operator evaluation)")
res_list = [self.operator(conc, source_rate if final_step else 0.0)]
StepResult.save(self.operator, [conc], res_list, [t, t],

39
openmc/deplete/chain.py
View file

@ -15,7 +15,7 @@ from numbers import Real, Integral
from warnings import warn
from openmc.checkvalue import check_type, check_greater_than
from openmc.data import gnd_name, zam, DataLibrary
from openmc.data import gnds_name, zam, DataLibrary
from openmc.exceptions import DataError
from .nuclide import FissionYieldDistribution
@ -135,14 +135,14 @@ def replace_missing(product, decay_data):
Parameters
----------
product : str
Name of product in GND format, e.g. 'Y86_m1'.
Name of product in GNDS format, e.g. 'Y86_m1'.
decay_data : dict
Dictionary of decay data
Returns
-------
product : str
Replacement for missing product in GND format.
Replacement for missing product in GNDS format.
"""
# Determine atomic number, mass number, and metastable state
@ -213,7 +213,7 @@ def replace_missing_fpy(actinide, fpy_data, decay_data):
# Check if metastable state has data (e.g., Am242m)
Z, A, m = zam(actinide)
if m == 0:
metastable = gnd_name(Z, A, 1)
metastable = gnds_name(Z, A, 1)
if metastable in fpy_data:
return metastable
@ -222,7 +222,7 @@ def replace_missing_fpy(actinide, fpy_data, decay_data):
while isotone in decay_data:
Z += 1
A += 1
isotone = gnd_name(Z, A, 0)
isotone = gnds_name(Z, A, 0)
if isotone in fpy_data:
return isotone
@ -231,7 +231,7 @@ def replace_missing_fpy(actinide, fpy_data, decay_data):
while isotone in decay_data:
Z -= 1
A -= 1
isotone = gnd_name(Z, A, 0)
isotone = gnds_name(Z, A, 0)
if isotone in fpy_data:
return isotone
@ -357,7 +357,7 @@ class Chain:
reactions = {}
for f in neutron_files:
evaluation = openmc.data.endf.Evaluation(f)
name = evaluation.gnd_name
name = evaluation.gnds_name
reactions[name] = {}
for mf, mt, nc, mod in evaluation.reaction_list:
if mf == 3:
@ -630,15 +630,27 @@ class Chain:
# Gain from radioactive decay
if nuc.n_decay_modes != 0:
for _, target, branching_ratio in nuc.decay_modes:
# Allow for total annihilation for debug purposes
if target is not None:
branch_val = branching_ratio * decay_constant
for decay_type, target, branching_ratio in nuc.decay_modes:
branch_val = branching_ratio * decay_constant
if branch_val != 0.0:
# Allow for total annihilation for debug purposes
if branch_val != 0.0:
if target is not None:
k = self.nuclide_dict[target]
matrix[k, i] += branch_val
# Produce alphas and protons from decay
if 'alpha' in decay_type:
k = self.nuclide_dict.get('He4')
if k is not None:
count = decay_type.count('alpha')
matrix[k, i] += count * branch_val
elif 'p' in decay_type:
k = self.nuclide_dict.get('H1')
if k is not None:
count = decay_type.count('p')
matrix[k, i] += count * branch_val
if nuc.name in rates.index_nuc:
# Extract all reactions for this nuclide in this cell
nuc_ind = rates.index_nuc[nuc.name]
@ -892,7 +904,7 @@ class Chain:
ground_target = grounds.get(parent_name)
if ground_target is None:
pz, pa, pm = zam(parent_name)
ground_target = gnd_name(pz, pa + 1, 0)
ground_target = gnds_name(pz, pa + 1, 0)
new_ratios[ground_target] = ground_br
parent.add_reaction(reaction, ground_target, rxn_Q, ground_br)
@ -1025,6 +1037,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

3
openmc/deplete/coupled_operator.py
View file

@ -460,7 +460,6 @@ class CoupledOperator(OpenMCOperator):
# Run OpenMC
openmc.lib.run()
openmc.lib.reset_timers()
# Extract results
rates = self._calculate_reaction_rates(source_rate)
@ -529,6 +528,8 @@ class CoupledOperator(OpenMCOperator):
"openmc_simulation_n{}.h5".format(step),
write_source=False)
openmc.lib.reset_timers()
def finalize(self):
"""Finalize a depletion simulation and release resources."""
if self.cleanup_when_done:

22
openmc/deplete/independent_operator.py
View file

@ -31,6 +31,9 @@ class IndependentOperator(OpenMCOperator):
passed to an integrator class, such as
:class:`openmc.deplete.CECMIntegrator`.
Note that passing an empty :class:`~openmc.deplete.MicroXS` instance to the
``micro_xs`` argument allows a decay-only calculation to be run.
.. versionadded:: 0.13.1
Parameters
@ -38,9 +41,11 @@ class IndependentOperator(OpenMCOperator):
materials : openmc.Materials
Materials to deplete.
micro_xs : MicroXS
One-group microscopic cross sections in [b] .
One-group microscopic cross sections in [b]. If the
:class:`~openmc.deplete.MicroXS` object is empty, a decay-only calculation will
be run.
chain_file : str
Path to the depletion chain XML file. Defaults to
Path to the depletion chain XML file. Defaults to
``openmc.config['chain_file']``.
keff : 2-tuple of float, optional
keff eigenvalue and uncertainty from transport calculation.
@ -150,7 +155,7 @@ class IndependentOperator(OpenMCOperator):
@classmethod
def from_nuclides(cls, volume, nuclides,
micro_xs,
chain_file,
chain_file=None,
nuc_units='atom/b-cm',
keff=None,
normalization_mode='fission-q',
@ -169,10 +174,13 @@ class IndependentOperator(OpenMCOperator):
Dictionary with nuclide names as keys and nuclide concentrations as
values.
micro_xs : MicroXS
One-group microscopic cross sections.
chain_file : str
Path to the depletion chain XML file.
nuc_units : {'atom/cm3', 'atom/b-cm'}
One-group microscopic cross sections in [b]. If the
:class:`~openmc.deplete.MicroXS` object is empty, a decay-only calculation
will be run.
chain_file : str, optional
Path to the depletion chain XML file. Defaults to
``openmc.config['chain_file']``.
nuc_units : {'atom/cm3', 'atom/b-cm'}, optional
Units for nuclide concentration.
keff : 2-tuple of float, optional
keff eigenvalue and uncertainty from transport calculation.

2
openmc/deplete/nuclide.py
View file

@ -82,7 +82,7 @@ class Nuclide:
Parameters
----------
name : str, optional
GND name of this nuclide, e.g. ``"He4"``, ``"Am242_m1"``
GNDS name of this nuclide, e.g. ``"He4"``, ``"Am242_m1"``
Attributes
----------

6
openmc/deplete/openmc_operator.py
View file

@ -541,10 +541,10 @@ class OpenMCOperator(TransportOperator):
self._normalization_helper.update(
tally_rates[:, fission_ind])
# Divide by total number and store
rates[i] = self._rate_helper.divide_by_adens(number)
# Divide by total number of atoms and store
rates[i] = self._rate_helper.divide_by_atoms(number)
# Scale reaction rates to obtain units of reactions/sec
# Scale reaction rates to obtain units of (reactions/sec)/atom
rates *= self._normalization_helper.factor(source_rate)
# Store new fission yields on the chain

14
openmc/deplete/results.py
View file

@ -1,7 +1,8 @@
import numbers
import bisect
import math
from typing import Iterable, Optional, Tuple, Union
import typing # required to prevent typing.Union namespace overwriting Union
from typing import Iterable, Optional, Tuple
from warnings import warn
import h5py
@ -95,7 +96,7 @@ class Results(list):
def get_atoms(
self,
mat: Union[Material, str],
mat: typing.Union[Material, str],
nuc: str,
nuc_units: str = "atoms",
time_units: str = "s"
@ -165,7 +166,7 @@ class Results(list):
def get_reaction_rate(
self,
mat: Union[Material, str],
mat: typing.Union[Material, str],
nuc: str,
rx: str
) -> Tuple[np.ndarray, np.ndarray]:
@ -374,7 +375,8 @@ class Results(list):
def export_to_materials(
self,
burnup_index: int,
nuc_with_data: Optional[Iterable[str]] = None
nuc_with_data: Optional[Iterable[str]] = None,
path: PathLike = 'materials.xml'
) -> Materials:
"""Return openmc.Materials object based on results at a given step
@ -393,6 +395,8 @@ class Results(list):
If not provided, nuclides from the cross_sections element of
materials.xml will be used. If that element is not present,
nuclides from openmc.config['cross_sections'] will be used.
path : PathLike
Path to materials XML file to read. Defaults to 'materials.xml'.
Returns
-------
@ -406,7 +410,7 @@ class Results(list):
# updated. If for some reason you have modified OpenMC to produce
# new materials as depletion takes place, this method will not
# work as expected and leave out that material.
mat_file = Materials.from_xml("materials.xml")
mat_file = Materials.from_xml(path)
# Only nuclides with valid transport data will be written to
# the new materials XML file. The precedence of nuclides to select

28
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
----------
@ -318,10 +318,11 @@ class StepResult:
chunks=(1, 1, n_mats, n_nuc_number),
dtype='float64')
handle.create_dataset("reaction rates", (1, n_stages, n_mats, n_nuc_rxn, n_rxn),
maxshape=(None, n_stages, n_mats, n_nuc_rxn, n_rxn),
chunks=(1, 1, n_mats, n_nuc_rxn, n_rxn),
dtype='float64')
if n_nuc_rxn > 0 and n_rxn > 0:
handle.create_dataset("reaction rates", (1, n_stages, n_mats, n_nuc_rxn, n_rxn),
maxshape=(None, n_stages, n_mats, n_nuc_rxn, n_rxn),
chunks=(1, 1, n_mats, n_nuc_rxn, n_rxn),
dtype='float64')
handle.create_dataset("eigenvalues", (1, n_stages, 2),
maxshape=(None, n_stages, 2), dtype='float64')
@ -358,7 +359,9 @@ class StepResult:
# Grab handles
number_dset = handle["/number"]
rxn_dset = handle["/reaction rates"]
has_reactions = ("reaction rates" in handle)
if has_reactions:
rxn_dset = handle["/reaction rates"]
eigenvalues_dset = handle["/eigenvalues"]
time_dset = handle["/time"]
source_rate_dset = handle["/source_rate"]
@ -375,9 +378,10 @@ class StepResult:
number_shape[0] = new_shape
number_dset.resize(number_shape)
rxn_shape = list(rxn_dset.shape)
rxn_shape[0] = new_shape
rxn_dset.resize(rxn_shape)
if has_reactions:
rxn_shape = list(rxn_dset.shape)
rxn_shape[0] = new_shape
rxn_dset.resize(rxn_shape)
eigenvalues_shape = list(eigenvalues_dset.shape)
eigenvalues_shape[0] = new_shape
@ -407,7 +411,8 @@ class StepResult:
high = max(inds)
for i in range(n_stages):
number_dset[index, i, low:high+1] = self.data[i]
rxn_dset[index, i, low:high+1] = self.rates[i]
if has_reactions:
rxn_dset[index, i, low:high+1] = self.rates[i]
if comm.rank == 0:
eigenvalues_dset[index, i] = self.k[i]
if comm.rank == 0:
@ -483,7 +488,8 @@ class StepResult:
for i in range(results.n_stages):
rate = ReactionRates(results.index_mat, rxn_nuc_to_ind, rxn_to_ind, True)
rate[:] = handle["/reaction rates"][step, i, :, :, :]
if "reaction rates" in handle:
rate[:] = handle["/reaction rates"][step, i, :, :, :]
results.rates.append(rate)
return results

69
openmc/executor.py
View file

@ -9,7 +9,7 @@ from .plots import _get_plot_image
def _process_CLI_arguments(volume=False, geometry_debug=False, particles=None,
plot=False, restart_file=None, threads=None,
tracks=False, event_based=None,
openmc_exec='openmc', mpi_args=None):
openmc_exec='openmc', mpi_args=None, path_input=None):
"""Converts user-readable flags in to command-line arguments to be run with
the OpenMC executable via subprocess.
@ -41,7 +41,10 @@ def _process_CLI_arguments(volume=False, geometry_debug=False, particles=None,
Path to OpenMC executable. Defaults to 'openmc'.
mpi_args : list of str, optional
MPI execute command and any additional MPI arguments to pass,
e.g. ['mpiexec', '-n', '8'].
e.g., ['mpiexec', '-n', '8'].
path_input : str or Pathlike
Path to a single XML file or a directory containing XML files for the
OpenMC executable to read.
.. versionadded:: 0.13.0
@ -82,6 +85,9 @@ def _process_CLI_arguments(volume=False, geometry_debug=False, particles=None,
if mpi_args is not None:
args = mpi_args + args
if path_input is not None:
args += [path_input]
return args
@ -118,7 +124,7 @@ def _run(args, output, cwd):
raise RuntimeError(error_msg)
def plot_geometry(output=True, openmc_exec='openmc', cwd='.'):
def plot_geometry(output=True, openmc_exec='openmc', cwd='.', path_input=None):
"""Run OpenMC in plotting mode
Parameters
@ -129,6 +135,11 @@ def plot_geometry(output=True, openmc_exec='openmc', cwd='.'):
Path to OpenMC executable
cwd : str, optional
Path to working directory to run in
path_input : str
Path to a single XML file or a directory containing XML files for the
OpenMC executable to read.
.. versionadded:: 0.13.3
Raises
------
@ -136,10 +147,13 @@ def plot_geometry(output=True, openmc_exec='openmc', cwd='.'):
If the `openmc` executable returns a non-zero status
"""
_run([openmc_exec, '-p'], output, cwd)
args = [openmc_exec, '-p']
if path_input is not None:
args += [path_input]
_run(args, output, cwd)
def plot_inline(plots, openmc_exec='openmc', cwd='.'):
def plot_inline(plots, openmc_exec='openmc', cwd='.', path_input=None):
"""Display plots inline in a Jupyter notebook.
.. versionchanged:: 0.13.0
@ -155,6 +169,11 @@ def plot_inline(plots, openmc_exec='openmc', cwd='.'):
Path to OpenMC executable
cwd : str, optional
Path to working directory to run in
path_input : str
Path to a single XML file or a directory containing XML files for the
OpenMC executable to read.
.. versionadded:: 0.13.3
Raises
------
@ -171,7 +190,7 @@ def plot_inline(plots, openmc_exec='openmc', cwd='.'):
openmc.Plots(plots).export_to_xml(cwd)
# Run OpenMC in geometry plotting mode
plot_geometry(False, openmc_exec, cwd)
plot_geometry(False, openmc_exec, cwd, path_input)
if plots is not None:
images = [_get_plot_image(p, cwd) for p in plots]
@ -179,7 +198,8 @@ def plot_inline(plots, openmc_exec='openmc', cwd='.'):
def calculate_volumes(threads=None, output=True, cwd='.',
openmc_exec='openmc', mpi_args=None):
openmc_exec='openmc', mpi_args=None,
path_input=None):
"""Run stochastic volume calculations in OpenMC.
This function runs OpenMC in stochastic volume calculation mode. To specify
@ -206,10 +226,14 @@ def calculate_volumes(threads=None, output=True, cwd='.',
Path to OpenMC executable. Defaults to 'openmc'.
mpi_args : list of str, optional
MPI execute command and any additional MPI arguments to pass,
e.g. ['mpiexec', '-n', '8'].
e.g., ['mpiexec', '-n', '8'].
cwd : str, optional
Path to working directory to run in. Defaults to the current working
directory.
path_input : str or Pathlike
Path to a single XML file or a directory containing XML files for the
OpenMC executable to read.
Raises
------
@ -223,14 +247,16 @@ def calculate_volumes(threads=None, output=True, cwd='.',
"""
args = _process_CLI_arguments(volume=True, threads=threads,
openmc_exec=openmc_exec, mpi_args=mpi_args)
openmc_exec=openmc_exec, mpi_args=mpi_args,
path_input=path_input)
_run(args, output, cwd)
def run(particles=None, threads=None, geometry_debug=False,
restart_file=None, tracks=False, output=True, cwd='.',
openmc_exec='openmc', mpi_args=None, event_based=False):
openmc_exec='openmc', mpi_args=None, event_based=False,
path_input=None):
"""Run an OpenMC simulation.
Parameters
@ -239,17 +265,17 @@ def run(particles=None, threads=None, geometry_debug=False,
Number of particles to simulate per generation.
threads : int, optional
Number of OpenMP threads. If OpenMC is compiled with OpenMP threading
enabled, the default is implementation-dependent but is usually equal
to the number of hardware threads available (or a value set by the
enabled, the default is implementation-dependent but is usually equal to
the number of hardware threads available (or a value set by the
:envvar:`OMP_NUM_THREADS` environment variable).
geometry_debug : bool, optional
Turn on geometry debugging during simulation. Defaults to False.
restart_file : str, optional
Path to restart file to use
tracks : bool, optional
Enables the writing of particles tracks. The number of particle
tracks written to tracks.h5 is limited to 1000 unless
Settings.max_tracks is set. Defaults to False.
Enables the writing of particles tracks. The number of particle tracks
written to tracks.h5 is limited to 1000 unless Settings.max_tracks is
set. Defaults to False.
output : bool
Capture OpenMC output from standard out
cwd : str, optional
@ -258,13 +284,19 @@ def run(particles=None, threads=None, geometry_debug=False,
openmc_exec : str, optional
Path to OpenMC executable. Defaults to 'openmc'.
mpi_args : list of str, optional
MPI execute command and any additional MPI arguments to pass,
e.g. ['mpiexec', '-n', '8'].
MPI execute command and any additional MPI arguments to pass, e.g.,
['mpiexec', '-n', '8'].
event_based : bool, optional
Turns on event-based parallelism, instead of default history-based
.. versionadded:: 0.12
path_input : str or Pathlike
Path to a single XML file or a directory containing XML files for the
OpenMC executable to read.
.. versionadded:: 0.13.3
Raises
------
RuntimeError
@ -275,6 +307,7 @@ def run(particles=None, threads=None, geometry_debug=False,
args = _process_CLI_arguments(
volume=False, geometry_debug=geometry_debug, particles=particles,
restart_file=restart_file, threads=threads, tracks=tracks,
event_based=event_based, openmc_exec=openmc_exec, mpi_args=mpi_args)
event_based=event_based, openmc_exec=openmc_exec, mpi_args=mpi_args,
path_input=path_input)
_run(args, output, cwd)

12
openmc/filter.py
View file

@ -533,7 +533,7 @@ class CellFromFilter(WithIDFilter):
expected_type = Cell
class CellbornFilter(WithIDFilter):
class CellBornFilter(WithIDFilter):
"""Bins tally events based on which cell the particle was born in.
Parameters
@ -557,6 +557,14 @@ class CellbornFilter(WithIDFilter):
expected_type = Cell
# Temporary alias for CellbornFilter
def CellbornFilter(*args, **kwargs):
warnings.warn('The name of "CellbornFilter" has changed to '
'"CellBornFilter". "CellbornFilter" will be '
'removed in the future.', FutureWarning)
return CellBornFilter(*args, **kwargs)
class CellInstanceFilter(Filter):
"""Bins tally events based on which cell instance a particle is in.
@ -2095,7 +2103,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):

242
openmc/geometry.py
View file

@ -1,12 +1,15 @@
import os
import typing
from collections import OrderedDict, defaultdict
from collections.abc import Iterable
from copy import deepcopy
from pathlib import Path
from xml.etree import ElementTree as ET
import warnings
import openmc
import openmc._xml as xml
from .checkvalue import check_type
from .checkvalue import check_type, check_less_than, check_greater_than, PathLike
class Geometry:
@ -25,12 +28,19 @@ class Geometry:
bounding_box : 2-tuple of numpy.array
Lower-left and upper-right coordinates of an axis-aligned bounding box
of the universe.
merge_surfaces : bool
Whether to remove redundant surfaces when the geometry is exported.
surface_precision : int
Number of decimal places to round to for comparing the coefficients of
surfaces for considering them topologically equivalent.
"""
def __init__(self, root=None):
self._root_universe = None
self._offsets = {}
self.merge_surfaces = False
self.surface_precision = 10
if root is not None:
if isinstance(root, openmc.UniverseBase):
self.root_universe = root
@ -48,11 +58,31 @@ class Geometry:
def bounding_box(self):
return self.root_universe.bounding_box
@property
def merge_surfaces(self):
return self._merge_surfaces
@property
def surface_precision(self):
return self._surface_precision
@root_universe.setter
def root_universe(self, root_universe):
check_type('root universe', root_universe, openmc.UniverseBase)
self._root_universe = root_universe
@merge_surfaces.setter
def merge_surfaces(self, merge_surfaces):
check_type('merge surfaces', merge_surfaces, bool)
self._merge_surfaces = merge_surfaces
@surface_precision.setter
def surface_precision(self, surface_precision):
check_type('surface precision', surface_precision, int)
check_less_than('surface_precision', surface_precision, 16)
check_greater_than('surface_precision', surface_precision, 0)
self._surface_precision = surface_precision
def add_volume_information(self, volume_calc):
"""Add volume information from a stochastic volume calculation.
@ -75,6 +105,39 @@ class Geometry:
if universe.id in volume_calc.volumes:
universe.add_volume_information(volume_calc)
def to_xml_element(self, remove_surfs=False):
"""Creates a 'geometry' element to be written to an XML file.
Parameters
----------
remove_surfs : bool
Whether or not to remove redundant surfaces from the geometry when
exporting
"""
# Find and remove redundant surfaces from the geometry
if remove_surfs:
warnings.warn("remove_surfs kwarg will be deprecated soon, please "
"set the Geometry.merge_surfaces attribute instead.")
self.merge_surfaces = True
if self.merge_surfaces:
self.remove_redundant_surfaces()
# Create XML representation
element = ET.Element("geometry")
self.root_universe.create_xml_subelement(element, memo=set())
# Sort the elements in the file
element[:] = sorted(element, key=lambda x: (
x.tag, int(x.get('id'))))
# Clean the indentation in the file to be user-readable
xml.clean_indentation(element)
xml.reorder_attributes(element) # TODO: Remove when support is Python 3.8+
return element
def export_to_xml(self, path='geometry.xml', remove_surfs=False):
"""Export geometry to an XML file.
@ -89,20 +152,7 @@ class Geometry:
.. versionadded:: 0.12
"""
# Find and remove redundant surfaces from the geometry
if remove_surfs:
self.remove_redundant_surfaces()
# Create XML representation
root_element = ET.Element("geometry")
self.root_universe.create_xml_subelement(root_element, memo=set())
# Sort the elements in the file
root_element[:] = sorted(root_element, key=lambda x: (
x.tag, int(x.get('id'))))
# Clean the indentation in the file to be user-readable
xml.clean_indentation(root_element)
root_element = self.to_xml_element(remove_surfs)
# Check if path is a directory
p = Path(path)
@ -110,18 +160,17 @@ class Geometry:
p /= 'geometry.xml'
# Write the XML Tree to the geometry.xml file
xml.reorder_attributes(root_element) # TODO: Remove when support is Python 3.8+
tree = ET.ElementTree(root_element)
tree.write(str(p), xml_declaration=True, encoding='utf-8')
@classmethod
def from_xml(cls, path='geometry.xml', materials=None):
"""Generate geometry from XML file
def from_xml_element(cls, elem, materials=None):
"""Generate geometry from an XML element
Parameters
----------
path : str, optional
Path to geometry XML file
elem : xml.etree.ElementTree.Element
XML element
materials : openmc.Materials or None
Materials used to assign to cells. If None, an attempt is made to
generate it from the materials.xml file.
@ -132,6 +181,11 @@ class Geometry:
Geometry object
"""
mats = dict()
if materials is not None:
mats.update({str(m.id): m for m in materials})
mats['void'] = None
# Helper function for keeping a cache of Universe instances
universes = {}
def get_universe(univ_id):
@ -140,13 +194,10 @@ class Geometry:
universes[univ_id] = univ
return universes[univ_id]
tree = ET.parse(path)
root = tree.getroot()
# Get surfaces
surfaces = {}
periodic = {}
for surface in root.findall('surface'):
for surface in elem.findall('surface'):
s = openmc.Surface.from_xml_element(surface)
surfaces[s.id] = s
@ -160,24 +211,24 @@ class Geometry:
surfaces[s1].periodic_surface = surfaces[s2]
# Add any DAGMC universes
for elem in root.findall('dagmc_universe'):
dag_univ = openmc.DAGMCUniverse.from_xml_element(elem)
for e in elem.findall('dagmc_universe'):
dag_univ = openmc.DAGMCUniverse.from_xml_element(e)
universes[dag_univ.id] = dag_univ
# Dictionary that maps each universe to a list of cells/lattices that
# contain it (needed to determine which universe is the root)
# contain it (needed to determine which universe is the elem)
child_of = defaultdict(list)
for elem in root.findall('lattice'):
lat = openmc.RectLattice.from_xml_element(elem, get_universe)
for e in elem.findall('lattice'):
lat = openmc.RectLattice.from_xml_element(e, get_universe)
universes[lat.id] = lat
if lat.outer is not None:
child_of[lat.outer].append(lat)
for u in lat.universes.ravel():
child_of[u].append(lat)
for elem in root.findall('hex_lattice'):
lat = openmc.HexLattice.from_xml_element(elem, get_universe)
for e in elem.findall('hex_lattice'):
lat = openmc.HexLattice.from_xml_element(e, get_universe)
universes[lat.id] = lat
if lat.outer is not None:
child_of[lat.outer].append(lat)
@ -191,15 +242,8 @@ class Geometry:
for u in ring:
child_of[u].append(lat)
# Create dictionary to easily look up materials
if materials is None:
filename = Path(path).parent / 'materials.xml'
materials = openmc.Materials.from_xml(str(filename))
mats = {str(m.id): m for m in materials}
mats['void'] = None
for elem in root.findall('cell'):
c = openmc.Cell.from_xml_element(elem, surfaces, mats, get_universe)
for e in elem.findall('cell'):
c = openmc.Cell.from_xml_element(e, surfaces, mats, get_universe)
if c.fill_type in ('universe', 'lattice'):
child_of[c.fill].append(c)
@ -211,6 +255,41 @@ class Geometry:
else:
raise ValueError('Error determining root universe.')
@classmethod
def from_xml(
cls,
path: PathLike = 'geometry.xml',
materials: typing.Optional[typing.Union[PathLike, 'openmc.Materials']] = 'materials.xml'
):
"""Generate geometry from XML file
Parameters
----------
path : PathLike, optional
Path to geometry XML file
materials : openmc.Materials or PathLike
Materials used to assign to cells. If PathLike, an attempt is made
to generate materials from the provided xml file.
Returns
-------
openmc.Geometry
Geometry object
"""
# Using str and os.Pathlike here to avoid error when using just the imported PathLike
# TypeError: Subscripted generics cannot be used with class and instance checks
check_type('materials', materials, (str, os.PathLike, openmc.Materials))
if isinstance(materials, (str, os.PathLike)):
materials = openmc.Materials.from_xml(materials)
tree = ET.parse(path)
root = tree.getroot()
return cls.from_xml_element(root, materials)
def find(self, point):
"""Find cells/universes/lattices which contain a given point
@ -396,28 +475,6 @@ class Geometry:
surfaces = cell.region.get_surfaces(surfaces)
return surfaces
def get_redundant_surfaces(self):
"""Return all of the topologically redundant surface IDs
.. versionadded:: 0.12
Returns
-------
dict
Dictionary whose keys are the ID of a redundant surface and whose
values are the topologically equivalent :class:`openmc.Surface`
that should replace it.
"""
tally = defaultdict(list)
for surf in self.get_all_surfaces().values():
coeffs = tuple(surf._coefficients[k] for k in surf._coeff_keys)
key = (surf._type,) + coeffs
tally[key].append(surf)
return {replace.id: keep
for keep, *redundant in tally.values()
for replace in redundant}
def _get_domains_by_name(self, name, case_sensitive, matching, domain_type):
if not case_sensitive:
name = name.lower()
@ -477,6 +534,27 @@ class Geometry:
"""
return self._get_domains_by_name(name, case_sensitive, matching, 'cell')
def get_surfaces_by_name(self, name, case_sensitive=False, matching=False):
"""Return a list of surfaces with matching names.
Parameters
----------
name : str
The name to search match
case_sensitive : bool
Whether to distinguish upper and lower case letters in each
surface's name (default is False)
matching : bool
Whether the names must match completely (default is False)
Returns
-------
list of openmc.Surface
Surfaces matching the queried name
"""
return self._get_domains_by_name(name, case_sensitive, matching, 'surface')
def get_cells_by_fill_name(self, name, case_sensitive=False, matching=False):
"""Return a list of cells with fills with matching names.
@ -565,16 +643,42 @@ class Geometry:
return self._get_domains_by_name(name, case_sensitive, matching, 'lattice')
def remove_redundant_surfaces(self):
"""Remove redundant surfaces from the geometry"""
"""Remove and return all of the redundant surfaces.
Uses surface_precision attribute of Geometry instance for rounding and
comparing surface coefficients.
.. versionadded:: 0.12
Returns
-------
redundant_surfaces
Dictionary whose keys are the ID of a redundant surface and whose
values are the topologically equivalent :class:`openmc.Surface`
that should replace it.
"""
# Get redundant surfaces
redundant_surfaces = self.get_redundant_surfaces()
redundancies = defaultdict(list)
for surf in self.get_all_surfaces().values():
coeffs = tuple(round(surf._coefficients[k],
self.surface_precision)
for k in surf._coeff_keys)
key = (surf._type,) + coeffs
redundancies[key].append(surf)
# Iterate through all cells contained in the geometry
for cell in self.get_all_cells().values():
# Recursively remove redundant surfaces from regions
if cell.region:
cell.region.remove_redundant_surfaces(redundant_surfaces)
redundant_surfaces = {replace.id: keep
for keep, *redundant in redundancies.values()
for replace in redundant}
if redundant_surfaces:
# Iterate through all cells contained in the geometry
for cell in self.get_all_cells().values():
# Recursively remove redundant surfaces from regions
if cell.region:
cell.region.remove_redundant_surfaces(redundant_surfaces)
return redundant_surfaces
def determine_paths(self, instances_only=False):
"""Determine paths through CSG tree for cells and materials.

6
openmc/lib/__init__.py
View file

@ -42,12 +42,18 @@ else:
def _dagmc_enabled():
return c_bool.in_dll(_dll, "DAGMC_ENABLED").value
def _ncrystal_enabled():
return c_bool.in_dll(_dll, "NCRYSTAL_ENABLED").value
def _coord_levels():
return c_int.in_dll(_dll, "n_coord_levels").value
def _libmesh_enabled():
return c_bool.in_dll(_dll, "LIBMESH_ENABLED").value
def _mcpl_enabled():
return c_bool.in_dll(_dll, "MCPL_ENABLED").value
from .error import *
from .core import *
from .nuclide import *

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()

267
openmc/material.py
View file

@ -6,7 +6,7 @@ from pathlib import Path
import re
import typing # imported separately as py3.8 requires typing.Iterable
import warnings
from typing import Optional, Union
from typing import Optional
from xml.etree import ElementTree as ET
import h5py
@ -92,12 +92,17 @@ class Material(IDManagerMixin):
fissionable_mass : float
Mass of fissionable nuclides in the material in [g]. Requires that the
:attr:`volume` attribute is set.
decay_photon_energy : openmc.stats.Univariate
decay_photon_energy : openmc.stats.Univariate or None
Energy distribution of photons emitted from decay of unstable nuclides
within the material. The integral of this distribution is the total
intensity of the photon source in [decay/sec].
within the material, or None if no photon source exists. The integral of
this distribution is the total intensity of the photon source in
[decay/sec].
.. versionadded:: 0.14.0
.. versionadded:: 0.13.2
ncrystal_cfg : str
NCrystal configuration string
.. versionadded:: 0.13.3
"""
@ -117,6 +122,7 @@ class Material(IDManagerMixin):
self._volume = None
self._atoms = {}
self._isotropic = []
self._ncrystal_cfg = None
# A list of tuples (nuclide, percent, percent type)
self._nuclides = []
@ -139,6 +145,9 @@ class Material(IDManagerMixin):
string += '{: <16}\n'.format('\tS(a,b) Tables')
if self._ncrystal_cfg:
string += '{: <16}=\t{}\n'.format('\tNCrystal conf', self._ncrystal_cfg)
for sab in self._sab:
string += '{: <16}=\t{}\n'.format('\tS(a,b)', sab)
@ -218,6 +227,10 @@ class Material(IDManagerMixin):
def volume(self):
return self._volume
@property
def ncrystal_cfg(self):
return self._ncrystal_cfg
@name.setter
def name(self, name: Optional[str]):
if name is not None:
@ -264,7 +277,7 @@ class Material(IDManagerMixin):
return density*self.volume
@property
def decay_photon_energy(self) -> Univariate:
def decay_photon_energy(self) -> Optional[Univariate]:
atoms = self.get_nuclide_atoms()
dists = []
probs = []
@ -273,7 +286,7 @@ class Material(IDManagerMixin):
if source_per_atom is not None:
dists.append(source_per_atom)
probs.append(num_atoms)
return openmc.data.combine_distributions(dists, probs)
return openmc.data.combine_distributions(dists, probs) if dists else None
@classmethod
def from_hdf5(cls, group: h5py.Group):
@ -330,6 +343,66 @@ class Material(IDManagerMixin):
return material
@classmethod
def from_ncrystal(cls, cfg, **kwargs):
"""Create material from NCrystal configuration string.
Density, temperature, and material composition, and (ultimately) thermal
neutron scattering will be automatically be provided by NCrystal based
on this string. The name and material_id parameters are simply passed on
to the Material constructor.
.. versionadded:: 0.13.3
Parameters
----------
cfg : str
NCrystal configuration string
**kwargs
Keyword arguments passed to :class:`openmc.Material`
Returns
-------
openmc.Material
Material instance
"""
import NCrystal
nc_mat = NCrystal.createInfo(cfg)
def openmc_natabund(Z):
#nc_mat.getFlattenedComposition might need natural abundancies.
#This call-back function is used so NCrystal can flatten composition
#using OpenMC's natural abundancies. In practice this function will
#only get invoked in the unlikely case where a material is specified
#by referring both to natural elements and specific isotopes of the
#same element.
elem_name = openmc.data.ATOMIC_SYMBOL[Z]
return [
(int(iso_name[len(elem_name):]), abund)
for iso_name, abund in openmc.data.isotopes(elem_name)
]
flat_compos = nc_mat.getFlattenedComposition(
preferNaturalElements=True, naturalAbundProvider=openmc_natabund)
# Create the Material
material = cls(temperature=nc_mat.getTemperature(), **kwargs)
for Z, A_vals in flat_compos:
elemname = openmc.data.ATOMIC_SYMBOL[Z]
for A, frac in A_vals:
if A:
material.add_nuclide(f'{elemname}{A}', frac)
else:
material.add_element(elemname, frac)
material.set_density('g/cm3', nc_mat.getDensity())
material._ncrystal_cfg = NCrystal.normaliseCfg(cfg)
return material
def add_volume_information(self, volume_calc):
"""Add volume information to a material.
@ -345,10 +418,10 @@ class Material(IDManagerMixin):
self._atoms = volume_calc.atoms[self.id]
else:
raise ValueError('No volume information found for material ID={}.'
.format(self.id))
.format(self.id))
else:
raise ValueError('No volume information found for material ID={}.'
.format(self.id))
.format(self.id))
def set_density(self, units: str, density: Optional[float] = None):
"""Set the density of the material
@ -404,6 +477,9 @@ class Material(IDManagerMixin):
'macroscopic data-set has already been added'.format(self._id)
raise ValueError(msg)
if self._ncrystal_cfg is not None:
raise ValueError("Cannot add nuclides to NCrystal material")
# If nuclide name doesn't look valid, give a warning
try:
Z, _, _ = openmc.data.zam(nuclide)
@ -457,11 +533,10 @@ class Material(IDManagerMixin):
params['percent_type'] = percent_type
## check if nuclide
# check if nuclide
if not component.isalpha():
self.add_nuclide(component, **params)
else: # is element
kwargs = params
else:
self.add_element(component, **params)
def remove_nuclide(self, nuclide: str):
@ -608,6 +683,9 @@ class Material(IDManagerMixin):
raise ValueError("Element name should be given by the "
"element's symbol or name, e.g., 'Zr', 'zirconium'")
if self._ncrystal_cfg is not None:
raise ValueError("Cannot add elements to NCrystal material")
# Allow for element identifier to be given as a symbol or name
if len(element) > 2:
el = element.lower()
@ -824,6 +902,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
@ -876,6 +956,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
@ -963,7 +1045,7 @@ class Material(IDManagerMixin):
Returns
-------
Union[dict, float]
typing.Union[dict, float]
If by_nuclide is True then a dictionary whose keys are nuclide
names and values are activity is returned. Otherwise the activity
of the material is returned as a float.
@ -986,6 +1068,49 @@ class Material(IDManagerMixin):
return activity if by_nuclide else sum(activity.values())
def get_decay_heat(self, units: str = 'W', by_nuclide: bool = False):
"""Returns the decay heat of the material or for each nuclide in the
material in units of [W], [W/g] or [W/cm3].
.. versionadded:: 0.13.3
Parameters
----------
units : {'W', 'W/g', 'W/cm3'}
Specifies the units of decay heat to return. Options include total
heat [W], specific [W/g] or volumetric heat [W/cm3].
Default is total heat [W].
by_nuclide : bool
Specifies if the decay heat should be returned for the material as a
whole or per nuclide. Default is False.
Returns
-------
Union[dict, float]
If `by_nuclide` is True then a dictionary whose keys are nuclide
names and values are decay heat is returned. Otherwise the decay heat
of the material is returned as a float.
"""
cv.check_value('units', units, {'W', 'W/g', 'W/cm3'})
cv.check_type('by_nuclide', by_nuclide, bool)
if units == 'W':
multiplier = self.volume
elif units == 'W/cm3':
multiplier = 1
elif units == 'W/g':
multiplier = 1.0 / self.get_mass_density()
decayheat = {}
for nuclide, atoms_per_bcm in self.get_nuclide_atom_densities().items():
decay_erg = openmc.data.decay_energy(nuclide)
inv_seconds = openmc.data.decay_constant(nuclide)
decay_erg *= openmc.data.JOULE_PER_EV
decayheat[nuclide] = inv_seconds * decay_erg * 1e24 * atoms_per_bcm * multiplier
return decayheat if by_nuclide else sum(decayheat.values())
def get_nuclide_atoms(self):
"""Return number of atoms of each nuclide in the material
@ -1134,6 +1259,14 @@ class Material(IDManagerMixin):
if self._volume:
element.set("volume", str(self._volume))
if self._ncrystal_cfg:
if self._sab:
raise ValueError("NCrystal materials are not compatible with S(a,b).")
if self._macroscopic is not None:
raise ValueError("NCrystal materials are not compatible with macroscopic cross sections.")
element.set("cfg", str(self._ncrystal_cfg))
# Create temperature XML subelement
if self.temperature is not None:
element.set("temperature", str(self.temperature))
@ -1401,6 +1534,56 @@ class Materials(cv.CheckedList):
for material in self:
material.make_isotropic_in_lab()
def _write_xml(self, file, header=True, level=0, spaces_per_level=2, trailing_indent=True):
"""Writes XML content of the materials to an open file handle.
Parameters
----------
file : IOTextWrapper
Open file handle to write content into.
header : bool
Whether or not to write the XML header
level : int
Indentation level of materials element
spaces_per_level : int
Number of spaces per indentation
trailing_indentation : bool
Whether or not to write a trailing indentation for the materials element
"""
indentation = level*spaces_per_level*' '
# Write the header and the opening tag for the root element.
if header:
file.write("<?xml version='1.0' encoding='utf-8'?>\n")
file.write(indentation+'<materials>\n')
# Write the <cross_sections> element.
if self.cross_sections is not None:
element = ET.Element('cross_sections')
element.text = str(self.cross_sections)
clean_indentation(element, level=level+1)
element.tail = element.tail.strip(' ')
file.write((level+1)*spaces_per_level*' ')
reorder_attributes(element) # TODO: Remove when support is Python 3.8+
ET.ElementTree(element).write(file, encoding='unicode')
# Write the <material> elements.
for material in sorted(self, key=lambda x: x.id):
element = material.to_xml_element()
clean_indentation(element, level=level+1)
element.tail = element.tail.strip(' ')
file.write((level+1)*spaces_per_level*' ')
reorder_attributes(element) # TODO: Remove when support is Python 3.8+
ET.ElementTree(element).write(file, encoding='unicode')
# Write the closing tag for the root element.
file.write(indentation+'</materials>\n')
# Write a trailing indentation for the next element
# at this level if needed
if trailing_indent:
file.write(indentation)
def export_to_xml(self, path: PathLike = 'materials.xml'):
"""Export material collection to an XML file.
@ -1420,32 +1603,34 @@ class Materials(cv.CheckedList):
# one go.
with open(str(p), 'w', encoding='utf-8',
errors='xmlcharrefreplace') as fh:
self._write_xml(fh)
# Write the header and the opening tag for the root element.
fh.write("<?xml version='1.0' encoding='utf-8'?>\n")
fh.write('<materials>\n')
@classmethod
def from_xml_element(cls, elem):
"""Generate materials collection from XML file
# Write the <cross_sections> element.
if self.cross_sections is not None:
element = ET.Element('cross_sections')
element.text = str(self.cross_sections)
clean_indentation(element, level=1)
element.tail = element.tail.strip(' ')
fh.write(' ')
reorder_attributes(element) # TODO: Remove when support is Python 3.8+
ET.ElementTree(element).write(fh, encoding='unicode')
Parameters
----------
elem : xml.etree.ElementTree.Element
XML element
# Write the <material> elements.
for material in sorted(self, key=lambda x: x.id):
element = material.to_xml_element()
clean_indentation(element, level=1)
element.tail = element.tail.strip(' ')
fh.write(' ')
reorder_attributes(element) # TODO: Remove when support is Python 3.8+
ET.ElementTree(element).write(fh, encoding='unicode')
Returns
-------
openmc.Materials
Materials collection
# Write the closing tag for the root element.
fh.write('</materials>\n')
"""
# Generate each material
materials = cls()
for material in elem.findall('material'):
materials.append(Material.from_xml_element(material))
# Check for cross sections settings
xs = elem.find('cross_sections')
if xs is not None:
materials.cross_sections = xs.text
return materials
@classmethod
def from_xml(cls, path: PathLike = 'materials.xml'):
@ -1465,14 +1650,4 @@ class Materials(cv.CheckedList):
tree = ET.parse(path)
root = tree.getroot()
# Generate each material
materials = cls()
for material in root.findall('material'):
materials.append(Material.from_xml_element(material))
# Check for cross sections settings
xs = tree.find('cross_sections')
if xs is not None:
materials.cross_sections = xs.text
return materials
return cls.from_xml_element(root)

Some files were not shown because too many files have changed in this diff Show more