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

Merge branch 'develop' of https://github.com/openmc-dev/openmc into cam_WP7_task4

Browse source
This commit is contained in:
Mikolaj Adam Kowalski 2020-03-11 16:04:54 +00:00
commit c2a769f1f2
23 changed files with 669 additions and 758 deletions
Download patch file Download diff file Expand all files Collapse all files

1
CMakeLists.txt
View file

@ -165,7 +165,6 @@ endif()
add_subdirectory(vendor/xtl)
set(xtl_DIR ${CMAKE_CURRENT_BINARY_DIR}/vendor/xtl)
add_subdirectory(vendor/xtensor)
target_link_libraries(xtensor INTERFACE xtl)
#===============================================================================
# GSL header-only library

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

@ -40,6 +40,12 @@ of an element, you specify the element itself. For example,
mat.add_element('C', 1.0)
This method can also accept case-insensitive element names such as
::
mat.add_element('aluminium', 1.0)
Internally, OpenMC stores data on the atomic masses and natural abundances of
all known isotopes and then uses this data to determine what isotopes should be
added to the material. When the material is later exported to XML for use by the
@ -55,6 +61,20 @@ following would add 3.2% enriched uranium to a material::
In addition to U235 and U238, concentrations of U234 and U236 will be present
and are determined through a correlation based on measured data.
It is also possible to perform enrichment of any element that is composed
of two naturally-occurring isotopes (e.g., Li or B) in terms of atomic percent.
To invoke this, provide the additional argument `enrichment_target` to
:meth:`Material.add_element`. For example the following would enrich B10
to 30ao%::
mat.add_element('B', 1.0, enrichment=30.0, enrichment_target='B10')
In order to enrich an isotope in terms of mass percent (wo%), provide the extra
argument `enrichment_type`. For example the following would enrich Li6 to 15wo%::
mat.add_element('Li', 1.0, enrichment=15.0, enrichment_target='Li6',
enrichment_type='wo')
Often, cross section libraries don't actually have all naturally-occurring
isotopes for a given element. For example, in ENDF/B-VII.1, cross section
evaluations are given for O16 and O17 but not for O18. If OpenMC is aware of
@ -135,6 +155,33 @@ attribute, e.g.,
:attr:`Material.temperature` or :attr:`Cell.temperature`
attributes, respectively.
-----------------
Material Mixtures
-----------------
In OpenMC it is possible to mix any number of materials to create a new material
with the correct nuclide composition and density. The
:meth:`Material.mix_materials` method takes a list of materials and
a list of their mixing fractions. Mixing fractions can be provided as atomic
fractions, weight fractions, or volume fractions. The fraction type
can be specified by passing 'ao', 'wo', or 'vo' as the third argument, respectively.
For example, assuming the required materials have already been defined, a MOX
material with 3% plutonium oxide by weight could be created using the following:
::
mox = openmc.Material.mix_materials([uo2, puo2], [0.97, 0.03], 'wo')
It should be noted that, if mixing fractions are specifed as atomic or weight
fractions, the supplied fractions should sum to one. If the fractions are specified
as volume fractions, and the sum of the fractions is less than one, then the remaining
fraction is set as void material.
.. warning:: Materials with :math:`S(\alpha,\beta)` thermal scattering data
cannot be used in :meth:`Material.mix_materials`. However, thermal
scattering data can be added to a material created by
:meth:`Material.mix_materials`.
--------------------
Material Collections
--------------------

217
examples/jupyter/pincell.ipynb
View file

File diff suppressed because one or more lines are too long

6
include/openmc/geometry_aux.h
View file

@ -86,11 +86,15 @@ void count_cell_instances(int32_t univ_indx);
//! Recursively search through universes and count universe instances.
//! \param search_univ The index of the universe to begin searching from.
//! \param target_univ_id The ID of the universe to be counted.
//! \param univ_count_memo Memoized counts that make this function faster for
//! large systems. The first call to this function for each target_univ_id
//! should start with an empty memo.
//! \return The number of instances of target_univ_id in the geometry tree under
//! search_univ.
//==============================================================================
int count_universe_instances(int32_t search_univ, int32_t target_univ_id);
int count_universe_instances(int32_t search_univ, int32_t target_univ_id,
std::unordered_map<int32_t, int32_t>& univ_count_memo);
//==============================================================================
//! Build a character array representing the path to a distribcell instance.

3
include/openmc/lattice.h
View file

@ -77,7 +77,8 @@ public:
{offsets_.resize(n_maps * universes_.size(), C_NONE);}
//! Populate the distribcell offset tables.
int32_t fill_offset_table(int32_t offset, int32_t target_univ_id, int map);
int32_t fill_offset_table(int32_t offset, int32_t target_univ_id, int map,
std::unordered_map<int32_t, int32_t>& univ_count_memo);
//! \brief Check lattice indices.
//! \param i_xyz[3] The indices for a lattice tile.

5
include/openmc/mgxs.h
View file

@ -40,8 +40,8 @@ class Mgxs {
xt::xtensor<double, 1> kTs; // temperature in eV (k * T)
AngleDistributionType scatter_format; // flag for if this is legendre, histogram, or tabular
int num_delayed_groups; // number of delayed neutron groups
int num_groups; // number of energy groups
int num_delayed_groups; // number of delayed neutron groups
std::vector<XsData> xs; // Cross section data
// MGXS Incoming Flux Angular grid information
bool is_isotropic; // used to skip search for angle indices if isotropic
@ -185,6 +185,9 @@ class Mgxs {
//! @param u Incoming particle direction.
void
set_angle_index(Direction u);
//! \brief Provide const access to list of XsData held by this
const std::vector<XsData>& get_xsdata() const { return xs; }
};
} // namespace openmc

4
include/openmc/scattdata.h
View file

@ -43,8 +43,8 @@ class ScattData {
double_2dvec energy; // Normalized p0 matrix for sampling Eout
double_2dvec mult; // nu-scatter multiplication (nu-scatt/scatt)
double_3dvec dist; // Angular distribution
xt::xtensor<double, 1> gmin; // minimum outgoing group
xt::xtensor<double, 1> gmax; // maximum outgoing group
xt::xtensor<int, 1> gmin; // minimum outgoing group
xt::xtensor<int, 1> gmax; // maximum outgoing group
xt::xtensor<double, 1> scattxs; // Isotropic Sigma_{s,g_{in}}
//! \brief Calculates the value of normalized f(mu).

44
openmc/data/data.py
View file

@ -10,6 +10,7 @@ from warnings import warn
# pp. 293-306 (2013). The "representative isotopic abundance" values from
# column 9 are used except where an interval is given, in which case the
# "best measurement" is used.
# Note that the abundances are given as atomic fractions!
NATURAL_ABUNDANCE = {
'H1': 0.99984426, 'H2': 0.00015574, 'He3': 0.000002,
'He4': 0.999998, 'Li6': 0.07589, 'Li7': 0.92411,
@ -110,6 +111,49 @@ NATURAL_ABUNDANCE = {
'U238': 0.992742
}
# Dictionary to give element symbols from IUPAC names
# (and some common mispellings)
ELEMENT_SYMBOL = {'neutron': 'n', 'hydrogen': 'H', 'helium': 'He',
'lithium': 'Li', 'beryllium': 'Be', 'boron': 'B',
'carbon': 'C', 'nitrogen': 'N', 'oxygen': 'O', 'fluorine': 'F',
'neon': 'Ne', 'sodium': 'Na', 'magnesium': 'Mg',
'aluminium': 'Al', 'aluminum': 'Al', 'silicon': 'Si',
'phosphorus': 'P', 'sulfur': 'S', 'sulphur': 'S',
'chlorine': 'Cl', 'argon': 'Ar', 'potassium': 'K',
'calcium': 'Ca', 'scandium': 'Sc', 'titanium': 'Ti',
'vanadium': 'V', 'chromium': 'Cr', 'manganese': 'Mn',
'iron': 'Fe', 'cobalt': 'Co', 'nickel': 'Ni', 'copper': 'Cu',
'zinc': 'Zn', 'gallium': 'Ga', 'germanium': 'Ge',
'arsenic': 'As', 'selenium': 'Se', 'bromine': 'Br',
'krypton': 'Kr', 'rubidium': 'Rb', 'strontium': 'Sr',
'yttrium': 'Y', 'zirconium': 'Zr', 'niobium': 'Nb',
'molybdenum': 'Mo', 'technetium': 'Tc', 'ruthenium': 'Ru',
'rhodium': 'Rh', 'palladium': 'Pd', 'silver': 'Ag',
'cadmium': 'Cd', 'indium': 'In', 'tin': 'Sn', 'antimony': 'Sb',
'tellurium': 'Te', 'iodine': 'I', 'xenon': 'Xe',
'caesium': 'Cs', 'cesium': 'Cs', 'barium': 'Ba',
'lanthanum': 'La', 'cerium': 'Ce', 'praseodymium': 'Pr',
'neodymium': 'Nd', 'promethium': 'Pm', 'samarium': 'Sm',
'europium': 'Eu', 'gadolinium': 'Gd', 'terbium': 'Tb',
'dysprosium': 'Dy', 'holmium': 'Ho', 'erbium': 'Er',
'thulium': 'Tm', 'ytterbium': 'Yb', 'lutetium': 'Lu',
'hafnium': 'Hf', 'tantalum': 'Ta', 'tungsten': 'W',
'wolfram': 'W', 'rhenium': 'Re', 'osmium': 'Os',
'iridium': 'Ir', 'platinum': 'Pt', 'gold': 'Au',
'mercury': 'Hg', 'thallium': 'Tl', 'lead': 'Pb',
'bismuth': 'Bi', 'polonium': 'Po', 'astatine': 'At',
'radon': 'Rn', 'francium': 'Fr', 'radium': 'Ra',
'actinium': 'Ac', 'thorium': 'Th', 'protactinium': 'Pa',
'uranium': 'U', 'neptunium': 'Np', 'plutonium': 'Pu',
'americium': 'Am', 'curium': 'Cm', 'berkelium': 'Bk',
'californium': 'Cf', 'einsteinium': 'Es', 'fermium': 'Fm',
'mendelevium': 'Md', 'nobelium': 'No', 'lawrencium': 'Lr',
'rutherfordium': 'Rf', 'dubnium': 'Db', 'seaborgium': 'Sg',
'bohrium': 'Bh', 'hassium': 'Hs', 'meitnerium': 'Mt',
'darmstadtium': 'Ds', 'roentgenium': 'Rg', 'copernicium': 'Cn',
'nihonium': 'Nh', 'flerovium': 'Fl', 'moscovium': 'Mc',
'livermorium': 'Lv', 'tennessine': 'Ts', 'oganesson': 'Og'}
ATOMIC_SYMBOL = {0: 'n', 1: 'H', 2: 'He', 3: 'Li', 4: 'Be', 5: 'B', 6: 'C',
7: 'N', 8: 'O', 9: 'F', 10: 'Ne', 11: 'Na', 12: 'Mg', 13: 'Al',
14: 'Si', 15: 'P', 16: 'S', 17: 'Cl', 18: 'Ar', 19: 'K',

135
openmc/element.py
View file

@ -4,6 +4,7 @@ import os
from xml.etree import ElementTree as ET
import openmc.checkvalue as cv
from numbers import Real
from openmc.data import NATURAL_ABUNDANCE, atomic_mass
@ -35,6 +36,7 @@ class Element(str):
return self
def expand(self, percent, percent_type, enrichment=None,
enrichment_target=None, enrichment_type=None,
cross_sections=None):
"""Expand natural element into its naturally-occurring isotopes.
@ -52,9 +54,15 @@ class Element(str):
percent_type : {'ao', 'wo'}
'ao' for atom percent and 'wo' for weight percent
enrichment : float, optional
Enrichment for U235 in weight percent. For example, input 4.95 for
4.95 weight percent enriched U. Default is None
(natural composition).
Enrichment of an enrichment_taget nuclide in percent (ao or wo).
If enrichment_taget is not supplied then it is enrichment for U235
in weight percent. For example, input 4.95 for 4.95 weight percent
enriched U. Default is None (natural composition).
enrichment_target: str, optional
Single nuclide name to enrich from a natural composition (e.g., 'O16')
enrichment_type: {'ao', 'wo'}, optional
'ao' for enrichment as atom percent and 'wo' for weight percent.
Default is: 'ao' for two-isotope enrichment; 'wo' for U enrichment
cross_sections : str, optional
Location of cross_sections.xml file. Default is None.
@ -65,16 +73,47 @@ class Element(str):
is a tuple consisting of a nuclide string, the atom/weight percent,
and the string 'ao' or 'wo'.
Raises
------
ValueError
No data is available for any of natural isotopes of the element
ValueError
If only some natural isotopes are available in the cross-section data
library and the element is not O, W, or Ta
ValueError
If a non-naturally-occurring isotope is requested
ValueError
If enrichment is requested of an element with more than two
naturally-occurring isotopes.
ValueError
If enrichment procedure for Uranium is used when element is not
Uranium.
ValueError
Uranium enrichment is requested with enrichment_type=='ao'
Notes
-----
When the `enrichment` argument is specified, a correlation from
`ORNL/CSD/TM-244 <https://doi.org/10.2172/5561567>`_ is used to
calculate the weight fractions of U234, U235, U236, and U238. Namely,
the weight fraction of U234 and U236 are taken to be 0.89% and 0.46%,
respectively, of the U235 weight fraction. The remainder of the isotopic
weight is assigned to U238.
respectively, of the U235 weight fraction. The remainder of the
isotopic weight is assigned to U238.
When the `enrichment` argument is specified with `enrichment_target`, a
general enrichment procedure is used for elements composed of exactly
two naturally-occurring isotopes. `enrichment` is interpreted as atom
percent by default but can be controlled by the `enrichment_type`
argument.
"""
# Check input
if enrichment_type is not None:
cv.check_value('enrichment_type', enrichment_type, {'ao', 'wo'})
if enrichment is not None:
cv.check_less_than('enrichment', enrichment, 100.0, equality=True)
cv.check_greater_than('enrichment', enrichment, 0., equality=True)
# Get the nuclides present in nature
natural_nuclides = set()
@ -110,8 +149,8 @@ class Element(str):
mutual_nuclides = sorted(list(mutual_nuclides))
absent_nuclides = sorted(list(absent_nuclides))
# If all natural nuclides are present in the library, expand element
# using all natural nuclides
# If all natural nuclides are present in the library,
# expand element using all natural nuclides
if len(absent_nuclides) == 0:
for nuclide in mutual_nuclides:
abundances[nuclide] = NATURAL_ABUNDANCE[nuclide]
@ -164,7 +203,20 @@ class Element(str):
abundances[nuclide] = NATURAL_ABUNDANCE[nuclide]
# Modify mole fractions if enrichment provided
if enrichment is not None:
# Old treatment for Uranium
if enrichment is not None and enrichment_target is None:
# Check that the element is Uranium
if self.name != 'U':
msg = ('Enrichment procedure for Uranium was requested, '
'but the isotope is {} not U'.format(self))
raise ValueError(msg)
# Check that enrichment_type is not 'ao'
if enrichment_type == 'ao':
msg = ('Enrichment procedure for Uranium requires that '
'enrichment value is provided as wo%.')
raise ValueError(msg)
# Calculate the mass fractions of isotopes
abundances['U234'] = 0.0089 * enrichment
@ -181,6 +233,73 @@ class Element(str):
for nuclide in abundances.keys():
abundances[nuclide] /= sum_abundances
# Modify mole fractions if enrichment provided
# New treatment for arbitrary element
elif enrichment is not None and enrichment_target is not None:
# Provide more informative error message for U235
if enrichment_target == 'U235':
msg = ("There is a special procedure for enrichment of U235 "
"in U. To invoke it, the arguments 'enrichment_target'"
"and 'enrichment_type' should be omitted. Provide "
"a value only for 'enrichment' in weight percent.")
raise ValueError(msg)
# Check if it is two-isotope mixture
if len(abundances) != 2:
msg = ('Element {} does not consist of two naturally-occurring '
'isotopes. Please enter isotopic abundances manually.'
.format(self))
raise ValueError(msg)
# Check if the target nuclide is present in the mixture
if enrichment_target not in abundances:
msg = ('The target nuclide {} is not one of the naturally-occurring '
'isotopes ({})'.format(enrichment_target, list(abundances)))
raise ValueError(msg)
# If weight percent enrichment is requested convert to mass fractions
if enrichment_type == 'wo':
# Convert the atomic abundances to weight fractions
# Compute the element atomic mass
element_am = sum(atomic_mass(nuc)*abundances[nuc] for nuc in abundances)
# Convert Molar Fractions to mass fractions
for nuclide in abundances:
abundances[nuclide] *= atomic_mass(nuclide) / element_am
# Normalize to one
sum_abundances = sum(abundances.values())
for nuclide in abundances:
abundances[nuclide] /= sum_abundances
# Enrich the mixture
# The procedure is more generic that it needs to be. It allows
# to enrich mixtures of more then 2 isotopes, keeping the ratios
# of non-enriched nuclides the same as in natural composition
# Get fraction of non-enriched isotopes in nat. composition
non_enriched = 1.0 - abundances[enrichment_target]
tail_fraction = 1.0 - enrichment / 100.0
# Enrich all nuclides
# Do bogus operation for enrichment target but overwrite immediatly
# to avoid if statement in the loop
for nuclide, fraction in abundances.items():
abundances[nuclide] = tail_fraction * fraction / non_enriched
abundances[enrichment_target] = enrichment / 100.0
# Convert back to atomic fractions if requested
if enrichment_type == 'wo':
# Convert the mass fractions to mole fractions
for nuclide in abundances:
abundances[nuclide] /= atomic_mass(nuclide)
# Normalize the mole fractions to one
sum_abundances = sum(abundances.values())
for nuclide in abundances:
abundances[nuclide] /= sum_abundances
# Compute the ratio of the nuclide atomic masses to the element
# atomic mass
if percent_type == 'wo':

55
openmc/material.py
View file

@ -499,34 +499,62 @@ class Material(IDManagerMixin):
if macroscopic == self._macroscopic:
self._macroscopic = None
def add_element(self, element, percent, percent_type='ao', enrichment=None):
def add_element(self, element, percent, percent_type='ao', enrichment=None,
enrichment_target=None, enrichment_type=None):
"""Add a natural element to the material
Parameters
----------
element : str
Element to add, e.g., 'Zr'
Element to add, e.g., 'Zr' or 'Zirconium'
percent : float
Atom or weight percent
percent_type : {'ao', 'wo'}, optional
'ao' for atom percent and 'wo' for weight percent. Defaults to atom
percent.
enrichment : float, optional
Enrichment for U235 in weight percent. For example, input 4.95 for
4.95 weight percent enriched U. Default is None
(natural composition).
Enrichment of an enrichment_taget nuclide in percent (ao or wo).
If enrichment_taget is not supplied then it is enrichment for U235
in weight percent. For example, input 4.95 for 4.95 weight percent
enriched U.
Default is None (natural composition).
enrichment_target: str, optional
Single nuclide name to enrich from a natural composition (e.g., 'O16')
enrichment_type: {'ao', 'wo'}, optional
'ao' for enrichment as atom percent and 'wo' for weight percent.
Default is: 'ao' for two-isotope enrichment; 'wo' for U enrichment
Notes
-----
General enrichment procedure is allowed only for elements composed of
two isotopes. If `enrichment_target` is given without `enrichment`
natural composition is added to the material.
"""
cv.check_type('nuclide', element, str)
cv.check_type('percent', percent, Real)
cv.check_value('percent type', percent_type, {'ao', 'wo'})
# Make sure element name is just that
if not element.isalpha():
raise ValueError("Element name should be given by the "
"element's symbol or name, e.g., 'Zr', 'zirconium'")
# Allow for element identifier to be given as a symbol or name
if len(element) > 2:
el = element.lower()
element = openmc.data.ELEMENT_SYMBOL.get(el)
if element is None:
msg = 'Element name "{}" not recognised'.format(el)
raise ValueError(msg)
if self._macroscopic is not None:
msg = 'Unable to add an Element to Material ID="{}" as a ' \
'macroscopic data-set has already been added'.format(self._id)
raise ValueError(msg)
if enrichment is not None:
if enrichment is not None and enrichment_target is None:
if not isinstance(enrichment, Real):
msg = 'Unable to add an Element to Material ID="{}" with a ' \
'non-floating point enrichment value "{}"'\
@ -551,14 +579,13 @@ class Material(IDManagerMixin):
format(enrichment, self._id)
warnings.warn(msg)
# Make sure element name is just that
if not element.isalpha():
raise ValueError("Element name should be given by the "
"element's symbol, e.g., 'Zr'")
# Add naturally-occuring isotopes
element = openmc.Element(element)
for nuclide in element.expand(percent, percent_type, enrichment):
for nuclide in element.expand(percent,
percent_type,
enrichment,
enrichment_target,
enrichment_type):
self.add_nuclide(*nuclide)
def add_s_alpha_beta(self, name, fraction=1.0):
@ -927,7 +954,7 @@ class Material(IDManagerMixin):
Fractions of each material to be combined
percent_type : {'ao', 'wo', 'vo'}
Type of percentage, must be one of 'ao', 'wo', or 'vo', to signify atom
percent (molar percent), weight percent, or volume percent,
percent (molar percent), weight percent, or volume percent,
optional. Defaults to 'ao'
name : str
The name for the new material, optional. Defaults to concatenated
@ -996,7 +1023,7 @@ class Material(IDManagerMixin):
zip(materials, fracs)])
new_mat = openmc.Material(name=name)
# Compute atom fractions of nuclides and add them to the new material
# Compute atom fractions of nuclides and add them to the new material
tot_nuclides_per_cc = np.sum([dens for dens in nuclides_per_cc.values()])
for nuc, atom_dens in nuclides_per_cc.items():
new_mat.add_nuclide(nuc, atom_dens/tot_nuclides_per_cc, 'ao')

752
openmc/surface.py
View file

@ -27,6 +27,35 @@ will not accept positional parameters for superclass arguments.\
"""
class SurfaceCoefficient:
"""Descriptor class for surface coefficients.
Parameters
-----------
value : float or str
Value of the coefficient (float) or the name of the coefficient that
it is equivalent to (str).
"""
def __init__(self, value):
self.value = value
def __get__(self, instance, owner=None):
if instance is None:
return self
else:
if isinstance(self.value, str):
return instance._coefficients[self.value]
else:
return self.value
def __set__(self, instance, value):
if isinstance(self.value, Real):
raise AttributeError('This coefficient is read-only')
check_type('{} coefficient'.format(self.value), value, Real)
instance._coefficients[self.value] = value
def _future_kwargs_warning_helper(cls, *args, **kwargs):
# Warn if Surface parameters are passed by position, not by keyword
argsdict = dict(zip(('boundary_type', 'name', 'surface_id'), args))
@ -484,7 +513,7 @@ class PlaneMixin(metaclass=ABCMeta):
if np.any(np.isclose(np.abs(nhat), 1., rtol=0., atol=self._atol)):
sign = nhat.sum()
a, b, c, d = self._get_base_coeffs()
vals = [d/val if not np.isclose(val, 0., rtol=0., atol=self._atol)
vals = [d/val if not np.isclose(val, 0., rtol=0., atol=self._atol)
else np.nan for val in (a, b, c)]
if side == '-':
if sign > 0:
@ -674,41 +703,10 @@ class Plane(PlaneMixin, Surface):
return True
return NotImplemented
@property
def a(self):
return self.coefficients['a']
@property
def b(self):
return self.coefficients['b']
@property
def c(self):
return self.coefficients['c']
@property
def d(self):
return self.coefficients['d']
@a.setter
def a(self, a):
check_type('A coefficient', a, Real)
self._coefficients['a'] = a
@b.setter
def b(self, b):
check_type('B coefficient', b, Real)
self._coefficients['b'] = b
@c.setter
def c(self, c):
check_type('C coefficient', c, Real)
self._coefficients['c'] = c
@d.setter
def d(self, d):
check_type('D coefficient', d, Real)
self._coefficients['d'] = d
a = SurfaceCoefficient('a')
b = SurfaceCoefficient('b')
c = SurfaceCoefficient('c')
d = SurfaceCoefficient('d')
@classmethod
def from_points(cls, p1, p2, p3, **kwargs):
@ -792,30 +790,11 @@ class XPlane(PlaneMixin, Surface):
super().__init__(**kwargs)
self.x0 = x0
@property
def x0(self):
return self.coefficients['x0']
@property
def a(self):
return 1.
@property
def b(self):
return 0.
@property
def c(self):
return 0.
@property
def d(self):
return self.x0
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
x0 = SurfaceCoefficient('x0')
a = SurfaceCoefficient(1.)
b = SurfaceCoefficient(0.)
c = SurfaceCoefficient(0.)
d = x0
def evaluate(self, point):
return point[0] - self.x0
@ -870,30 +849,11 @@ class YPlane(PlaneMixin, Surface):
super().__init__(**kwargs)
self.y0 = y0
@property
def y0(self):
return self.coefficients['y0']
@property
def a(self):
return 0.
@property
def b(self):
return 1.
@property
def c(self):
return 0.
@property
def d(self):
return self.y0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
y0 = SurfaceCoefficient('y0')
a = SurfaceCoefficient(0.)
b = SurfaceCoefficient(1.)
c = SurfaceCoefficient(0.)
d = y0
def evaluate(self, point):
return point[1] - self.y0
@ -948,30 +908,11 @@ class ZPlane(PlaneMixin, Surface):
super().__init__(**kwargs)
self.z0 = z0
@property
def z0(self):
return self.coefficients['z0']
@property
def a(self):
return 0.
@property
def b(self):
return 0.
@property
def c(self):
return 1.
@property
def d(self):
return self.z0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
z0 = SurfaceCoefficient('z0')
a = SurfaceCoefficient(0.)
b = SurfaceCoefficient(0.)
c = SurfaceCoefficient(1.)
d = z0
def evaluate(self, point):
return point[2] - self.z0
@ -1095,7 +1036,7 @@ class QuadricMixin(metaclass=ABCMeta):
return surf
def rotate(self, rotation, pivot=(0., 0., 0.), order='xyz', inplace=False):
# Get pivot and rotation matrix
# Get pivot and rotation matrix
pivot = np.asarray(pivot)
rotation = np.asarray(rotation, dtype=float)
@ -1225,68 +1166,13 @@ class Cylinder(QuadricMixin, Surface):
return True
return NotImplemented
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r(self):
return self.coefficients['r']
@property
def dx(self):
return self.coefficients['dx']
@property
def dy(self):
return self.coefficients['dy']
@property
def dz(self):
return self.coefficients['dz']
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r.setter
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
@dx.setter
def dx(self, dx):
check_type('dx coefficient', dx, Real)
self._coefficients['dx'] = dx
@dy.setter
def dy(self, dy):
check_type('dy coefficient', dy, Real)
self._coefficients['dy'] = dy
@dz.setter
def dz(self, dz):
check_type('dz coefficient', dz, Real)
self._coefficients['dz'] = dz
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r = SurfaceCoefficient('r')
dx = SurfaceCoefficient('dx')
dy = SurfaceCoefficient('dy')
dz = SurfaceCoefficient('dz')
def bounding_box(self, side):
if side == '-':
@ -1305,7 +1191,7 @@ class Cylinder(QuadricMixin, Surface):
# Get x, y, z coordinates of two points
x1, y1, z1 = self._origin
x2, y2, z2 = self._origin + self._axis
r = self.r
r = self.r
# Define intermediate terms
dx = x2 - x1
@ -1375,7 +1261,7 @@ class Cylinder(QuadricMixin, Surface):
with catch_warnings():
simplefilter('ignore', IDWarning)
kwargs = {'boundary_type': self.boundary_type, 'name': self.name,
'surface_id': self.id}
'surface_id': self.id}
quad_rep = Quadric(*self._get_base_coeffs(), **kwargs)
return quad_rep.to_xml_element()
@ -1440,48 +1326,13 @@ class XCylinder(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (y0, z0, r)):
setattr(self, key, val)
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r(self):
return self.coefficients['r']
@property
def x0(self):
return 0.
@property
def dx(self):
return 1.
@property
def dy(self):
return 0.
@property
def dz(self):
return 0.
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r.setter
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
x0 = SurfaceCoefficient(0.)
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r = SurfaceCoefficient('r')
dx = SurfaceCoefficient(1.)
dy = SurfaceCoefficient(0.)
dz = SurfaceCoefficient(0.)
def _get_base_coeffs(self):
y0, z0, r = self.y0, self.z0, self.r
@ -1566,48 +1417,13 @@ class YCylinder(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, z0, r)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r(self):
return self.coefficients['r']
@property
def y0(self):
return 0.
@property
def dx(self):
return 0.
@property
def dy(self):
return 1.
@property
def dz(self):
return 0.
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r.setter
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient(0.)
z0 = SurfaceCoefficient('z0')
r = SurfaceCoefficient('r')
dx = SurfaceCoefficient(0.)
dy = SurfaceCoefficient(1.)
dz = SurfaceCoefficient(0.)
def _get_base_coeffs(self):
x0, z0, r = self.x0, self.z0, self.r
@ -1692,48 +1508,13 @@ class ZCylinder(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, y0, r)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def r(self):
return self.coefficients['r']
@property
def z0(self):
return 0.
@property
def dx(self):
return 0.
@property
def dy(self):
return 0.
@property
def dz(self):
return 1.
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@r.setter
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient(0.)
r = SurfaceCoefficient('r')
dx = SurfaceCoefficient(0.)
dy = SurfaceCoefficient(0.)
dz = SurfaceCoefficient(1.)
def _get_base_coeffs(self):
x0, y0, r = self.x0, self.y0, self.r
@ -1820,41 +1601,10 @@ class Sphere(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, y0, z0, r)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r(self):
return self.coefficients['r']
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r.setter
def r(self, r):
check_type('r coefficient', r, Real)
self._coefficients['r'] = r
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r = SurfaceCoefficient('r')
def _get_base_coeffs(self):
x0, y0, z0, r = self.x0, self.y0, self.z0, self.r
@ -1966,68 +1716,13 @@ class Cone(QuadricMixin, Surface):
return True
return NotImplemented
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r2(self):
return self.coefficients['r2']
@property
def dx(self):
return self.coefficients['dx']
@property
def dy(self):
return self.coefficients['dy']
@property
def dz(self):
return self.coefficients['dz']
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r2.setter
def r2(self, r2):
check_type('r^2 coefficient', r2, Real)
self._coefficients['r2'] = r2
@dx.setter
def dx(self, dx):
check_type('dx coefficient', dx, Real)
self._coefficients['dx'] = dx
@dy.setter
def dy(self, dy):
check_type('dy coefficient', dy, Real)
self._coefficients['dy'] = dy
@dz.setter
def dz(self, dz):
check_type('dz coefficient', dz, Real)
self._coefficients['dz'] = dz
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r2 = SurfaceCoefficient('r2')
dx = SurfaceCoefficient('dx')
dy = SurfaceCoefficient('dy')
dz = SurfaceCoefficient('dz')
def _get_base_coeffs(self):
# The equation for a general cone with vertex at point p = (x0, y0, z0)
@ -2074,7 +1769,7 @@ class Cone(QuadricMixin, Surface):
with catch_warnings():
simplefilter('ignore', IDWarning)
kwargs = {'boundary_type': self.boundary_type, 'name': self.name,
'surface_id': self.id}
'surface_id': self.id}
quad_rep = Quadric(*self._get_base_coeffs(), **kwargs)
return quad_rep.to_xml_element()
@ -2142,53 +1837,13 @@ class XCone(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, y0, z0, r2)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r2(self):
return self.coefficients['r2']
@property
def dx(self):
return 1.
@property
def dy(self):
return 0.
@property
def dz(self):
return 0.
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r2.setter
def r2(self, r2):
check_type('r^2 coefficient', r2, Real)
self._coefficients['r2'] = r2
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r2 = SurfaceCoefficient('r2')
dx = SurfaceCoefficient(1.)
dy = SurfaceCoefficient(0.)
dz = SurfaceCoefficient(0.)
def _get_base_coeffs(self):
x0, y0, z0, r2 = self.x0, self.y0, self.z0, self.r2
@ -2271,53 +1926,13 @@ class YCone(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, y0, z0, r2)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r2(self):
return self.coefficients['r2']
@property
def dx(self):
return 0.
@property
def dy(self):
return 1.
@property
def dz(self):
return 0.
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r2.setter
def r2(self, r2):
check_type('r^2 coefficient', r2, Real)
self._coefficients['r2'] = r2
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r2 = SurfaceCoefficient('r2')
dx = SurfaceCoefficient(0.)
dy = SurfaceCoefficient(1.)
dz = SurfaceCoefficient(0.)
def _get_base_coeffs(self):
x0, y0, z0, r2 = self.x0, self.y0, self.z0, self.r2
@ -2400,53 +2015,13 @@ class ZCone(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (x0, y0, z0, r2)):
setattr(self, key, val)
@property
def x0(self):
return self.coefficients['x0']
@property
def y0(self):
return self.coefficients['y0']
@property
def z0(self):
return self.coefficients['z0']
@property
def r2(self):
return self.coefficients['r2']
@property
def dx(self):
return 0.
@property
def dy(self):
return 0.
@property
def dz(self):
return 1.
@x0.setter
def x0(self, x0):
check_type('x0 coefficient', x0, Real)
self._coefficients['x0'] = x0
@y0.setter
def y0(self, y0):
check_type('y0 coefficient', y0, Real)
self._coefficients['y0'] = y0
@z0.setter
def z0(self, z0):
check_type('z0 coefficient', z0, Real)
self._coefficients['z0'] = z0
@r2.setter
def r2(self, r2):
check_type('r^2 coefficient', r2, Real)
self._coefficients['r2'] = r2
x0 = SurfaceCoefficient('x0')
y0 = SurfaceCoefficient('y0')
z0 = SurfaceCoefficient('z0')
r2 = SurfaceCoefficient('r2')
dx = SurfaceCoefficient(0.)
dy = SurfaceCoefficient(0.)
dz = SurfaceCoefficient(1.)
def _get_base_coeffs(self):
x0, y0, z0, r2 = self.x0, self.y0, self.z0, self.r2
@ -2513,95 +2088,16 @@ class Quadric(QuadricMixin, Surface):
for key, val in zip(self._coeff_keys, (a, b, c, d, e, f, g, h, j, k)):
setattr(self, key, val)
@property
def a(self):
return self.coefficients['a']
@property
def b(self):
return self.coefficients['b']
@property
def c(self):
return self.coefficients['c']
@property
def d(self):
return self.coefficients['d']
@property
def e(self):
return self.coefficients['e']
@property
def f(self):
return self.coefficients['f']
@property
def g(self):
return self.coefficients['g']
@property
def h(self):
return self.coefficients['h']
@property
def j(self):
return self.coefficients['j']
@property
def k(self):
return self.coefficients['k']
@a.setter
def a(self, a):
check_type('a coefficient', a, Real)
self._coefficients['a'] = a
@b.setter
def b(self, b):
check_type('b coefficient', b, Real)
self._coefficients['b'] = b
@c.setter
def c(self, c):
check_type('c coefficient', c, Real)
self._coefficients['c'] = c
@d.setter
def d(self, d):
check_type('d coefficient', d, Real)
self._coefficients['d'] = d
@e.setter
def e(self, e):
check_type('e coefficient', e, Real)
self._coefficients['e'] = e
@f.setter
def f(self, f):
check_type('f coefficient', f, Real)
self._coefficients['f'] = f
@g.setter
def g(self, g):
check_type('g coefficient', g, Real)
self._coefficients['g'] = g
@h.setter
def h(self, h):
check_type('h coefficient', h, Real)
self._coefficients['h'] = h
@j.setter
def j(self, j):
check_type('j coefficient', j, Real)
self._coefficients['j'] = j
@k.setter
def k(self, k):
check_type('k coefficient', k, Real)
self._coefficients['k'] = k
a = SurfaceCoefficient('a')
b = SurfaceCoefficient('b')
c = SurfaceCoefficient('c')
d = SurfaceCoefficient('d')
e = SurfaceCoefficient('e')
f = SurfaceCoefficient('f')
g = SurfaceCoefficient('g')
h = SurfaceCoefficient('h')
j = SurfaceCoefficient('j')
k = SurfaceCoefficient('k')
def _get_base_coeffs(self):
return tuple(getattr(self, c) for c in self._coeff_keys)
@ -2734,7 +2230,7 @@ class Halfspace(Region):
Parameters
----------
redundant_surfaces : dict
Dictionary mapping redundant surface IDs to surface IDs for the
Dictionary mapping redundant surface IDs to surface IDs for the
:class:`openmc.Surface` instances that should replace them.
"""

29
src/geometry_aux.cpp
View file

@ -385,8 +385,10 @@ prepare_distribcell()
}
// Fill the cell and lattice offset tables.
#pragma omp parallel for
for (int map = 0; map < target_univ_ids.size(); map++) {
auto target_univ_id = target_univ_ids[map];
std::unordered_map<int32_t, int32_t> univ_count_memo;
for (const auto& univ : model::universes) {
int32_t offset = 0;
for (int32_t cell_indx : univ->cells_) {
@ -395,11 +397,13 @@ prepare_distribcell()
if (c.type_ == Fill::UNIVERSE) {
c.offset_[map] = offset;
int32_t search_univ = c.fill_;
offset += count_universe_instances(search_univ, target_univ_id);
offset += count_universe_instances(search_univ, target_univ_id,
univ_count_memo);
} else if (c.type_ == Fill::LATTICE) {
Lattice& lat = *model::lattices[c.fill_];
offset = lat.fill_offset_table(offset, target_univ_id, map);
offset = lat.fill_offset_table(offset, target_univ_id, map,
univ_count_memo);
}
}
}
@ -411,7 +415,6 @@ prepare_distribcell()
void
count_cell_instances(int32_t univ_indx)
{
const auto univ_counts = model::universe_cell_counts.find(univ_indx);
if (univ_counts != model::universe_cell_counts.end()) {
for (const auto& it : univ_counts->second) {
@ -442,30 +445,42 @@ count_cell_instances(int32_t univ_indx)
//==============================================================================
int
count_universe_instances(int32_t search_univ, int32_t target_univ_id)
count_universe_instances(int32_t search_univ, int32_t target_univ_id,
std::unordered_map<int32_t, int32_t>& univ_count_memo)
{
// If this is the target, it can't contain itself.
// If this is the target, it can't contain itself.
if (model::universes[search_univ]->id_ == target_univ_id) {
return 1;
}
// If we have already counted the number of instances, reuse that value.
auto search = univ_count_memo.find(search_univ);
if (search != univ_count_memo.end()) {
return search->second;
}
int count {0};
for (int32_t cell_indx : model::universes[search_univ]->cells_) {
Cell& c = *model::cells[cell_indx];
if (c.type_ == Fill::UNIVERSE) {
int32_t next_univ = c.fill_;
count += count_universe_instances(next_univ, target_univ_id);
count += count_universe_instances(next_univ, target_univ_id,
univ_count_memo);
} else if (c.type_ == Fill::LATTICE) {
Lattice& lat = *model::lattices[c.fill_];
for (auto it = lat.begin(); it != lat.end(); ++it) {
int32_t next_univ = *it;
count += count_universe_instances(next_univ, target_univ_id);
count += count_universe_instances(next_univ, target_univ_id,
univ_count_memo);
}
}
}
// Remember the number of instances in this universe.
univ_count_memo[search_univ] = count;
return count;
}

5
src/lattice.cpp
View file

@ -93,11 +93,12 @@ Lattice::adjust_indices()
//==============================================================================
int32_t
Lattice::fill_offset_table(int32_t offset, int32_t target_univ_id, int map)
Lattice::fill_offset_table(int32_t offset, int32_t target_univ_id, int map,
std::unordered_map<int32_t, int32_t>& univ_count_memo)
{
for (LatticeIter it = begin(); it != end(); ++it) {
offsets_[map * universes_.size() + it.indx_] = offset;
offset += count_universe_instances(*it, target_univ_id);
offset += count_universe_instances(*it, target_univ_id, univ_count_memo);
}
return offset;
}

2
src/main.cpp
View file

@ -43,6 +43,8 @@ int main(int argc, char* argv[]) {
case RunMode::VOLUME:
err = openmc_calculate_volumes();
break;
default:
break;
}
if (err) fatal_error(openmc_err_msg);

2
src/particle.cpp
View file

@ -675,6 +675,8 @@ Particle::write_restart() const
case RunMode::PARTICLE:
write_dataset(file_id, "run_mode", "particle restart");
break;
default:
break;
}
write_dataset(file_id, "id", id_);
write_dataset(file_id, "type", static_cast<int>(type_));

2
src/scattdata.cpp
View file

@ -227,6 +227,8 @@ ScattData::get_xs(MgxsType xstype, int gin, const int* gout, const double* mu)
fatal_error("Invalid call to get_xs");
}
break;
default:
break;
}
return val;
}

2
src/state_point.cpp
View file

@ -87,6 +87,8 @@ openmc_statepoint_write(const char* filename, bool* write_source)
case RunMode::EIGENVALUE:
write_dataset(file_id, "run_mode", "eigenvalue");
break;
default:
break;
}
write_attribute(file_id, "photon_transport", settings::photon_transport);
write_dataset(file_id, "n_particles", settings::n_particles);

4
src/tallies/tally_scoring.cpp
View file

@ -27,7 +27,7 @@ namespace openmc {
//==============================================================================
FilterBinIter::FilterBinIter(const Tally& tally, Particle* p)
: tally_{tally}, filter_matches_{p->filter_matches_}
: filter_matches_{p->filter_matches_}, tally_{tally}
{
// Find all valid bins in each relevant filter if they have not already been
// found for this event.
@ -57,7 +57,7 @@ FilterBinIter::FilterBinIter(const Tally& tally, Particle* p)
FilterBinIter::FilterBinIter(const Tally& tally, bool end,
std::vector<FilterMatch>* particle_filter_matches)
: tally_{tally}, filter_matches_{*particle_filter_matches}
: filter_matches_{*particle_filter_matches}, tally_{tally}
{
// Handle the special case for an iterator that points to the end.
if (end) {

4
src/volume_calc.cpp
View file

@ -295,6 +295,8 @@ std::vector<VolumeCalculation::Result> VolumeCalculation::execute() const
case TriggerMetric::variance:
val = result.volume[1] * result.volume[1];
break;
default:
break;
}
// update max if entry is valid
if (val > 0.0) { trigger_val = std::max(trigger_val, val); }
@ -374,6 +376,8 @@ void VolumeCalculation::to_hdf5(const std::string& filename,
case TriggerMetric::relative_error:
trigger_str = "rel_err";
break;
default:
break;
}
write_attribute(file_id, "trigger_type", trigger_str);
} else {

81
tests/unit_tests/test_element.py Normal file
View file

@ -0,0 +1,81 @@
import openmc
from pytest import approx, raises
from openmc.data import NATURAL_ABUNDANCE, atomic_mass
def test_expand_no_enrichment():
""" Expand Li in natural compositions"""
lithium = openmc.Element('Li')
# Verify the expansion into ATOMIC fraction against natural composition
for isotope in lithium.expand(100.0, 'ao'):
assert isotope[1] == approx(NATURAL_ABUNDANCE[isotope[0]] * 100.0)
# Verify the expansion into WEIGHT fraction against natural composition
natural = {'Li6': NATURAL_ABUNDANCE['Li6'] * atomic_mass('Li6'),
'Li7': NATURAL_ABUNDANCE['Li7'] * atomic_mass('Li7')}
li_am = sum(natural.values())
for key in natural:
natural[key] /= li_am
for isotope in lithium.expand(100.0, 'wo'):
assert isotope[1] == approx(natural[isotope[0]] * 100.0)
def test_expand_enrichment():
""" Expand and verify enrichment of Li """
lithium = openmc.Element('Li')
# Verify the enrichment by atoms
ref = {'Li6': 75.0, 'Li7': 25.0}
for isotope in lithium.expand(100.0, 'ao', 25.0, 'Li7', 'ao'):
assert isotope[1] == approx(ref[isotope[0]])
# Verify the enrichment by weight
for isotope in lithium.expand(100.0, 'wo', 25.0, 'Li7', 'wo'):
assert isotope[1] == approx(ref[isotope[0]])
def test_expand_exceptions():
""" Test that correct exceptions are raised for invalid input """
# 1 Isotope Element
with raises(ValueError):
element = openmc.Element('Be')
element.expand(70.0, 'ao', 4.0, 'Be9')
# 3 Isotope Element
with raises(ValueError):
element = openmc.Element('Cr')
element.expand(70.0, 'ao', 4.0, 'Cr52')
# Non-present Enrichment Target
with raises(ValueError):
element = openmc.Element('H')
element.expand(70.0, 'ao', 4.0, 'H4')
# Enrichment Procedure for Uranium if not Uranium
with raises(ValueError):
element = openmc.Element('Li')
element.expand(70.0, 'ao', 4.0)
# Missing Enrichment Target
with raises(ValueError):
element = openmc.Element('Li')
element.expand(70.0, 'ao', 4.0, enrichment_type='ao')
# Invalid Enrichment Type Entry
with raises(ValueError):
element = openmc.Element('Li')
element.expand(70.0, 'ao', 4.0, 'Li7', 'Grand Moff Tarkin')
# Trying to enrich Uranium
with raises(ValueError):
element = openmc.Element('U')
element.expand(70.0, 'ao', 4.0, 'U235', 'wo')
# Trying to enrich Uranium with wrong enrichment_target
with raises(ValueError):
element = openmc.Element('U')
element.expand(70.0, 'ao', 4.0, enrichment_type='ao')

23
tests/unit_tests/test_material.py
View file

@ -29,11 +29,34 @@ def test_elements():
m = openmc.Material()
m.add_element('Zr', 1.0)
m.add_element('U', 1.0, enrichment=4.5)
m.add_element('Li', 1.0, enrichment=60.0, enrichment_target='Li7')
m.add_element('H', 1.0, enrichment=50.0, enrichment_target='H2',
enrichment_type='wo')
with pytest.raises(ValueError):
m.add_element('U', 1.0, enrichment=100.0)
with pytest.raises(ValueError):
m.add_element('Pu', 1.0, enrichment=3.0)
with pytest.raises(ValueError):
m.add_element('U', 1.0, enrichment=70.0, enrichment_target='U235')
with pytest.raises(ValueError):
m.add_element('He', 1.0, enrichment=17.0, enrichment_target='He6')
def test_elements_by_name():
"""Test adding elements by name"""
m = openmc.Material()
m.add_element('woLfrAm', 1.0)
with pytest.raises(ValueError):
m.add_element('uranum', 1.0)
m.add_element('uRaNiUm', 1.0)
m.add_element('Aluminium', 1.0)
a = openmc.Material()
b = openmc.Material()
c = openmc.Material()
a.add_element('sulfur', 1.0)
b.add_element('SulPhUR', 1.0)
c.add_element('S', 1.0)
assert a._nuclides == b._nuclides
assert b._nuclides == c._nuclides
def test_density():
m = openmc.Material()

2
vendor/xtensor vendored

@ -1 +1 @@
Subproject commit ef091807f7ed0e5ba7e251a6c46f4af7bba79e2e
Subproject commit 31acec1e90bbea6d4bc17af0710a123bd5da6689

2
vendor/xtl vendored

@ -1 +1 @@
Subproject commit f5d13e6c4f856becc178939365fcdcf9a657ffb5
Subproject commit 0024346605bd92bcc4009caad7f4be88687e063a