crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 21:55:41 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Update TRISO tests

Browse source
This commit is contained in:
amandalund 2018-11-14 10:54:24 -06:00
parent c443208e21
commit 66fe2185f3
5 changed files with 121 additions and 130 deletions
Download patch file Download diff file Expand all files Collapse all files

26
openmc/model/triso.py
View file

@ -882,21 +882,21 @@ def _create_container(region, sphere_radius):
if p1.surface.x0 > p2.surface.x0:
p1, p2 = p2, p1
length = p2.surface.x0 - p1.surface.x0
center = (length/2, cyl.surface.y0, cyl.surface.z0)
center = (p1.surface.x0 + length/2, cyl.surface.y0, cyl.surface.z0)
# Calculate the parameters for a cylinder along the y-axis
elif axis == 'y':
if p1.surface.y0 > p2.surface.y0:
p1, p2 = p2, p1
length = p2.surface.y0 - p1.surface.y0
center = (cyl.surface.x0, length/2, cyl.surface.z0)
center = (cyl.surface.x0, p1.surface.y0 + length/2, cyl.surface.z0)
# Calculate the parameters for a cylinder along the z-axis
else:
if p1.surface.z0 > p2.surface.z0:
p1, p2 = p2, p1
length = p2.surface.z0 - p1.surface.z0
center = (cyl.surface.x0, cyl.surface.y0, length/2)
center = (cyl.surface.x0, cyl.surface.y0, p1.surface.z0 + length/2)
# Make sure the half-spaces are on the correct side of the surfaces
if cyl.side != '-' or p1.side != '+' or p2.side != '-':
@ -1424,23 +1424,3 @@ def pack_spheres(radius, region, pf=None, num_spheres=None, initial_pf=0.3,
_close_random_pack(domain, spheres, contraction_rate)
return spheres + domain.center
def create_trisos(outer_radius, fill, centers):
"""Create TRISO particles at the given coordinates.
Parameters
----------
outer_radius : float
Outer radius of TRISO particle
fill : openmc.Universe
Universe which contains all layers of the TRISO particle
centers : Iterable of float
Cartesian coordinates of the centers of the TRISO particles in cm
Returns
-------
list of openmc.model.TRISO
List of TRISO particles in the domain.
"""
return [TRISO(outer_radius, fill, c) for c in centers]

14
tests/regression_tests/triso/inputs_true.dat
View file

@ -5,7 +5,7 @@
<cell id="15" material="15" region="10 -11" universe="9" />
<cell id="16" material="16" region="11 -12" universe="9" />
<cell id="17" material="17" region="12" universe="9" />
<cell fill="10" id="134" region="122 -123 124 -125 126 -127" universe="0" />
<cell fill="10" id="18" region="13 -14 15 -16 17 -18" universe="0" />
<cell fill="9" id="135" region="-128" translation="-0.0014062011817633224 -0.015257747167296776 -0.09184766966024316" universe="35" />
<cell fill="9" id="136" region="-129" translation="0.10909902639777394 -0.004179177767328068 -0.046240930185910245" universe="23" />
<cell fill="9" id="137" region="-130" translation="0.13878336652768952 -0.06900710536877097 -0.03848821497485455" universe="18" />
@ -220,12 +220,12 @@
<surface coeffs="0.0 0.0 0.0 0.03125" id="10" type="sphere" />
<surface coeffs="0.0 0.0 0.0 0.03475" id="11" type="sphere" />
<surface coeffs="0.0 0.0 0.0 0.03825" id="12" type="sphere" />
<surface boundary="reflective" coeffs="-0.5" id="122" type="x-plane" />
<surface boundary="reflective" coeffs="0.5" id="123" type="x-plane" />
<surface boundary="reflective" coeffs="-0.5" id="124" type="y-plane" />
<surface boundary="reflective" coeffs="0.5" id="125" type="y-plane" />
<surface boundary="reflective" coeffs="-0.5" id="126" type="z-plane" />
<surface boundary="reflective" coeffs="0.5" id="127" type="z-plane" />
<surface boundary="reflective" coeffs="-0.5" id="13" type="x-plane" />
<surface boundary="reflective" coeffs="0.5" id="14" type="x-plane" />
<surface boundary="reflective" coeffs="-0.5" id="15" type="y-plane" />
<surface boundary="reflective" coeffs="0.5" id="16" type="y-plane" />
<surface boundary="reflective" coeffs="-0.5" id="17" type="z-plane" />
<surface boundary="reflective" coeffs="0.5" id="18" type="z-plane" />
<surface coeffs="-0.0014062011817633224 -0.015257747167296776 -0.09184766966024316 0.04225" id="128" type="sphere" />
<surface coeffs="0.10909902639777394 -0.004179177767328068 -0.046240930185910245 0.04225" id="129" type="sphere" />
<surface coeffs="0.13878336652768952 -0.06900710536877097 -0.03848821497485455 0.04225" id="130" type="sphere" />

2
tests/regression_tests/triso/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
1.683226E+00 7.383559E-02
1.683227E+00 7.383566E-02

16
tests/regression_tests/triso/test.py
View file

@ -53,19 +53,21 @@ class TRISOTestHarness(PyAPITestHarness):
c5 = openmc.Cell(fill=opyc, region=+spheres[3])
inner_univ = openmc.Universe(cells=[c1, c2, c3, c4, c5])
outer_radius = 422.5*1e-4
trisos = openmc.model.pack_trisos(
radius=outer_radius, fill=inner_univ, domain_shape='cube',
domain_length=1., domain_center=(0., 0., 0.), n_particles=100)
# Define box to contain lattice
# Define box to contain lattice and to pack TRISO particles in
min_x = openmc.XPlane(x0=-0.5, boundary_type='reflective')
max_x = openmc.XPlane(x0=0.5, boundary_type='reflective')
min_y = openmc.YPlane(y0=-0.5, boundary_type='reflective')
max_y = openmc.YPlane(y0=0.5, boundary_type='reflective')
min_z = openmc.ZPlane(z0=-0.5, boundary_type='reflective')
max_z = openmc.ZPlane(z0=0.5, boundary_type='reflective')
box = openmc.Cell(region=+min_x & -max_x & +min_y & -max_y & +min_z & -max_z)
box_region = +min_x & -max_x & +min_y & -max_y & +min_z & -max_z
box = openmc.Cell(region=box_region)
outer_radius = 422.5*1e-4
centers = openmc.model.pack_spheres(radius=outer_radius,
region=box_region, num_spheres=100)
trisos = [openmc.model.TRISO(outer_radius, inner_univ, c)
for c in centers]
# Create lattice
ll, ur = box.region.bounding_box

193
tests/unit_tests/test_model_triso.py
View file

@ -10,19 +10,54 @@ import pytest
import scipy.spatial
_RADIUS = 0.1
_PACKING_FRACTION = 0.35
_RADIUS = 4.25e-2
domain_params = [
{'shape': 'cube', 'length': 0.75, 'radius': 0., 'volume': 0.75**3},
{'shape': 'cylinder', 'length': 0.5, 'radius': 0.5, 'volume': 0.5*pi*0.5**2},
{'shape': 'sphere', 'length': 0., 'radius': 0.5, 'volume': 4/3*pi*0.5**3}
]
_SHAPES = ['rectangular_prism', 'cylinder', 'sphere', 'spherical_shell']
_VOLUMES = [1.**3, 1.*pi*1.**2, 4/3*pi*1.**3, 4/3*pi*(1.**3 - 0.5**3)]
@pytest.fixture(scope='module', params=domain_params,
ids=['cube', 'cylinder', 'sphere'])
def domain(request):
return request.param
@pytest.fixture(scope='module')
def container(request):
return request.getfixturevalue(request.param)
@pytest.fixture(scope='module')
def centers(request):
container = request.getfixturevalue(request.param)
return openmc.model.pack_spheres(radius=_RADIUS, region=container,
pf=_PACKING_FRACTION, initial_pf=0.2)
@pytest.fixture(scope='module')
def rectangular_prism():
min_x = openmc.XPlane(x0=-0.5)
max_x = openmc.XPlane(x0=0.5)
min_y = openmc.YPlane(y0=-0.5)
max_y = openmc.YPlane(y0=0.5)
min_z = openmc.ZPlane(z0=-0.5)
max_z = openmc.ZPlane(z0=0.5)
return +min_x & -max_x & +min_y & -max_y & +min_z & -max_z
@pytest.fixture(scope='module')
def cylinder():
cylinder = openmc.ZCylinder(R=1.)
min_z = openmc.ZPlane(z0=-0.5)
max_z = openmc.ZPlane(z0=0.5)
return +min_z & -max_z & -cylinder
@pytest.fixture(scope='module')
def sphere():
sphere = openmc.Sphere(R=1.)
return -sphere
@pytest.fixture(scope='module')
def spherical_shell():
sphere = openmc.Sphere(R=1.)
inner_sphere = openmc.Sphere(R=0.5)
return -sphere & +inner_sphere
@pytest.fixture(scope='module')
@ -33,73 +68,68 @@ def triso_universe():
return univ
@pytest.fixture(scope='module')
def trisos(domain, triso_universe):
trisos = openmc.model.pack_trisos(
radius=_RADIUS,
fill=triso_universe,
domain_shape=domain['shape'],
domain_length=domain['length'],
domain_radius=domain['radius'],
domain_center=(0., 0., 0.),
initial_packing_fraction=0.2,
packing_fraction=_PACKING_FRACTION
)
return trisos
def test_overlap(trisos):
"""Check that no TRISO particles overlap."""
centers = [t.center for t in trisos]
@pytest.mark.parametrize('centers', _SHAPES, indirect=True)
def test_overlap(centers):
"""Check that none of the spheres in the packed configuration overlap."""
# Create KD tree for quick nearest neighbor search
tree = scipy.spatial.cKDTree(centers)
# Find distance to nearest neighbor for all particles
# Find distance to nearest neighbor for all spheres
d = tree.query(centers, k=2)[0]
# Get the smallest distance between any two particles
# Get the smallest distance between any two spheres
d_min = min(d[:, 1])
assert d_min > 2*_RADIUS or d_min == pytest.approx(2*_RADIUS)
def test_contained(trisos, domain):
"""Make sure all particles are entirely contained within the domain."""
if domain['shape'] == 'cube':
x = max(np.hstack([abs(t.center) for t in trisos])) + _RADIUS
assert x < 0.5*domain['length'] or x == pytest.approx(0.5*domain['length'])
@pytest.mark.parametrize('centers,shape', zip(_SHAPES, _SHAPES),
indirect=['centers'])
def test_contained(centers, shape):
"""Make sure all spheres are entirely contained within the domain."""
if shape == 'rectangular_prism':
x = np.amax(abs(centers)) + _RADIUS
assert x < 0.5 or x == pytest.approx(0.5)
elif domain['shape'] == 'cylinder':
r = max([norm(t.center[0:2]) for t in trisos]) + _RADIUS
z = max([abs(t.center[2]) for t in trisos]) + _RADIUS
assert r < domain['radius'] or r == pytest.approx(domain['radius'])
assert z < 0.5*domain['length'] or z == pytest.approx(0.5*domain['length'])
elif shape == 'cylinder':
r = max(np.linalg.norm(centers[:,0:2], axis=1)) + _RADIUS
z = max(abs(centers[:,2])) + _RADIUS
assert r < 1. or r == pytest.approx(1.)
assert z < 0.5 or z == pytest.approx(0.5)
elif domain['shape'] == 'sphere':
r = max([norm(t.center) for t in trisos]) + _RADIUS
assert r < domain['radius'] or r == pytest.approx(domain['radius'])
elif shape == 'sphere':
r = max(np.linalg.norm(centers, axis=1)) + _RADIUS
assert r < 1. or r == pytest.approx(1.)
elif shape == 'spherical_shell':
d = np.linalg.norm(centers, axis=1)
r_max = max(d) + _RADIUS
r_min = min(d) - _RADIUS
assert r_max < 1. or r_max == pytest.approx(1.)
assert r_min > 0.5 or r_min == pytest.approx(0.5)
def test_packing_fraction(trisos, domain):
@pytest.mark.parametrize('centers,volume', zip(_SHAPES, _VOLUMES),
indirect=['centers'])
def test_packing_fraction(centers, volume):
"""Check that the actual PF is close to the requested PF."""
pf = len(trisos)*4/3*pi*_RADIUS**3/domain['volume']
pf = len(centers) * 4/3 * pi *_RADIUS**3 / volume
assert pf == pytest.approx(_PACKING_FRACTION, rel=1e-2)
def test_n_particles(triso_universe):
"""Check that the function returns the correct number of particles"""
trisos = openmc.model.pack_trisos(
radius=_RADIUS, fill=triso_universe, domain_shape='cube',
domain_length=1.0, n_particles=800
@pytest.mark.parametrize('container', _SHAPES, indirect=True)
def test_num_spheres(container):
"""Check that the function returns the correct number of spheres"""
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=container, num_spheres=50
)
assert len(trisos) == 800
assert len(centers) == 50
def test_triso_lattice(triso_universe):
trisos = openmc.model.pack_trisos(
radius=_RADIUS, fill=triso_universe, domain_shape='cube',
domain_length=1.0, domain_center=(0., 0., 0.), packing_fraction=0.2
def test_triso_lattice(triso_universe, rectangular_prism):
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=rectangular_prism, pf=0.2
)
trisos = [openmc.model.TRISO(_RADIUS, triso_universe, c) for c in centers]
lower_left = np.array((-.5, -.5, -.5))
upper_right = np.array((.5, .5, .5))
@ -112,50 +142,29 @@ def test_triso_lattice(triso_universe):
)
def test_domain_input(triso_universe):
# Invalid domain shape
def test_container_input(triso_universe):
# Invalid container shape
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, n_particles=100,
domain_shape='circle'
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=+openmc.Sphere(R=1.), num_spheres=100
)
# Don't specify domain length on a cube
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, n_particles=100,
domain_shape='cube'
)
# Don't specify domain radius on a sphere
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, n_particles=100,
domain_shape='sphere'
)
def test_packing_fraction_input(triso_universe):
# Provide neither packing fraction nor number of particles
def test_packing_fraction_input(sphere):
# Provide neither packing fraction nor number of spheres
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, domain_shape='cube',
domain_length=10
)
# Provide both packing fraction and number of particles
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, domain_shape='cube',
domain_length=10, n_particles=100, packing_fraction=0.2
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=sphere
)
# Specify a packing fraction that is too high for CRP
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, domain_shape='cube',
domain_length=10, packing_fraction=1
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=sphere, pf=1.
)
# Specify a packing fraction that is too high for RSP
with pytest.raises(ValueError):
trisos = openmc.model.pack_trisos(
radius=1, fill=triso_universe, domain_shape='cube',
domain_length=10, packing_fraction=0.5,
initial_packing_fraction=0.4
centers = openmc.model.pack_spheres(
radius=_RADIUS, region=sphere, pf=0.5, initial_pf=0.4
)