mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-28 06:05:58 -04:00
Merge branch 'develop' into adding_typehints_to_match_docs_material.py
This commit is contained in:
Jonathan Shimwell
GitHub
commit
5985adbccb
18 changed files with 493 additions and 102 deletions
|
|
@ -15,7 +15,7 @@
|
|||
|
||||
# sudo docker run image_name:tag_name or ID with no tag sudo docker run ID number
|
||||
|
||||
FROM debian:bullseye-slim
|
||||
FROM debian:bullseye-slim AS dependencies
|
||||
|
||||
# By default this Dockerfile builds OpenMC without DAGMC and LIBMESH support
|
||||
ARG build_dagmc=off
|
||||
|
|
@ -172,6 +172,8 @@ RUN if [ "$build_libmesh" = "on" ]; then \
|
|||
&& rm -rf ${LIBMESH_INSTALL_DIR}/build ${LIBMESH_INSTALL_DIR}/libmesh ; \
|
||||
fi
|
||||
|
||||
FROM dependencies AS build
|
||||
|
||||
# clone and install openmc
|
||||
RUN mkdir -p ${HOME}/OpenMC && cd ${HOME}/OpenMC \
|
||||
&& git clone --shallow-submodules --recurse-submodules --single-branch -b ${openmc_branch} --depth=1 ${OPENMC_REPO} \
|
||||
|
|
@ -207,5 +209,7 @@ RUN mkdir -p ${HOME}/OpenMC && cd ${HOME}/OpenMC \
|
|||
&& cd ../openmc && pip install .[test,depletion-mpi] \
|
||||
&& python -c "import openmc"
|
||||
|
||||
FROM build AS release
|
||||
|
||||
# Download cross sections (NNDC and WMP) and ENDF data needed by test suite
|
||||
RUN ${HOME}/OpenMC/openmc/tools/ci/download-xs.sh
|
||||
|
|
|
|||
|
|
@ -61,6 +61,7 @@ Core Functions
|
|||
|
||||
atomic_mass
|
||||
atomic_weight
|
||||
decay_constant
|
||||
dose_coefficients
|
||||
gnd_name
|
||||
half_life
|
||||
|
|
|
|||
|
|
@ -36,7 +36,6 @@ public:
|
|||
int threshold;
|
||||
double n_electrons;
|
||||
double binding_energy;
|
||||
xt::xtensor<double, 1> cross_section;
|
||||
vector<Transition> transitions;
|
||||
};
|
||||
|
||||
|
|
@ -80,8 +79,11 @@ public:
|
|||
Tabulated1D coherent_anomalous_real_;
|
||||
Tabulated1D coherent_anomalous_imag_;
|
||||
|
||||
// Photoionization and atomic relaxation data
|
||||
// Photoionization and atomic relaxation data. Subshell cross sections are
|
||||
// stored separately to improve memory access pattern when calculating the
|
||||
// total cross section
|
||||
vector<ElectronSubshell> shells_;
|
||||
xt::xtensor<double, 2> cross_sections_;
|
||||
|
||||
// Compton profile data
|
||||
xt::xtensor<double, 2> profile_pdf_;
|
||||
|
|
|
|||
|
|
@ -214,8 +214,8 @@ class Cell(IDManagerMixin):
|
|||
self._atoms = self._fill.get_nuclide_atom_densities()
|
||||
|
||||
# Convert to total number of atoms
|
||||
for key, nuclide in self._atoms.items():
|
||||
atom = nuclide[1] * self._volume * 1.0e+24
|
||||
for key, atom_per_bcm in self._atoms.items():
|
||||
atom = atom_per_bcm * self._volume * 1.0e+24
|
||||
self._atoms[key] = atom
|
||||
|
||||
elif self.fill_type == 'distribmat':
|
||||
|
|
@ -224,12 +224,12 @@ class Cell(IDManagerMixin):
|
|||
partial_volume = self.volume / len(self.fill)
|
||||
self._atoms = OrderedDict()
|
||||
for mat in self.fill:
|
||||
for key, nuclide in mat.get_nuclide_atom_densities().items():
|
||||
for key, atom_per_bcm in mat.get_nuclide_atom_densities().items():
|
||||
# To account for overlap of nuclides between distribmat
|
||||
# we need to append new atoms to any existing value
|
||||
# hence it is necessary to ask for default.
|
||||
atom = self._atoms.setdefault(key, 0)
|
||||
atom += nuclide[1] * partial_volume * 1.0e+24
|
||||
atom += atom_per_bcm * partial_volume * 1.0e+24
|
||||
self._atoms[key] = atom
|
||||
|
||||
else:
|
||||
|
|
|
|||
|
|
@ -3,7 +3,7 @@ import json
|
|||
import os
|
||||
import re
|
||||
from pathlib import Path
|
||||
from math import sqrt
|
||||
from math import sqrt, log
|
||||
from warnings import warn
|
||||
|
||||
# Isotopic abundances from Meija J, Coplen T B, et al, "Isotopic compositions
|
||||
|
|
@ -202,7 +202,7 @@ _GND_NAME_RE = re.compile(r'([A-Zn][a-z]*)(\d+)((?:_[em]\d+)?)')
|
|||
|
||||
# Used in half_life function as a cache
|
||||
_HALF_LIFE = {}
|
||||
|
||||
_LOG_TWO = log(2.0)
|
||||
|
||||
def atomic_mass(isotope):
|
||||
"""Return atomic mass of isotope in atomic mass units.
|
||||
|
|
@ -283,6 +283,8 @@ def half_life(isotope):
|
|||
Half-life values are from the `ENDF/B-VIII.0 decay sublibrary
|
||||
<https://www.nndc.bnl.gov/endf-b8.0/download.html>`_.
|
||||
|
||||
.. versionadded:: 0.13.1
|
||||
|
||||
Parameters
|
||||
----------
|
||||
isotope : str
|
||||
|
|
@ -302,6 +304,35 @@ def half_life(isotope):
|
|||
|
||||
return _HALF_LIFE.get(isotope.lower())
|
||||
|
||||
|
||||
def decay_constant(isotope):
|
||||
"""Return decay constant of isotope in [s^-1]
|
||||
|
||||
Decay constants are based on half-life values from the
|
||||
:func:`~openmc.data.half_life` function. When the isotope is stable, a decay
|
||||
constant of zero is returned.
|
||||
|
||||
.. versionadded:: 0.13.1
|
||||
|
||||
Parameters
|
||||
----------
|
||||
isotope : str
|
||||
Name of isotope, e.g., 'Pu239'
|
||||
|
||||
Returns
|
||||
-------
|
||||
float
|
||||
Decay constant of isotope in [s^-1]
|
||||
|
||||
See also
|
||||
--------
|
||||
openmc.data.half_life
|
||||
|
||||
"""
|
||||
t = half_life(isotope)
|
||||
return _LOG_TWO / t if t else 0.0
|
||||
|
||||
|
||||
def water_density(temperature, pressure=0.1013):
|
||||
"""Return the density of liquid water at a given temperature and pressure.
|
||||
|
||||
|
|
|
|||
|
|
@ -465,11 +465,10 @@ class Decay(EqualityMixin):
|
|||
|
||||
@property
|
||||
def decay_constant(self):
|
||||
if hasattr(self.half_life, 'n'):
|
||||
return log(2.)/self.half_life
|
||||
else:
|
||||
mu, sigma = self.half_life
|
||||
return ufloat(log(2.)/mu, log(2.)/mu**2*sigma)
|
||||
if self.half_life.n == 0.0:
|
||||
name = self.nuclide['name']
|
||||
raise ValueError(f"{name} is listed as unstable but has a zero half-life.")
|
||||
return log(2.)/self.half_life
|
||||
|
||||
@property
|
||||
def decay_energy(self):
|
||||
|
|
|
|||
|
|
@ -528,9 +528,9 @@ class Operator(TransportOperator):
|
|||
"""
|
||||
mat_id = str(mat.id)
|
||||
|
||||
for nuclide, density in mat.get_nuclide_atom_densities().values():
|
||||
number = density * 1.0e24
|
||||
self.number.set_atom_density(mat_id, nuclide, number)
|
||||
for nuclide, atom_per_bcm in mat.get_nuclide_atom_densities().items():
|
||||
atom_per_cc = atom_per_bcm * 1.0e24
|
||||
self.number.set_atom_density(mat_id, nuclide, atom_per_cc)
|
||||
|
||||
def _set_number_from_results(self, mat, prev_res):
|
||||
"""Extracts material nuclides and number densities.
|
||||
|
|
@ -556,16 +556,15 @@ class Operator(TransportOperator):
|
|||
# Merge lists of nuclides, with the same order for every calculation
|
||||
geom_nuc_densities.update(depl_nuc)
|
||||
|
||||
for nuclide in geom_nuc_densities.keys():
|
||||
for nuclide, atom_per_bcm in geom_nuc_densities.items():
|
||||
if nuclide in depl_nuc:
|
||||
concentration = prev_res.get_atoms(mat_id, nuclide)[1][-1]
|
||||
volume = prev_res[-1].volume[mat_id]
|
||||
number = concentration / volume
|
||||
atom_per_cc = concentration / volume
|
||||
else:
|
||||
density = geom_nuc_densities[nuclide][1]
|
||||
number = density * 1.0e24
|
||||
atom_per_cc = atom_per_bcm * 1.0e24
|
||||
|
||||
self.number.set_atom_density(mat_id, nuclide, number)
|
||||
self.number.set_atom_density(mat_id, nuclide, atom_per_cc)
|
||||
|
||||
def initial_condition(self):
|
||||
"""Performs final setup and returns initial condition.
|
||||
|
|
|
|||
|
|
@ -4,7 +4,6 @@ from copy import deepcopy
|
|||
from numbers import Real
|
||||
from pathlib import Path
|
||||
import os
|
||||
import math
|
||||
import re
|
||||
import typing # imported separately as py3.8 requires typing.Iterable
|
||||
import warnings
|
||||
|
|
@ -90,6 +89,9 @@ class Material(IDManagerMixin):
|
|||
fissionable_mass : float
|
||||
Mass of fissionable nuclides in the material in [g]. Requires that the
|
||||
:attr:`volume` attribute is set.
|
||||
activity : float
|
||||
Activity of the material in [Bq]. Requires that the :attr:`volume`
|
||||
attribute is set.
|
||||
|
||||
"""
|
||||
|
||||
|
|
@ -146,20 +148,10 @@ class Material(IDManagerMixin):
|
|||
|
||||
return string
|
||||
|
||||
|
||||
@property
|
||||
def activity(self):
|
||||
"""Returns the total activity of the material in Becquerels."""
|
||||
|
||||
atoms_per_barn_cm = self.get_nuclide_atom_densities()
|
||||
total_activity = 0
|
||||
for key, value in atoms_per_barn_cm.items():
|
||||
half_life = openmc.data.half_life(key)
|
||||
if half_life:
|
||||
total_activity += value[1] / half_life
|
||||
total_activity *= math.log(2) * 1e24 * self.volume
|
||||
|
||||
return total_activity
|
||||
return sum(self.get_nuclide_activity().values())
|
||||
|
||||
@property
|
||||
def name(self):
|
||||
|
|
@ -263,10 +255,10 @@ class Material(IDManagerMixin):
|
|||
if self.volume is None:
|
||||
raise ValueError("Volume must be set in order to determine mass.")
|
||||
density = 0.0
|
||||
for nuc, atoms_per_cc in self.get_nuclide_atom_densities().values():
|
||||
for nuc, atoms_per_bcm in self.get_nuclide_atom_densities().items():
|
||||
Z = openmc.data.zam(nuc)[0]
|
||||
if Z >= 90:
|
||||
density += 1e24 * atoms_per_cc * openmc.data.atomic_mass(nuc) \
|
||||
density += 1e24 * atoms_per_bcm * openmc.data.atomic_mass(nuc) \
|
||||
/ openmc.data.AVOGADRO
|
||||
return density*self.volume
|
||||
|
||||
|
|
@ -794,11 +786,15 @@ class Material(IDManagerMixin):
|
|||
"""Returns all nuclides in the material and their atomic densities in
|
||||
units of atom/b-cm
|
||||
|
||||
.. versionchanged:: 0.13.1
|
||||
The values in the dictionary were changed from a tuple containing
|
||||
the nuclide name and the density to just the density.
|
||||
|
||||
Returns
|
||||
-------
|
||||
nuclides : dict
|
||||
Dictionary whose keys are nuclide names and values are tuples of
|
||||
(nuclide, density in atom/b-cm)
|
||||
Dictionary whose keys are nuclide names and values are densities in
|
||||
[atom/b-cm]
|
||||
|
||||
"""
|
||||
|
||||
|
|
@ -858,10 +854,46 @@ class Material(IDManagerMixin):
|
|||
|
||||
nuclides = OrderedDict()
|
||||
for n, nuc in enumerate(nucs):
|
||||
nuclides[nuc] = (nuc, nuc_densities[n])
|
||||
nuclides[nuc] = nuc_densities[n]
|
||||
|
||||
return nuclides
|
||||
|
||||
def get_nuclide_activity(self):
|
||||
"""Return activity in [Bq] for each nuclide in the material
|
||||
|
||||
.. versionadded:: 0.13.1
|
||||
|
||||
Returns
|
||||
-------
|
||||
dict
|
||||
Dictionary whose keys are nuclide names and values are activity in
|
||||
[Bq].
|
||||
"""
|
||||
activity = {}
|
||||
for nuclide, atoms in self.get_nuclide_atoms().items():
|
||||
inv_seconds = openmc.data.decay_constant(nuclide)
|
||||
activity[nuclide] = inv_seconds * atoms
|
||||
return activity
|
||||
|
||||
def get_nuclide_atoms(self):
|
||||
"""Return number of atoms of each nuclide in the material
|
||||
|
||||
.. versionadded:: 0.13.1
|
||||
|
||||
Returns
|
||||
-------
|
||||
dict
|
||||
Dictionary whose keys are nuclide names and values are number of
|
||||
atoms present in the material.
|
||||
|
||||
"""
|
||||
if self.volume is None:
|
||||
raise ValueError("Volume must be set in order to determine atoms.")
|
||||
atoms = {}
|
||||
for nuclide, atom_per_bcm in self.get_nuclide_atom_densities().items():
|
||||
atoms[nuclide] = 1.0e24 * atom_per_bcm * self.volume
|
||||
return atoms
|
||||
|
||||
def get_mass_density(self, nuclide: Optional[str] = None):
|
||||
"""Return mass density of one or all nuclides
|
||||
|
||||
|
|
@ -878,9 +910,9 @@ class Material(IDManagerMixin):
|
|||
|
||||
"""
|
||||
mass_density = 0.0
|
||||
for nuc, atoms_per_cc in self.get_nuclide_atom_densities().values():
|
||||
for nuc, atoms_per_bcm in self.get_nuclide_atom_densities().items():
|
||||
if nuclide is None or nuclide == nuc:
|
||||
density_i = 1e24 * atoms_per_cc * openmc.data.atomic_mass(nuc) \
|
||||
density_i = 1e24 * atoms_per_bcm * openmc.data.atomic_mass(nuc) \
|
||||
/ openmc.data.AVOGADRO
|
||||
mass_density += density_i
|
||||
return mass_density
|
||||
|
|
@ -1101,7 +1133,7 @@ class Material(IDManagerMixin):
|
|||
nuclides_per_cc = defaultdict(float)
|
||||
mass_per_cc = defaultdict(float)
|
||||
for mat, wgt in zip(materials, wgts):
|
||||
for nuc, atoms_per_bcm in mat.get_nuclide_atom_densities().values():
|
||||
for nuc, atoms_per_bcm in mat.get_nuclide_atom_densities().items():
|
||||
nuc_per_cc = wgt*1.e24*atoms_per_bcm
|
||||
nuclides_per_cc[nuc] += nuc_per_cc
|
||||
mass_per_cc[nuc] += nuc_per_cc*openmc.data.atomic_mass(nuc) / \
|
||||
|
|
|
|||
242
openmc/mesh.py
242
openmc/mesh.py
|
|
@ -193,6 +193,81 @@ class StructuredMesh(MeshBase):
|
|||
s1 = (slice(1, None),)*ndim + (slice(None),)
|
||||
return (vertices[s0] + vertices[s1]) / 2
|
||||
|
||||
@property
|
||||
def num_mesh_cells(self):
|
||||
return np.prod(self.dimension)
|
||||
|
||||
def write_data_to_vtk(self, points, filename, datasets, volume_normalization=True):
|
||||
"""Creates a VTK object of the mesh
|
||||
|
||||
Parameters
|
||||
----------
|
||||
points : list or np.array
|
||||
List of (X,Y,Z) tuples.
|
||||
filename : str
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
and values are the data sets.
|
||||
volume_normalization : bool, optional
|
||||
Whether or not to normalize the data by
|
||||
the volume of the mesh elements.
|
||||
|
||||
Raises
|
||||
------
|
||||
RuntimeError
|
||||
When the size of a dataset doesn't match the number of cells
|
||||
|
||||
Returns
|
||||
-------
|
||||
vtk.vtkStructuredGrid
|
||||
the VTK object
|
||||
"""
|
||||
|
||||
import vtk
|
||||
from vtk.util import numpy_support as nps
|
||||
|
||||
# check that the data sets are appropriately sized
|
||||
errmsg = "The size of the dataset {} should be equal to the number of cells"
|
||||
for label, dataset in datasets.items():
|
||||
if isinstance(dataset, np.ndarray):
|
||||
if not dataset.size == self.dimension[0] * self.dimension[1]* self.dimension[2]:
|
||||
raise RuntimeError(errmsg.format(label))
|
||||
else:
|
||||
if len(dataset) == self.dimension[0] * self.dimension[1]* self.dimension[2]:
|
||||
raise RuntimeError(errmsg.format(label))
|
||||
cv.check_type('label', label, str)
|
||||
|
||||
vtk_grid = vtk.vtkStructuredGrid()
|
||||
|
||||
vtk_grid.SetDimensions(*self.dimension)
|
||||
|
||||
vtkPts = vtk.vtkPoints()
|
||||
vtkPts.SetData(nps.numpy_to_vtk(points, deep=True))
|
||||
vtk_grid.SetPoints(vtkPts)
|
||||
|
||||
# create VTK arrays for each of
|
||||
# the data sets
|
||||
for label, dataset in datasets.items():
|
||||
dataset = np.asarray(dataset).flatten()
|
||||
|
||||
if volume_normalization:
|
||||
dataset /= self.volumes.flatten()
|
||||
|
||||
dataset_array = vtk.vtkDoubleArray()
|
||||
dataset_array.SetName(label)
|
||||
dataset_array.SetArray(nps.numpy_to_vtk(dataset),
|
||||
dataset.size,
|
||||
True)
|
||||
vtk_grid.GetCellData().AddArray(dataset_array)
|
||||
|
||||
# write the .vtk file
|
||||
writer = vtk.vtkStructuredGridWriter()
|
||||
writer.SetFileName(str(filename))
|
||||
writer.SetInputData(vtk_grid)
|
||||
writer.Write()
|
||||
|
||||
return vtk_grid
|
||||
|
||||
class RegularMesh(StructuredMesh):
|
||||
"""A regular Cartesian mesh in one, two, or three dimensions
|
||||
|
|
@ -273,10 +348,6 @@ class RegularMesh(StructuredMesh):
|
|||
dims = self._dimension
|
||||
return [(u - l) / d for u, l, d in zip(us, ls, dims)]
|
||||
|
||||
@property
|
||||
def num_mesh_cells(self):
|
||||
return np.prod(self._dimension)
|
||||
|
||||
@property
|
||||
def volumes(self):
|
||||
"""Return Volumes for every mesh cell
|
||||
|
|
@ -621,7 +692,42 @@ class RegularMesh(StructuredMesh):
|
|||
root_cell.fill = lattice
|
||||
|
||||
return root_cell, cells
|
||||
|
||||
def write_data_to_vtk(self, filename, datasets, volume_normalization=True):
|
||||
"""Creates a VTK object of the mesh
|
||||
|
||||
Parameters
|
||||
----------
|
||||
filename : str or pathlib.Path
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
and values are the data sets.
|
||||
volume_normalization : bool, optional
|
||||
Whether or not to normalize the data by
|
||||
the volume of the mesh elements.
|
||||
Defaults to True.
|
||||
|
||||
Returns
|
||||
-------
|
||||
vtk.vtkStructuredGrid
|
||||
the VTK object
|
||||
"""
|
||||
|
||||
ll, ur = self.lower_left, self.upper_right
|
||||
x_vals = np.linspace(ll[0], ur[0], num=self.dimension[0] + 1)
|
||||
y_vals = np.linspace(ll[1], ur[1], num=self.dimension[1] + 1)
|
||||
z_vals = np.linspace(ll[2], ur[2], num=self.dimension[2] + 1)
|
||||
|
||||
# create points
|
||||
pts_cartesian = np.array([[x, y, z] for z in z_vals for y in y_vals for x in x_vals])
|
||||
|
||||
return super().write_data_to_vtk(
|
||||
points=pts_cartesian,
|
||||
filename=filename,
|
||||
datasets=datasets,
|
||||
volume_normalization=volume_normalization
|
||||
)
|
||||
|
||||
def Mesh(*args, **kwargs):
|
||||
warnings.warn("Mesh has been renamed RegularMesh. Future versions of "
|
||||
|
|
@ -821,6 +927,36 @@ class RectilinearMesh(StructuredMesh):
|
|||
|
||||
return element
|
||||
|
||||
def write_data_to_vtk(self, filename, datasets, volume_normalization=True):
|
||||
"""Creates a VTK object of the mesh
|
||||
|
||||
Parameters
|
||||
----------
|
||||
filename : str or pathlib.Path
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
and values are the data sets.
|
||||
volume_normalization : bool, optional
|
||||
Whether or not to normalize the data by
|
||||
the volume of the mesh elements.
|
||||
Defaults to True.
|
||||
|
||||
Returns
|
||||
-------
|
||||
vtk.vtkStructuredGrid
|
||||
the VTK object
|
||||
"""
|
||||
# create points
|
||||
pts_cartesian = np.array([[x, y, z] for z in self.z_grid for y in self.y_grid for x in self.x_grid])
|
||||
|
||||
return super().write_data_to_vtk(
|
||||
points=pts_cartesian,
|
||||
filename=filename,
|
||||
datasets=datasets,
|
||||
volume_normalization=volume_normalization
|
||||
)
|
||||
|
||||
|
||||
class CylindricalMesh(StructuredMesh):
|
||||
"""A 3D cylindrical mesh
|
||||
|
|
@ -1012,6 +1148,46 @@ class CylindricalMesh(StructuredMesh):
|
|||
|
||||
return np.multiply.outer(np.outer(V_r, V_p), V_z)
|
||||
|
||||
def write_data_to_vtk(self, filename, datasets, volume_normalization=True):
|
||||
"""Creates a VTK object of the mesh
|
||||
|
||||
Parameters
|
||||
----------
|
||||
filename : str or pathlib.Path
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
and values are the data sets.
|
||||
volume_normalization : bool, optional
|
||||
Whether or not to normalize the data by
|
||||
the volume of the mesh elements.
|
||||
Defaults to True.
|
||||
|
||||
Returns
|
||||
-------
|
||||
vtk.vtkStructuredGrid
|
||||
the VTK object
|
||||
"""
|
||||
# create points
|
||||
pts_cylindrical = np.array(
|
||||
[
|
||||
[r, phi, z]
|
||||
for z in self.z_grid
|
||||
for phi in self.phi_grid
|
||||
for r in self.r_grid
|
||||
]
|
||||
)
|
||||
pts_cartesian = np.copy(pts_cylindrical)
|
||||
r, phi = pts_cylindrical[:, 0], pts_cylindrical[:, 1]
|
||||
pts_cartesian[:, 0] = r * np.cos(phi)
|
||||
pts_cartesian[:, 1] = r * np.sin(phi)
|
||||
|
||||
return super().write_data_to_vtk(
|
||||
points=pts_cartesian,
|
||||
filename=filename,
|
||||
datasets=datasets,
|
||||
volume_normalization=volume_normalization
|
||||
)
|
||||
|
||||
class SphericalMesh(StructuredMesh):
|
||||
"""A 3D spherical mesh
|
||||
|
|
@ -1204,6 +1380,49 @@ class SphericalMesh(StructuredMesh):
|
|||
|
||||
return np.multiply.outer(np.outer(V_r, V_t), V_p)
|
||||
|
||||
def write_data_to_vtk(self, filename, datasets, volume_normalization=True):
|
||||
"""Creates a VTK object of the mesh
|
||||
|
||||
Parameters
|
||||
----------
|
||||
filename : str or pathlib.Path
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
and values are the data sets.
|
||||
volume_normalization : bool, optional
|
||||
Whether or not to normalize the data by
|
||||
the volume of the mesh elements.
|
||||
Defaults to True.
|
||||
|
||||
Returns
|
||||
-------
|
||||
vtk.vtkStructuredGrid
|
||||
the VTK object
|
||||
"""
|
||||
|
||||
# create points
|
||||
pts_spherical = np.array(
|
||||
[
|
||||
[r, theta, phi]
|
||||
for phi in self.phi_grid
|
||||
for theta in self.theta_grid
|
||||
for r in self.r_grid
|
||||
]
|
||||
)
|
||||
pts_cartesian = np.copy(pts_spherical)
|
||||
r, theta, phi = pts_spherical[:, 0], pts_spherical[:, 1], pts_spherical[:, 2]
|
||||
pts_cartesian[:, 0] = r * np.sin(phi) * np.cos(theta)
|
||||
pts_cartesian[:, 1] = r * np.sin(phi) * np.sin(theta)
|
||||
pts_cartesian[:, 2] = r * np.cos(phi)
|
||||
|
||||
return super().write_data_to_vtk(
|
||||
points=pts_cartesian,
|
||||
filename=filename,
|
||||
datasets=datasets,
|
||||
volume_normalization=volume_normalization
|
||||
)
|
||||
|
||||
|
||||
class UnstructuredMesh(MeshBase):
|
||||
"""A 3D unstructured mesh
|
||||
|
|
@ -1328,6 +1547,7 @@ class UnstructuredMesh(MeshBase):
|
|||
"been loaded from a statepoint file.")
|
||||
return len(self._centroids)
|
||||
|
||||
|
||||
@centroids.setter
|
||||
def centroids(self, centroids):
|
||||
cv.check_type("Unstructured mesh centroids", centroids,
|
||||
|
|
@ -1373,7 +1593,7 @@ class UnstructuredMesh(MeshBase):
|
|||
|
||||
Parameters
|
||||
----------
|
||||
filename : str
|
||||
filename : str or pathlib.Path
|
||||
Name of the VTK file to write.
|
||||
datasets : dict
|
||||
Dictionary whose keys are the data labels
|
||||
|
|
@ -1381,6 +1601,11 @@ class UnstructuredMesh(MeshBase):
|
|||
volume_normalization : bool
|
||||
Whether or not to normalize the data by the
|
||||
volume of the mesh elements
|
||||
|
||||
Raises
|
||||
------
|
||||
RuntimeError
|
||||
when the size of a dataset doesn't match the number of cells
|
||||
"""
|
||||
|
||||
import vtk
|
||||
|
|
@ -1397,11 +1622,14 @@ class UnstructuredMesh(MeshBase):
|
|||
" mesh information from a statepoint file.")
|
||||
|
||||
# check that the data sets are appropriately sized
|
||||
errmsg = "The size of the dataset {} should be equal to the number of cells"
|
||||
for label, dataset in datasets.items():
|
||||
if isinstance(dataset, np.ndarray):
|
||||
assert dataset.size == self.n_elements
|
||||
if not dataset.size == self.dimension[0] * self.dimension[1]* self.dimension[2]:
|
||||
raise RuntimeError(errmsg.format(label))
|
||||
else:
|
||||
assert len(dataset) == self.n_elements
|
||||
if len(dataset) == self.dimension[0] * self.dimension[1]* self.dimension[2]:
|
||||
raise RuntimeError(errmsg.format(label))
|
||||
cv.check_type('label', label, str)
|
||||
|
||||
# create data arrays for the cells/points
|
||||
|
|
|
|||
|
|
@ -542,15 +542,11 @@ def _calculate_cexs_elem_mat(this, types, temperature=294.,
|
|||
|
||||
if isinstance(this, openmc.Material):
|
||||
# Expand elements in to nuclides with atomic densities
|
||||
nuclides = this.get_nuclide_atom_densities()
|
||||
# For ease of processing split out the nuclide and its fraction
|
||||
nuc_fractions = {nuclide[1][0]: nuclide[1][1]
|
||||
for nuclide in nuclides.items()}
|
||||
nuc_fractions = this.get_nuclide_atom_densities()
|
||||
# Create a dict of [nuclide name] = nuclide object to carry forward
|
||||
# with a common nuclides format between openmc.Material and
|
||||
# openmc.Element objects
|
||||
nuclides = {nuclide[1][0]: nuclide[1][0]
|
||||
for nuclide in nuclides.items()}
|
||||
nuclides = {nuclide: nuclide for nuclide in nuc_fractions}
|
||||
else:
|
||||
# Expand elements in to nuclides with atomic densities
|
||||
nuclides = this.expand(1., 'ao', enrichment=enrichment,
|
||||
|
|
@ -885,13 +881,13 @@ def _calculate_mgxs_elem_mat(this, types, library, orders=None,
|
|||
# Check to see if we have nuclides/elements or a macroscopic object
|
||||
if this._macroscopic is not None:
|
||||
# We have macroscopics
|
||||
nuclides = {this._macroscopic: (this._macroscopic, this.density)}
|
||||
nuclides = {this._macroscopic: this.density}
|
||||
else:
|
||||
# Expand elements in to nuclides with atomic densities
|
||||
nuclides = this.get_nuclide_atom_densities()
|
||||
|
||||
# For ease of processing split out nuc and nuc_density
|
||||
nuc_fraction = [nuclide[1][1] for nuclide in nuclides.items()]
|
||||
nuc_fraction = list(nuclides.values())
|
||||
else:
|
||||
T = temperature
|
||||
# Expand elements in to nuclides with atomic densities
|
||||
|
|
|
|||
|
|
@ -14,7 +14,9 @@
|
|||
#include "openmc/settings.h"
|
||||
|
||||
#include "xtensor/xbuilder.hpp"
|
||||
#include "xtensor/xmath.hpp"
|
||||
#include "xtensor/xoperation.hpp"
|
||||
#include "xtensor/xslice.hpp"
|
||||
#include "xtensor/xview.hpp"
|
||||
|
||||
#include <cmath>
|
||||
|
|
@ -44,6 +46,8 @@ vector<unique_ptr<PhotonInteraction>> elements;
|
|||
|
||||
PhotonInteraction::PhotonInteraction(hid_t group)
|
||||
{
|
||||
using namespace xt::placeholders;
|
||||
|
||||
// Set index of element in global vector
|
||||
index_ = data::elements.size();
|
||||
|
||||
|
|
@ -125,6 +129,7 @@ PhotonInteraction::PhotonInteraction(hid_t group)
|
|||
}
|
||||
|
||||
shells_.resize(n_shell);
|
||||
cross_sections_ = xt::zeros<double>({energy_.size(), n_shell});
|
||||
|
||||
// Create mapping from designator to index
|
||||
std::unordered_map<int, int> shell_map;
|
||||
|
|
@ -155,13 +160,15 @@ PhotonInteraction::PhotonInteraction(hid_t group)
|
|||
read_attribute(tgroup, "num_electrons", shell.n_electrons);
|
||||
|
||||
// Read subshell cross section
|
||||
xt::xtensor<double, 1> xs;
|
||||
dset = open_dataset(tgroup, "xs");
|
||||
read_attribute(dset, "threshold_idx", shell.threshold);
|
||||
close_dataset(dset);
|
||||
read_dataset(tgroup, "xs", shell.cross_section);
|
||||
read_dataset(tgroup, "xs", xs);
|
||||
|
||||
auto& xs = shell.cross_section;
|
||||
xs = xt::where(xs > 0.0, xt::log(xs), -500.0);
|
||||
auto cross_section =
|
||||
xt::view(cross_sections_, xt::range(shell.threshold, _), i);
|
||||
cross_section = xt::where(xs > 0, xt::log(xs), 0);
|
||||
|
||||
if (object_exists(tgroup, "transitions")) {
|
||||
// Determine dimensions of transitions
|
||||
|
|
@ -566,18 +573,13 @@ void PhotonInteraction::calculate_xs(Particle& p) const
|
|||
|
||||
// Calculate microscopic photoelectric cross section
|
||||
xs.photoelectric = 0.0;
|
||||
for (const auto& shell : shells_) {
|
||||
// Check threshold of reaction
|
||||
int i_start = shell.threshold;
|
||||
if (i_grid < i_start)
|
||||
continue;
|
||||
const auto& xs_lower = xt::row(cross_sections_, i_grid);
|
||||
const auto& xs_upper = xt::row(cross_sections_, i_grid + 1);
|
||||
|
||||
// Evaluation subshell photoionization cross section
|
||||
xs.photoelectric +=
|
||||
std::exp(shell.cross_section(i_grid - i_start) +
|
||||
f * (shell.cross_section(i_grid + 1 - i_start) -
|
||||
shell.cross_section(i_grid - i_start)));
|
||||
}
|
||||
for (int i = 0; i < xs_upper.size(); ++i)
|
||||
if (xs_lower(i) != 0)
|
||||
xs.photoelectric +=
|
||||
std::exp(xs_lower(i) + f * (xs_upper(i) - xs_lower(i)));
|
||||
|
||||
// Calculate microscopic pair production cross section
|
||||
xs.pair_production = std::exp(
|
||||
|
|
|
|||
|
|
@ -29,6 +29,7 @@
|
|||
|
||||
#include <algorithm> // for max, min, max_element
|
||||
#include <cmath> // for sqrt, exp, log, abs, copysign
|
||||
#include <xtensor/xview.hpp>
|
||||
|
||||
namespace openmc {
|
||||
|
||||
|
|
@ -344,25 +345,26 @@ void sample_photon_reaction(Particle& p)
|
|||
|
||||
// Photoelectric effect
|
||||
double prob_after = prob + micro.photoelectric;
|
||||
|
||||
if (prob_after > cutoff) {
|
||||
// Get grid index, interpolation factor, and bounding subshell
|
||||
// cross sections
|
||||
int i_grid = micro.index_grid;
|
||||
double f = micro.interp_factor;
|
||||
const auto& xs_lower = xt::row(element.cross_sections_, i_grid);
|
||||
const auto& xs_upper = xt::row(element.cross_sections_, i_grid + 1);
|
||||
|
||||
for (int i_shell = 0; i_shell < element.shells_.size(); ++i_shell) {
|
||||
const auto& shell {element.shells_[i_shell]};
|
||||
|
||||
// Get grid index and interpolation factor
|
||||
int i_grid = micro.index_grid;
|
||||
double f = micro.interp_factor;
|
||||
|
||||
// Check threshold of reaction
|
||||
int i_start = shell.threshold;
|
||||
if (i_grid < i_start)
|
||||
if (xs_lower(i_shell) == 0)
|
||||
continue;
|
||||
|
||||
// Evaluation subshell photoionization cross section
|
||||
double xs = std::exp(shell.cross_section(i_grid - i_start) +
|
||||
f * (shell.cross_section(i_grid + 1 - i_start) -
|
||||
shell.cross_section(i_grid - i_start)));
|
||||
// Evaluation subshell photoionization cross section
|
||||
prob += std::exp(
|
||||
xs_lower(i_shell) + f * (xs_upper(i_shell) - xs_lower(i_shell)));
|
||||
|
||||
prob += xs;
|
||||
if (prob > cutoff) {
|
||||
double E_electron = p.E() - shell.binding_energy;
|
||||
|
||||
|
|
|
|||
|
|
@ -1,6 +1,6 @@
|
|||
#!/usr/bin/env python
|
||||
|
||||
from collections.abc import Mapping
|
||||
from math import log
|
||||
import os
|
||||
from pathlib import Path
|
||||
|
||||
|
|
@ -126,3 +126,14 @@ def test_zam():
|
|||
assert openmc.data.zam('Am242_m10') == (95, 242, 10)
|
||||
with pytest.raises(ValueError):
|
||||
openmc.data.zam('garbage')
|
||||
|
||||
|
||||
def test_half_life():
|
||||
assert openmc.data.half_life('H2') is None
|
||||
assert openmc.data.half_life('U235') == pytest.approx(2.22102e16)
|
||||
assert openmc.data.half_life('Am242') == pytest.approx(57672.0)
|
||||
assert openmc.data.half_life('Am242_m1') == pytest.approx(4449622000.0)
|
||||
assert openmc.data.decay_constant('H2') == 0.0
|
||||
assert openmc.data.decay_constant('U235') == pytest.approx(log(2.0)/2.22102e16)
|
||||
assert openmc.data.decay_constant('Am242') == pytest.approx(log(2.0)/57672.0)
|
||||
assert openmc.data.decay_constant('Am242_m1') == pytest.approx(log(2.0)/4449622000.0)
|
||||
|
|
|
|||
|
|
@ -92,7 +92,7 @@ def test_activation(run_in_tmpdir, model, reaction_rate_mode, reaction_rate_opts
|
|||
|
||||
w = model.geometry.get_materials_by_name('tungsten')[0]
|
||||
atom_densities = w.get_nuclide_atom_densities()
|
||||
atom_per_cc = 1e24 * atom_densities['W186'][1] # Density in atom/cm^3
|
||||
atom_per_cc = 1e24 * atom_densities['W186'] # Density in atom/cm^3
|
||||
n0 = atom_per_cc * w.volume # Absolute number of atoms
|
||||
|
||||
# Pick a random irradiation time and then determine necessary source rate to
|
||||
|
|
|
|||
|
|
@ -100,7 +100,7 @@ def test_add_elements_by_formula():
|
|||
m.add_elements_from_formula('Li4SiO4')
|
||||
# checking the ratio of elements is 4:1:4 for Li:Si:O
|
||||
elem = defaultdict(float)
|
||||
for nuclide, adens in m.get_nuclide_atom_densities().values():
|
||||
for nuclide, adens in m.get_nuclide_atom_densities().items():
|
||||
if nuclide.startswith("Li"):
|
||||
elem["Li"] += adens
|
||||
if nuclide.startswith("Si"):
|
||||
|
|
@ -117,7 +117,7 @@ def test_add_elements_by_formula():
|
|||
'O16': 0.443386, 'O17': 0.000168}
|
||||
nuc_dens = m.get_nuclide_atom_densities()
|
||||
for nuclide in ref_dens:
|
||||
assert nuc_dens[nuclide][1] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
assert nuc_dens[nuclide] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
|
||||
# testing the correct nuclides are added to the Material when enriched
|
||||
m = openmc.Material()
|
||||
|
|
@ -129,7 +129,7 @@ def test_add_elements_by_formula():
|
|||
'O16': 0.443386, 'O17': 0.000168}
|
||||
nuc_dens = m.get_nuclide_atom_densities()
|
||||
for nuclide in ref_dens:
|
||||
assert nuc_dens[nuclide][1] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
assert nuc_dens[nuclide] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
|
||||
# testing the use of brackets
|
||||
m = openmc.Material()
|
||||
|
|
@ -137,7 +137,7 @@ def test_add_elements_by_formula():
|
|||
|
||||
# checking the ratio of elements is 2:2:6 for Mg:N:O
|
||||
elem = defaultdict(float)
|
||||
for nuclide, adens in m.get_nuclide_atom_densities().values():
|
||||
for nuclide, adens in m.get_nuclide_atom_densities().items():
|
||||
if nuclide.startswith("Mg"):
|
||||
elem["Mg"] += adens
|
||||
if nuclide.startswith("N"):
|
||||
|
|
@ -155,7 +155,7 @@ def test_add_elements_by_formula():
|
|||
'O16': 0.599772, 'O17': 0.000227}
|
||||
nuc_dens = m.get_nuclide_atom_densities()
|
||||
for nuclide in ref_dens:
|
||||
assert nuc_dens[nuclide][1] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
assert nuc_dens[nuclide] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
|
||||
# testing non integer multiplier results in a value error
|
||||
m = openmc.Material()
|
||||
|
|
@ -289,12 +289,21 @@ def test_get_nuclide_densities(uo2):
|
|||
|
||||
|
||||
def test_get_nuclide_atom_densities(uo2):
|
||||
nucs = uo2.get_nuclide_atom_densities()
|
||||
for nuc, density in nucs.values():
|
||||
for nuc, density in uo2.get_nuclide_atom_densities().items():
|
||||
assert nuc in ('U235', 'O16')
|
||||
assert density > 0
|
||||
|
||||
|
||||
def test_get_nuclide_atoms():
|
||||
mat = openmc.Material()
|
||||
mat.add_nuclide('Li6', 1.0)
|
||||
mat.set_density('atom/cm3', 3.26e20)
|
||||
mat.volume = 100.0
|
||||
|
||||
atoms = mat.get_nuclide_atoms()
|
||||
assert atoms['Li6'] == pytest.approx(mat.density * mat.volume)
|
||||
|
||||
|
||||
def test_mass():
|
||||
m = openmc.Material()
|
||||
m.add_nuclide('Zr90', 1.0, 'wo')
|
||||
|
|
@ -341,7 +350,7 @@ def test_borated_water():
|
|||
'O16':2.4672e-02}
|
||||
nuc_dens = m.get_nuclide_atom_densities()
|
||||
for nuclide in ref_dens:
|
||||
assert nuc_dens[nuclide][1] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
assert nuc_dens[nuclide] == pytest.approx(ref_dens[nuclide], 1e-2)
|
||||
assert m.id == 50
|
||||
|
||||
# Test the Celsius conversion.
|
||||
|
|
|
|||
76
tests/unit_tests/test_mesh_to_vtk.py
Normal file
76
tests/unit_tests/test_mesh_to_vtk.py
Normal file
|
|
@ -0,0 +1,76 @@
|
|||
import numpy as np
|
||||
from pathlib import Path
|
||||
import pytest
|
||||
|
||||
vtk = pytest.importorskip("vtk")
|
||||
from vtk.util import numpy_support as nps
|
||||
|
||||
import openmc
|
||||
|
||||
|
||||
regular_mesh = openmc.RegularMesh()
|
||||
regular_mesh.lower_left = (0, 0, 0)
|
||||
regular_mesh.upper_right = (1, 1, 1)
|
||||
regular_mesh.dimension = [30, 20, 10]
|
||||
|
||||
rectilinear_mesh = openmc.RectilinearMesh()
|
||||
rectilinear_mesh.x_grid = np.linspace(1, 2, num=30)
|
||||
rectilinear_mesh.y_grid = np.linspace(1, 2, num=30)
|
||||
rectilinear_mesh.z_grid = np.linspace(1, 2, num=30)
|
||||
|
||||
cylinder_mesh = openmc.CylindricalMesh()
|
||||
cylinder_mesh.r_grid = np.linspace(1, 2, num=30)
|
||||
cylinder_mesh.phi_grid = np.linspace(0, np.pi, num=50)
|
||||
cylinder_mesh.z_grid = np.linspace(0, 1, num=30)
|
||||
|
||||
spherical_mesh = openmc.SphericalMesh()
|
||||
spherical_mesh.r_grid = np.linspace(1, 2, num=30)
|
||||
spherical_mesh.phi_grid = np.linspace(0, np.pi, num=50)
|
||||
spherical_mesh.theta_grid = np.linspace(0, np.pi / 2, num=30)
|
||||
|
||||
@pytest.mark.parametrize("mesh", [cylinder_mesh, regular_mesh, rectilinear_mesh, spherical_mesh])
|
||||
def test_write_data_to_vtk(mesh, tmpdir):
|
||||
# BUILD
|
||||
filename = Path(tmpdir) / "out.vtk"
|
||||
|
||||
data = np.random.random(mesh.num_mesh_cells)
|
||||
|
||||
# RUN
|
||||
mesh.write_data_to_vtk(filename=filename, datasets={"label1": data, "label2": data})
|
||||
|
||||
# TEST
|
||||
assert filename.is_file()
|
||||
|
||||
# read file
|
||||
reader = vtk.vtkStructuredGridReader()
|
||||
reader.SetFileName(filename)
|
||||
reader.Update()
|
||||
|
||||
# check name of datasets
|
||||
vtk_grid = reader.GetOutput()
|
||||
array1 = vtk_grid.GetCellData().GetArray(0)
|
||||
array2 = vtk_grid.GetCellData().GetArray(1)
|
||||
|
||||
assert array1.GetName() == "label1"
|
||||
assert array2.GetName() == "label2"
|
||||
|
||||
# check size of datasets
|
||||
assert nps.vtk_to_numpy(array1).size == data.size
|
||||
assert nps.vtk_to_numpy(array2).size == data.size
|
||||
|
||||
@pytest.mark.parametrize("mesh", [cylinder_mesh, regular_mesh, rectilinear_mesh, spherical_mesh])
|
||||
def test_write_data_to_vtk_size_mismatch(mesh):
|
||||
"""Checks that an error is raised when the size of the dataset
|
||||
doesn't match the mesh number of cells
|
||||
|
||||
Parameters
|
||||
----------
|
||||
mesh : openmc.StructuredMesh
|
||||
The mesh to test
|
||||
"""
|
||||
right_size = mesh.num_mesh_cells
|
||||
data = np.random.random(right_size + 1)
|
||||
|
||||
expected_error_msg = "The size of the dataset label should be equal to the number of cells"
|
||||
with pytest.raises(RuntimeError, match=expected_error_msg):
|
||||
mesh.write_data_to_vtk(filename="out.vtk", datasets={"label": data})
|
||||
|
|
@ -391,16 +391,14 @@ def test_py_lib_attributes(run_in_tmpdir, pin_model_attributes, mpi_intracomm):
|
|||
assert openmc.lib.materials[1].get_density('atom/b-cm') == \
|
||||
pytest.approx(0.06891296988603757, abs=1e-13)
|
||||
mat_a_dens = np.sum(
|
||||
[v[1] for v in test_model.materials[0].
|
||||
get_nuclide_atom_densities().values()])
|
||||
list(test_model.materials[0].get_nuclide_atom_densities().values()))
|
||||
assert mat_a_dens == pytest.approx(0.06891296988603757, abs=1e-8)
|
||||
# Change the density
|
||||
test_model.update_densities(['UO2'], 2.)
|
||||
assert openmc.lib.materials[1].get_density('atom/b-cm') == \
|
||||
pytest.approx(2., abs=1e-13)
|
||||
mat_a_dens = np.sum(
|
||||
[v[1] for v in test_model.materials[0].
|
||||
get_nuclide_atom_densities().values()])
|
||||
list(test_model.materials[0].get_nuclide_atom_densities().values()))
|
||||
assert mat_a_dens == pytest.approx(2., abs=1e-8)
|
||||
|
||||
# Now lets do the cell temperature updates.
|
||||
|
|
@ -441,7 +439,7 @@ def test_deplete(run_in_tmpdir, pin_model_attributes, mpi_intracomm):
|
|||
test_model = openmc.Model(geom, mats, settings, tals, plots)
|
||||
|
||||
initial_mat = mats[0].clone()
|
||||
initial_u = initial_mat.get_nuclide_atom_densities()['U235'][1]
|
||||
initial_u = initial_mat.get_nuclide_atom_densities()['U235']
|
||||
|
||||
# Note that the chain file includes only U-235 fission to a stable Xe136 w/
|
||||
# a yield of 100%. Thus all the U235 we lose becomes Xe136
|
||||
|
|
@ -453,8 +451,8 @@ def test_deplete(run_in_tmpdir, pin_model_attributes, mpi_intracomm):
|
|||
operator_kwargs=op_kwargs,
|
||||
power=1., output=False)
|
||||
# Get the new Xe136 and U235 atom densities
|
||||
after_xe = mats[0].get_nuclide_atom_densities()['Xe136'][1]
|
||||
after_u = mats[0].get_nuclide_atom_densities()['U235'][1]
|
||||
after_xe = mats[0].get_nuclide_atom_densities()['Xe136']
|
||||
after_u = mats[0].get_nuclide_atom_densities()['U235']
|
||||
assert after_xe + after_u == pytest.approx(initial_u, abs=1e-15)
|
||||
assert test_model.is_initialized is False
|
||||
|
||||
|
|
@ -462,7 +460,7 @@ def test_deplete(run_in_tmpdir, pin_model_attributes, mpi_intracomm):
|
|||
mats[0].nuclides.clear()
|
||||
densities = initial_mat.get_nuclide_atom_densities()
|
||||
tot_density = 0.
|
||||
for nuc, density in densities.values():
|
||||
for nuc, density in densities.items():
|
||||
mats[0].add_nuclide(nuc, density)
|
||||
tot_density += density
|
||||
mats[0].set_density('atom/b-cm', tot_density)
|
||||
|
|
@ -473,8 +471,8 @@ def test_deplete(run_in_tmpdir, pin_model_attributes, mpi_intracomm):
|
|||
operator_kwargs=op_kwargs,
|
||||
power=1., output=False)
|
||||
# Get the new Xe136 and U235 atom densities
|
||||
after_lib_xe = mats[0].get_nuclide_atom_densities()['Xe136'][1]
|
||||
after_lib_u = mats[0].get_nuclide_atom_densities()['U235'][1]
|
||||
after_lib_xe = mats[0].get_nuclide_atom_densities()['Xe136']
|
||||
after_lib_u = mats[0].get_nuclide_atom_densities()['U235']
|
||||
assert after_lib_xe + after_lib_u == pytest.approx(initial_u, abs=1e-15)
|
||||
assert test_model.is_initialized is True
|
||||
|
||||
|
|
|
|||
|
|
@ -1,4 +1,5 @@
|
|||
import numpy as np
|
||||
|
||||
import openmc
|
||||
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue