Skip to content

Commit

Permalink
add optimized model alongside serpent-idential model
Browse files Browse the repository at this point in the history
  • Loading branch information
@yardasol
yardasol committed Oct 20, 2022
1 parent 375a100 commit 7afaf2f
Show file tree
Hide file tree
Showing 9 changed files with 3,424 additions and 134 deletions.
70 changes: 28 additions & 42 deletions examples/msbr/_control_rods.py
View file
Original file line number Diff line number Diff line change
@@ -1,55 +1,41 @@
import openmc
import numpy as np

def shared_cr_geometry():
fuel_hole = openmc.ZCylinder(r=5.08, name='cr_fuel_hole')
elem_bound = openmc.rectangular_prism(7.62*2, 7.62*2) # Pin outer boundary

gr_sq_neg = openmc.rectangular_prism(7.23646*2, 7.23646*2, corner_radius=0.99) # Graphite square
# params for
r_d = 1.16
e_d = 2 * r_d / np.sqrt(3)
r_dt = 0.8
r_c = 0.18
l1 = 5.8801
l2 = 6.505
l3 = 8.03646
ul = openmc.model.hexagonal_prism(origin=(-l1, l2), edge_length=e_d, orientation='x', corner_radius=r_c)
br = openmc.model.hexagonal_prism(origin=(l1, -l2), edge_length=e_d, orientation='x',corner_radius=r_c)
lb = openmc.model.hexagonal_prism(origin=(-l2, -l1), edge_length=e_d, orientation='y',corner_radius=r_c)
ru = openmc.model.hexagonal_prism(origin=(l2, l1), edge_length=e_d, orientation='y',corner_radius=r_c)
ul_t = openmc.ZCylinder(-l1, -l3, r_dt, name='corner_ul_tip')
br_t = openmc.ZCylinder(l1, l3, r_dt, name='corner_br_tip')
ru_t = openmc.ZCylinder(-l3, l1, r_dt, name='corner_ru_tip')
lb_t = openmc.ZCylinder(l3, -l1, r_dt, name='corner_lb_tip')

gr_corners = elem_bound & (ul | br | lb | ru |
-ul_t | -br_t | -ru_t | -lb_t)

inter_elem_channel = (~gr_sq_neg & elem_bound &
~ul & ~br & ~lb & ~ru &
+ul_t & +ru_t & +ru_t & +lb_t)


return fuel_hole, gr_sq_neg, gr_corners, inter_elem_channel

def control_rod(fuel_hole, gr_sq_neg, gr_corners, inter_elem_channel, fuel, moder):
def control_rod(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, fuel_hole, fuel, moder):
s1 = openmc.ZCylinder(r=4.7625, name='control_rod')

c1 = openmc.Cell(fill=moder, region=-s1, name='control_rod')
c2 = openmc.Cell(fill=fuel, region=(+s1 & -fuel_hole), name='cr_fuel_inner')
moderator_block = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg | gr_corners), name='cr_moderator')
fuel_outer = openmc.Cell(fill=fuel, region=inter_elem_channel, name='cr_fuel_outer')
c3 = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg & inter_elem_channel),
name='cr_moderator')
c4 = openmc.Cell(fill=fuel, region= (~gr_sq_neg & elem_bound & inter_elem_channel),
name='cr_fuel_outer')

#universe_id=3
cr = openmc.Universe(name='control_rod')
cr.add_cells([c1, c2, moderator_block, fuel_outer])
cr = openmc.Universe(name='control_rod',
cells=[c1, c2, c3, c4])

for i, (corner, name) in enumerate(gr_corners):
if i < 4:
cr.add_cell(openmc.Cell(fill=moder, region=~corner, name=f'cr_moderator_{name}'))
else:
cr.add_cell(openmc.Cell(fill=moder, region=+corner, name=f'cr_moderator_{name}'))

return cr

def control_rod_channel(fuel_hole, gr_sq_neg, gr_corners, inter_elem_channel, fuel, moder):
def control_rod_channel(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, fuel_hole, fuel, moder):
c1 = openmc.Cell(fill=fuel, region=(-fuel_hole), name='crc_fuel_inner')
moderator_block = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg | gr_corners), name='crc_moderator')
fuel_outer = openmc.Cell(fill=fuel, region=inter_elem_channel, name='crc_fuel_outer')
c2 = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg & inter_elem_channel), name='crc_moderator')
c3 = openmc.Cell(fill=fuel, region=(~gr_sq_neg & elem_bound & inter_elem_channel), name='crc_fuel_outer')

# universe_id=4
crc = openmc.Universe(name='control_rod_channel')
crc.add_cells([c1, moderator_block, fuel_outer])
crc = openmc.Universe(name='control_rod_channel',
cells=[c1, c2, c3])

for i, (corner, name) in enumerate(gr_corners):
if i < 4:
crc.add_cell(openmc.Cell(fill=moder, region=~corner, name=f'crc_moderator_{name}'))
else:
crc.add_cell(openmc.Cell(fill=moder, region=+corner, name=f'crc_moderator_{name}'))

return crc
106 changes: 106 additions & 0 deletions examples/msbr/_control_rods_optimized.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
import openmc
import numpy as np

def shared_cr_geometry():
fuel_hole = openmc.ZCylinder(r=5.08, name='cr_fuel_hole')
elem_bound = openmc.rectangular_prism(7.62*2, 7.62*2) # Pin outer boundary
eb_minx, eb_maxx, eb_miny, eb_maxy = list(elem_bound.get_surfaces().values())

gr_sq_neg = openmc.rectangular_prism(7.23646*2, 7.23646*2, corner_radius=0.99) # Graphite square
# params for
(gr_minx, gr_maxx, gr_miny, gr_maxy,
gr_cyl_lb, gr_cyl_minx, gr_cyl_miny,
gr_cyl_ul, gr_cyl_maxy, gr_cyl_br,
gr_cyl_maxx, gr_cyl_ru) = list(gr_sq_neg.get_surfaces().values())

r_d = 1.16
e_d = 2 * r_d / np.sqrt(3)
r_dt = 0.8
r_c = 0.18
l1 = 5.8801
l2 = 6.505
l3 = 8.03646
ul = openmc.model.hexagonal_prism(origin=(-l1, l2), edge_length=e_d, orientation='x', corner_radius=r_c)
br = openmc.model.hexagonal_prism(origin=(l1, -l2), edge_length=e_d, orientation='x',corner_radius=r_c)
lb = openmc.model.hexagonal_prism(origin=(-l2, -l1), edge_length=e_d, orientation='y',corner_radius=r_c)
ru = openmc.model.hexagonal_prism(origin=(l2, l1), edge_length=e_d, orientation='y',corner_radius=r_c)
ul_t = openmc.ZCylinder(-l1, -l3, r_dt, name='corner_ul_tip')
br_t = openmc.ZCylinder(l1, l3, r_dt, name='corner_br_tip')
ru_t = openmc.ZCylinder(-l3, l1, r_dt, name='corner_ru_tip')
lb_t = openmc.ZCylinder(l3, -l1, r_dt, name='corner_lb_tip')

gr_ul = ul & -eb_maxy & +gr_maxy
gr_ul_fill = ul & +gr_cyl_ul & -gr_maxy & -gr_cyl_minx & +gr_cyl_maxy
gr_br = br & +eb_miny & -gr_miny
gr_br_fill = br & +gr_cyl_br & +gr_miny & +gr_cyl_maxx & -gr_cyl_miny
gr_lb = lb & +eb_minx & -gr_minx
gr_lb_fill = lb & +gr_cyl_lb & +gr_minx & -gr_cyl_minx & -gr_cyl_miny
gr_ru = ru & -eb_maxx & +gr_maxx
gr_ru_fill = ru & +gr_cyl_ru & -gr_maxx & +gr_cyl_maxx & +gr_cyl_maxy

gr_ul_t = -ul_t & +eb_miny & -gr_miny
gr_br_t = -br_t & -eb_maxy & +gr_maxy
gr_ru_t = -ru_t & +eb_minx & -gr_minx
gr_lb_t = -lb_t & -eb_maxx & +gr_maxx

gr_corners = (gr_ul, gr_br, gr_lb, gr_ru, gr_ul_fill, gr_br_fill, gr_lb_fill, gr_ru_fill)
gr_t = (gr_ul_t, gr_br_t, gr_ru_t, gr_lb_t)

inter_elem_channel = (~gr_sq_neg & elem_bound &
~ul & ~br & ~lb & ~ru &
+ul_t & +br_t & +ru_t & +lb_t)


return fuel_hole, gr_sq_neg, gr_corners, gr_t, inter_elem_channel

def control_rod(fuel_hole, gr_sq_neg, gr_corners, gr_t, inter_elem_channel, fuel, moder):
s1 = openmc.ZCylinder(r=4.7625, name='control_rod')

c1 = openmc.Cell(fill=moder, region=-s1, name='control_rod')
c2 = openmc.Cell(fill=fuel, region=(+s1 & -fuel_hole), name='cr_fuel_inner')
c3 = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg), name='cr_moderator')
c4 = openmc.Cell(fill=fuel, region=inter_elem_channel, name='cr_fuel_outer')

c3_ul = openmc.Cell(fill=moder, region=gr_corners[0], name='cr_moderator_ul')
c3_br = openmc.Cell(fill=moder, region=gr_corners[1], name='cr_moderator_br')
c3_ru = openmc.Cell(fill=moder, region=gr_corners[2], name='cr_moderator_ru')
c3_lb = openmc.Cell(fill=moder, region=gr_corners[3], name='cr_moderator_lb')
c3_ulf = openmc.Cell(fill=moder, region=gr_corners[4], name='cr_moderator_ul_fill')
c3_brf = openmc.Cell(fill=moder, region=gr_corners[5], name='cr_moderator_br_fill')
c3_ruf = openmc.Cell(fill=moder, region=gr_corners[6], name='cr_moderator_ru_fill')
c3_lbf = openmc.Cell(fill=moder, region=gr_corners[7], name='cr_moderator_lb_fill')
c3_ul_t = openmc.Cell(fill=moder, region=gr_t[0], name='cr_moderator_ul_t')
c3_br_t = openmc.Cell(fill=moder, region=gr_t[1], name='cr_moderator_br_t')
c3_ru_t = openmc.Cell(fill=moder, region=gr_t[2], name='cr_moderator_ru_t')
c3_lb_t = openmc.Cell(fill=moder, region=gr_t[3], name='cr_moderator_lb_t')


#universe_id=3
cr = openmc.Universe(name='control_rod')
cr.add_cells([c1, c2, c3, c4, c3_ul, c3_br, c3_ru, c3_lb, c3_ulf, c3_brf, c3_ruf, c3_lbf, c3_ul_t, c3_br_t, c3_ru_t, c3_lb_t])

return cr

def control_rod_channel(fuel_hole, gr_sq_neg, gr_corners, gr_t, inter_elem_channel, fuel, moder):
c1 = openmc.Cell(fill=fuel, region=(-fuel_hole), name='crc_fuel_inner')
c2 = openmc.Cell(fill=moder, region=(+fuel_hole & gr_sq_neg), name='crc_moderator')
c3 = openmc.Cell(fill=fuel, region=inter_elem_channel, name='crc_fuel_outer')

c2_ul = openmc.Cell(fill=moder, region=gr_corners[0], name='crc_moderator_ul')
c2_br = openmc.Cell(fill=moder, region=gr_corners[1], name='crc_moderator_br')
c2_ru = openmc.Cell(fill=moder, region=gr_corners[2], name='crc_moderator_ru')
c2_lb = openmc.Cell(fill=moder, region=gr_corners[3], name='crc_moderator_lb')
c2_ulf = openmc.Cell(fill=moder, region=gr_corners[4], name='crc_moderator_ul_fill')
c2_brf = openmc.Cell(fill=moder, region=gr_corners[5], name='crc_moderator_br_fill')
c2_ruf = openmc.Cell(fill=moder, region=gr_corners[6], name='crc_moderator_ru_fill')
c2_lbf = openmc.Cell(fill=moder, region=gr_corners[7], name='crc_moderator_lb_fill')
c2_ul_t = openmc.Cell(fill=moder, region=gr_t[0], name='crc_moderator_ul_t')
c2_br_t = openmc.Cell(fill=moder, region=gr_t[1], name='crc_moderator_br_t')
c2_ru_t = openmc.Cell(fill=moder, region=gr_t[2], name='crc_moderator_ru_t')
c2_lb_t = openmc.Cell(fill=moder, region=gr_t[3], name='crc_moderator_lb_t')

# universe_id=4
crc = openmc.Universe(name='control_rod_channel')
crc.add_cells([c1, c2, c3, c2_ul, c2_br, c2_ru, c2_lb, c2_ulf, c2_brf, c2_ruf, c2_lbf, c2_ul_t, c2_br_t, c2_ru_t, c2_lb_t])

return crc
87 changes: 36 additions & 51 deletions examples/msbr/_core_elements.py
View file
Original file line number Diff line number Diff line change
@@ -1,53 +1,33 @@
import openmc
import numpy as np

def _bound_zone_cells(cells_tuples, levels):
"""Helper function that moves Zone IA and Zone IIA cells to the
appropriate height."""
cell_list = []
n_levels = len(cells_tuples)
for i, cells in enumerate(cells_tuples):
lower_bound = levels[i]
upper_bound = levels[i+1]
if i == 0:
lower_bound = None
upper_bound = levels[i]
elif i == n_levels - 1:
lower_bound = levels[i-1]
upper_bound = None
else:
lower_bound = levels[i-1]
upper_bound = levels[i]
for j, cell in enumerate(cells):
cell.region = cell.region & +lower_bound & -upper_bound
if lower_bound is None:
cell.region = cell.region & -upper_bound
elif upper_bound is None:
cell.region = cell.region & +lower_bound
else:
cell.region = cell.region & +lower_bound & -upper_bound
cell_list.append(cell)
return cell_list

def shared_elem_geometry():
"""Surfaces and regions shared by Zone IA and Zone IIA elements.
Specs found in Robertson, 1971 Fig 3.4 (p. 17) and Fig 3.5 (p.18)"""

elem_bound = openmc.rectangular_prism(5.08*2, 5.08*2) # Pin outer boundary

gr_sq_neg = openmc.rectangular_prism(4.953*2, 4.953*2, corner_radius=0.46) # Graphite square

# params for main pin section for both I-A and II-A
r_d = 0.66802
l1 = 4.28498
l2 = 4.53898
l3 = 5.62102
ul = openmc.ZCylinder(-l1, l2, r_d, name='corner_ul')
br = openmc.ZCylinder(l1, -l2, r_d, name='corner_br')
lb = openmc.ZCylinder(-l2, -l1, r_d, name='corner_lb')
ru = openmc.ZCylinder(l2, l1, r_d, name='corner_ru')
ul_t = openmc.ZCylinder(-l1, -l3, r_d, name='corner_ul_tip')
br_t = openmc.ZCylinder(l1, l3, r_d, name='corner_br_tip')
ru_t = openmc.ZCylinder(-l3, l1, r_d, name='corner_ru_tip')
lb_t = openmc.ZCylinder(l3, -l1, r_d, name='corner_lb_tip')

gr_corners = elem_bound & (-ul | -br | -lb | -ru |
-ul_t | -br_t | -ru_t | -lb_t)

inter_elem_channel = (~gr_sq_neg & elem_bound &
+ul & +br & +lb & +ru &
+ul_t & +ru_t & +ru_t & +lb_t)
gr_round_4 = openmc.ZCylinder(r=2.2225, name='gr_round_4')

return elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4

def zoneIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, moder, fuel, hast):
"""Zone IA element. Specs found in Robertson, 1971 Fig 3.4 (p. 17)"""
elem_levels = [0.0, 22.86, 419.10, 438.15, 445.135, 449.58]
elem_levels = [22.86, 419.10, 438.15, 445.135]
level_bounds = []
for level in elem_levels:
level_bounds.append(openmc.ZPlane(z0=level))
Expand All @@ -57,7 +37,7 @@ def zoneIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, mo
h = 12.66
theta = np.arctan(4.953 / h)
r2 = (1 / np.cos(theta))**2 - 1
s3 = openmc.ZCone(z0=h + elem_levels[3], r2=r2, name='cone_i')
s3 = openmc.ZCone(z0=h + elem_levels[2], r2=r2, name='cone_i')

c1 = openmc.Cell(fill=fuel, region=(-s2), name='ia_fuel_inner_1')
c2 = openmc.Cell(fill=moder, region=(+s2 & -s1), name='ia_moderator_1')
Expand All @@ -67,11 +47,13 @@ def zoneIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, mo
# I-A main (lower 2)
c4 = c1.clone(clone_materials=False)
c4.name = 'ia_fuel_inner_main'
c5 = openmc.Cell(fill=moder, region=(+s2 &
gr_sq_neg |
gr_corners), name='ia_moderator_main')
c6 = openmc.Cell(fill=fuel, region=(inter_elem_channel), name='ia_fuel_outer_main')
iam = (c4, c5, c6)
c5 = openmc.Cell(fill=moder, region=(+s2 & gr_sq_neg & inter_elem_channel),
name='ia_moderator_main')
c6 = openmc.Cell(fill=fuel, region=(~gr_sq_neg & elem_bound & inter_elem_channel),
name='ia_fuel_outer_main')
iam = [c4, c5, c6]
for (corner, name) in gr_corners:
iam.append(openmc.Cell(fill=moder, region=+corner, name=f'ia_moderator_main_{name}'))

# I-A 2 (upper 1)
c7 = c1.clone(clone_materials=False)
Expand Down Expand Up @@ -103,7 +85,7 @@ def zoneIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, mo

def zoneIIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, moder, fuel):
"""Zone IIA element. Specs found in Robertson, 1971 Fig 3.5 (p. 18)"""
elem_levels = [0.0, 434.34, 436.88, 439.42, 441.96, 449.58]
elem_levels = [434.34, 436.88, 439.42, 441.96]
level_bounds = []
for level in elem_levels:
level_bounds.append(openmc.ZPlane(z0=level))
Expand All @@ -113,21 +95,24 @@ def zoneIIA(elem_bound, gr_corners, gr_sq_neg, inter_elem_channel, gr_round_4, m
h = 6.5
theta = np.arctan(3.65125 / h)
r2 = (1 / np.cos(theta))**2 - 1
s4 = openmc.ZCone(z0=h + elem_levels[3], r2=r2, name='cone_ii')
s4 = openmc.ZCone(z0=h + elem_levels[2], r2=r2, name='cone_ii')

# II-A main (lower 1)
c1 = openmc.Cell(fill=fuel, region=(-s1), name='iia_fuel_inner_main')
c2 = openmc.Cell(fill=moder, region=(+s1 & gr_sq_neg | gr_corners), name='iia_moderator_main')
c3 = openmc.Cell(fill=fuel, region=(inter_elem_channel), name='iia_fuel_outer_main')
c3.name = 'iia_fuel_outer_main'
iiam = (c1, c2, c3)
c2 = openmc.Cell(fill=moder, region=(+s1 & gr_sq_neg & inter_elem_channel), name='iia_moderator_main')
c3 = openmc.Cell(fill=fuel, region=(~gr_sq_neg & elem_bound & inter_elem_channel), name='iia_fuel_outer_main')
iiam = [c1, c2, c3]

# II-A 2 (upper 1)
c4 = openmc.Cell(fill=fuel, region=(-s2), name='iia_fuel_inner_2')
c5 = openmc.Cell(fill=moder, region=(+s2 & gr_sq_neg | gr_corners), name='iia_moderator_2')
c5 = openmc.Cell(fill=moder, region=(+s2 & gr_sq_neg & inter_elem_channel), name='iia_moderator_2')
c6 = c3.clone(clone_materials=False)
c6.name = 'iia_fuel_outer_2'
iia2 = (c4, c5, c6)
iia2 = [c4, c5, c6]

for (corner, name) in gr_corners:
iiam.append(openmc.Cell(fill=moder, region=+corner, name=f'iia_moderator_main_{name}'))
iia2.append(openmc.Cell(fill=moder, region=+corner, name=f'iia_moderator_2_{name}'))

# II-A 3 (upper 2)
c7 = c4.clone(clone_materials=False)
Expand Down
Loading

0 comments on commit 7afaf2f

Please sign in to comment.