removed element expansion in Fortran and added it to Python API

This commit is contained in:
Sam Shaner 2016-10-22 16:36:02 -04:00
parent 5708ff9d77
commit 716b6a92f5
27 changed files with 497 additions and 1081 deletions

View file

@ -1,12 +1,13 @@
import re
import sys
import os
from six import string_types
import openmc
from openmc.checkvalue import check_type, check_length
from openmc.data import NATURAL_ABUNDANCE
from xml.etree import ElementTree as ET
class Element(object):
"""A natural element used in a material via <element>. Internally, OpenMC
@ -90,9 +91,29 @@ class Element(object):
self._scattering = scattering
def expand(self):
def expand(self, percent, percent_type, enrichment=None,
cross_sections=None):
"""Expand natural element into its naturally-occurring isotopes.
An optional cross_sections argument or the OPENMC_CROSS_SECTIONS
environment variable is used to specify a cross_sections.xml file.
If the cross_sections.xml file is found, the element is expanded only
into the isotopes/nuclides present in cross_sections.xml. If no
cross_sections.xml file is found, the element is expanded based on its
naturally occurring isotopes.
Parameters
----------
percent : float
Atom or weight percent
percent_type : {'ao', 'wo'}
'ao' for atom percent and 'wo' for weight percent
enrichment : float, optional
Enrichment percent for U235 in U. Default is None
(natural composition).
cross_sections : str, optional
Location of cross_sections.xml file. Default is None.
Returns
-------
isotopes : list
@ -102,9 +123,150 @@ class Element(object):
"""
# Get the nuclides present in nature
natural_nuclides = set()
for nuclide in sorted(NATURAL_ABUNDANCE.keys()):
if re.match(r'{}\d+'.format(self.name), nuclide):
natural_nuclides.add(int(nuclide[len(self.name):]))
# Create lists to store the expanded nuclides and abundances
nuclides = []
abundances = []
# If cross_sections is None, get the cross sections from the
# OPENMC_CROSS_SECTIONS environment variable
if cross_sections is None:
cross_sections = os.environ.get('OPENMC_CROSS_SECTIONS')
# If a cross_sections library is present, check natural nuclides
# against the nuclides in the library
if cross_sections is not None:
library_nuclides = set()
tree = ET.parse(cross_sections)
root = tree.getroot()
for child in root:
nuclide = child.attrib['materials']
if re.match(r'{}\d+'.format(self.name), nuclide) and \
'_m' not in nuclide:
library_nuclides.add(int(nuclide[len(self.name):]))
# Get a set of the mutual and absent nuclides
mutual_nuclides = natural_nuclides.intersection(library_nuclides)
absent_nuclides = natural_nuclides.difference(mutual_nuclides)
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 len(absent_nuclides) == 0:
for nuclide, abundance in sorted(NATURAL_ABUNDANCE.items()):
if re.match(r'{}\d+'.format(self.name), nuclide):
nuclides.append(openmc.Nuclide(nuclide))
abundances.append(abundance)
# If no natural elements are present in the library, check if the
# 0 element is present. If so, set the abundance to 1 for this
# nuclide. Else, raise an error.
elif len(mutual_nuclides) == 0:
if 0 in library_nuclides:
nuclides.append(openmc.Nuclide(self.name + '0'))
abundances.append(1.0)
else:
msg = 'Unable to expand element {0} because the cross '\
'section library provided does not contain any of '\
'the natural isotopes for that element.'\
.format(self.name)
raise ValueError(msg)
# If some, but not all, natural nuclides are in the library, add
# the mutual nuclides. For the absent nuclides, increment the
# abundance of the nearest mutual nuclide with their abundances.
else:
# Add the mutual isotopes
nuclides_a = []
for nuclide, abundance in sorted(NATURAL_ABUNDANCE.items()):
if re.match(r'{}\d+'.format(self.name), nuclide) and \
int(nuclide[len(self.name):]) in mutual_nuclides:
nuclides.append(openmc.Nuclide(nuclide))
nuclides_a.append(int(nuclide[len(self.name):]))
abundances.append(abundance)
# Adjust the abundances for the absent nuclides
for nuclide, abundance in sorted(NATURAL_ABUNDANCE.items()):
if re.match(r'{}\d+'.format(self.name), nuclide) and \
int(nuclide[len(self.name):]) in absent_nuclides:
# Get index to the nearest nuclide
a = int(nuclide[len(self.name):])
i = min(list(range(len(nuclides_a))), key=lambda j: \
abs(nuclides_a[j] - a))
# Increment abundance of the nearest nuclide
abundances[i] += abundance
# If a cross_section library is not present, expand the element into
# its natural nuclides
else:
for nuclide, abundance in sorted(NATURAL_ABUNDANCE.items()):
if re.match(r'{}\d+'.format(self.name), nuclide):
nuclides.append(openmc.Nuclide(nuclide))
abundances.append(abundance)
# Create a list of atomic masses
n_nuclides = len(nuclides)
atomic_masses = []
for nuclide in nuclides:
atomic_masses.append(openmc.data.atomic_mass(nuclide.name))
# Modify mole fractions if enrichment provided
if enrichment is not None:
# Get the indices for the uranium nuclides
for i,nuc in enumerate(nuclides):
if nuc.name == 'U234':
u234 = i
elif nuc.name == 'U235':
u235 = i
elif nuc.name == 'U238':
u238 = i
# Calculate the mass fractions of isotopes
abundances[u234] = 0.008 * enrichment
abundances[u235] = enrichment
abundances[u238] = 1.0 - 1.008 * enrichment
# Convert the mass fractions to mole fractions
for i in range(n_nuclides):
abundances[i] /= atomic_masses[i]
# Normalize the mole fractions to one
sum_abundances = sum(abundances)
for i in range(n_nuclides):
abundances[i] /= sum_abundances
# Compute the ratio of the nuclide atomic massess to the element
# atomic mass
if percent_type == 'wo':
# Compute the element atomic mass
element_am = 0.
for i in range(n_nuclides):
element_am += atomic_masses[i] * abundances[i]
# Convert the molar fractions to mass fractions
for i in range(n_nuclides):
abundances[i] *= atomic_masses[i] / element_am
# Normalize the mass fractions to one
sum_abundances = sum(abundances)
for i in range(n_nuclides):
abundances[i] /= sum_abundances
isotopes = []
for isotope, abundance in sorted(NATURAL_ABUNDANCE.items()):
if re.match(r'{}\d+'.format(self.name), isotope):
nuc = openmc.Nuclide(isotope)
isotopes.append((nuc, abundance))
for nuclide, abundance in zip(nuclides, abundances):
pct = float('{:2.10f}'.format(percent*abundance))
isotopes.append((nuclide, pct, percent_type))
return isotopes

View file

@ -57,9 +57,9 @@ class Material(object):
'atom/b-cm', 'atom/cm3', 'sum', or 'macro'. The 'macro' unit only
applies in the case of a multi-group calculation.
elements : list of tuple
List in which each item is a 3-tuple consisting of an
:class:`openmc.Element` instance, the percent density, and the percent
type ('ao' or 'wo').
List in which each item is a 4-tuple consisting of an
:class:`openmc.Element` instance, the percent density, the percent
type ('ao' or 'wo'), and enrichment.
nuclides : list of tuple
List in which each item is a 3-tuple consisting of an
:class:`openmc.Nuclide` instance, the percent density, and the percent
@ -82,7 +82,7 @@ class Material(object):
# (only one is allowed, hence this is different than _nuclides, etc)
self._macroscopic = None
# A list of tuples (element, percent, percent type)
# A list of tuples (element, percent, percent type, enrichment)
self._elements = []
# If specified, a list of table names
@ -148,9 +148,13 @@ class Material(object):
string += '{0: <16}\n'.format('\tElements')
for element, percent, percent_type in self._elements:
for element, percent, percent_type, enr in self._elements:
string += '{0: <16}'.format('\t{0.name}'.format(element))
string += '=\t{0: <12} [{1}]\n'.format(percent, percent_type)
if enr is None:
string += '=\t{0: <12} [{1}]\n'.format(percent, percent_type)
else:
string += '=\t{0: <12} [{1}] @ {2} w/o enrichment\n'\
.format(percent, percent_type, enr)
return string
@ -401,8 +405,7 @@ class Material(object):
if macroscopic.name == self._macroscopic.name:
self._macroscopic = None
def add_element(self, element, percent, percent_type='ao', enrichment=None,
expand=False):
def add_element(self, element, percent, percent_type='ao', enrichment=None):
"""Add a natural element to the material
Parameters
@ -416,9 +419,6 @@ class Material(object):
percent.
enrichment : float, optional
Optional weight percent enrichment for uranium. Defaults to None.
expand : bool, optional
Whether to expand the natural element into its naturally-occurring
isotopes. Defaults to False.
"""
@ -470,68 +470,7 @@ class Material(object):
format(enrichment, self._id)
warnings.warn(msg)
if expand or enrichment is not None:
# Get the isotopes and their natural abundances
isotopes, abundances = zip(*element.expand())
abundances = list(abundances)
isotopes = list(isotopes)
n_isotopes = len(isotopes)
# Create a list of atomic masses
atomic_masses = []
for i,isotope in enumerate(isotopes):
atomic_masses.append(openmc.data.atomic_mass(isotope.name))
# Modify mole fractions if enrichment provided
if enrichment is not None:
# Get the indices for the uranium isotopes
for i,iso in enumerate(isotopes):
if iso.name == 'U234':
u234 = i
elif iso.name == 'U235':
u235 = i
elif iso.name == 'U238':
u238 = i
# Calculate the mass fractions of isotopes
abundances[u234] = 0.008 * enrichment
abundances[u235] = enrichment
abundances[u238] = 1.0 - 1.008 * enrichment
# Convert the mass fractions to mole fractions
for i in range(n_isotopes):
abundances[i] /= atomic_masses[i]
# Normalize the mole fractions to one
sum_abundances = sum(abundances)
for i in range(n_isotopes):
abundances[i] /= sum_abundances
# Compute the ratio of the isotope atomic massess to the element
# atomic mass
if percent_type == 'wo':
# Compute the element atomic mass
element_am = 0.
for i in range(n_isotopes):
element_am += atomic_masses[i] * abundances[i]
# Convert the molar fractions to mass fractions
for i in range(n_isotopes):
abundances[i] *= atomic_masses[i] / element_am
# Normalize the mass fractions to one
sum_abundances = sum(abundances)
for i in range(n_isotopes):
abundances[i] /= sum_abundances
for isotope, abundance in zip(isotopes, abundances):
self._nuclides.append((isotope, percent*abundance, percent_type))
else:
self._elements.append((element, percent, percent_type))
self._elements.append((element, percent, percent_type, enrichment))
def remove_element(self, element):
"""Remove a natural element from the material
@ -553,11 +492,6 @@ class Material(object):
if element == elm:
self._elements.remove(elm)
# If the Material contains Nuclides of the Element, delete them
for nuc in self._nuclides:
if nuc.get_element().name == element.name:
self._nuclides.remove(nuc)
def add_s_alpha_beta(self, name):
r"""Add an :math:`S(\alpha,\beta)` table to the material
@ -586,7 +520,6 @@ class Material(object):
self._sab.append(new_name)
def make_isotropic_in_lab(self):
for nuclide, percent, percent_type in self._nuclides:
nuclide.scattering = 'iso-in-lab'
@ -605,13 +538,15 @@ class Material(object):
nuclides = []
for nuclide, density, density_type in self._nuclides:
for nuclide, percent, percent_type in self._nuclides:
nuclides.append(nuclide.name)
for element, density, density_type in self._elements:
for ele, ele_pct, ele_pct_type, enr in self._elements:
# Expand natural element into isotopes
for isotope, abundance in element.expand():
nuclides.append(isotope.name)
isotopes = ele.expand(ele_pct, ele_pct_type, enr)
for iso, iso_pct, iso_pct_type in isotopes:
nuclides.append(iso.name)
return nuclides
@ -621,8 +556,8 @@ class Material(object):
Returns
-------
nuclides : dict
Dictionary whose keys are nuclide names and values are 2-tuples of
(nuclide, density)
Dictionary whose keys are nuclide names and values are 3-tuples of
(nuclide, density percent, density percent type)
"""
@ -631,10 +566,12 @@ class Material(object):
for nuclide, density, density_type in self._nuclides:
nuclides[nuclide.name] = (nuclide, density)
for element, density, density_type in self._elements:
for ele, ele_pct, ele_pct_type, enr in self._elements:
# Expand natural element into isotopes
for isotope, abundance in element.expand():
nuclides[isotope.name] = (isotope, density*abundance)
isotopes = ele.expand(ele_pct, ele_pct_type, enr)
for iso, iso_pct, iso_pct_type in isotopes:
nuclides[iso.name] = (iso, iso_pct, iso_pct_type)
return nuclides
@ -659,22 +596,20 @@ class Material(object):
return xml_element
def _get_element_xml(self, element, distrib=False):
xml_element = ET.Element("element")
xml_element.set("name", str(element[0].name))
def _get_element_xml(self, element, cross_sections, distrib=False):
if not distrib:
if element[2] == 'ao':
xml_element.set("ao", str(element[1]))
else:
xml_element.set("wo", str(element[1]))
# Get the nuclides in this element
nuclides = element[0].expand(element[1], element[2], element[3],
cross_sections)
if not element[0].scattering is None:
xml_element.set("scattering", element[0].scattering)
xml_elements = []
for nuclide in nuclides:
xml_elements.append(self._get_nuclide_xml(nuclide, distrib))
return xml_element
return xml_elements
def _get_nuclides_xml(self, nuclides, distrib=False):
xml_elements = []
for nuclide in nuclides:
@ -682,15 +617,19 @@ class Material(object):
return xml_elements
def _get_elements_xml(self, elements, distrib=False):
def _get_elements_xml(self, elements, cross_sections, distrib=False):
xml_elements = []
for element in elements:
xml_elements.append(self._get_element_xml(element, distrib))
nuclide_elements = self._get_element_xml(element, cross_sections,
distrib)
for nuclide_element in nuclide_elements:
xml_elements.append(nuclide_element)
return xml_elements
def get_material_xml(self):
def get_material_xml(self, cross_sections):
"""Return XML representation of the material
Returns
@ -726,7 +665,8 @@ class Material(object):
element.append(subelement)
# Create element XML subelements
subelements = self._get_elements_xml(self._elements)
subelements = self._get_elements_xml(self._elements,
cross_sections)
for subelement in subelements:
element.append(subelement)
else:
@ -740,12 +680,12 @@ class Material(object):
comps = []
allnucs = self._nuclides + self._elements
dist_per_type = allnucs[0][2]
for nuc, per, typ in allnucs:
if not typ == dist_per_type:
for nuc in allnucs:
if not nuc[2] == dist_per_type:
msg = 'All nuclides and elements in a distributed ' \
'material must have the same type, either ao or wo'
raise ValueError(msg)
comps.append(per)
comps.append(nuc[1])
if self._distrib_otf_file is None:
# Create values and units subelements
@ -760,12 +700,15 @@ class Material(object):
if self._macroscopic is None:
# Create nuclide XML subelements
subelements = self._get_nuclides_xml(self._nuclides, distrib=True)
subelements = self._get_nuclides_xml(self._nuclides,
distrib=True)
for subelement_nuc in subelements:
subelement.append(subelement_nuc)
# Create element XML subelements
subelements = self._get_elements_xml(self._elements, distrib=True)
subelements = self._get_elements_xml(self._elements,
cross_sections,
distrib=True)
for subsubelement in subelements:
subelement.append(subsubelement)
else:
@ -800,15 +743,81 @@ class Materials(cv.CheckedList):
----------
materials : Iterable of openmc.Material
Materials to add to the collection
cross_sections : str
Indicates the path to an XML cross section listing file (usually named
cross_sections.xml). If it is not set, the
:envvar:`OPENMC_CROSS_SECTIONS` environment variable will be used for
continuous-energy calculations and
:envvar:`OPENMC_MG_CROSS_SECTIONS` will be used for multi-group
calculations to find the path to the XML cross section file.
multipole_library : str
Indicates the path to a directory containing a windowed multipole
cross section library. If it is not set, the
:envvar:`OPENMC_MULTIPOLE_LIBRARY` environment variable will be used. A
multipole library is optional.
temperature : dict
Defines a default temperature and method for treating intermediate
temperatures at which nuclear data doesn't exist. Accepted keys are
'default', 'method', 'tolerance', and 'multipole'. The value for
'default' should be a float representing the default temperature in
Kelvin. The value for 'method' should be 'nearest' or 'interpolation'.
If the method is 'nearest', 'tolerance' indicates a range of temperature
within which cross sections may be used. 'multipole' is a boolean
indicating whether or not the windowed multipole method should be used
to evaluate resolved resonance cross sections.
"""
def __init__(self, materials=None):
super(Materials, self).__init__(Material, 'materials collection')
self._materials_file = ET.Element("materials")
self._cross_sections = None
self._multipole_library = None
self._temperature = {}
if materials is not None:
self += materials
@property
def cross_sections(self):
return self._cross_sections
@property
def multipole_library(self):
return self._multipole_library
@property
def temperature(self):
return self._temperature
@cross_sections.setter
def cross_sections(self, cross_sections):
cv.check_type('cross sections', cross_sections, string_types)
self._cross_sections = cross_sections
@multipole_library.setter
def multipole_library(self, multipole_library):
cv.check_type('cross sections', multipole_library, string_types)
self._multipole_library = multipole_library
@temperature.setter
def temperature(self, temperature):
cv.check_type('temperature settings', temperature, Mapping)
for key, value in temperature.items():
cv.check_value('temperature key', key,
['default', 'method', 'tolerance', 'multipole'])
if key == 'default':
cv.check_type('default temperature', value, Real)
elif key == 'method':
cv.check_value('temperature method', value,
['nearest', 'interpolation'])
elif key == 'tolerance':
cv.check_type('temperature tolerance', value, Real)
elif key == 'multipole':
cv.check_type('temperature multipole', value, bool)
self._temperature = temperature
def add_material(self, material):
"""Append material to collection
@ -891,9 +900,26 @@ class Materials(cv.CheckedList):
def _create_material_subelements(self):
for material in self:
xml_element = material.get_material_xml()
xml_element = material.get_material_xml(self.cross_sections)
self._materials_file.append(xml_element)
def _create_cross_sections_subelement(self):
if self._cross_sections is not None:
element = ET.SubElement(self._materials_file, "cross_sections")
element.text = str(self._cross_sections)
def _create_multipole_library_subelement(self):
if self._multipole_library is not None:
element = ET.SubElement(self._materials_file, "multipole_library")
element.text = str(self._multipole_library)
def _create_temperature_subelements(self):
if self.temperature:
for key, value in sorted(self.temperature.items()):
element = ET.SubElement(self._materials_file,
"temperature_{}".format(key))
element.text = str(value)
def export_to_xml(self, path='materials.xml'):
"""Export material collection to an XML file.
@ -908,6 +934,9 @@ class Materials(cv.CheckedList):
self._materials_file.clear()
self._create_material_subelements()
self._create_cross_sections_subelement()
self._create_multipole_library_subelement()
self._create_temperature_subelements()
# Clean the indentation in the file to be user-readable
sort_xml_elements(self._materials_file)

View file

@ -603,12 +603,19 @@ class Settings(object):
@cross_sections.setter
def cross_sections(self, cross_sections):
warnings.warn('Settings.cross_sections has been deprecated and will be '
'removed in a future version. Materials.cross_sections '
'should defined instead.', DeprecationWarning)
cv.check_type('cross sections', cross_sections, string_types)
self._cross_sections = cross_sections
@multipole_library.setter
def multipole_library(self, multipole_library):
cv.check_type('cross sections', multipole_library, string_types)
warnings.warn('Settings.multipole_library has been deprecated and will '
'be removed in a future version. '
'Materials.multipole_library should defined instead.',
DeprecationWarning)
cv.check_type('multipole library', multipole_library, string_types)
self._multipole_library = multipole_library
@ptables.setter
@ -713,6 +720,11 @@ class Settings(object):
@temperature.setter
def temperature(self, temperature):
warnings.warn('Settings.temperature has been deprecated and will '
'be removed in a future version. Materials.temperature '
'should defined instead.', DeprecationWarning)
cv.check_type('temperature settings', temperature, Mapping)
for key, value in temperature.items():
cv.check_value('temperature key', key,

View file

@ -257,25 +257,6 @@ module constants
! Maximum number of partial fission reactions
integer, parameter :: PARTIAL_FISSION_MAX = 4
! Major cross section libraries
integer, parameter :: &
ENDF_BVII0 = 1, &
ENDF_BVII1 = 2, &
JEFF_311 = 3, &
JEFF_312 = 4, &
JEFF_32 = 5, &
JENDL_32 = 6, &
JENDL_33 = 7, &
JENDL_40 = 8
! Number of natural elements and nuclides in the global natural_elements
! and natural_nuclides arrays. Note that some nuclides are repeated in order
! to make it easy to expand natural elements for different cross section
! libraries.
integer, parameter :: &
NUM_NATURAL_ELEMENTS = 84, &
NUM_NATURAL_NUCLIDES = 301
! Temperature treatment method
integer, parameter :: &
TEMPERATURE_NEAREST = 1, &

View file

@ -96,9 +96,6 @@ module global
! Unreoslved resonance probablity tables
logical :: urr_ptables_on = .true.
! What to assume for expanding natural elements
integer :: default_expand = ENDF_BVII1
! Default temperature and method for choosing temperatures
integer :: temperature_method = TEMPERATURE_NEAREST
logical :: temperature_multipole = .false.
@ -132,333 +129,6 @@ module global
! Number of points to use in the Legendre to tabular conversion
integer :: legendre_to_tabular_points = 33
! ============================================================================
! ELEMENT AND NUCLIDE RELATED VARIABLES
! List of all possible natural element expansion isotopes
character (len=5), dimension (NUM_NATURAL_NUCLIDES), parameter :: &
natural_nuclides = [character(len=5) :: &
'H1' , 'H2' , 'He3' , 'He4' , 'Li6' , 'Li7' , 'Be9' , 'B10' , &
'B11' , 'C0' , 'N14' , 'N15' , 'O16' , 'O17' , 'O18' , 'O16' , &
'O16' , 'O17' , 'F19' , 'Ne20' , 'Ne21' , 'Ne22' , 'Na23' , 'Mg24' , &
'Mg25' , 'Mg26' , 'Al27' , 'Si28' , 'Si29' , 'Si30' , 'P31' , 'S32' , &
'S33' , 'S34' , 'S35' , 'Cl35' , 'Cl37' , 'Ar36' , 'Ar38' , 'Ar40' , &
'K39' , 'K40' , 'K41' , 'Ca40' , 'Ca42' , 'Ca43' , 'Ca44' , 'Ca46' , &
'Ca48' , 'Sc45' , 'Ti46' , 'Ti47' , 'Ti48' , 'Ti49' , 'Ti50' , 'V0' , &
'V50' , 'V51' , 'Cr50' , 'Cr52' , 'Cr53' , 'Cr54' , 'Mn55' , 'Fe54' , &
'Fe56' , 'Fe57' , 'Fe58' , 'Co59' , 'Ni58' , 'Ni60' , 'Ni61' , 'Ni62' , &
'Ni64' , 'Cu63' , 'Cu65' , 'Zn0' , 'Zn64' , 'Zn66' , 'Zn67' , 'Zn68' , &
'Zn70' , 'Ga0' , 'Ga69' , 'Ga71' , 'Ge70' , 'Ge72' , 'Ge73' , 'Ge74' , &
'Ge76' , 'As75' , 'Se74' , 'Se76' , 'Se77' , 'Se78' , 'Se80' , 'Se82' , &
'Br79' , 'Br81' , 'Kr78' , 'Kr80' , 'Kr82' , 'Kr83' , 'Kr84' , 'Kr86' , &
'Rb85' , 'Rb87' , 'Sr84' , 'Sr86' , 'Sr87' , 'Sr88' , 'Y89' , 'Zr90' , &
'Zr91' , 'Zr92' , 'Zr94' , 'Zr96' , 'Nb93' , 'Mo92' , 'Mo94' , 'Mo95' , &
'Mo96' , 'Mo97' , 'Mo98' , 'Mo100', 'Ru96' , 'Ru98' , 'Ru99' , 'Ru100', &
'Ru101', 'Ru102', 'Ru104', 'Rh103', 'Pd102', 'Pd104', 'Pd105', 'Pd106', &
'Pd108', 'Pd110', 'Ag107', 'Ag109', 'Cd106', 'Cd108', 'Cd110', 'Cd111', &
'Cd112', 'Cd113', 'Cd114', 'Cd116', 'In113', 'In115', 'Sn112', 'Sn114', &
'Sn115', 'Sn116', 'Sn117', 'Sn118', 'Sn119', 'Sn120', 'Sn122', 'Sn124', &
'Sb121', 'Sb123', 'Te120', 'Te122', 'Te123', 'Te124', 'Te125', 'Te126', &
'Te128', 'Te130', 'I127' , 'Xe124', 'Xe126', 'Xe128', 'Xe129', 'Xe130', &
'Xe131', 'Xe132', 'Xe134', 'Xe136', 'Cs133', 'Ba130', 'Ba132', 'Ba134', &
'Ba135', 'Ba136', 'Ba137', 'Ba138', 'La138', 'La139', 'Ce136', 'Ce138', &
'Ce140', 'Ce142', 'Pr141', 'Nd142', 'Nd143', 'Nd144', 'Nd145', 'Nd146', &
'Nd148', 'Nd150', 'Sm144', 'Sm147', 'Sm148', 'Sm149', 'Sm150', 'Sm152', &
'Sm154', 'Eu151', 'Eu153', 'Gd152', 'Gd154', 'Gd155', 'Gd156', 'Gd157', &
'Gd158', 'Gd160', 'Tb159', 'Dy156', 'Dy158', 'Dy160', 'Dy161', 'Dy162', &
'Dy163', 'Dy164', 'Ho165', 'Er162', 'Er164', 'Er166', 'Er167', 'Er168', &
'Er170', 'Tm169', 'Yb168', 'Yb170', 'Yb171', 'Yb172', 'Yb173', 'Yb174', &
'Yb176', 'Lu175', 'Lu176', 'Hf174', 'Hf176', 'Hf177', 'Hf178', 'Hf179', &
'Hf180', 'Ta181', 'Ta180', 'Ta181', 'W182' , 'W183' , 'W184' , 'W186' , &
'W180' , 'W182' , 'W183' , 'W184' , 'W186' , 'Re185', 'Re187', 'Os0' , &
'Os184', 'Os186', 'Os187', 'Os188', 'Os189', 'Os190', 'Os192', 'Ir191', &
'Ir193', 'Pt0' , 'Pt190', 'Pt192', 'Pt194', 'Pt195', 'Pt196', 'Pt198', &
'Au197', 'Hg196', 'Hg198', 'Hg199', 'Hg200', 'Hg201', 'Hg202', 'Hg204', &
'Tl0' , 'Tl203', 'Tl205', 'Pb204', 'Pb206', 'Pb207', 'Pb208', 'Bi209', &
'Th232', 'Pa231', 'U234' , 'U235' , 'U238']
character (len=2), dimension (NUM_NATURAL_ELEMENTS), parameter :: &
natural_elements = [character(len=2) :: &
'h' , 'he', 'li', 'be', 'b' , 'c' , 'n' , 'o' , 'f' , 'ne', &
'na', 'mg', 'al', 'si', 'p' , 's' , 'cl', 'ar', 'k' , 'ca', &
'sc', 'ti', 'v' , 'cr', 'mn', 'fe', 'co', 'ni', 'cu', 'zn', &
'ga', 'ge', 'as', 'se', 'br', 'kr', 'rb', 'sr', 'y' , 'zr', &
'nb', 'mo', 'ru', 'rh', 'pd', 'ag', 'cd', 'in', 'sn', 'sb', &
'te', 'i' , 'xe', 'cs', 'ba', 'la', 'ce', 'pr', 'nd', 'sm', &
'eu', 'gd', 'tb', 'dy', 'ho', 'er', 'tm', 'yb', 'lu', 'hf', &
'ta', 'w' , 're', 'os', 'ir', 'pt', 'au', 'hg', 'tl', 'pb', &
'bi', 'th', 'pa', 'u']
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_ENDF_BVII0_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 17 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 253 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_ENDF_BVII0_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 19 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 57 , 63 , 64 , 68 , 69 , 74 , 76 , 77 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 257 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_ENDF_BVII1_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 17 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 57 , 59 , 63 , 64 , 68 , 69 , 74 , 77 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
251 , 257 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_ENDF_BVII1_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 19 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , 82 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
253 , 262 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_311_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 17 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , &
82 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 253 , 262 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_311_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 19 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 57 , 63 , 64 , 68 , 69 , 74 , 76 , 77 , &
83 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 257 , 264 , 265 , 274 , 275 , 282 , 289 , 290 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_312_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 17 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 57 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , &
82 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 253 , 262 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_312_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 19 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , 77 , &
83 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 257 , 264 , 265 , 274 , 275 , 282 , 289 , 290 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_32_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 13 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 77 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
251 , 257 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JEFF_32_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 16 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 57 , 63 , 64 , 68 , 69 , 74 , 76 , 82 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
253 , 262 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_32_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 16 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 77 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 253 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_32_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 17 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 57 , 63 , 64 , 68 , 69 , 74 , 76 , 82 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 257 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_33_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 16 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 77 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 253 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_33_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 17 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 57 , 63 , 64 , 68 , 69 , 74 , 76 , 82 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 257 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_40_start = (/ &
1 , 3 , 5 , 7 , 8 , 10 , 11 , 16 , 19 , 20 , &
23 , 24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , &
50 , 51 , 57 , 59 , 63 , 64 , 68 , 69 , 74 , 77 , &
83 , 85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , &
117 , 118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , &
163 , 171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , &
210 , 212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , &
250 , 257 , 262 , 265 , 272 , 275 , 281 , 282 , 290 , 292 , &
296 , 297 , 298 , 299 /)
integer, dimension(NUM_NATURAL_ELEMENTS) :: &
natural_nuclides_JENDL_40_end = (/ &
3 , 5 , 7 , 8 , 10 , 11 , 13 , 17 , 20 , 23 , &
24 , 27 , 28 , 31 , 32 , 36 , 38 , 41 , 44 , 50 , &
51 , 56 , 59 , 63 , 64 , 68 , 69 , 74 , 76 , 82 , &
85 , 90 , 90 , 97 , 99 , 105 , 107 , 111 , 112 , 117 , &
118 , 125 , 132 , 133 , 139 , 141 , 149 , 151 , 161 , 163 , &
171 , 172 , 181 , 182 , 189 , 191 , 195 , 196 , 203 , 210 , &
212 , 219 , 220 , 227 , 228 , 234 , 235 , 242 , 244 , 250 , &
251 , 262 , 264 , 272 , 274 , 281 , 282 , 289 , 292 , 296 , &
297 , 298 , 299 , 302 /)
real(8), dimension (NUM_NATURAL_NUCLIDES) :: natural_nuclides_mf = (/ &
0.999885_8 , 0.000115_8 , 0.00000134_8, 0.99999866_8, &
0.0759_8 , 0.9241_8 , ONE , 0.199_8 , &
0.801_8 , ONE , 0.99636_8 , 0.00364_8 , &
0.99757_8 , 0.00038_8 , 0.00205_8 , ONE , &
0.99962_8 , 0.00038_8 , ONE , 0.9048_8 , &
0.0027_8 , 0.0925_8 , ONE , 0.7899_8 , &
0.1000_8 , 0.1101_8 , ONE , 0.92223_8 , &
0.04685_8 , 0.03092_8 , ONE , 0.9499_8 , &
0.0075_8 , 0.0425_8 , 0.0001_8 , 0.7576_8 , &
0.2424_8 , 0.003336_8 , 0.000629_8 , 0.996035_8 , &
0.932581_8 , 0.000117_8 , 0.067302_8 , 0.96941_8 , &
0.00647_8 , 0.00135_8 , 0.02086_8 , 0.00004_8 , &
0.00187_8 , ONE , 0.0825_8 , 0.0744_8 , &
0.7372_8 , 0.0541_8 , 0.0518_8 , ONE , &
0.0025_8 , 0.9975_8 , 0.04345_8 , 0.83789_8 , &
0.09501_8 , 0.02365_8 , ONE , 0.05845_8 , &
0.91754_8 , 0.02119_8 , 0.00282_8 , ONE , &
0.68077_8 , 0.26223_8 , 0.011399_8 , 0.036346_8 , &
0.009255_8 , 0.6915_8 , 0.3085_8 , ONE , &
0.4917_8 , 0.2773_8 , 0.0404_8 , 0.1845_8 , &
0.0061_8 , ONE , 0.60108_8 , 0.39892_8 , &
0.2057_8 , 0.2745_8 , 0.0775_8 , 0.3650_8 , &
0.0773_8 , ONE , 0.0089_8 , 0.0937_8 , &
0.0763_8 , 0.2377_8 , 0.4961_8 , 0.0873_8 , &
0.5069_8 , 0.4931_8 , 0.00355_8 , 0.02286_8 , &
0.11593_8 , 0.11500_8 , 0.56987_8 , 0.17279_8 , &
0.7217_8 , 0.2783_8 , 0.0056_8 , 0.0986_8 , &
0.0700_8 , 0.8258_8 , ONE , 0.5145_8 , &
0.1122_8 , 0.1715_8 , 0.1738_8 , 0.0280_8 , &
ONE , 0.1453_8 , 0.0915_8 , 0.1584_8 , &
0.1667_8 , 0.0960_8 , 0.2439_8 , 0.0982_8 , &
0.0554_8 , 0.0187_8 , 0.1276_8 , 0.1260_8 , &
0.1706_8 , 0.3155_8 , 0.1862_8 , ONE , &
0.0102_8 , 0.1114_8 , 0.2233_8 , 0.2733_8 , &
0.2646_8 , 0.1172_8 , 0.51839_8 , 0.48161_8 , &
0.0125_8 , 0.0089_8 , 0.1249_8 , 0.1280_8 , &
0.2413_8 , 0.1222_8 , 0.2873_8 , 0.0749_8 , &
0.0429_8 , 0.9571_8 , 0.0097_8 , 0.0066_8 , &
0.0034_8 , 0.1454_8 , 0.0768_8 , 0.2422_8 , &
0.0859_8 , 0.3258_8 , 0.0463_8 , 0.0579_8 , &
0.5721_8 , 0.4279_8 , 0.0009_8 , 0.0255_8 , &
0.0089_8 , 0.0474_8 , 0.0707_8 , 0.1884_8 , &
0.3174_8 , 0.3408_8 , ONE , 0.000952_8 , &
0.000890_8 , 0.019102_8 , 0.264006_8 , 0.040710_8 , &
0.212324_8 , 0.269086_8 , 0.104357_8 , 0.088573_8 , &
ONE , 0.00106_8 , 0.00101_8 , 0.02417_8 , &
0.06592_8 , 0.07854_8 , 0.11232_8 , 0.71698_8 , &
0.0008881_8 , 0.9991119_8 , 0.00185_8 , 0.00251_8 , &
0.88450_8 , 0.11114_8 , ONE , 0.27152_8 , &
0.12174_8 , 0.23798_8 , 0.08293_8 , 0.17189_8 , &
0.05756_8 , 0.05638_8 , 0.0307_8 , 0.1499_8 , &
0.1124_8 , 0.1382_8 , 0.0738_8 , 0.2675_8 , &
0.2275_8 , 0.4781_8 , 0.5219_8 , 0.0020_8 , &
0.0218_8 , 0.1480_8 , 0.2047_8 , 0.1565_8 , &
0.2484_8 , 0.2186_8 , ONE , 0.00056_8 , &
0.00095_8 , 0.02329_8 , 0.18889_8 , 0.25475_8 , &
0.24896_8 , 0.28260_8 , ONE , 0.00139_8 , &
0.01601_8 , 0.33503_8 , 0.22869_8 , 0.26978_8 , &
0.14910_8 , ONE , 0.00123_8 , 0.02982_8 , &
0.1409_8 , 0.2168_8 , 0.16103_8 , 0.32026_8 , &
0.12996_8 , 0.97401_8 , 0.02599_8 , 0.0016_8 , &
0.0526_8 , 0.1860_8 , 0.2728_8 , 0.1362_8 , &
0.3508_8 , ONE , 0.0001201_8 , 0.9998799_8 , &
0.2662_8 , 0.1431_8 , 0.3064_8 , 0.2843_8 , &
0.0012_8 , 0.2650_8 , 0.1431_8 , 0.3064_8 , &
0.2843_8 , 0.3740_8 , 0.6260_8 , ONE , &
0.0002_8 , 0.0159_8 , 0.0196_8 , 0.1324_8 , &
0.1615_8 , 0.2626_8 , 0.4078_8 , 0.373_8 , &
0.627_8 , ONE , 0.00012_8 , 0.00782_8 , &
0.3286_8 , 0.3378_8 , 0.2521_8 , 0.07356_8 , &
ONE , 0.0015_8 , 0.0997_8 , 0.1687_8 , &
0.2310_8 , 0.1318_8 , 0.2986_8 , 0.0687_8 , &
ONE , 0.2952_8 , 0.7048_8 , 0.014_8 , &
0.241_8 , 0.221_8 , 0.524_8 , ONE , &
ONE , ONE , 0.000054_8 , 0.007204_8 , &
0.992742_8 /)
! ============================================================================
! TALLY-RELATED VARIABLES
@ -647,14 +317,14 @@ module global
logical :: restart_run = .false.
integer :: restart_batch
character(MAX_FILE_LEN) :: path_input ! Path to input file
character(MAX_FILE_LEN) :: path_cross_sections ! Path to cross_sections.xml
character(MAX_FILE_LEN) :: path_multipole ! Path to wmp library
character(MAX_FILE_LEN) :: path_source = '' ! Path to binary source
character(MAX_FILE_LEN) :: path_state_point ! Path to binary state point
character(MAX_FILE_LEN) :: path_source_point ! Path to binary source point
character(MAX_FILE_LEN) :: path_particle_restart ! Path to particle restart
character(MAX_FILE_LEN) :: path_output = '' ! Path to output directory
character(MAX_FILE_LEN) :: path_input ! Path to input file
character(MAX_FILE_LEN) :: path_cross_sections = '' ! Path to cross_sections.xml
character(MAX_FILE_LEN) :: path_multipole ! Path to wmp library
character(MAX_FILE_LEN) :: path_source = '' ! Path to binary source
character(MAX_FILE_LEN) :: path_state_point ! Path to binary state point
character(MAX_FILE_LEN) :: path_source_point ! Path to binary source point
character(MAX_FILE_LEN) :: path_particle_restart ! Path to particle restart
character(MAX_FILE_LEN) :: path_output = '' ! Path to output directory
! The verbosity controls how much information will be printed to the
! screen and in logs

View file

@ -81,7 +81,6 @@ contains
real(8), allocatable :: temp_real(:)
integer :: n_tracks
logical :: file_exists
character(MAX_FILE_LEN) :: env_variable
character(MAX_WORD_LEN) :: type
character(MAX_LINE_LEN) :: filename
type(Node), pointer :: doc => null()
@ -138,55 +137,24 @@ contains
end if
end if
! Find cross_sections.xml file -- the first place to look is the
! settings.xml file. If no file is found there, then we check the
! CROSS_SECTIONS environment variable
if (.not. check_for_node(doc, "cross_sections")) then
! No cross_sections.xml file specified in settings.xml, check
! environment variable
if (run_CE) then
call get_environment_variable("OPENMC_CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0) then
call get_environment_variable("CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0) then
call fatal_error("No cross_sections.xml file was specified in &
&settings.xml or in the OPENMC_CROSS_SECTIONS environment &
&variable. OpenMC needs such a file to identify where to &
&find ACE cross section libraries. Please consult the &
&user's guide at http://mit-crpg.github.io/openmc for &
&information on how to set up ACE cross section libraries.")
else
call warning("The CROSS_SECTIONS environment variable is &
&deprecated. Please update your environment to use &
&OPENMC_CROSS_SECTIONS instead.")
end if
end if
path_cross_sections = trim(env_variable)
else
call get_environment_variable("OPENMC_MG_CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0) then
call fatal_error("No mgxs.xml file was specified in &
&settings.xml or in the OPENMC_MG_CROSS_SECTIONS environment &
&variable. OpenMC needs such a file to identify where to &
&find ACE cross section libraries. Please consult the user's &
&guide at http://mit-crpg.github.io/openmc for information on &
&how to set up ACE cross section libraries.")
else
path_cross_sections = trim(env_variable)
end if
end if
else
! Look for deprecated cross_sections.xml file in settings.xml
if (check_for_node(doc, "cross_sections")) then
call warning("Setting cross_sections in settings.xml has been deprecated. &
&The cross_sections are now set in materials.xml and the &
&cross_sections input to materials.xml and OPENMC_CROSS_SECTIONS &
&environment variable will take precendent over setting &
&cross_sections in settings.xml.")
call get_node_value(doc, "cross_sections", path_cross_sections)
end if
! Find the windowed multipole library
! Look for deprecated windowed_multipole file in settings.xml
if (run_mode /= MODE_PLOTTING) then
if (.not. check_for_node(doc, "multipole_library")) then
! No library location specified in settings.xml, check
! environment variable
call get_environment_variable("OPENMC_MULTIPOLE_LIBRARY", env_variable)
path_multipole = trim(env_variable)
else
if (check_for_node(doc, "multipole_library")) then
call warning("Setting multipole_library in settings.xml has been &
&deprecated. The multipole_library is now set in materials.xml and&
& the multipole_library input to materials.xml and &
&OPENMC_MULTIPOLE_LIBRARY environment variable will take &
&precendent over setting multipole_library in settings.xml.")
call get_node_value(doc, "multipole_library", path_multipole)
end if
if (.not. ends_with(path_multipole, "/")) &
@ -1055,32 +1023,6 @@ contains
end if
end if
! Natural element expansion option
if (check_for_node(doc, "natural_elements")) then
call get_node_value(doc, "natural_elements", temp_str)
select case (to_lower(temp_str))
case ('endf/b-vii.0')
default_expand = ENDF_BVII0
case ('endf/b-vii.1')
default_expand = ENDF_BVII1
case ('jeff-3.1.1')
default_expand = JEFF_311
case ('jeff-3.1.2')
default_expand = JEFF_312
case ('jeff-3.2')
default_expand = JEFF_32
case ('jendl-3.2')
default_expand = JENDL_32
case ('jendl-3.3')
default_expand = JENDL_33
case ('jendl-4.0')
default_expand = JENDL_40
case default
call fatal_error("Unknown natural element expansion option: " &
// trim(temp_str))
end select
end if
call get_node_list(doc, "volume_calc", node_vol_list)
n = get_list_size(node_vol_list)
allocate(volume_calcs(n))
@ -1091,9 +1033,17 @@ contains
! Get temperature settings
if (check_for_node(doc, "temperature_default")) then
call warning("Setting temperature_default in settings.xml has been &
&deprecated. temperature_default is now set in materials.xml and &
&its input to materials.xml will take precendent &
&over setting it in settings.xml.")
call get_node_value(doc, "temperature_default", temperature_default)
end if
if (check_for_node(doc, "temperature_method")) then
call warning("Setting temperature_method in settings.xml has been &
&deprecated. temperature_method is now set in materials.xml and &
&its input to materials.xml will take precendent &
&over setting it in settings.xml.")
call get_node_value(doc, "temperature_method", temp_str)
select case (to_lower(temp_str))
case ('nearest')
@ -1105,9 +1055,17 @@ contains
end select
end if
if (check_for_node(doc, "temperature_tolerance")) then
call warning("Setting temperature_tolerance in settings.xml has been &
&deprecated. temperature_tolerance is now set in materials.xml and &
&its input to materials.xml will take precendent &
&over setting it in settings.xml.")
call get_node_value(doc, "temperature_tolerance", temperature_tolerance)
end if
if (check_for_node(doc, "temperature_multipole")) then
call warning("Setting temperature_multipole in settings.xml has been &
&deprecated. temperature_multipole is now set in materials.xml and &
&its input to materials.xml will take precendent &
&over setting it in settings.xml.")
call get_node_value(doc, "temperature_multipole", temp_str)
select case (to_lower(temp_str))
case ('true', '1')
@ -2095,7 +2053,115 @@ contains
type(Library), allocatable :: libraries(:)
type(VectorReal), allocatable :: nuc_temps(:) ! List of T to read for each nuclide
type(VectorReal), allocatable :: sab_temps(:) ! List of T to read for each S(a,b)
real(8), allocatable :: material_temps(:)
character(MAX_LINE_LEN) :: temp_str
real(8), allocatable :: material_temps(:)
logical :: file_exists
character(MAX_FILE_LEN) :: env_variable
character(MAX_LINE_LEN) :: filename
type(Node), pointer :: doc => null()
! Display output message
call write_message("Reading materials XML file...", 5)
! Check is materials.xml exists
filename = trim(path_input) // "materials.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
call fatal_error("Material XML file '" // trim(filename) // "' does not &
&exist!")
end if
! Parse materials.xml file
call open_xmldoc(doc, filename)
! Find cross_sections.xml file -- the first place to look is the
! materials.xml file. If no file is found there, then we check the
! CROSS_SECTIONS environment variable
if (.not. check_for_node(doc, "cross_sections")) then
! No cross_sections.xml file specified in settings.xml, check
! environment variable
if (run_CE) then
call get_environment_variable("OPENMC_CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0) then
call get_environment_variable("CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0 .and. path_cross_sections == '') then
call fatal_error("No cross_sections.xml file was specified in &
&materials.xml, settings.xml, or in the OPENMC_CROSS_SECTIONS&
& environment variable. OpenMC needs such a file to identify &
&where to find ACE cross section libraries. Please consult the&
& user's guide at http://mit-crpg.github.io/openmc for &
&information on how to set up ACE cross section libraries.")
else
call warning("The CROSS_SECTIONS environment variable is &
&deprecated. Please update your environment to use &
&OPENMC_CROSS_SECTIONS instead.")
end if
end if
path_cross_sections = trim(env_variable)
else
call get_environment_variable("OPENMC_MG_CROSS_SECTIONS", env_variable)
if (len_trim(env_variable) == 0 .and. path_cross_sections == '') then
call fatal_error("No mgxs.xml file was specified in &
&materials.xml or in the OPENMC_MG_CROSS_SECTIONS environment &
&variable. OpenMC needs such a file to identify where to &
&find ACE cross section libraries. Please consult the user's &
&guide at http://mit-crpg.github.io/openmc for information on &
&how to set up ACE cross section libraries.")
else
path_cross_sections = trim(env_variable)
end if
end if
else
call get_node_value(doc, "cross_sections", path_cross_sections)
end if
! Find the windowed multipole library
if (run_mode /= MODE_PLOTTING) then
if (.not. check_for_node(doc, "multipole_library")) then
! No library location specified in materials.xml, check
! environment variable
call get_environment_variable("OPENMC_MULTIPOLE_LIBRARY", env_variable)
path_multipole = trim(env_variable)
else
call get_node_value(doc, "multipole_library", path_multipole)
end if
if (.not. ends_with(path_multipole, "/")) &
path_multipole = trim(path_multipole) // "/"
end if
! Get temperature settings
if (check_for_node(doc, "temperature_default")) then
call get_node_value(doc, "temperature_default", temperature_default)
end if
if (check_for_node(doc, "temperature_method")) then
call get_node_value(doc, "temperature_method", temp_str)
select case (to_lower(temp_str))
case ('nearest')
temperature_method = TEMPERATURE_NEAREST
case ('interpolation')
temperature_method = TEMPERATURE_INTERPOLATION
case default
call fatal_error("Unknown temperature method: " // trim(temp_str))
end select
end if
if (check_for_node(doc, "temperature_tolerance")) then
call get_node_value(doc, "temperature_tolerance", temperature_tolerance)
end if
if (check_for_node(doc, "temperature_multipole")) then
call get_node_value(doc, "temperature_multipole", temp_str)
select case (to_lower(temp_str))
case ('true', '1')
temperature_multipole = .true.
case ('false', '0')
temperature_multipole = .false.
case default
call fatal_error("Unrecognized value for <use_windowed_multipole> in &
&materials.xml")
end select
end if
! Close materials XML file
call close_xmldoc(doc)
if (run_CE) then
call read_ce_cross_sections_xml(libraries)
@ -2153,10 +2219,8 @@ contains
integer :: i ! loop index for materials
integer :: j ! loop index for nuclides
integer :: k ! loop index for elements
integer :: n ! number of nuclides
integer :: n_sab ! number of sab tables for a material
integer :: n_nuc_ele ! number of nuclides in an element
integer :: i_library ! index in libraries array
integer :: index_nuclide ! index in nuclides
integer :: index_sab ! index in sab_tables
@ -2170,17 +2234,12 @@ contains
type(Node), pointer :: doc => null()
type(Node), pointer :: node_mat => null()
type(Node), pointer :: node_nuc => null()
type(Node), pointer :: node_ele => null()
type(Node), pointer :: node_sab => null()
type(NodeList), pointer :: node_mat_list => null()
type(NodeList), pointer :: node_nuc_list => null()
type(NodeList), pointer :: node_macro_list => null()
type(NodeList), pointer :: node_ele_list => null()
type(NodeList), pointer :: node_sab_list => null()
! Display output message
call write_message("Reading materials XML file...", 5)
! Check is materials.xml exists
filename = trim(path_input) // "materials.xml"
inquire(FILE=filename, EXIST=file_exists)
@ -2240,10 +2299,9 @@ contains
! Check to ensure material has at least one nuclide
if (.not. check_for_node(node_mat, "nuclide") .and. &
.not. check_for_node(node_mat, "element") .and. &
.not. check_for_node(node_mat, "macroscopic")) then
call fatal_error("No macroscopic data, nuclides or natural elements &
&specified on material " // trim(to_str(mat % id)))
call fatal_error("No macroscopic data or nuclides specified on &
&material " // trim(to_str(mat % id)))
end if
! Create list of macroscopic x/s based on those specified, just treat
@ -2279,7 +2337,7 @@ contains
! Get pointer list of XML <nuclide>
call get_node_list(node_mat, "nuclide", node_nuc_list)
! Create list of nuclides based on those specified plus natural elements
! Create list of nuclides based on those specified
INDIVIDUAL_NUCLIDES: do j = 1, get_list_size(node_nuc_list)
! Combine nuclide identifier and cross section and copy into names
call get_list_item(node_nuc_list, j, node_nuc)
@ -2317,66 +2375,6 @@ contains
end do INDIVIDUAL_NUCLIDES
end if
! =======================================================================
! READ AND PARSE <element> TAGS
! Get pointer list of XML <element>
call get_node_list(node_mat, "element", node_ele_list)
NATURAL_ELEMENTS: do j = 1, get_list_size(node_ele_list)
call get_list_item(node_ele_list, j, node_ele)
! Check for empty name on natural element
if (.not. check_for_node(node_ele, "name")) then
call fatal_error("No name specified on nuclide in material " &
// trim(to_str(mat % id)))
end if
call get_node_value(node_ele, "name", name)
! Check if no atom/weight percents were specified or if both atom and
! weight percents were specified
if (.not. check_for_node(node_ele, "ao") .and. &
.not. check_for_node(node_ele, "wo")) then
call fatal_error("No atom or weight percent specified for element " &
// trim(name))
elseif (check_for_node(node_ele, "ao") .and. &
check_for_node(node_ele, "wo")) then
call fatal_error("Cannot specify both atom and weight percents for &
&element: " // trim(name))
end if
! Get current number of nuclides
n_nuc_ele = names % size()
! Expand element into naturally-occurring nuclides
call expand_natural_element_names(name, names)
! Compute number of new nuclides from the natural element expansion
n_nuc_ele = names % size() - n_nuc_ele
! Check enforced isotropic lab scattering
if (run_CE) then
if (check_for_node(node_ele, "scattering")) then
call get_node_value(node_ele, "scattering", temp_str)
else
temp_str = "data"
end if
! Set ace or iso-in-lab scattering for each nuclide in element
do k = 1, n_nuc_ele
if (adjustl(to_lower(temp_str)) == "iso-in-lab") then
call list_iso_lab % push_back(1)
else if (adjustl(to_lower(temp_str)) == "data") then
call list_iso_lab % push_back(0)
else
call fatal_error("Scattering must be isotropic in lab or follow&
& the ACE file data")
end if
end do
end if
end do NATURAL_ELEMENTS
! ========================================================================
! COPY NUCLIDES TO ARRAYS IN MATERIAL
@ -4588,205 +4586,6 @@ contains
end subroutine read_mg_cross_sections_header
!===============================================================================
! EXPAND_NATURAL_ELEMENT_NAMES converts natural elements specified using an
! <element> tag within a material and adds the names to the input names array.
! In some cases, modifications have been made to work with ENDF/B-VII.1 where
! evaluations of particular nuclides don't exist.
!===============================================================================
subroutine expand_natural_element_names(name, names)
character(*), intent(in) :: name
type(VectorChar), intent(inout) :: names
integer :: i
character(2) :: element_name
integer :: natural_elements_loc
integer :: nuclide_start
integer :: nuclide_end
! Convert the element name to lower case
element_name = to_lower(name(1:2))
! Find location of element name in natural_elements
i = -1
do i = 1, NUM_NATURAL_ELEMENTS
if (natural_elements(i) == element_name) then
natural_elements_loc = i
end if
end do
! Issue error if element not found
if (i == -1) then
call fatal_error("Cannot expand element: " // name)
end if
! Get start and end locations in natural_nuclides arrays
select case (default_expand)
case (ENDF_BVII0)
nuclide_start = natural_nuclides_ENDF_BVII0_start(natural_elements_loc)
nuclide_end = natural_nuclides_ENDF_BVII0_end(natural_elements_loc)
case (ENDF_BVII1)
nuclide_start = natural_nuclides_ENDF_BVII1_start(natural_elements_loc)
nuclide_end = natural_nuclides_ENDF_BVII1_end(natural_elements_loc)
case (JEFF_311)
nuclide_start = natural_nuclides_JEFF_311_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_311_end(natural_elements_loc)
case (JEFF_312)
nuclide_start = natural_nuclides_JEFF_312_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_312_end(natural_elements_loc)
case (JEFF_32)
nuclide_start = natural_nuclides_JEFF_32_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_32_end(natural_elements_loc)
case (JENDL_32)
nuclide_start = natural_nuclides_JENDL_32_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_32_end(natural_elements_loc)
case (JENDL_33)
nuclide_start = natural_nuclides_JENDL_33_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_33_end(natural_elements_loc)
case (JENDL_40)
nuclide_start = natural_nuclides_JENDL_40_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_40_end(natural_elements_loc)
end select
! Add the nuclide names to the names array
do i = nuclide_start, nuclide_end-1
call names % push_back(natural_nuclides(i))
end do
end subroutine expand_natural_element_names
!===============================================================================
! EXPAND_NATURAL_ELEMENT_DENSITIES converts natural elements specified using an
! <element> tag within a material into individual nuclides based on IUPAC
! Isotopic Compositions of the Elements 2009 (doi:10.1351/PAC-REP-10-06-02). In
! some cases, modifications have been made to work with ENDF/B-VII.1 where
! evaluations of particular nuclides don't exist.
!===============================================================================
subroutine expand_natural_element_densities(name, expand_by, density, &
densities)
character(*), intent(in) :: name ! element name
character(*), intent(in) :: expand_by ! "ao" or "wo"
real(8), intent(in) :: density ! value for "ao" or "wo"
type(VectorReal), intent(inout) :: densities ! nuclide densities vector
integer :: i ! iterators
integer :: n_nuclides ! number of nuclides in the element
integer :: natural_elements_loc ! location in global array
integer :: nuclide_start ! ending nuclide in global array
integer :: nuclide_end ! starting nuclide in global array
character(2) :: element_name ! element atomic symbol
real(8) :: element_awr ! element atomic weight ratio
real(8), allocatable :: awr(:) ! nuclide atomic weight ratios
real(8), allocatable :: mf(:) ! nuclide mole or mass fractions
element_name = to_lower(name(1:2))
! Find location of element name in natural_elements
i = -1
do i = 1, NUM_NATURAL_ELEMENTS
if (natural_elements(i) == element_name) then
natural_elements_loc = i
end if
end do
! Issue error if element not found
if (i == -1) then
call fatal_error("Cannot expand element: " // name)
end if
! Get start and end locations in natural_nuclides arrays
select case (default_expand)
case (ENDF_BVII0)
nuclide_start = natural_nuclides_ENDF_BVII0_start(natural_elements_loc)
nuclide_end = natural_nuclides_ENDF_BVII0_end(natural_elements_loc)
case (ENDF_BVII1)
nuclide_start = natural_nuclides_ENDF_BVII1_start(natural_elements_loc)
nuclide_end = natural_nuclides_ENDF_BVII1_end(natural_elements_loc)
case (JEFF_311)
nuclide_start = natural_nuclides_JEFF_311_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_311_end(natural_elements_loc)
case (JEFF_312)
nuclide_start = natural_nuclides_JEFF_312_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_312_end(natural_elements_loc)
case (JEFF_32)
nuclide_start = natural_nuclides_JEFF_32_start(natural_elements_loc)
nuclide_end = natural_nuclides_JEFF_32_end(natural_elements_loc)
case (JENDL_32)
nuclide_start = natural_nuclides_JENDL_32_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_32_end(natural_elements_loc)
case (JENDL_33)
nuclide_start = natural_nuclides_JENDL_33_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_33_end(natural_elements_loc)
case (JENDL_40)
nuclide_start = natural_nuclides_JENDL_40_start(natural_elements_loc)
nuclide_end = natural_nuclides_JENDL_40_end(natural_elements_loc)
end select
! Set the number of nuclides in this element
n_nuclides = nuclide_end - nuclide_start
! Allocate and initialize array for mole fractions
allocate(mf(n_nuclides))
! Add the nuclide names to the names array
do i = nuclide_start, nuclide_end-1
mf(i - nuclide_start + 1) = natural_nuclides_mf(i)
end do
! Compute the ratio of the nuclide atomic weights to the element atomic
! weight
if (expand_by == "wo") then
! Allocate array for the atomic weight ratios
allocate(awr(n_nuclides))
! Get the atomic weight ratios
do i = nuclide_start, nuclide_end-1
awr(i - nuclide_start + 1) = nuclides(nuclide_dict % &
get_key(to_lower(natural_nuclides(i)))) % awr
end do
! Compute the element awr
element_awr = sum(awr * mf)
! Convert the mole fractions to mass fractions
mf = mf * awr / element_awr
! Normalize the mass fractions to ONE
mf = mf / sum(mf)
! Deallocate array for the atomic weight ratios
deallocate(awr)
end if
! Add the densities to the master array
do i = 1, n_nuclides
call densities % push_back(density * mf(i))
end do
! Deallocate array for mole or mass fractions
deallocate(mf)
end subroutine expand_natural_element_densities
!===============================================================================
! GENERATE_RPN implements the shunting-yard algorithm to generate a Reverse
! Polish notation (RPN) expression for the region specification of a cell given
@ -4894,11 +4693,9 @@ contains
type(Node), pointer :: node_mat => null()
type(Node), pointer :: node_dens => null()
type(Node), pointer :: node_nuc => null()
type(Node), pointer :: node_ele => null()
type(NodeList), pointer :: node_mat_list => null()
type(NodeList), pointer :: node_nuc_list => null()
type(NodeList), pointer :: node_macro_list => null()
type(NodeList), pointer :: node_ele_list => null()
! Display output message
call write_message("Reading material densities from XML file...", 5)
@ -5006,7 +4803,7 @@ contains
! Get pointer list of XML <nuclide>
call get_node_list(node_mat, "nuclide", node_nuc_list)
! Create list of nuclides based on those specified plus natural elements
! Create list of nuclides based on those specified
INDIVIDUAL_NUCLIDES: do j = 1, get_list_size(node_nuc_list)
! Combine nuclide identifier and cross section and copy into names
@ -5043,42 +4840,6 @@ contains
end do INDIVIDUAL_NUCLIDES
end if
! =======================================================================
! READ AND PARSE <element> TAGS
! Get pointer list of XML <element>
call get_node_list(node_mat, "element", node_ele_list)
NATURAL_ELEMENTS: do j = 1, get_list_size(node_ele_list)
call get_list_item(node_ele_list, j, node_ele)
call get_node_value(node_ele, "name", name)
! Check if no atom/weight percents were specified or if both atom and
! weight percents were specified
if (.not. check_for_node(node_ele, "ao") .and. &
.not. check_for_node(node_ele, "wo")) then
call fatal_error("No atom or weight percent specified for element " &
// trim(name))
elseif (check_for_node(node_ele, "ao") .and. &
check_for_node(node_ele, "wo")) then
call fatal_error("Cannot specify both atom and weight percents for &
&element: " // trim(name))
end if
! Expand element into naturally-occurring nuclides
if (check_for_node(node_ele, "ao")) then
call get_node_value(node_ele, "ao", temp_dble)
call expand_natural_element_densities(name, "ao", temp_dble, &
densities)
else
call get_node_value(node_ele, "wo", temp_dble)
call expand_natural_element_densities(name, "wo", -temp_dble, &
densities)
end if
end do NATURAL_ELEMENTS
! ========================================================================
! SET MATERIAL ATOM DENSITIES

View file

@ -757,6 +757,7 @@ class MGInputSet(InputSet):
# Define the materials file
self.xs_data = xs
self.materials += mats
self.materials.cross_sections = "../1d_mgxs.h5"
# Define surfaces.
# Assembly/Problem Boundary
@ -793,7 +794,6 @@ class MGInputSet(InputSet):
self.settings.source = Source(space=Box([0.0, 0.0, 0.0],
[10.0, 10.0, 5.]))
self.settings.energy_mode = "multi-group"
self.settings.cross_sections = "../1d_mgxs.h5"
def build_defualt_plots(self):
plot = openmc.Plot()

View file

@ -11,11 +11,11 @@
</material>
<material id="3">
<density value="2.0" units="g/cc" />
<element name="Zr" ao="1.0" />
<nuclide name="Zr90" ao="1.0" />
</material>
<material id="4">
<density value="0.1" units="g/cc" />
<element name="N" ao="1.0" />
<nuclide name="N14" ao="1.0" />
</material>
</materials>

View file

@ -1,11 +1,11 @@
k-combined:
2.531110E-01 3.041974E-03
2.511523E-01 2.296777E-03
tally 1:
2.594626E+00
1.346701E+00
2.683653E+00
1.440725E+00
9.933862E-01
1.977011E-01
1.112289E-01
2.476655E-03
2.578905E+00
1.331155E+00
2.688693E+00
1.447930E+00
9.863774E-01
1.948549E-01
1.123802E-01
2.527538E-03

View file

@ -1 +1 @@
0d54ed29b2a348854d99d36f602b188766465b41e70536647578c61d076909a2351b5c4cafbbdc62dd00f36bf2685f231fae7c75d1289befddd2da79042e31a7
c472877dd96f51dcfd7b45ec1c230c3847b26eb1f0a97751bb175b70eafff2ad8ee97e3ed8691a3bc8bf8f27c074a9d34021322988185bb04c3e55958f40589b

View file

@ -1,5 +1,5 @@
k-combined:
8.617147E-01 2.915803E-02
8.781218E-01 2.426380E-02
tally 1:
7.879002E+01
1.242562E+03
8.505928E+01
1.447208E+03

View file

@ -77,7 +77,6 @@ class ElementWOTestHarness(PyAPITestHarness):
# Instantiate a Geometry, register the root Universe, and export to XML
self._input_set.geometry.root_universe = root
mat_filter = openmc.MaterialFilter((fuel.id,))
flux_tally = openmc.Tally()
flux_tally.filters = [mat_filter]

View file

@ -1 +1 @@
47d753380f995651e6cec8e3648ff9e4ac254fc80919014124f4f6674afc95262006eda3821b4d6dcf09e29f28f5e21840c9b76b95c1c7f8ae522faac2386090
6f8d495cd7537f8ca29355664b31d429fe13f656a868852769b27ef9566e7bb683973322e41a92bab65beb1ba6cf477698b30a7ae6bc17d2e34cb44751e8c81d

View file

@ -1,5 +1,5 @@
k-combined:
1.385412E+00 1.198466E-01
1.508419E+00 3.763047E-02
tally 1:
6.197393E+01
7.688717E+02
6.274668E+01
7.902869E+02

View file

@ -73,7 +73,6 @@ class EnrichmentTestHarness(PyAPITestHarness):
# Instantiate a Geometry, register the root Universe, and export to XML
self._input_set.geometry.root_universe = root
mat_filter = openmc.MaterialFilter((fuel.id,))
flux_tally = openmc.Tally()
flux_tally.filters = [mat_filter]

View file

@ -1 +1 @@
115218221500e60ea43145875324a9a6800366fcc97357ef4cd8182f9253e114cbe1c7481dd3bdd8d7dce5983e2cb757ae0a2e05ede236d94cba7909703ba7f7
12abd69924751d1cf2a296ce799bfb2a6a4e744ffc49eef5a5cee8b764254caa3783c406cff71238c5608b454f451953479fb4caad45460a013225ab8de2271e

View file

@ -1 +1 @@
139d760cd83eab001ed3fd52d7146fb5f30b7021b26ada6a1e425f5446185ac6fa61ac2ac58aa9e895981e10540160b7ae7428c833024679ef1ca762715e2d51
1f82e622b36562a0a14536fe24c95ee9c901c428f7f082f245afae5c699ccd647aaf31a49790822ceb7ec136bdb75f8e25c3677642f73702ca2e30b058f6b440

View file

@ -1 +1 @@
6d03b988671fe21dfe2a908685536edb3059e5bdda7974ddfa65c7f63cecb7f6b3fce49bc1533ef2787e284121779653cb14a66cfd709258e4b07d8bd7571843
ec7f0697d22668dcc4cdaf134c58b9f84137730d52cb44264ee9205369abcf995b2f66847d620715cd609bd2f5f865294374978f8fb89f182962e40801171cba

View file

@ -1 +1 @@
8c496e13f9456bedcc0789ffdca81913596cca3c7867717ced38c3063b6cdc8949bd204c98723b72c1b456e3fc85fc90e6710482b0c74e5e29283902a5d9ed50
54cd9db30bc8e671591d76e8da2b5dead54eb4cead1eef09a005e96ca2c0fd707872e0cc21f31d198915b156605de78322fbb7957698ea8febceb23f048e44ad

View file

@ -60,8 +60,8 @@ class MGXSTestHarness(PyAPITestHarness):
self._input_set.mgxs_file, self._input_set.materials, \
self._input_set.geometry = self.mgxs_lib.create_mg_mode()
# Modify settings so we can run in MG mode
self._input_set.settings.cross_sections = './mgxs.h5'
# Modify materials and settings so we can run in MG mode
self._input_set.materials.cross_sections = './mgxs.h5'
self._input_set.settings.energy_mode = 'multi-group'
# Write modified input files

View file

@ -1,94 +0,0 @@
<?xml version="1.0"?>
<geometry>
<!-- pu-met-fast-021 case 2 -->
<surface id="1" type="z-plane" coeffs="-17.51" boundary="vacuum" />
<surface id="2" type="z-plane" coeffs="-2.57" />
<!-- Bottom Plutonium Core -->
<surface id="3" type="z-plane" coeffs="-2.55" />
<surface id="4" type="z-plane" coeffs="-2.37" />
<surface id="5" type="z-plane" coeffs="-2.1" />
<surface id="6" type="z-plane" coeffs="-2.06" />
<surface id="7" type="z-plane" coeffs="-1.61" />
<surface id="8" type="z-plane" coeffs="-1.57" />
<surface id="9" type="z-plane" coeffs="-1.12" />
<surface id="10" type="z-plane" coeffs="-1.08" />
<surface id="11" type="z-plane" coeffs="-0.63" />
<surface id="12" type="z-plane" coeffs="-0.59" />
<surface id="13" type="z-plane" coeffs="-0.14" />
<!-- Middle Gap -->
<surface id="14" type="z-plane" coeffs="-0.12" />
<surface id="15" type="z-plane" coeffs="0.12" />
<!-- Top Plutonium Core -->
<surface id="16" type="z-plane" coeffs="0.14" />
<surface id="17" type="z-plane" coeffs="0.59" />
<surface id="18" type="z-plane" coeffs="0.63" />
<surface id="19" type="z-plane" coeffs="1.08" />
<surface id="20" type="z-plane" coeffs="1.12" />
<surface id="21" type="z-plane" coeffs="1.57" />
<surface id="22" type="z-plane" coeffs="1.61" />
<surface id="23" type="z-plane" coeffs="2.06" />
<surface id="24" type="z-plane" coeffs="2.1" />
<surface id="25" type="z-plane" coeffs="2.37" />
<surface id="26" type="z-plane" coeffs="2.55" />
<surface id="27" type="z-plane" coeffs="2.57" />
<!-- Top Beryllium Reflector -->
<surface id="28" type="z-plane" coeffs="2.59" />
<surface id="29" type="z-plane" coeffs="17.51" boundary="vacuum" />
<!-- Cylinders -->
<surface id="30" type="z-cylinder" coeffs="0 0 9.995" boundary="vacuum" />
<surface id="31" type="z-cylinder" coeffs="0 0 6.263" />
<surface id="32" type="z-cylinder" coeffs="0 0 6.063" />
<surface id="33" type="z-cylinder" coeffs="0 0 5.995" />
<!-- Beryllium Reflector -->
<cell id="1" material="2" region=" 1 -2 -30" />
<cell id="2" material="2" region="28 -29 -30" />
<!-- Duralumin base and top -->
<cell id="3" material="4" region=" 2 -4 -30 31" />
<cell id="4" material="4" region="25 -28 -30 31" />
<!-- Duralumin top cylindrical shell -->
<cell id="6" material="5" region="15 -28 -31 32" />
<!-- Duralumin bottom cylindrical shell -->
<cell id="7" material="6" region=" 2 -14 -31 32" />
<!-- Steel Cover -->
<cell id="8" material="3" region=" 2 -3 -33" />
<cell id="9" material="3" region=" 5 -6 -33" />
<cell id="10" material="3" region=" 7 -8 -33" />
<cell id="11" material="3" region=" 9 -10 -33" />
<cell id="12" material="3" region="11 -12 -33" />
<cell id="13" material="3" region="13 -14 -33" />
<cell id="14" material="3" region="15 -16 -33" />
<cell id="15" material="3" region="17 -18 -33" />
<cell id="16" material="3" region="19 -20 -33" />
<cell id="17" material="3" region="21 -22 -33" />
<cell id="18" material="3" region="23 -24 -33" />
<cell id="19" material="3" region="26 -27 -33" />
<cell id="20" material="3" region=" 2 -14 33 -32" />
<cell id="21" material="3" region="15 -27 33 -32" />
<!-- Plutonium -->
<cell id="22" material="1" region=" 3 -5 -33" />
<cell id="23" material="1" region=" 6 -7 -33" />
<cell id="24" material="1" region=" 8 -9 -33" />
<cell id="25" material="1" region="10 -11 -33" />
<cell id="26" material="1" region="12 -13 -33" />
<cell id="27" material="1" region="16 -17 -33" />
<cell id="28" material="1" region="18 -19 -33" />
<cell id="29" material="1" region="20 -21 -33" />
<cell id="30" material="1" region="22 -23 -33" />
<cell id="31" material="1" region="24 -26 -33" />
<!-- Void -->
<cell id="32" material="void" region="14 -15 -31" />
<cell id="33" material="void" region=" 4 -25 31 -30" />
<cell id="34" material="void" region="27 -28 -30" />
</geometry>

View file

@ -1,62 +0,0 @@
<?xml version="1.0"?>
<materials>
<!-- pu-met-fast-021 case 2 -->
<!-- Plutonium -->
<material id="1">
<density value="18.5345" units="g/cm3" />
<nuclide name="Pu239" ao="4.4422e-02" />
<nuclide name="Pu240" ao="2.1326e-03" />
<nuclide name="Pu241" ao="9.2538e-05" />
<element name="C" ao="1.9515e-04" />
<element name="Fe" ao="8.1943e-05" />
</material>
<!-- Berylium-oxide Reflector -->
<material id="2">
<density value="2.86750" units="g/cm3" />
<nuclide name="Be9" ao="6.9041e-02" />
<element name="O" ao="6.9041e-02" />
</material>
<!-- Steel Cover -->
<material id="3">
<density value="6.9323" units="g/cm3" />
<element name="Fe" ao="5.1280e-02" />
<element name="C" ao="3.4757e-04" />
<element name="Si" ao="8.9185e-04" />
<element name="Ti" ao="6.1034e-04" />
<element name="Cr" ao="1.4452e-02" />
<element name="Mn" ao="1.5198e-03" />
<element name="Ni" ao="7.1131e-03" />
</material>
<!-- Duralumin base and top -->
<material id="4">
<density value="2.78" units="g/cm3" />
<element name="Al" ao="5.8077e-02" />
<element name="Mg" ao="1.0332e-03" />
<element name="Mn" ao="1.8284e-04" />
<element name="Cu" ao="1.1329e-03" />
</material>
<!-- Duralumin top cylindrical shell -->
<material id="5">
<density value="0.417" units="g/cm3" />
<element name="Al" ao="8.7115e-03" />
<element name="Mg" ao="1.5498e-04" />
<element name="Mn" ao="2.7426e-05" />
<element name="Cu" ao="1.6993e-04" />
</material>
<!-- Duralumin bottom cylindrical shell -->
<material id="6">
<density value="1.8155" units="g/cm3" />
<element name="Al" ao="3.7928e-02" />
<element name="Mg" ao="6.7475e-04" />
<element name="Mn" ao="1.1941e-04" />
<element name="Cu" ao="7.3982e-04" />
</material>
</materials>

View file

@ -1,2 +0,0 @@
k-combined:
1.034427E+00 1.583807E-02

View file

@ -1,22 +0,0 @@
<?xml version="1.0"?>
<settings>
<!-- pu-met-fast-021 case 2 -->
<eigenvalue>
<batches>10</batches>
<inactive>5</inactive>
<particles>400</particles>
</eigenvalue>
<natural_elements>endf/b-vii.1</natural_elements>
<source>
<space>
<type>box</type>
<parameters>0 0 0 1 1 1</parameters>
</space>
</source>
</settings>

View file

@ -1,11 +0,0 @@
#!/usr/bin/env python
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness
if __name__ == '__main__':
harness = TestHarness('statepoint.10.*')
harness.main()

View file

@ -1 +1 @@
b22973093e2b0690b30fb1262a11e27004555b796c446d256cb58a1d7329888ab60c1b93729b3d74bb015b98aa434416daa2dc2aee526eb8df0e9052911f94b4
3fdba48bb553aa75723e9ab226aeb873d4ebfc5934081263cf27063645e63036ec1088f1977fa38f45a95ed793bdb7b2a0e7a1f948b2a5b5b366110df72e8f0d

View file

@ -1,6 +0,0 @@
<?xml version="1.0"?>
<plots>
<plot id="1" type="slice" basis="xy" color="material"
origin="0.0 0.0 0.0" width="1.0 1.0" pixels="400 400">
</plot>
</plots>