mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-28 14:15:42 -04:00
Merge pull request #1955 from paulromano/from-xml-complete
Add missing from_xml capabilities
This commit is contained in:
Patrick Shriwise
GitHub
commit
6450d669cb
22 changed files with 714 additions and 141 deletions
102
openmc/filter.py
102
openmc/filter.py
|
|
@ -16,6 +16,7 @@ from .material import Material
|
|||
from .mixin import IDManagerMixin
|
||||
from .surface import Surface
|
||||
from .universe import UniverseBase
|
||||
from ._xml import get_text
|
||||
|
||||
|
||||
_FILTER_TYPES = (
|
||||
|
|
@ -231,9 +232,43 @@ class Filter(IDManagerMixin, metaclass=FilterMeta):
|
|||
|
||||
subelement = ET.SubElement(element, 'bins')
|
||||
subelement.text = ' '.join(str(b) for b in self.bins)
|
||||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
"""Generate a filter from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
**kwargs
|
||||
Keyword arguments (e.g., mesh information)
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Filter
|
||||
Filter object
|
||||
|
||||
"""
|
||||
filter_type = elem.get('type')
|
||||
|
||||
# If the filter type matches this class's short_name, then
|
||||
# there is no overriden from_xml_element method
|
||||
if filter_type == cls.short_name.lower():
|
||||
# Get bins from element -- the default here works for any filters
|
||||
# that just store a list of bins that can be represented as integers
|
||||
filter_id = int(elem.get('id'))
|
||||
bins = [int(x) for x in get_text(elem, 'bins').split()]
|
||||
return cls(bins, filter_id=filter_id)
|
||||
|
||||
# Search through all subclasses and find the one matching the HDF5
|
||||
# 'type'. Call that class's from_hdf5 method
|
||||
for subclass in cls._recursive_subclasses():
|
||||
if filter_type == subclass.short_name.lower():
|
||||
return subclass.from_xml_element(elem, **kwargs)
|
||||
|
||||
|
||||
def can_merge(self, other):
|
||||
"""Determine if filter can be merged with another.
|
||||
|
||||
|
|
@ -622,6 +657,13 @@ class CellInstanceFilter(Filter):
|
|||
subelement.text = ' '.join(str(i) for i in self.bins.ravel())
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
bins = [int(x) for x in get_text(elem, 'bins').split()]
|
||||
cell_instances = list(zip(bins[::2], bins[1::2]))
|
||||
return cls(cell_instances, filter_id=filter_id)
|
||||
|
||||
|
||||
class SurfaceFilter(WithIDFilter):
|
||||
"""Filters particles by surface crossing
|
||||
|
|
@ -661,8 +703,8 @@ class ParticleFilter(Filter):
|
|||
|
||||
Attributes
|
||||
----------
|
||||
bins : Iterable of Integral
|
||||
The Particles to tally
|
||||
bins : iterable of str
|
||||
The particles to tally
|
||||
id : int
|
||||
Unique identifier for the filter
|
||||
num_bins : Integral
|
||||
|
|
@ -698,6 +740,12 @@ class ParticleFilter(Filter):
|
|||
filter_id = int(group.name.split('/')[-1].lstrip('filter '))
|
||||
return cls(particles, filter_id=filter_id)
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
bins = get_text(elem, 'bins').split()
|
||||
return cls(bins, filter_id=filter_id)
|
||||
|
||||
|
||||
class MeshFilter(Filter):
|
||||
"""Bins tally event locations onto a regular, rectangular mesh.
|
||||
|
|
@ -877,6 +925,18 @@ class MeshFilter(Filter):
|
|||
element.set('translation', ' '.join(map(str, self.translation)))
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
mesh_id = int(get_text(elem, 'bins'))
|
||||
mesh_obj = kwargs['meshes'][mesh_id]
|
||||
filter_id = int(elem.get('id'))
|
||||
out = cls(mesh_obj, filter_id=filter_id)
|
||||
|
||||
translation = elem.get('translation')
|
||||
if translation:
|
||||
out.translation = [float(x) for x in translation.split()]
|
||||
return out
|
||||
|
||||
|
||||
class MeshSurfaceFilter(MeshFilter):
|
||||
"""Filter events by surface crossings on a regular, rectangular mesh.
|
||||
|
|
@ -1019,35 +1079,12 @@ class CollisionFilter(Filter):
|
|||
self.bins = np.asarray(bins)
|
||||
self.id = filter_id
|
||||
|
||||
def __repr__(self):
|
||||
string = type(self).__name__ + '\n'
|
||||
string += '{: <16}=\t{}\n'.format('\tValues', self.bins)
|
||||
string += '{: <16}=\t{}\n'.format('\tID', self.id)
|
||||
return string
|
||||
|
||||
@Filter.bins.setter
|
||||
def bins(self, bins):
|
||||
Filter.bins.__set__(self, np.asarray(bins))
|
||||
|
||||
def check_bins(self, bins):
|
||||
for x in bins:
|
||||
# Values should be integers
|
||||
cv.check_type('filter value', x, Integral)
|
||||
cv.check_greater_than('filter value', x, 0, equality=True)
|
||||
|
||||
def to_xml_element(self):
|
||||
"""Return XML Element representing the Filter.
|
||||
|
||||
Returns
|
||||
-------
|
||||
element : xml.etree.ElementTree.Element
|
||||
XML element containing filter data
|
||||
|
||||
"""
|
||||
element = super().to_xml_element()
|
||||
element[0].text = ' '.join(str(x) for x in self.bins)
|
||||
return element
|
||||
|
||||
|
||||
class RealFilter(Filter):
|
||||
"""Tally modifier that describes phase-space and other characteristics
|
||||
|
|
@ -1236,6 +1273,12 @@ class RealFilter(Filter):
|
|||
element[0].text = ' '.join(str(x) for x in self.values)
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
bins = [float(x) for x in get_text(elem, 'bins').split()]
|
||||
return cls(bins, filter_id=filter_id)
|
||||
|
||||
|
||||
class EnergyFilter(RealFilter):
|
||||
"""Bins tally events based on incident particle energy.
|
||||
|
|
@ -1969,6 +2012,13 @@ class EnergyFunctionFilter(Filter):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
energy = [float(x) for x in get_text(elem, 'energy').split()]
|
||||
y = [float(x) for x in get_text(elem, 'y').split()]
|
||||
return cls(energy, y, filter_id=filter_id)
|
||||
|
||||
def can_merge(self, other):
|
||||
return False
|
||||
|
||||
|
|
|
|||
|
|
@ -46,6 +46,12 @@ class ExpansionFilter(Filter):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
order = int(elem.find('order').text)
|
||||
return cls(order, filter_id=filter_id)
|
||||
|
||||
|
||||
class LegendreFilter(ExpansionFilter):
|
||||
r"""Score Legendre expansion moments up to specified order.
|
||||
|
|
@ -226,6 +232,15 @@ class SpatialLegendreFilter(ExpansionFilter):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
order = int(elem.find('order').text)
|
||||
axis = elem.find('axis').text
|
||||
minimum = float(elem.find('min').text)
|
||||
maximum = float(elem.find('max').text)
|
||||
return cls(order, axis, minimum, maximum, filter_id=filter_id)
|
||||
|
||||
|
||||
class SphericalHarmonicsFilter(ExpansionFilter):
|
||||
r"""Score spherical harmonic expansion moments up to specified order.
|
||||
|
|
@ -316,6 +331,14 @@ class SphericalHarmonicsFilter(ExpansionFilter):
|
|||
element.set('cosine', self.cosine)
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
order = int(elem.find('order').text)
|
||||
filter = cls(order, filter_id=filter_id)
|
||||
filter.cosine = elem.get('cosine')
|
||||
return filter
|
||||
|
||||
|
||||
class ZernikeFilter(ExpansionFilter):
|
||||
r"""Score Zernike expansion moments in space up to specified order.
|
||||
|
|
@ -358,7 +381,7 @@ class ZernikeFilter(ExpansionFilter):
|
|||
x-coordinate of center of circle for normalization
|
||||
y : float
|
||||
y-coordinate of center of circle for normalization
|
||||
r : int or None
|
||||
r : float
|
||||
Radius of circle for normalization
|
||||
|
||||
Attributes
|
||||
|
|
@ -369,7 +392,7 @@ class ZernikeFilter(ExpansionFilter):
|
|||
x-coordinate of center of circle for normalization
|
||||
y : float
|
||||
y-coordinate of center of circle for normalization
|
||||
r : int or None
|
||||
r : float
|
||||
Radius of circle for normalization
|
||||
id : int
|
||||
Unique identifier for the filter
|
||||
|
|
@ -464,6 +487,15 @@ class ZernikeFilter(ExpansionFilter):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
filter_id = int(elem.get('id'))
|
||||
order = int(elem.find('order').text)
|
||||
x = float(elem.find('x').text)
|
||||
y = float(elem.find('y').text)
|
||||
r = float(elem.find('r').text)
|
||||
return cls(order, x, y, r, filter_id=filter_id)
|
||||
|
||||
|
||||
class ZernikeRadialFilter(ZernikeFilter):
|
||||
r"""Score the :math:`m = 0` (radial variation only) Zernike moments up to
|
||||
|
|
@ -499,7 +531,7 @@ class ZernikeRadialFilter(ZernikeFilter):
|
|||
x-coordinate of center of circle for normalization
|
||||
y : float
|
||||
y-coordinate of center of circle for normalization
|
||||
r : int or None
|
||||
r : float
|
||||
Radius of circle for normalization
|
||||
|
||||
Attributes
|
||||
|
|
@ -510,7 +542,7 @@ class ZernikeRadialFilter(ZernikeFilter):
|
|||
x-coordinate of center of circle for normalization
|
||||
y : float
|
||||
y-coordinate of center of circle for normalization
|
||||
r : int or None
|
||||
r : float
|
||||
Radius of circle for normalization
|
||||
id : int
|
||||
Unique identifier for the filter
|
||||
|
|
|
|||
|
|
@ -285,6 +285,18 @@ class RectilinearMesh(Mesh):
|
|||
return (lower_left, upper_right, dimension, width)
|
||||
|
||||
def set_grid(self, x_grid, y_grid, z_grid):
|
||||
"""Set grid values
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x_grid : iterable of float
|
||||
Mesh boundary points along the x-axis
|
||||
y_grid : iterable of float
|
||||
Mesh boundary points along the y-axis
|
||||
z_grid : iterable of float
|
||||
Mesh boundary points along the z-axis
|
||||
|
||||
"""
|
||||
nx = len(x_grid)
|
||||
x_grid = (c_double*nx)(*x_grid)
|
||||
ny = len(y_grid)
|
||||
|
|
@ -371,15 +383,27 @@ class CylindricalMesh(Mesh):
|
|||
|
||||
return (lower_left, upper_right, dimension, width)
|
||||
|
||||
def set_grid(self, x_grid, y_grid, z_grid):
|
||||
nx = len(x_grid)
|
||||
x_grid = (c_double*nx)(*x_grid)
|
||||
ny = len(y_grid)
|
||||
y_grid = (c_double*ny)(*y_grid)
|
||||
def set_grid(self, r_grid, phi_grid, z_grid):
|
||||
"""Set grid values
|
||||
|
||||
Parameters
|
||||
----------
|
||||
r_grid : iterable of float
|
||||
Mesh boundary points along the r-axis
|
||||
phi_grid : Iterable of float
|
||||
Mesh boundary points along the phi-axis
|
||||
z_grid : Iterable of float
|
||||
Mesh boundary points along the z-axis
|
||||
|
||||
"""
|
||||
nr = len(r_grid)
|
||||
r_grid = (c_double*nr)(*r_grid)
|
||||
nphi = len(phi_grid)
|
||||
phi_grid = (c_double*nphi)(*phi_grid)
|
||||
nz = len(z_grid)
|
||||
z_grid = (c_double*nz)(*z_grid)
|
||||
_dll.openmc_cylindrical_mesh_set_grid(self._index, x_grid, nx, y_grid,
|
||||
ny, z_grid, nz)
|
||||
_dll.openmc_cylindrical_mesh_set_grid(self._index, r_grid, nr, phi_grid,
|
||||
nphi, z_grid, nz)
|
||||
|
||||
class SphericalMesh(Mesh):
|
||||
"""SphericalMesh stored internally.
|
||||
|
|
@ -457,15 +481,27 @@ class SphericalMesh(Mesh):
|
|||
|
||||
return (lower_left, upper_right, dimension, width)
|
||||
|
||||
def set_grid(self, x_grid, y_grid, z_grid):
|
||||
nx = len(x_grid)
|
||||
x_grid = (c_double*nx)(*x_grid)
|
||||
ny = len(y_grid)
|
||||
y_grid = (c_double*ny)(*y_grid)
|
||||
nz = len(z_grid)
|
||||
z_grid = (c_double*nz)(*z_grid)
|
||||
_dll.openmc_spherical_mesh_set_grid(self._index, x_grid, nx, y_grid,
|
||||
ny, z_grid, nz)
|
||||
def set_grid(self, r_grid, theta_grid, phi_grid):
|
||||
"""Set grid values
|
||||
|
||||
Parameters
|
||||
----------
|
||||
r_grid : iterable of float
|
||||
Mesh boundary points along the r-axis
|
||||
theta_grid : Iterable of float
|
||||
Mesh boundary points along the theta-axis
|
||||
phi_grid : Iterable of float
|
||||
Mesh boundary points along the phi-axis
|
||||
|
||||
"""
|
||||
nr = len(r_grid)
|
||||
r_grid = (c_double*nr)(*r_grid)
|
||||
ntheta = len(theta_grid)
|
||||
theta_grid = (c_double*ntheta)(*theta_grid)
|
||||
nphi = len(phi_grid)
|
||||
phi_grid = (c_double*nphi)(*phi_grid)
|
||||
_dll.openmc_spherical_mesh_set_grid(self._index, r_grid, nr, theta_grid,
|
||||
ntheta, phi_grid, nphi)
|
||||
|
||||
|
||||
_MESH_TYPE_MAP = {
|
||||
|
|
|
|||
103
openmc/mesh.py
103
openmc/mesh.py
|
|
@ -1,5 +1,6 @@
|
|||
from abc import ABC
|
||||
from collections.abc import Iterable
|
||||
from math import pi
|
||||
from numbers import Real, Integral
|
||||
import warnings
|
||||
from xml.etree import ElementTree as ET
|
||||
|
|
@ -89,7 +90,7 @@ class MeshBase(IDManagerMixin, ABC):
|
|||
raise ValueError('Unrecognized mesh type: "' + mesh_type + '"')
|
||||
|
||||
@classmethod
|
||||
def from_xml(cls, elem):
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generates a mesh from an XML element
|
||||
|
||||
Parameters
|
||||
|
|
@ -109,6 +110,10 @@ class MeshBase(IDManagerMixin, ABC):
|
|||
return RegularMesh.from_xml_element(elem)
|
||||
elif mesh_type == 'rectilinear':
|
||||
return RectilinearMesh.from_xml_element(elem)
|
||||
elif mesh_type == 'cylindrical':
|
||||
return CylindricalMesh.from_xml_element(elem)
|
||||
elif mesh_type == 'spherical':
|
||||
return SphericalMesh.from_xml_element(elem)
|
||||
elif mesh_type == 'unstructured':
|
||||
return UnstructuredMesh.from_xml_element(elem)
|
||||
else:
|
||||
|
|
@ -697,11 +702,11 @@ class CylindricalMesh(MeshBase):
|
|||
n_dimension : int
|
||||
Number of mesh dimensions (always 3 for a CylindricalMesh).
|
||||
r_grid : Iterable of float
|
||||
Mesh boundary points along the r-axis.
|
||||
Mesh boundary points along the r-axis.
|
||||
Requirement is r >= 0.
|
||||
phi_grid : Iterable of float
|
||||
Mesh boundary points along the phi-axis.
|
||||
The default value is [0, 360], i.e. the full phi range.
|
||||
Mesh boundary points along the phi-axis.
|
||||
The default value is [0, 2π], i.e. the full phi range.
|
||||
z_grid : Iterable of float
|
||||
Mesh boundary points along the z-axis.
|
||||
indices : Iterable of tuple
|
||||
|
|
@ -714,7 +719,7 @@ class CylindricalMesh(MeshBase):
|
|||
super().__init__(mesh_id, name)
|
||||
|
||||
self._r_grid = None
|
||||
self._phi_grid = [0, 360]
|
||||
self._phi_grid = [0.0, 2*pi]
|
||||
self._z_grid = None
|
||||
|
||||
@property
|
||||
|
|
@ -757,12 +762,12 @@ class CylindricalMesh(MeshBase):
|
|||
@phi_grid.setter
|
||||
def phi_grid(self, grid):
|
||||
cv.check_type('mesh phi_grid', grid, Iterable, Real)
|
||||
self._phi_grid = np.array(grid)
|
||||
self._phi_grid = np.asarray(grid)
|
||||
|
||||
@z_grid.setter
|
||||
def z_grid(self, grid):
|
||||
cv.check_type('mesh z_grid', grid, Iterable, Real)
|
||||
self._z_grid = np.array(grid)
|
||||
self._z_grid = np.asarray(grid)
|
||||
|
||||
def __repr__(self):
|
||||
fmt = '{0: <16}{1}{2}\n'
|
||||
|
|
@ -792,7 +797,7 @@ class CylindricalMesh(MeshBase):
|
|||
# Read and assign mesh properties
|
||||
mesh = cls(mesh_id)
|
||||
mesh.r_grid = group['r_grid'][()]
|
||||
mesh.phi_grid = 180 / np.pi * group['phi_grid'][()]
|
||||
mesh.phi_grid = group['phi_grid'][()]
|
||||
mesh.z_grid = group['z_grid'][()]
|
||||
|
||||
return mesh
|
||||
|
|
@ -822,6 +827,29 @@ class CylindricalMesh(MeshBase):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate a cylindrical mesh from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.CylindricalMesh
|
||||
Cylindrical mesh object
|
||||
|
||||
"""
|
||||
|
||||
mesh_id = int(get_text(elem, 'id'))
|
||||
mesh = cls(mesh_id)
|
||||
mesh.r_grid = [float(x) for x in get_text(elem, "r_grid").split()]
|
||||
mesh.phi_grid = [float(x) for x in get_text(elem, "phi_grid").split()]
|
||||
mesh.z_grid = [float(x) for x in get_text(elem, "z_grid").split()]
|
||||
return mesh
|
||||
|
||||
def calc_mesh_volumes(self):
|
||||
"""Return Volumes for every mesh cell
|
||||
|
||||
|
|
@ -832,9 +860,9 @@ class CylindricalMesh(MeshBase):
|
|||
|
||||
"""
|
||||
|
||||
V_r = np.diff(np.array(self.r_grid)**2 / 2)
|
||||
V_p = np.diff(np.array(self.phi_grid) * np.pi / 180.0)
|
||||
V_z = np.diff(np.array(self.z_grid))
|
||||
V_r = np.diff(np.asarray(self.r_grid)**2 / 2)
|
||||
V_p = np.diff(self.phi_grid)
|
||||
V_z = np.diff(self.z_grid)
|
||||
|
||||
return np.multiply.outer(np.outer(V_r, V_p), V_z)
|
||||
|
||||
|
|
@ -864,10 +892,10 @@ class SphericalMesh(MeshBase):
|
|||
Requirement is r >= 0.
|
||||
theta_grid : Iterable of float
|
||||
Mesh boundary points along the theta-axis in degrees.
|
||||
The default value is [0, 180], i.e. the full theta range.
|
||||
The default value is [0, π], i.e. the full theta range.
|
||||
phi_grid : Iterable of float
|
||||
Mesh boundary points along the phi-axis in degrees.
|
||||
The default value is [0, 360], i.e. the full phi range.
|
||||
The default value is [0, 2π], i.e. the full phi range.
|
||||
indices : Iterable of tuple
|
||||
An iterable of mesh indices for each mesh element, e.g. [(1, 1, 1),
|
||||
(2, 1, 1), ...]
|
||||
|
|
@ -878,8 +906,8 @@ class SphericalMesh(MeshBase):
|
|||
super().__init__(mesh_id, name)
|
||||
|
||||
self._r_grid = None
|
||||
self._theta_grid = [0, 180]
|
||||
self._phi_grid = [0, 360]
|
||||
self._theta_grid = [0, pi]
|
||||
self._phi_grid = [0, 2*pi]
|
||||
|
||||
@property
|
||||
def dimension(self):
|
||||
|
|
@ -921,12 +949,12 @@ class SphericalMesh(MeshBase):
|
|||
@theta_grid.setter
|
||||
def theta_grid(self, grid):
|
||||
cv.check_type('mesh theta_grid', grid, Iterable, Real)
|
||||
self._theta_grid = np.array(grid)
|
||||
self._theta_grid = np.asarray(grid)
|
||||
|
||||
@phi_grid.setter
|
||||
def phi_grid(self, grid):
|
||||
cv.check_type('mesh phi_grid', grid, Iterable, Real)
|
||||
self._phi_grid = np.array(grid)
|
||||
self._phi_grid = np.asarray(grid)
|
||||
|
||||
def __repr__(self):
|
||||
fmt = '{0: <16}{1}{2}\n'
|
||||
|
|
@ -956,8 +984,8 @@ class SphericalMesh(MeshBase):
|
|||
# Read and assign mesh properties
|
||||
mesh = cls(mesh_id)
|
||||
mesh.r_grid = group['r_grid'][()]
|
||||
mesh.theta_grid = 180 / np.pi * group['theta_grid'][()]
|
||||
mesh.phi_grid = 180 / np.pi * group['phi_grid'][()]
|
||||
mesh.theta_grid = group['theta_grid'][()]
|
||||
mesh.phi_grid = group['phi_grid'][()]
|
||||
|
||||
return mesh
|
||||
|
||||
|
|
@ -986,6 +1014,29 @@ class SphericalMesh(MeshBase):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate a spherical mesh from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.SphericalMesh
|
||||
Spherical mesh object
|
||||
|
||||
"""
|
||||
|
||||
mesh_id = int(get_text(elem, 'id'))
|
||||
mesh = cls(mesh_id)
|
||||
mesh.r_grid = [float(x) for x in get_text(elem, "r_grid").split()]
|
||||
mesh.theta_grid = [float(x) for x in get_text(elem, "theta_grid").split()]
|
||||
mesh.phi_grid = [float(x) for x in get_text(elem, "phi_grid").split()]
|
||||
return mesh
|
||||
|
||||
def calc_mesh_volumes(self):
|
||||
"""Return Volumes for every mesh cell
|
||||
|
||||
|
|
@ -996,15 +1047,13 @@ class SphericalMesh(MeshBase):
|
|||
|
||||
"""
|
||||
|
||||
V_r = np.diff(np.array(self.r_grid)**3 / 3)
|
||||
V_t = np.diff(-np.cos(np.pi * np.array(self.theta_grid) / 180.0))
|
||||
V_p = np.diff(np.array(self.phi_grid) * np.pi / 180)
|
||||
V_r = np.diff(np.asarray(self.r_grid)**3 / 3)
|
||||
V_t = np.diff(-np.cos(self.theta_grid))
|
||||
V_p = np.diff(self.phi_grid)
|
||||
|
||||
return np.multiply.outer(np.outer(V_r, V_t), V_p)
|
||||
|
||||
|
||||
|
||||
|
||||
class UnstructuredMesh(MeshBase):
|
||||
"""A 3D unstructured mesh
|
||||
|
||||
|
|
@ -1017,6 +1066,8 @@ class UnstructuredMesh(MeshBase):
|
|||
----------
|
||||
filename : str
|
||||
Location of the unstructured mesh file
|
||||
library : {'moab', 'libmesh'}
|
||||
Mesh library used for the unstructured mesh tally
|
||||
mesh_id : int
|
||||
Unique identifier for the mesh
|
||||
name : str
|
||||
|
|
@ -1034,7 +1085,7 @@ class UnstructuredMesh(MeshBase):
|
|||
Name of the file containing the unstructured mesh
|
||||
length_multiplier: float
|
||||
Multiplicative factor to apply to mesh coordinates
|
||||
library : str
|
||||
library : {'moab', 'libmesh'}
|
||||
Mesh library used for the unstructured mesh tally
|
||||
output : bool
|
||||
Indicates whether or not automatic tally output should
|
||||
|
|
@ -1291,4 +1342,4 @@ class UnstructuredMesh(MeshBase):
|
|||
library = get_text(elem, 'library')
|
||||
length_multiplier = float(get_text(elem, 'length_multiplier', 1.0))
|
||||
|
||||
return cls(filename, library, mesh_id, '', length_multiplier)
|
||||
return cls(filename, library, mesh_id, '', length_multiplier)
|
||||
|
|
|
|||
|
|
@ -5,7 +5,6 @@ import os
|
|||
from pathlib import Path
|
||||
from numbers import Integral
|
||||
from tempfile import NamedTemporaryFile
|
||||
import time
|
||||
|
||||
import h5py
|
||||
|
||||
|
|
@ -204,10 +203,9 @@ class Model:
|
|||
|
||||
@classmethod
|
||||
def from_xml(cls, geometry='geometry.xml', materials='materials.xml',
|
||||
settings='settings.xml'):
|
||||
settings='settings.xml', tallies='tallies.xml',
|
||||
plots='plots.xml'):
|
||||
"""Create model from existing XML files
|
||||
When initializing this way, the user must manually load plots and
|
||||
tallies.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
|
|
@ -217,6 +215,14 @@ class Model:
|
|||
Path to materials.xml file
|
||||
settings : str
|
||||
Path to settings.xml file
|
||||
tallies : str
|
||||
Path to tallies.xml file
|
||||
|
||||
.. versionadded:: 0.13.0
|
||||
plots : str
|
||||
Path to plots.xml file
|
||||
|
||||
.. versionadded:: 0.13.0
|
||||
|
||||
Returns
|
||||
-------
|
||||
|
|
@ -227,7 +233,9 @@ class Model:
|
|||
materials = openmc.Materials.from_xml(materials)
|
||||
geometry = openmc.Geometry.from_xml(geometry, materials)
|
||||
settings = openmc.Settings.from_xml(settings)
|
||||
return cls(geometry, materials, settings)
|
||||
tallies = openmc.Tallies.from_xml(tallies) if Path(tallies).exists() else None
|
||||
plots = openmc.Plots.from_xml(plots) if Path(plots).exists() else None
|
||||
return cls(geometry, materials, settings, tallies, plots)
|
||||
|
||||
def init_lib(self, threads=None, geometry_debug=False, restart_file=None,
|
||||
tracks=False, output=True, event_based=None, intracomm=None):
|
||||
|
|
|
|||
103
openmc/plots.py
103
openmc/plots.py
|
|
@ -684,6 +684,85 @@ class Plot(IDManagerMixin):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate plot object from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Plot
|
||||
Plot object
|
||||
|
||||
"""
|
||||
plot_id = int(elem.get("id"))
|
||||
plot = cls(plot_id)
|
||||
if "filename" in elem.keys():
|
||||
plot.filename = elem.get("filename")
|
||||
plot.color_by = elem.get("color_by")
|
||||
plot.type = elem.get("type")
|
||||
plot.basis = elem.get("basis")
|
||||
|
||||
# Helper function to get a tuple of values
|
||||
def get_tuple(elem, name, dtype=int):
|
||||
subelem = elem.find(name)
|
||||
if subelem is not None:
|
||||
return tuple([dtype(x) for x in subelem.text.split()])
|
||||
|
||||
plot.origin = get_tuple(elem, "origin", float)
|
||||
plot.width = get_tuple(elem, "width", float)
|
||||
plot.pixels = get_tuple(elem, "pixels")
|
||||
plot._background = get_tuple(elem, "background")
|
||||
|
||||
# Set plot colors
|
||||
colors = {}
|
||||
for color_elem in elem.findall("color"):
|
||||
uid = color_elem.get("id")
|
||||
colors[uid] = tuple([int(x) for x in color_elem.get("rgb").split()])
|
||||
# TODO: set colors (needs geometry information)
|
||||
|
||||
# Set masking information
|
||||
mask_elem = elem.find("mask")
|
||||
if mask_elem is not None:
|
||||
mask_components = [int(x) for x in mask_elem.get("components").split()]
|
||||
# TODO: set mask components (needs geometry information)
|
||||
background = mask_elem.get("background")
|
||||
if background is not None:
|
||||
plot.mask_background = tuple([int(x) for x in background.split()])
|
||||
|
||||
# show overlaps
|
||||
overlap_elem = elem.find("show_overlaps")
|
||||
if overlap_elem is not None:
|
||||
plot.show_overlaps = (overlap_elem.text in ('true', '1'))
|
||||
overlap_color = get_tuple(elem, "overlap_color")
|
||||
if overlap_color is not None:
|
||||
plot.overlap_color = overlap_color
|
||||
|
||||
# Set universe level
|
||||
level = elem.find("level")
|
||||
if level is not None:
|
||||
plot.level = int(level.text)
|
||||
|
||||
# Set meshlines
|
||||
mesh_elem = elem.find("meshlines")
|
||||
if mesh_elem is not None:
|
||||
meshlines = {'type': mesh_elem.get('meshtype')}
|
||||
if 'id' in mesh_elem.keys():
|
||||
meshlines['id'] = int(mesh_elem.get('id'))
|
||||
if 'linewidth' in mesh_elem.keys():
|
||||
meshlines['linewidth'] = int(mesh_elem.get('linewidth'))
|
||||
if 'color' in mesh_elem.keys():
|
||||
meshlines['color'] = tuple(
|
||||
[int(x) for x in mesh_elem.get('color').split()]
|
||||
)
|
||||
plot.meshlines = meshlines
|
||||
|
||||
return plot
|
||||
|
||||
def to_ipython_image(self, openmc_exec='openmc', cwd='.'):
|
||||
"""Render plot as an image
|
||||
|
||||
|
|
@ -848,3 +927,27 @@ class Plots(cv.CheckedList):
|
|||
reorder_attributes(self._plots_file) # TODO: Remove when support is Python 3.8+
|
||||
tree = ET.ElementTree(self._plots_file)
|
||||
tree.write(str(p), xml_declaration=True, encoding='utf-8')
|
||||
|
||||
@classmethod
|
||||
def from_xml(cls, path='plots.xml'):
|
||||
"""Generate plots collection from XML file
|
||||
|
||||
Parameters
|
||||
----------
|
||||
path : str, optional
|
||||
Path to plots XML file
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Plots
|
||||
Plots collection
|
||||
|
||||
"""
|
||||
tree = ET.parse(path)
|
||||
root = tree.getroot()
|
||||
|
||||
# Generate each plot
|
||||
plots = cls()
|
||||
for elem in root.findall('plot'):
|
||||
plots.append(Plot.from_xml_element(elem))
|
||||
return plots
|
||||
|
|
|
|||
|
|
@ -1247,6 +1247,12 @@ class Settings:
|
|||
for elem in root.findall('source'):
|
||||
self.source.append(Source.from_xml_element(elem))
|
||||
|
||||
def _volume_calcs_from_xml_element(self, root):
|
||||
volume_elems = root.findall("volume_calc")
|
||||
if volume_elems:
|
||||
self.volume_calculations = [VolumeCalculation.from_xml_element(elem)
|
||||
for elem in volume_elems]
|
||||
|
||||
def _output_from_xml_element(self, root):
|
||||
elem = root.find('output')
|
||||
if elem is not None:
|
||||
|
|
@ -1585,6 +1591,7 @@ class Settings:
|
|||
settings._generations_per_batch_from_xml_element(root)
|
||||
settings._keff_trigger_from_xml_element(root)
|
||||
settings._source_from_xml_element(root)
|
||||
settings._volume_calcs_from_xml_element(root)
|
||||
settings._output_from_xml_element(root)
|
||||
settings._statepoint_from_xml_element(root)
|
||||
settings._sourcepoint_from_xml_element(root)
|
||||
|
|
|
|||
|
|
@ -492,6 +492,7 @@ class SphericalIndependent(Spatial):
|
|||
origin = [float(x) for x in elem.get('origin').split()]
|
||||
return cls(r, theta, phi, origin=origin)
|
||||
|
||||
|
||||
class CylindricalIndependent(Spatial):
|
||||
r"""Spatial distribution represented in cylindrical coordinates.
|
||||
|
||||
|
|
|
|||
|
|
@ -10,12 +10,14 @@ from xml.etree import ElementTree as ET
|
|||
import h5py
|
||||
import numpy as np
|
||||
import pandas as pd
|
||||
from scipy.misc import derivative
|
||||
import scipy.sparse as sps
|
||||
|
||||
import openmc
|
||||
import openmc.checkvalue as cv
|
||||
from ._xml import clean_indentation, reorder_attributes
|
||||
from ._xml import clean_indentation, reorder_attributes, get_text
|
||||
from .mixin import IDManagerMixin
|
||||
from .mesh import MeshBase
|
||||
|
||||
|
||||
# The tally arithmetic product types. The tensor product performs the full
|
||||
|
|
@ -844,13 +846,8 @@ class Tally(IDManagerMixin):
|
|||
'not contain any scores'
|
||||
raise ValueError(msg)
|
||||
|
||||
else:
|
||||
scores = ''
|
||||
for score in self.scores:
|
||||
scores += f'{score} '
|
||||
|
||||
subelement = ET.SubElement(element, "scores")
|
||||
subelement.text = scores.rstrip(' ')
|
||||
subelement = ET.SubElement(element, "scores")
|
||||
subelement.text = ' '.join(str(x) for x in self.scores)
|
||||
|
||||
# Tally estimator type
|
||||
if self.estimator is not None:
|
||||
|
|
@ -859,7 +856,7 @@ class Tally(IDManagerMixin):
|
|||
|
||||
# Optional Triggers
|
||||
for trigger in self.triggers:
|
||||
trigger.get_trigger_xml(element)
|
||||
element.append(trigger.to_xml_element())
|
||||
|
||||
# Optional derivatives
|
||||
if self.derivative is not None:
|
||||
|
|
@ -868,6 +865,60 @@ class Tally(IDManagerMixin):
|
|||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem, **kwargs):
|
||||
"""Generate tally object from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Tally
|
||||
Tally object
|
||||
|
||||
"""
|
||||
tally_id = int(elem.get('id'))
|
||||
name = elem.get('name', '')
|
||||
tally = cls(tally_id=tally_id, name=name)
|
||||
|
||||
# Read filters
|
||||
filters_elem = elem.find('filters')
|
||||
if filters_elem is not None:
|
||||
filter_ids = [int(x) for x in filters_elem.text.split()]
|
||||
tally.filters = [kwargs['filters'][uid] for uid in filter_ids]
|
||||
|
||||
# Read nuclides
|
||||
nuclides_elem = elem.find('nuclides')
|
||||
if nuclides_elem is not None:
|
||||
tally.nuclides = nuclides_elem.text.split()
|
||||
|
||||
# Read scores
|
||||
scores_elem = elem.find('scores')
|
||||
if scores_elem is not None:
|
||||
tally.scores = scores_elem.text.split()
|
||||
|
||||
# Set estimator
|
||||
estimator_elem = elem.find('estimator')
|
||||
if estimator_elem is not None:
|
||||
tally.estimator = estimator_elem.text
|
||||
|
||||
# Read triggers
|
||||
tally.triggers = [
|
||||
openmc.Trigger.from_xml_element(trigger_elem)
|
||||
for trigger_elem in elem.findall('trigger')
|
||||
]
|
||||
|
||||
# Read tally derivative
|
||||
deriv_elem = elem.find('derivative')
|
||||
if deriv_elem is not None:
|
||||
deriv_id = int(deriv_elem.text)
|
||||
tally.derivative = kwargs['derivatives'][deriv_id]
|
||||
|
||||
return tally
|
||||
|
||||
def contains_filter(self, filter_type):
|
||||
"""Looks for a filter in the tally that matches a specified type
|
||||
|
||||
|
|
@ -3143,3 +3194,49 @@ class Tallies(cv.CheckedList):
|
|||
reorder_attributes(root_element) # TODO: Remove when support is Python 3.8+
|
||||
tree = ET.ElementTree(root_element)
|
||||
tree.write(str(p), xml_declaration=True, encoding='utf-8')
|
||||
|
||||
@classmethod
|
||||
def from_xml(cls, path='tallies.xml'):
|
||||
"""Generate tallies from XML file
|
||||
|
||||
Parameters
|
||||
----------
|
||||
path : str, optional
|
||||
Path to tallies XML file
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Tallies
|
||||
Tallies object
|
||||
|
||||
"""
|
||||
tree = ET.parse(path)
|
||||
root = tree.getroot()
|
||||
|
||||
# Read mesh elements
|
||||
meshes = {}
|
||||
for elem in root.findall('mesh'):
|
||||
mesh = MeshBase.from_xml_element(elem)
|
||||
meshes[mesh.id] = mesh
|
||||
|
||||
# Read filter elements
|
||||
filters = {}
|
||||
for elem in root.findall('filter'):
|
||||
filter = openmc.Filter.from_xml_element(elem, meshes=meshes)
|
||||
filters[filter.id] = filter
|
||||
|
||||
# Read derivative elements
|
||||
derivatives = {}
|
||||
for elem in root.findall('derivative'):
|
||||
deriv = openmc.TallyDerivative.from_xml_element(elem)
|
||||
derivatives[deriv.id] = deriv
|
||||
|
||||
# Read tally elements
|
||||
tallies = []
|
||||
for elem in root.findall('tally'):
|
||||
tally = openmc.Tally.from_xml_element(
|
||||
elem, filters=filters, derivatives=derivatives
|
||||
)
|
||||
tallies.append(tally)
|
||||
|
||||
return cls(tallies)
|
||||
|
|
|
|||
|
|
@ -105,3 +105,24 @@ class TallyDerivative(EqualityMixin, IDManagerMixin):
|
|||
if self.variable == 'nuclide_density':
|
||||
element.set("nuclide", self.nuclide)
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate tally derivative from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.TallyDerivative
|
||||
Tally derivative object
|
||||
|
||||
"""
|
||||
derivative_id = int(elem.get("id"))
|
||||
variable = elem.get("variable")
|
||||
material = int(elem.get("material"))
|
||||
nuclide = elem.get("nuclide") if variable == "nuclide_density" else None
|
||||
return cls(derivative_id, variable, material, nuclide)
|
||||
|
|
|
|||
|
|
@ -74,7 +74,7 @@ class Trigger(EqualityMixin):
|
|||
if score not in self._scores:
|
||||
self._scores.append(score)
|
||||
|
||||
def get_trigger_xml(self, element):
|
||||
def to_xml_element(self):
|
||||
"""Return XML representation of the trigger
|
||||
|
||||
Returns
|
||||
|
|
@ -84,8 +84,36 @@ class Trigger(EqualityMixin):
|
|||
|
||||
"""
|
||||
|
||||
subelement = ET.SubElement(element, "trigger")
|
||||
subelement.set("type", self._trigger_type)
|
||||
subelement.set("threshold", str(self._threshold))
|
||||
element = ET.Element("trigger")
|
||||
element.set("type", self._trigger_type)
|
||||
element.set("threshold", str(self._threshold))
|
||||
if len(self._scores) != 0:
|
||||
subelement.set("scores", ' '.join(map(str, self._scores)))
|
||||
element.set("scores", ' '.join(self._scores))
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate trigger object from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.Trigger
|
||||
Trigger object
|
||||
|
||||
"""
|
||||
# Generate trigger object
|
||||
trigger_type = elem.get("type")
|
||||
threshold = float(elem.get("threshold"))
|
||||
trigger = cls(trigger_type, threshold)
|
||||
|
||||
# Add scores if present
|
||||
scores = elem.get("scores")
|
||||
if scores is not None:
|
||||
trigger.scores = scores.split()
|
||||
|
||||
return trigger
|
||||
|
|
|
|||
|
|
@ -11,6 +11,7 @@ from uncertainties import ufloat
|
|||
|
||||
import openmc
|
||||
import openmc.checkvalue as cv
|
||||
from openmc._xml import get_text
|
||||
|
||||
_VERSION_VOLUME = 1
|
||||
|
||||
|
|
@ -352,3 +353,49 @@ class VolumeCalculation:
|
|||
trigger_elem.set("type", self.trigger_type)
|
||||
trigger_elem.set("threshold", str(self.threshold))
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate volume calculation object from an XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.VolumeCalculation
|
||||
Volume calculation object
|
||||
|
||||
"""
|
||||
domain_type = get_text(elem, "domain_type")
|
||||
domain_ids = get_text(elem, "domain_ids").split()
|
||||
ids = [int(x) for x in domain_ids]
|
||||
samples = int(get_text(elem, "samples"))
|
||||
lower_left = get_text(elem, "lower_left").split()
|
||||
lower_left = tuple([float(x) for x in lower_left])
|
||||
upper_right = get_text(elem, "upper_right").split()
|
||||
upper_right = tuple([float(x) for x in upper_right])
|
||||
|
||||
# Instantiate some throw-away domains that are used by the constructor
|
||||
# to assign IDs
|
||||
with warnings.catch_warnings():
|
||||
warnings.simplefilter('ignore', openmc.IDWarning)
|
||||
if domain_type == 'cell':
|
||||
domains = [openmc.Cell(uid) for uid in ids]
|
||||
elif domain_type == 'material':
|
||||
domains = [openmc.Material(uid) for uid in ids]
|
||||
elif domain_type == 'universe':
|
||||
domains = [openmc.Universe(uid) for uid in ids]
|
||||
|
||||
vol = cls(domains, samples, lower_left, upper_right)
|
||||
|
||||
# Check for trigger
|
||||
trigger_elem = elem.find("threshold")
|
||||
if trigger_elem is not None:
|
||||
trigger_type = get_text(trigger_elem, "type")
|
||||
threshold = float(get_text(trigger_elem, "threshold"))
|
||||
vol.set_trigger(threshold, trigger_type)
|
||||
|
||||
return vol
|
||||
|
|
|
|||
|
|
@ -270,7 +270,7 @@ class WeightWindows(IDManagerMixin):
|
|||
path = f"./mesh[@id='{mesh_id}']"
|
||||
mesh_elem = root.find(path)
|
||||
if mesh_elem is not None:
|
||||
mesh = MeshBase.from_xml(mesh_elem)
|
||||
mesh = MeshBase.from_xml_element(mesh_elem)
|
||||
|
||||
# Read all other parameters
|
||||
lower_ww_bounds = [float(l) for l in get_text(elem, 'lower_ww_bounds').split()]
|
||||
|
|
|
|||
33
src/mesh.cpp
33
src/mesh.cpp
|
|
@ -1001,7 +1001,7 @@ double CylindricalMesh::find_r_crossing(
|
|||
double CylindricalMesh::find_phi_crossing(
|
||||
const Position& r, const Direction& u, double l, int shell) const
|
||||
{
|
||||
// Phi grid is [0, 360], thus there is no real surface to cross
|
||||
// Phi grid is [0, 2π], thus there is no real surface to cross
|
||||
if (full_phi_ && (shape_[1] == 1))
|
||||
return INFTY;
|
||||
|
||||
|
|
@ -1104,18 +1104,14 @@ int CylindricalMesh::set_grid()
|
|||
"cylindrical meshes must start at phi >= 0.");
|
||||
return OPENMC_E_INVALID_ARGUMENT;
|
||||
}
|
||||
if (grid_[1].back() > 360) {
|
||||
if (grid_[1].back() > 2.0 * PI) {
|
||||
set_errmsg("phi-grids for "
|
||||
"cylindrical meshes must end with theta <= 360 degree.");
|
||||
"cylindrical meshes must end with theta <= 2*pi.");
|
||||
|
||||
return OPENMC_E_INVALID_ARGUMENT;
|
||||
}
|
||||
|
||||
full_phi_ = (grid_[1].front() == 0.0) && (grid_[1].back() == 360.0);
|
||||
|
||||
// Transform phi-grid from degrees to radians
|
||||
std::transform(grid_[1].begin(), grid_[1].end(), grid_[1].begin(),
|
||||
[](double d) { return M_PI * d / 180.0; });
|
||||
full_phi_ = (grid_[1].front() == 0.0) && (grid_[1].back() == 2.0 * PI);
|
||||
|
||||
lower_left_ = {grid_[0].front(), grid_[1].front(), grid_[2].front()};
|
||||
upper_right_ = {grid_[0].back(), grid_[1].back(), grid_[2].back()};
|
||||
|
|
@ -1225,7 +1221,7 @@ double SphericalMesh::find_r_crossing(
|
|||
double SphericalMesh::find_theta_crossing(
|
||||
const Position& r, const Direction& u, double l, int shell) const
|
||||
{
|
||||
// Theta grid is [0, 180], thus there is no real surface to cross
|
||||
// Theta grid is [0, π], thus there is no real surface to cross
|
||||
if (full_theta_ && (shape_[1] == 1))
|
||||
return INFTY;
|
||||
|
||||
|
|
@ -1287,7 +1283,7 @@ double SphericalMesh::find_theta_crossing(
|
|||
double SphericalMesh::find_phi_crossing(
|
||||
const Position& r, const Direction& u, double l, int shell) const
|
||||
{
|
||||
// Phi grid is [0, 360], thus there is no real surface to cross
|
||||
// Phi grid is [0, 2π], thus there is no real surface to cross
|
||||
if (full_phi_ && (shape_[2] == 1))
|
||||
return INFTY;
|
||||
|
||||
|
|
@ -1368,25 +1364,20 @@ int SphericalMesh::set_grid()
|
|||
return OPENMC_E_INVALID_ARGUMENT;
|
||||
}
|
||||
}
|
||||
if (grid_[1].back() > 180) {
|
||||
if (grid_[1].back() > PI) {
|
||||
set_errmsg("theta-grids for "
|
||||
"spherical meshes must end with theta <= 180 degree.");
|
||||
"spherical meshes must end with theta <= pi.");
|
||||
|
||||
return OPENMC_E_INVALID_ARGUMENT;
|
||||
}
|
||||
if (grid_[2].back() > 360) {
|
||||
if (grid_[2].back() > 2 * PI) {
|
||||
set_errmsg("phi-grids for "
|
||||
"spherical meshes must end with phi <= 180 degree.");
|
||||
"spherical meshes must end with phi <= 2*pi.");
|
||||
return OPENMC_E_INVALID_ARGUMENT;
|
||||
}
|
||||
|
||||
full_theta_ = (grid_[1].front() == 0.0) && (grid_[1].back() == 180.0);
|
||||
full_phi_ = (grid_[2].front() == 0.0) && (grid_[2].back() == 360.0);
|
||||
|
||||
// Transform theta- and phi-grid from degrees to radians
|
||||
for (int i = 1; i < 3; i++)
|
||||
std::transform(grid_[i].begin(), grid_[i].end(), grid_[i].begin(),
|
||||
[](double d) { return M_PI * d / 180.0; });
|
||||
full_theta_ = (grid_[1].front() == 0.0) && (grid_[1].back() == PI);
|
||||
full_phi_ = (grid_[2].front() == 0.0) && (grid_[2].back() == 2 * PI);
|
||||
|
||||
lower_left_ = {grid_[0].front(), grid_[1].front(), grid_[2].front()};
|
||||
upper_right_ = {grid_[0].back(), grid_[1].back(), grid_[2].back()};
|
||||
|
|
|
|||
|
|
@ -55,13 +55,13 @@
|
|||
</mesh>
|
||||
<mesh id="5" type="cylindrical">
|
||||
<r_grid>0.0 0.4411764705882353 0.8823529411764706 1.3235294117647058 1.7647058823529411 2.2058823529411766 2.6470588235294117 3.0882352941176467 3.5294117647058822 3.9705882352941178 4.411764705882353 4.852941176470588 5.294117647058823 5.735294117647059 6.1764705882352935 6.617647058823529 7.0588235294117645 7.5</r_grid>
|
||||
<phi_grid>0.0 20.0 40.0 60.0 80.0 100.0 120.0 140.0 160.0 180.0 200.0 220.0 240.0 260.0 280.0 300.0 320.0 340.0 360.0</phi_grid>
|
||||
<phi_grid>0.0 0.3490658503988659 0.6981317007977318 1.0471975511965976 1.3962634015954636 1.7453292519943295 2.0943951023931953 2.443460952792061 2.792526803190927 3.141592653589793 3.490658503988659 3.839724354387525 4.1887902047863905 4.537856055185257 4.886921905584122 5.235987755982989 5.585053606381854 5.93411945678072 6.283185307179586</phi_grid>
|
||||
<z_grid>-7.5 -6.5625 -5.625 -4.6875 -3.75 -2.8125 -1.875 -0.9375 0.0 0.9375 1.875 2.8125 3.75 4.6875 5.625 6.5625 7.5</z_grid>
|
||||
</mesh>
|
||||
<mesh id="6" type="spherical">
|
||||
<r_grid>0.0 0.4411764705882353 0.8823529411764706 1.3235294117647058 1.7647058823529411 2.2058823529411766 2.6470588235294117 3.0882352941176467 3.5294117647058822 3.9705882352941178 4.411764705882353 4.852941176470588 5.294117647058823 5.735294117647059 6.1764705882352935 6.617647058823529 7.0588235294117645 7.5</r_grid>
|
||||
<theta_grid>0.0 22.5 45.0 67.5 90.0 112.5 135.0 157.5 180.0</theta_grid>
|
||||
<phi_grid>0.0 20.0 40.0 60.0 80.0 100.0 120.0 140.0 160.0 180.0 200.0 220.0 240.0 260.0 280.0 300.0 320.0 340.0 360.0</phi_grid>
|
||||
<theta_grid>0.0 0.39269908169872414 0.7853981633974483 1.1780972450961724 1.5707963267948966 1.9634954084936207 2.356194490192345 2.748893571891069 3.141592653589793</theta_grid>
|
||||
<phi_grid>0.0 0.3490658503988659 0.6981317007977318 1.0471975511965976 1.3962634015954636 1.7453292519943295 2.0943951023931953 2.443460952792061 2.792526803190927 3.141592653589793 3.490658503988659 3.839724354387525 4.1887902047863905 4.537856055185257 4.886921905584122 5.235987755982989 5.585053606381854 5.93411945678072 6.283185307179586</phi_grid>
|
||||
</mesh>
|
||||
<filter id="1" type="mesh">
|
||||
<bins>1</bins>
|
||||
|
|
|
|||
|
|
@ -1,4 +1,5 @@
|
|||
import numpy as np
|
||||
from math import pi
|
||||
|
||||
import openmc
|
||||
import pytest
|
||||
|
|
@ -53,13 +54,13 @@ def model():
|
|||
|
||||
cyl_mesh = openmc.CylindricalMesh()
|
||||
cyl_mesh.r_grid = np.linspace(0, 7.5, 18)
|
||||
cyl_mesh.phi_grid = np.linspace(0, 360, 19)
|
||||
cyl_mesh.phi_grid = np.linspace(0, 2*pi, 19)
|
||||
cyl_mesh.z_grid = np.linspace(-7.5, 7.5, 17)
|
||||
|
||||
sph_mesh = openmc.SphericalMesh()
|
||||
sph_mesh.r_grid = np.linspace(0, 7.5, 18)
|
||||
sph_mesh.theta_grid = np.linspace(0, 180, 9)
|
||||
sph_mesh.phi_grid = np.linspace(0, 360, 19)
|
||||
sph_mesh.theta_grid = np.linspace(0, pi, 9)
|
||||
sph_mesh.phi_grid = np.linspace(0, 2*pi, 19)
|
||||
|
||||
# Create filters
|
||||
reg_filters = [
|
||||
|
|
|
|||
|
|
@ -1,5 +1,4 @@
|
|||
from math import sqrt, pi
|
||||
|
||||
import numpy as np
|
||||
import openmc
|
||||
from pytest import fixture, approx
|
||||
|
||||
|
|
@ -38,6 +37,11 @@ def test_cell_instance():
|
|||
assert all(x == c1.id for x in bins[:6:2])
|
||||
assert all(x == c2.id for x in bins[6::2])
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
assert new_f.id == f.id
|
||||
assert np.all(new_f.bins == f.bins)
|
||||
|
||||
# get_pandas_dataframe()
|
||||
df = f.get_pandas_dataframe(f.num_bins, 1)
|
||||
cells = df['cellinstance', 'cell']
|
||||
|
|
@ -61,6 +65,11 @@ def test_collision():
|
|||
assert elem.tag == 'filter'
|
||||
assert elem.attrib['type'] == 'collision'
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
assert new_f.id == f.id
|
||||
assert np.all(new_f.bins == f.bins)
|
||||
|
||||
|
||||
def test_legendre():
|
||||
n = 5
|
||||
|
|
@ -79,6 +88,11 @@ def test_legendre():
|
|||
assert elem.attrib['type'] == 'legendre'
|
||||
assert elem.find('order').text == str(n)
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
assert new_f.id == f.id
|
||||
assert new_f.bins, f.bins
|
||||
|
||||
|
||||
def test_spatial_legendre():
|
||||
n = 5
|
||||
|
|
@ -102,6 +116,12 @@ def test_spatial_legendre():
|
|||
assert elem.find('order').text == str(n)
|
||||
assert elem.find('axis').text == str(axis)
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
assert new_f.id == f.id
|
||||
assert new_f.order == f.order
|
||||
assert new_f.axis == f.axis
|
||||
|
||||
|
||||
def test_spherical_harmonics():
|
||||
n = 3
|
||||
|
|
@ -122,6 +142,12 @@ def test_spherical_harmonics():
|
|||
assert elem.attrib['cosine'] == f.cosine
|
||||
assert elem.find('order').text == str(n)
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
assert new_f.id == f.id
|
||||
assert new_f.order == f.order
|
||||
assert new_f.cosine == f.cosine
|
||||
|
||||
|
||||
def test_zernike():
|
||||
n = 4
|
||||
|
|
@ -140,6 +166,12 @@ def test_zernike():
|
|||
assert elem.attrib['type'] == 'zernike'
|
||||
assert elem.find('order').text == str(n)
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
for attr in ('id', 'order', 'x', 'y', 'r'):
|
||||
assert getattr(new_f, attr) == getattr(f, attr)
|
||||
|
||||
|
||||
def test_zernike_radial():
|
||||
n = 4
|
||||
f = openmc.ZernikeRadialFilter(n, 0., 0., 1.)
|
||||
|
|
@ -157,6 +189,11 @@ def test_zernike_radial():
|
|||
assert elem.attrib['type'] == 'zernikeradial'
|
||||
assert elem.find('order').text == str(n)
|
||||
|
||||
# from_xml_element()
|
||||
new_f = openmc.Filter.from_xml_element(elem)
|
||||
for attr in ('id', 'order', 'x', 'y', 'r'):
|
||||
assert getattr(new_f, attr) == getattr(f, attr)
|
||||
|
||||
|
||||
def test_first_moment(run_in_tmpdir, box_model):
|
||||
plain_tally = openmc.Tally()
|
||||
|
|
|
|||
|
|
@ -583,17 +583,17 @@ def test_rectilinear_mesh(lib_init):
|
|||
def test_cylindrical_mesh(lib_init):
|
||||
deg2rad = lambda deg: deg*np.pi/180
|
||||
mesh = openmc.lib.CylindricalMesh()
|
||||
x_grid = [0., 5., 10.]
|
||||
y_grid = [0., 10., 20.]
|
||||
r_grid = [0., 5., 10.]
|
||||
phi_grid = np.radians([0., 10., 20.])
|
||||
z_grid = [10., 20., 30.]
|
||||
mesh.set_grid(x_grid, y_grid, z_grid)
|
||||
mesh.set_grid(r_grid, phi_grid, z_grid)
|
||||
assert np.all(mesh.lower_left == (0., 0., 10.))
|
||||
assert np.all(mesh.upper_right == (10., deg2rad(20.), 30.))
|
||||
assert np.all(mesh.dimension == (2, 2, 2))
|
||||
for i, diff_x in enumerate(np.diff(x_grid)):
|
||||
for j, diff_y in enumerate(np.diff(y_grid)):
|
||||
for k, diff_z in enumerate(np.diff(z_grid)):
|
||||
assert np.all(mesh.width[i, j, k, :] == (5, deg2rad(10), 10))
|
||||
for i, _ in enumerate(np.diff(r_grid)):
|
||||
for j, _ in enumerate(np.diff(phi_grid)):
|
||||
for k, _ in enumerate(np.diff(z_grid)):
|
||||
assert np.allclose(mesh.width[i, j, k, :], (5, deg2rad(10), 10))
|
||||
|
||||
with pytest.raises(exc.AllocationError):
|
||||
mesh2 = openmc.lib.CylindricalMesh(mesh.id)
|
||||
|
|
@ -614,17 +614,17 @@ def test_cylindrical_mesh(lib_init):
|
|||
def test_spherical_mesh(lib_init):
|
||||
deg2rad = lambda deg: deg*np.pi/180
|
||||
mesh = openmc.lib.SphericalMesh()
|
||||
x_grid = [0., 5., 10.]
|
||||
y_grid = [0., 10., 20.]
|
||||
z_grid = [10., 20., 30.]
|
||||
mesh.set_grid(x_grid, y_grid, z_grid)
|
||||
r_grid = [0., 5., 10.]
|
||||
theta_grid = np.radians([0., 10., 20.])
|
||||
phi_grid = np.radians([10., 20., 30.])
|
||||
mesh.set_grid(r_grid, theta_grid, phi_grid)
|
||||
assert np.all(mesh.lower_left == (0., 0., deg2rad(10.)))
|
||||
assert np.all(mesh.upper_right == (10., deg2rad(20.), deg2rad(30.)))
|
||||
assert np.all(mesh.dimension == (2, 2, 2))
|
||||
for i, diff_x in enumerate(np.diff(x_grid)):
|
||||
for j, diff_y in enumerate(np.diff(y_grid)):
|
||||
for k, diff_z in enumerate(np.diff(z_grid)):
|
||||
assert np.all(abs(mesh.width[i, j, k, :] - (5, deg2rad(10), deg2rad(10))) < 1e-16)
|
||||
for i, _ in enumerate(np.diff(r_grid)):
|
||||
for j, _ in enumerate(np.diff(theta_grid)):
|
||||
for k, _ in enumerate(np.diff(phi_grid)):
|
||||
assert np.allclose(mesh.width[i, j, k, :], (5, deg2rad(10), deg2rad(10)))
|
||||
|
||||
with pytest.raises(exc.AllocationError):
|
||||
mesh2 = openmc.lib.SphericalMesh(mesh.id)
|
||||
|
|
|
|||
|
|
@ -195,8 +195,8 @@ def test_from_xml(run_in_tmpdir, pin_model_attributes):
|
|||
keys = sorted(k for k in settings.__dict__.keys() if k not in no_test)
|
||||
for ref_k in keys:
|
||||
assert test_model.settings.__dict__[ref_k] == settings.__dict__[ref_k]
|
||||
assert len(test_model.tallies) == 0
|
||||
assert len(test_model.plots) == 0
|
||||
assert len(test_model.tallies) == 1
|
||||
assert len(test_model.plots) == 2
|
||||
assert test_model._materials_by_id == \
|
||||
{1: test_model.materials[0], 2: test_model.materials[1],
|
||||
3: test_model.materials[2]}
|
||||
|
|
|
|||
|
|
@ -1,3 +1,4 @@
|
|||
import numpy as np
|
||||
import openmc
|
||||
import openmc.examples
|
||||
import pytest
|
||||
|
|
@ -12,8 +13,8 @@ def myplot():
|
|||
plot.filename = 'myplot'
|
||||
plot.type = 'slice'
|
||||
plot.basis = 'yz'
|
||||
plot.background = (0, 0, 0)
|
||||
plot.background = 'black'
|
||||
plot.background = (0, 0, 0)
|
||||
|
||||
plot.color_by = 'material'
|
||||
m1, m2 = openmc.Material(), openmc.Material()
|
||||
|
|
@ -21,8 +22,8 @@ def myplot():
|
|||
plot.colors = {m1: 'green', m2: 'blue'}
|
||||
|
||||
plot.mask_components = [openmc.Material()]
|
||||
plot.mask_background = (255, 255, 255)
|
||||
plot.mask_background = 'white'
|
||||
plot.mask_background = (255, 255, 255)
|
||||
|
||||
plot.overlap_color = (255, 211, 0)
|
||||
plot.overlap_color = 'yellow'
|
||||
|
|
@ -70,7 +71,7 @@ def test_highlight_domains():
|
|||
plots.highlight_domains(model.geometry, mats)
|
||||
|
||||
|
||||
def test_to_xml_element(myplot):
|
||||
def test_xml_element(myplot):
|
||||
elem = myplot.to_xml_element()
|
||||
assert 'id' in elem.attrib
|
||||
assert 'color_by' in elem.attrib
|
||||
|
|
@ -80,10 +81,19 @@ def test_to_xml_element(myplot):
|
|||
assert elem.find('pixels') is not None
|
||||
assert elem.find('background').text == '0 0 0'
|
||||
|
||||
newplot = openmc.Plot.from_xml_element(elem)
|
||||
attributes = ('id', 'color_by', 'filename', 'type', 'basis', 'level',
|
||||
'meshlines', 'show_overlaps', 'origin', 'width', 'pixels',
|
||||
'background', 'mask_background')
|
||||
for attr in attributes:
|
||||
assert getattr(newplot, attr) == getattr(myplot, attr), attr
|
||||
|
||||
|
||||
def test_plots(run_in_tmpdir):
|
||||
p1 = openmc.Plot(name='plot1')
|
||||
p1.origin = (5., 5., 5.)
|
||||
p2 = openmc.Plot(name='plot2')
|
||||
p2.origin = (-3., -3., -3.)
|
||||
plots = openmc.Plots([p1, p2])
|
||||
assert len(plots) == 2
|
||||
|
||||
|
|
@ -92,3 +102,9 @@ def test_plots(run_in_tmpdir):
|
|||
assert len(plots) == 3
|
||||
|
||||
plots.export_to_xml()
|
||||
|
||||
# from_xml
|
||||
new_plots = openmc.Plots.from_xml()
|
||||
assert len(plots)
|
||||
assert plots[0].origin == p1.origin
|
||||
assert plots[1].origin == p2.origin
|
||||
|
|
|
|||
|
|
@ -112,3 +112,10 @@ def test_export_to_xml(run_in_tmpdir):
|
|||
assert not s.photon_transport
|
||||
assert s.electron_treatment == 'led'
|
||||
assert s.write_initial_source == True
|
||||
assert len(s.volume_calculations) == 1
|
||||
vol = s.volume_calculations[0]
|
||||
assert vol.domain_type == 'cell'
|
||||
assert len(vol.ids) == 1
|
||||
assert vol.samples == 1000
|
||||
assert vol.lower_left == (-10., -10., -10.)
|
||||
assert vol.upper_right == (10., 10., 10.)
|
||||
|
|
|
|||
40
tests/unit_tests/test_tallies.py
Normal file
40
tests/unit_tests/test_tallies.py
Normal file
|
|
@ -0,0 +1,40 @@
|
|||
import numpy as np
|
||||
import openmc
|
||||
|
||||
|
||||
def test_xml_roundtrip(run_in_tmpdir):
|
||||
# Create a tally with all possible gizmos
|
||||
mesh = openmc.RegularMesh()
|
||||
mesh.lower_left = (-10., -10., -10.)
|
||||
mesh.upper_right = (10., 10., 10.,)
|
||||
mesh.dimension = (5, 5, 5)
|
||||
mesh_filter = openmc.MeshFilter(mesh)
|
||||
tally = openmc.Tally()
|
||||
tally.filters = [mesh_filter]
|
||||
tally.nuclides = ['U235', 'I135', 'Li6']
|
||||
tally.scores = ['total', 'fission', 'heating']
|
||||
tally.derivative = openmc.TallyDerivative(
|
||||
variable='nuclide_density', material=1, nuclide='Li6'
|
||||
)
|
||||
tally.triggers = [openmc.Trigger('rel_err', 0.025)]
|
||||
tally.triggers[0].scores = ['total', 'fission']
|
||||
tallies = openmc.Tallies([tally])
|
||||
|
||||
# Roundtrip through XML and make sure we get what we started with
|
||||
tallies.export_to_xml()
|
||||
new_tallies = openmc.Tallies.from_xml()
|
||||
assert len(new_tallies) == 1
|
||||
new_tally = new_tallies[0]
|
||||
assert new_tally.id == tally.id
|
||||
assert len(new_tally.filters) == 1
|
||||
assert isinstance(new_tally.filters[0], openmc.MeshFilter)
|
||||
assert np.allclose(new_tally.filters[0].mesh.lower_left, mesh.lower_left)
|
||||
assert new_tally.nuclides == tally.nuclides
|
||||
assert new_tally.scores == tally.scores
|
||||
assert new_tally.derivative.variable == tally.derivative.variable
|
||||
assert new_tally.derivative.material == tally.derivative.material
|
||||
assert new_tally.derivative.nuclide == tally.derivative.nuclide
|
||||
assert len(new_tally.triggers) == 1
|
||||
assert new_tally.triggers[0].trigger_type == tally.triggers[0].trigger_type
|
||||
assert new_tally.triggers[0].threshold == tally.triggers[0].threshold
|
||||
assert new_tally.triggers[0].scores == tally.triggers[0].scores
|
||||
Loading…
Add table
Add a link
Reference in a new issue