From 19579740738ed99afa453735111086ffba8a3643 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 09:06:41 -0600 Subject: [PATCH 1/8] Don't deepcopy surfaces PR #1393 in OpenMC changed how memoization during XML writing is done. Specifically, the check for surfaces no longer looks at the IDs already written. --- smr/smr/surfaces.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/smr/smr/surfaces.py b/smr/smr/surfaces.py index 5a53ef2..d97133f 100644 --- a/smr/smr/surfaces.py +++ b/smr/smr/surfaces.py @@ -166,8 +166,8 @@ surfs['BA IR 7'] = openmc.ZCylinder( r=burn_abs_r7, name='BA IR 7') surfs['BA IR 8'] = openmc.ZCylinder( r=burn_abs_r8, name='BA IR 8') -surfs['IT IR'] = copy.deepcopy(surfs['BA IR 5']) -surfs['IT OR'] = copy.deepcopy(surfs['BA IR 6']) +surfs['IT IR'] = surfs['BA IR 5'] +surfs['IT OR'] = surfs['BA IR 6'] # Rectangular prisms for grid spacers surfs['rod grid box'] = \ @@ -184,8 +184,8 @@ surfs['bot support plate'] = openmc.ZPlane( surfs['top support plate'] = openmc.ZPlane( z0=top_support_plate, name='top support plate') surfs['bottom FR'] = openmc.ZPlane(z0=bottom_fuel_rod, name='bottom FR') -surfs['top lower nozzle'] = copy.deepcopy(surfs['bottom FR']) -surfs['bot lower nozzle'] = copy.deepcopy(surfs['top support plate']) +surfs['top lower nozzle'] = surfs['bottom FR'] +surfs['bot lower nozzle'] = surfs['top support plate'] # axial surfaces surfs['bot active core'] = openmc.ZPlane( @@ -193,7 +193,7 @@ surfs['bot active core'] = openmc.ZPlane( surfs['top active core'] = openmc.ZPlane( z0=top_active_core, name='top active core') -surfs['top lower thimble'] = copy.deepcopy(surfs['bot active core']) +surfs['top lower thimble'] = surfs['bot active core'] surfs['BA bot'] = openmc.ZPlane( z0=bot_burn_abs, name='bottom of BA') From 6ae045801edf6ebbcd5f301f97985e06c15483c5 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 07:09:56 -0600 Subject: [PATCH 2/8] Remove default temperatures on materials --- smr/smr/materials.py | 7 ------- 1 file changed, 7 deletions(-) diff --git a/smr/smr/materials.py b/smr/smr/materials.py index 989bda5..72a39d6 100644 --- a/smr/smr/materials.py +++ b/smr/smr/materials.py @@ -192,7 +192,6 @@ aB_bsg = aB10_bsg + aB11_bsg # Create borosilicate glass material mats['BSG'] = openmc.Material(name='Borosilicate Glass') -mats['BSG'].temperature = 300 mats['BSG'].set_density('g/cc', 2.26) mats['BSG'].add_element('O', aO_bsg, 'ao') mats['BSG'].add_element('Si', aSi_bsg, 'ao') @@ -205,7 +204,6 @@ mats['BSG'].add_nuclide('B11', aB11_bsg, 'ao') # Create 1.6% enriched UO2 fuel material mat = openmc.Material(name='1.6% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.31341) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=1.61006) @@ -213,7 +211,6 @@ mats['UO2 1.6 fresh'] = mat # Create 2.4% enriched UO2 fuel material mat = openmc.Material(name='2.4% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.29748) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=2.39993) @@ -221,7 +218,6 @@ mats['UO2 2.4 fresh'] = mat # Create 3.1% enriched UO2 fuel material mat = openmc.Material(name='3.1% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.30166) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=3.10221) @@ -229,7 +225,6 @@ mats['UO2 3.1 fresh'] = mat # Depleted versions of 1.6%, 2.4%, 3.1% fuel mat = openmc.Material(name='2.4% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.29748) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=2.39993) @@ -238,7 +233,6 @@ for nuc in _DEPLETION_NUCLIDES: mats['UO2 2.4 depleted'] = mat mat = openmc.Material(name='1.6% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.31341) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=1.61006) @@ -247,7 +241,6 @@ for nuc in _DEPLETION_NUCLIDES: mats['UO2 1.6 depleted'] = mat mat = openmc.Material(name='3.1% Enr. UO2 Fuel') -mat.temperature = 300 mat.set_density('g/cc', 10.30166) mat.add_element('O', 2., 'ao') mat.add_element('U', 1., 'ao', enrichment=3.10221) From 314331f27541c5efde8c2283365e59a6ae7435cb Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 07:13:15 -0600 Subject: [PATCH 3/8] Don't override empty guide tube universe (without spacers) --- smr/smr/assemblies.py | 2 +- smr/smr/pins.py | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/smr/smr/assemblies.py b/smr/smr/assemblies.py index dd71d97..aeba677 100644 --- a/smr/smr/assemblies.py +++ b/smr/smr/assemblies.py @@ -144,7 +144,7 @@ def assembly_universes(num_rings, num_axial, depleted): univs = {} # commonly needed universes - gtu = pins['GT empty'] + gtu = pins['GT empty stack'] gti = pins['GT empty instr'] bas = pins['BA stack'] ins = pins['IT stack'] diff --git a/smr/smr/pins.py b/smr/smr/pins.py index 676fc8a..e71adf5 100644 --- a/smr/smr/pins.py +++ b/smr/smr/pins.py @@ -230,8 +230,8 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): surfs['top upper nozzle'] ] - univs['GT empty'] = make_stack( - 'GT empty', surfaces=stack_surfs, + univs['GT empty stack'] = make_stack( + 'GT empty stack', surfaces=stack_surfs, universes=[univs['water pin'], univs['water pin'], univs['water pin'], From 5481118862b3cbf04f0d41228a344b00432847b3 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 07:38:41 -0600 Subject: [PATCH 4/8] Allow configurable fuel in ring radii --- smr/build-assembly.py | 5 +++-- smr/build-core-fresh.py | 5 +++-- smr/smr/assemblies.py | 9 +++++---- smr/smr/core.py | 9 +++++---- smr/smr/pins.py | 33 +++++++++++++++++++++------------ 5 files changed, 37 insertions(+), 24 deletions(-) diff --git a/smr/build-assembly.py b/smr/build-assembly.py index 2970847..e399b45 100644 --- a/smr/build-assembly.py +++ b/smr/build-assembly.py @@ -9,7 +9,7 @@ from tqdm import tqdm import openmc from smr.materials import materials -from smr.surfaces import surfs, lattice_pitch, bottom_fuel_stack, top_active_core +from smr.surfaces import surfs, lattice_pitch, bottom_fuel_stack, top_active_core, pellet_OR from smr.assemblies import assembly_universes from smr.plots import assembly_plots from smr import inlet_temperature @@ -41,7 +41,8 @@ else: directory.mkdir(exist_ok=True) # Define geometry with a single assembly -assembly = assembly_universes(args.rings, args.axial, args.depleted) +ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +assembly = assembly_universes(ring_radii, args.axial, args.depleted) lattice_sides = openmc.model.get_rectangular_prism(lattice_pitch, lattice_pitch, boundary_type='reflective') main_cell = openmc.Cell( diff --git a/smr/build-core-fresh.py b/smr/build-core-fresh.py index 6908b23..caba492 100644 --- a/smr/build-core-fresh.py +++ b/smr/build-core-fresh.py @@ -11,7 +11,7 @@ import numpy as np import openmc from smr.materials import materials from smr.plots import core_plots -from smr.surfaces import lattice_pitch, bottom_fuel_stack, top_active_core +from smr.surfaces import lattice_pitch, bottom_fuel_stack, top_active_core, pellet_OR from smr.core import core_geometry from smr import inlet_temperature @@ -48,7 +48,8 @@ else: directory = args.output_dir directory.mkdir(exist_ok=True) -geometry = core_geometry(args.rings, args.axial, args.depleted) +ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +geometry = core_geometry(ring_radii, args.axial, args.depleted) #### "Differentiate" the geometry if using distribmats if args.tallies == 'mat': diff --git a/smr/smr/assemblies.py b/smr/smr/assemblies.py index aeba677..2fffb51 100644 --- a/smr/smr/assemblies.py +++ b/smr/smr/assemblies.py @@ -120,13 +120,14 @@ def make_assembly(name, universes): return universe -def assembly_universes(num_rings, num_axial, depleted): +def assembly_universes(ring_radii, num_axial, depleted): """Generate universes for SMR fuel assemblies. Parameters ---------- - num_rings : int - Number of annual regions in fuel + ring_radii : iterable of float + Radii of rings in fuel (note that this doesn't need to include the + full fuel pin radius) num_axial : int Number of axial subdivisions in fuel depleted : bool @@ -138,7 +139,7 @@ def assembly_universes(num_rings, num_axial, depleted): Dictionary mapping a universe name to a openmc.Universe object """ - pins = pin_universes(num_rings, num_axial, depleted) + pins = pin_universes(ring_radii, num_axial, depleted) # Create dictionary to store assembly universes univs = {} diff --git a/smr/smr/core.py b/smr/smr/core.py index aef9444..c150e06 100644 --- a/smr/smr/core.py +++ b/smr/smr/core.py @@ -10,13 +10,14 @@ from .reflector import reflector_universes from .assemblies import assembly_universes -def core_geometry(num_rings, num_axial, depleted): +def core_geometry(ring_radii, num_axial, depleted): """Generate full core SMR geometry. Parameters ---------- - num_rings : int - Number of annual regions in fuel + ring_radii : iterable of float + Radii of rings in fuel (note that this doesn't need to include the + full fuel pin radius) num_axial : int Number of axial subdivisions in fuel depleted : bool @@ -28,7 +29,7 @@ def core_geometry(num_rings, num_axial, depleted): SMR full core geometry """ - assembly = assembly_universes(num_rings, num_axial, depleted) + assembly = assembly_universes(ring_radii, num_axial, depleted) reflector = reflector_universes() # Construct main core lattice diff --git a/smr/smr/pins.py b/smr/smr/pins.py index e71adf5..37de78d 100644 --- a/smr/smr/pins.py +++ b/smr/smr/pins.py @@ -142,13 +142,14 @@ def make_pin_stack(name, zsurfaces, universes, boundary, fuel_fill): return universe -def pin_universes(num_rings=10, num_axial=196, depleted=False): +def pin_universes(ring_radii=None, num_axial=196, depleted=False): """Generate universes for SMR fuel pins. Parameters ---------- - num_rings : int - Number of annual regions in fuel + ring_radii : iterable of float + Radii of rings in fuel (note that this doesn't need to include the + full fuel pin radius) num_axial : int Number of axial subdivisions in fuel depleted : bool @@ -706,12 +707,11 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): axial_splits = np.linspace(bottom_fuel_stack, top_active_core, num_axial + 1)[1:-1] axial_surfs = [openmc.ZPlane(z0=z) for z in axial_splits] - if num_rings > 1: + if ring_radii is not None: # Get z-cylinder surfaces for each ring rings = [] - for i in range(1, num_rings): - R = sqrt(i*pellet_OR**2/num_rings) - cyl = openmc.ZCylinder(r=R, name='fuel ring {}'.format(i)) + for i, r in enumerate(ring_radii): + cyl = openmc.ZCylinder(r=r, name='fuel ring {}'.format(i)) rings.append(cyl) def subdivided_fuel(fill): @@ -719,14 +719,14 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): uo2_cells = [] if num_axial > 1: for axial_region in subdivide(axial_surfs): - if num_rings > 1: + if ring_radii is not None: for ring_region in subdivide(rings): cell = openmc.Cell(fill=fill, region=axial_region & ring_region) uo2_cells.append(cell) else: uo2_cells.append(openmc.Cell(fill=fill, region=axial_region)) else: - if num_rings > 1: + if ring_radii is not None: for ring_region in subdivide(rings): cell = openmc.Cell(fill=fill, region=ring_region) uo2_cells.append(cell) @@ -737,7 +737,7 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): # If rings/axial segments are present, create a universe for the subdivided # fuel. Otherwise just use a plain material. - if num_rings > 1 or num_axial > 1: + if ring_radii is not None or num_axial > 1: fuel_fill = subdivided_fuel(mats['UO2 1.6 {}'.format(fuel)]) else: fuel_fill = mats['UO2 1.6 {}'.format(fuel)] @@ -829,7 +829,7 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): # If rings/axial segments are present, create a universe for the subdivided # fuel. Otherwise just use a plain material. - if num_rings > 1 or num_axial > 1: + if ring_radii is not None or num_axial > 1: fuel_fill = subdivided_fuel(mats['UO2 2.4 {}'.format(fuel)]) else: fuel_fill = mats['UO2 2.4 {}'.format(fuel)] @@ -875,11 +875,20 @@ def pin_universes(num_rings=10, num_axial=196, depleted=False): # If rings/axial segments are present, create a universe for the subdivided # fuel. Otherwise just use a plain material. - if num_rings > 1 or num_axial > 1: + if ring_radii is not None or num_axial > 1: fuel_fill = subdivided_fuel(mats['UO2 3.1 {}'.format(fuel)]) else: fuel_fill = mats['UO2 3.1 {}'.format(fuel)] + if num_axial > 1: + water_cells = [] + for i, r in enumerate(subdivide(axial_surfs)): + cell = openmc.Cell(fill=mats['H2O'], region=r, name=f'Water ({i})') + water_cells.append(cell) + water_fill = openmc.Universe(cells=water_cells) + else: + water_fill = mats['H2O'] + # Stack all axial pieces of 3.1% enriched fuel pin cell univs['Fuel pin (3.1%) stack'] = make_pin_stack( From 0e409cba870c2d97c035b1b78385578cdbc3931f Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 07:40:58 -0600 Subject: [PATCH 5/8] Add script for short single assembly ExaSMR coupled runs --- smr/build-assembly-short.py | 198 ++++++++++++++++++++++++++++++++++++ smr/smr/pins.py | 18 ++++ smr/smr/surfaces.py | 6 +- 3 files changed, 220 insertions(+), 2 deletions(-) create mode 100644 smr/build-assembly-short.py diff --git a/smr/build-assembly-short.py b/smr/build-assembly-short.py new file mode 100644 index 0000000..64999fc --- /dev/null +++ b/smr/build-assembly-short.py @@ -0,0 +1,198 @@ +#!/usr/bin/env python3 + +import argparse +import copy +from math import pi, isclose +from pathlib import Path + +import numpy as np +from tqdm import tqdm +import openmc + +from smr.materials import materials, mats +from smr.surfaces import surfs, lattice_pitch, pin_pitch, bottom_fuel_stack, \ + top_active_core, pellet_OR, clad_OR, clad_IR, guide_tube_IR, guide_tube_OR +from smr.pins import pin_universes + + +# Define command-line options +parser = argparse.ArgumentParser() +parser.add_argument('-m', '--multipole', action='store_true', + help='Whether to use multipole cross sections') +parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='mat', + help='Whether to use distribmats or distribcells for tallies') +parser.add_argument('-a', '--axial', type=int, default=10, + help='Number of axial subdivisions in fuel') +parser.add_argument('-d', '--depleted', action='store_true', + help='Whether UO2 compositions should represent depleted fuel') +parser.add_argument('-o', '--output-dir', type=Path, default=None) +args = parser.parse_args() + +# Make directory for inputs +if args.output_dir is None: + if args.depleted: + directory = Path('assembly-short-depleted') + else: + directory = Path('assembly-short-fresh') +else: + directory = args.output_dir +directory.mkdir(exist_ok=True) + +rings = [0.1*pin_pitch, 0.2*pin_pitch] + +# Define the NumPy array indices for assembly locations where there +# may be CR guide tubes, instrument tubes and burnable absorbers +nonfuel_y = np.array([2,2,2,3,3,5,5,5,5,5,8,8,8,8,8,11,11,11,11,11,13,13,14,14,14]) +nonfuel_x = np.array([5,8,11,3,13,2,5,8,11,14,2,5,8,11,14,2,5,8,11,14,3,13,5,8,11]) + +# NO BURNABLE ABSORBERS +pins = pin_universes(rings, args.axial, args.depleted) +gtu = pins['GT empty'] +#gti = pins['GT empty instr'] +universes = np.empty((17,17), dtype=openmc.Universe) +universes[:,:] = pins['Fuel pin (3.1%) no grid'] +universes[nonfuel_y, nonfuel_x] = [ gtu, gtu, gtu, + gtu, gtu, + gtu, gtu, gtu, gtu, gtu, + gtu, gtu, gtu, gtu, gtu, + gtu, gtu, gtu, gtu, gtu, + gtu, gtu, + gtu, gtu, gtu ] + +# Instantiate the lattice +lattice = openmc.RectLattice(name='Pin lattice') +lattice.lower_left = (-17.*pin_pitch/2., -17.*pin_pitch/2.) +lattice.pitch = (pin_pitch, pin_pitch) +lattice.universes = universes + +# Add lattice to bounding cell +root_universe = openmc.Universe(name='Root universe') +cell = openmc.Cell(name='Lattice cell') +cell.fill = lattice +z_bounds = +surfs['bot active core'] & -surfs['top active core'] +cell.region = surfs['lat grid box inner'] & z_bounds +root_universe.add_cell(cell) + +# Apply reflective boundaries +surfs['bot active core'].boundary_type = 'reflective' +surfs['top active core'].boundary_type = 'reflective' +for halfspace in surfs['lat grid box inner']: + halfspace.surface.boundary_type = 'reflective' + +# Define geometry with a single assembly +geometry = openmc.Geometry(root_universe) + + +def clone(material): + """Perform copy of material but share nuclide densities""" + shared_mat = copy.copy(material) + shared_mat.id = None + return shared_mat + + +#### "Differentiate" the geometry if using distribmats +h = 10.0*pin_pitch / args.axial +if args.tallies == 'mat': + # Count the number of instances for each cell and material + geometry.determine_paths(instances_only=True) + + for cell in tqdm(geometry.get_all_material_cells().values(), + desc='Differentiating materials'): + if cell.fill in materials: + # Fill cell with list of "differentiated" materials + cell.fill = [clone(cell.fill) for i in range(cell.num_instances)] + + # Determine volume of each fuel material + if 'UO2 Fuel' in cell.fill[0].name: + lower_left, _ = cell.region.bounding_box + if isclose(lower_left[0], rings[0]): + ri, ro = 0.0, rings[0] + elif isclose(lower_left[0], rings[1]): + ri, ro = rings[0], rings[1] + else: + ri, ro = rings[1], pellet_OR + for mat in cell.fill: + mat.volume = pi * (ro*ro - ri*ri) * h + elif cell.fill[0].name == 'Borated Water': + for mat in cell.fill: + mat.volume = pin_pitch**2 - pi*clad_OR**2 * h + elif cell.fill[0].name == 'Helium': + for mat in cell.fill: + mat.volume = pi * (clad_IR**2 - pellet_OR**2) * h + elif cell.fill[0].name == 'M5': + for mat in cell.fill: + mat.volume = pi * (clad_OR**2 - clad_IR**2) * h + elif cell.fill[0].name == 'Zircaloy-4': + for mat in cell.fill: + mat.volume = pi * (guide_tube_OR**2 - guide_tube_IR**2) * h + +#### Create OpenMC "materials.xml" file +print('Getting materials...') +all_materials = geometry.get_all_materials() +print('Creating materials collection...') +materials = openmc.Materials(all_materials.values()) +print('Exporting materials to XML...') +materials.export_to_xml(str(directory / 'materials.xml')) + + +#### Create OpenMC "geometry.xml" file +geometry.export_to_xml(str(directory / 'geometry.xml')) + + +#### Create OpenMC "settings.xml" file + +# Construct uniform initial source distribution over fissionable zones +lower_left = (-lattice_pitch/2, -lattice_pitch/2, bottom_fuel_stack) +upper_right = (lattice_pitch/2, lattice_pitch/2, top_active_core) +source = openmc.source.Source(space=openmc.stats.Box(lower_left, upper_right)) +source.space.only_fissionable = True + +settings = openmc.Settings() +settings.batches = 200 +settings.inactive = 100 +settings.particles = 10000 +settings.output = {'tallies': False, 'summary': False} +settings.source = source +settings.sourcepoint = {'write': False} + +if args.multipole: + settings.temperature = { + 'multipole': True, + 'tolerance': 1000, + 'default': 531.5, + 'method': 'interpolation', + 'range': (500.0, 1300.0) + } + +settings.export_to_xml(str(directory / 'settings.xml')) + + +#### Create OpenMC "tallies.xml" file +tallies = openmc.Tallies() + +# Extract all fuel materials +materials = geometry.get_materials_by_name(name='Fuel', matching=False) + +# If using distribcells, create distribcell tally needed for depletion +if args.tallies == 'cell': + # Extract all cells filled by a fuel material + fuel_cells = [] + for cell in geometry.get_all_cells().values(): + if cell.fill in materials: + tally = openmc.Tally(name='depletion tally') + tally.scores = ['(n,p)', '(n,a)', '(n,gamma)', + 'fission', '(n,2n)', '(n,3n)', '(n,4n)'] + tally.nuclides = cell.fill.get_nuclides() + tally.filters.append(openmc.DistribcellFilter([cell])) + tallies.append(tally) + +# If using distribmats, create material tally needed for depletion +elif args.tallies == 'mat': + tally = openmc.Tally(name='depletion tally') + tally.scores = ['(n,p)', '(n,a)', '(n,gamma)', + 'fission', '(n,2n)', '(n,3n)', '(n,4n)'] + tally.nuclides = materials[0].get_nuclides() + tally.filters = [openmc.MaterialFilter(materials)] + tallies.append(tally) + +tallies.export_to_xml(str(directory / 'tallies.xml')) diff --git a/smr/smr/pins.py b/smr/smr/pins.py index 37de78d..d804951 100644 --- a/smr/smr/pins.py +++ b/smr/smr/pins.py @@ -762,6 +762,12 @@ def pin_universes(ring_radii=None, num_axial=196, depleted=False): materials=outside_pin_mats, grid='intermediate') + univs['Fuel pin (1.6%) no grid'] = make_pin( + 'Pin no grid', + surfaces=[surfs['pellet OR']] + outside_pin_surfaces, + materials=[fuel_fill] + outside_pin_mats + ) + # Stack all axial pieces of 1.6% enriched fuel pin cell within_fuel_surfs = [ @@ -834,6 +840,12 @@ def pin_universes(ring_radii=None, num_axial=196, depleted=False): else: fuel_fill = mats['UO2 2.4 {}'.format(fuel)] + univs['Fuel pin (2.4%) no grid'] = make_pin( + 'Pin no grid', + surfaces=[surfs['pellet OR']] + outside_pin_surfaces, + materials=[fuel_fill] + outside_pin_mats + ) + # Stack all axial pieces of 2.4% enriched fuel pin cell univs['Fuel pin (2.4%) stack'] = make_pin_stack( @@ -889,6 +901,12 @@ def pin_universes(ring_radii=None, num_axial=196, depleted=False): else: water_fill = mats['H2O'] + univs['Fuel pin (3.1%) no grid'] = make_pin( + 'Pin no grid', + surfaces=[surfs['pellet OR']] + outside_pin_surfaces, + materials=[fuel_fill, mats['He'], mats['M5'], water_fill] + ) + # Stack all axial pieces of 3.1% enriched fuel pin cell univs['Fuel pin (3.1%) stack'] = make_pin_stack( diff --git a/smr/smr/surfaces.py b/smr/smr/surfaces.py index d97133f..96d9678 100644 --- a/smr/smr/surfaces.py +++ b/smr/smr/surfaces.py @@ -37,7 +37,7 @@ pellet_OR = 0.3195*INCHES/2 # ML17013A274, Table 4.1-2 pellet_length = 0.4*INCHES # ML17013A274, Table 4.1-2 clad_IR = 0.326*INCHES/2 # ML17013A274, Table 4.1-2 clad_OR = 0.374*INCHES/2 # ML17013A274, Table 4.1-2 -active_fuel_length = 78.74*INCHES # ML17013A274, Figure 4.2-10 +#active_fuel_length = 78.74*INCHES # ML17013A274, Figure 4.2-10 plenum_length = 5.311*INCHES # ML17013A274, Figure 4.2-10 fuel_rod_length = 85.00*INCHES # ML17013A274, Table 4.1-2 lower_end_cap_length = 0.575*INCHES # ML17007A001, Table 3-2 @@ -75,6 +75,8 @@ grid_strap_side = 21.47270 top_nozzle_height = 3.551*INCHES # ML17013A274, Figure 4.2-2 top_nozzle_width = 8.406*INCHES # ML17013A274, Figure 4.2-2 +active_fuel_length = 10.0*pin_pitch + # core radial parameters core_barrel_IR = 74*INCHES/2 # ML17013A274, Table 4.1-2 core_barrel_OR = 78*INCHES/2 # ML17013A274, Table 4.1-2 @@ -83,7 +85,7 @@ rpv_IR = 96.5*INCHES/2 # ML17013A274, Table 5.3-1 rpv_OR = 105*INCHES/2 # ML17013A274, Table 5.3-1 # axial parameters -reference_z = 0.0 +reference_z = -36.007 lowest_extent = reference_z bottom_support_plate = lowest_extent + 20.000 top_support_plate = bottom_support_plate + 5.000 From 2a24952857e29e6f218fb2c048b94937c86a0614 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Wed, 24 Jul 2019 22:55:08 -0500 Subject: [PATCH 6/8] Add script to build long assembly --- smr/build-assembly-long.py | 203 +++++++++++++++++++++++++++++++++++++ smr/smr/surfaces.py | 6 +- 2 files changed, 205 insertions(+), 4 deletions(-) create mode 100644 smr/build-assembly-long.py diff --git a/smr/build-assembly-long.py b/smr/build-assembly-long.py new file mode 100644 index 0000000..75373b9 --- /dev/null +++ b/smr/build-assembly-long.py @@ -0,0 +1,203 @@ +#!/usr/bin/env python3 + +import argparse +import copy +from math import pi, isclose +from pathlib import Path + +import numpy as np +from tqdm import tqdm +import openmc + +from smr.materials import materials, mats +from smr.surfaces import surfs, lattice_pitch, pin_pitch, bottom_fuel_stack, \ + top_active_core, pellet_OR, clad_OR, clad_IR, guide_tube_IR, guide_tube_OR, \ + active_fuel_length +from smr.pins import pin_universes, make_stack + + +# Define command-line options +parser = argparse.ArgumentParser() +parser.add_argument('-m', '--multipole', action='store_true', + help='Whether to use multipole cross sections') +parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='mat', + help='Whether to use distribmats or distribcells for tallies') +parser.add_argument('-a', '--axial', type=int, default=92, + help='Number of axial subdivisions in fuel') +parser.add_argument('-d', '--depleted', action='store_true', + help='Whether UO2 compositions should represent depleted fuel') +parser.add_argument('-o', '--output-dir', type=Path, default=None) +args = parser.parse_args() + +# Make directory for inputs +if args.output_dir is None: + if args.depleted: + directory = Path('assembly-long-depleted') + else: + directory = Path('assembly-long-fresh') +else: + directory = args.output_dir +directory.mkdir(exist_ok=True) + +rings = [0.1*pin_pitch, 0.2*pin_pitch] + +assembly_long_surfs = [ + surfs['bottom FR'], + surfs['bot active core'], + surfs['top active core'], + surfs['top pin plenum'], + surfs['top FR'], + surfs['bot upper nozzle'], + surfs['top upper nozzle'] +] + +univs = pin_universes(rings, args.axial, args.depleted) +fuel_univ = make_stack( + 'Fuel (3.1%) stack no grid', + surfaces=assembly_long_surfs, + universes=[ + univs['water pin'], + univs['end plug'], + univs['Fuel pin (3.1%) no grid'], + univs['pin plenum'], + univs['end plug'], + univs['water pin'] + ] +) + +# Define the NumPy array indices for assembly locations where there +# may be CR guide tubes, instrument tubes and burnable absorbers +nonfuel_y = np.array([2,2,2,3,3,5,5,5,5,5,8,8,8,8,8,11,11,11,11,11,13,13,14,14,14]) +nonfuel_x = np.array([5,8,11,3,13,2,5,8,11,14,2,5,8,11,14,2,5,8,11,14,3,13,5,8,11]) + +universes = np.full((17,17), fuel_univ) +universes[nonfuel_y, nonfuel_x] = univs['GT empty'] + +# Instantiate the lattice +lattice = openmc.RectLattice(name='Pin lattice') +lattice.lower_left = (-17.*pin_pitch/2., -17.*pin_pitch/2.) +lattice.pitch = (pin_pitch, pin_pitch) +lattice.universes = universes + +# Add lattice to bounding cell +root_universe = openmc.Universe(name='Root universe') +cell = openmc.Cell(name='Lattice cell') +cell.fill = lattice +z_bounds = +surfs['bottom FR'] & -surfs['top FR'] +cell.region = surfs['lat grid box inner'] & z_bounds +root_universe.add_cell(cell) + +# Apply reflective boundaries on sides and vacuum on bottom/top +surfs['bottom FR'].boundary_type = 'vacuum' +surfs['top FR'].boundary_type = 'vacuum' +for halfspace in surfs['lat grid box inner']: + halfspace.surface.boundary_type = 'reflective' + +# Define geometry with a single assembly +geometry = openmc.Geometry(root_universe) + + +def clone(material): + """Perform copy of material but share nuclide densities""" + shared_mat = copy.copy(material) + shared_mat.id = None + return shared_mat + + +#### "Differentiate" the geometry if using distribmats +h = active_fuel_length / args.axial +if args.tallies == 'mat': + # Count the number of instances for each cell and material + geometry.determine_paths(instances_only=True) + + for cell in tqdm(geometry.get_all_material_cells().values(), + desc='Differentiating materials'): + if cell.fill in materials: + # Fill cell with list of "differentiated" materials + cell.fill = [clone(cell.fill) for i in range(cell.num_instances)] + + # Determine volume of each fuel material + if 'UO2 Fuel' in cell.fill[0].name: + upper_right = cell.region.bounding_box[1] + if isclose(upper_right[0], rings[0]): + ri, ro = 0.0, rings[0] + elif isclose(upper_right[0], rings[1]): + ri, ro = rings[0], rings[1] + else: + ri, ro = rings[1], pellet_OR + for mat in cell.fill: + mat.volume = pi * (ro*ro - ri*ri) * h + else: + for mat in cell.fill: + mat.volume = 1.0 + +#### Create OpenMC "materials.xml" file +print('Getting materials...') +all_materials = geometry.get_all_materials() +print('Creating materials collection...') +materials = openmc.Materials(all_materials.values()) +print('Exporting materials to XML...') +materials.export_to_xml(str(directory / 'materials.xml')) + + +#### Create OpenMC "geometry.xml" file +geometry.export_to_xml(str(directory / 'geometry.xml')) + + +#### Create OpenMC "settings.xml" file + +# Construct uniform initial source distribution over fissionable zones +lower_left = (-lattice_pitch/2, -lattice_pitch/2, bottom_fuel_stack) +upper_right = (lattice_pitch/2, lattice_pitch/2, top_active_core) +source = openmc.source.Source(space=openmc.stats.Box(lower_left, upper_right)) +source.space.only_fissionable = True + +settings = openmc.Settings() +settings.batches = 200 +settings.inactive = 100 +settings.particles = 10000 +settings.output = {'tallies': False, 'summary': False} +settings.source = source +settings.sourcepoint = {'write': False} + +if args.multipole: + settings.temperature = { + 'multipole': True, + 'tolerance': 1000, + 'default': 531.5, + 'method': 'interpolation', + 'range': (500.0, 1300.0) + } + +settings.export_to_xml(str(directory / 'settings.xml')) + + +#### Create OpenMC "tallies.xml" file +tallies = openmc.Tallies() + +# Extract all fuel materials +materials = geometry.get_materials_by_name(name='Fuel', matching=False) + +# If using distribcells, create distribcell tally needed for depletion +if args.tallies == 'cell': + # Extract all cells filled by a fuel material + fuel_cells = [] + for cell in geometry.get_all_cells().values(): + if cell.fill in materials: + tally = openmc.Tally(name='depletion tally') + tally.scores = ['(n,p)', '(n,a)', '(n,gamma)', + 'fission', '(n,2n)', '(n,3n)', '(n,4n)'] + tally.nuclides = cell.fill.get_nuclides() + tally.filters.append(openmc.DistribcellFilter([cell])) + tallies.append(tally) + +# If using distribmats, create material tally needed for depletion +elif args.tallies == 'mat': + tally = openmc.Tally(name='depletion tally') + tally.scores = ['(n,p)', '(n,a)', '(n,gamma)', + 'fission', '(n,2n)', '(n,3n)', '(n,4n)'] + tally.nuclides = materials[0].get_nuclides() + tally.filters = [openmc.MaterialFilter(materials)] + tallies.append(tally) + +tallies.export_to_xml(str(directory / 'tallies.xml')) diff --git a/smr/smr/surfaces.py b/smr/smr/surfaces.py index 96d9678..0adf0d8 100644 --- a/smr/smr/surfaces.py +++ b/smr/smr/surfaces.py @@ -37,7 +37,7 @@ pellet_OR = 0.3195*INCHES/2 # ML17013A274, Table 4.1-2 pellet_length = 0.4*INCHES # ML17013A274, Table 4.1-2 clad_IR = 0.326*INCHES/2 # ML17013A274, Table 4.1-2 clad_OR = 0.374*INCHES/2 # ML17013A274, Table 4.1-2 -#active_fuel_length = 78.74*INCHES # ML17013A274, Figure 4.2-10 +active_fuel_length = 78.74*INCHES # ML17013A274, Figure 4.2-10 plenum_length = 5.311*INCHES # ML17013A274, Figure 4.2-10 fuel_rod_length = 85.00*INCHES # ML17013A274, Table 4.1-2 lower_end_cap_length = 0.575*INCHES # ML17007A001, Table 3-2 @@ -75,8 +75,6 @@ grid_strap_side = 21.47270 top_nozzle_height = 3.551*INCHES # ML17013A274, Figure 4.2-2 top_nozzle_width = 8.406*INCHES # ML17013A274, Figure 4.2-2 -active_fuel_length = 10.0*pin_pitch - # core radial parameters core_barrel_IR = 74*INCHES/2 # ML17013A274, Table 4.1-2 core_barrel_OR = 78*INCHES/2 # ML17013A274, Table 4.1-2 @@ -85,7 +83,7 @@ rpv_IR = 96.5*INCHES/2 # ML17013A274, Table 5.3-1 rpv_OR = 105*INCHES/2 # ML17013A274, Table 5.3-1 # axial parameters -reference_z = -36.007 +reference_z = -36.6205 lowest_extent = reference_z bottom_support_plate = lowest_extent + 20.000 top_support_plate = bottom_support_plate + 5.000 From 6448eaff093f95e75cf9c27c9a8d775264b89447 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 11:33:13 -0600 Subject: [PATCH 7/8] Use multipole cross sections by default --- smr/build-assembly-long.py | 7 +++++-- smr/build-assembly-short.py | 7 +++++-- smr/build-assembly.py | 7 +++++-- smr/build-core-fresh.py | 7 +++++-- 4 files changed, 20 insertions(+), 8 deletions(-) diff --git a/smr/build-assembly-long.py b/smr/build-assembly-long.py index 75373b9..cca7640 100644 --- a/smr/build-assembly-long.py +++ b/smr/build-assembly-long.py @@ -18,8 +18,10 @@ from smr.pins import pin_universes, make_stack # Define command-line options parser = argparse.ArgumentParser() -parser.add_argument('-m', '--multipole', action='store_true', - help='Whether to use multipole cross sections') +parser.add_argument('--multipole', action='store_true', + help='Use multipole cross sections') +parser.add_argument('--no-multipole', action='store_false', + help='Do not use multipole cross sections') parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='mat', help='Whether to use distribmats or distribcells for tallies') parser.add_argument('-a', '--axial', type=int, default=92, @@ -27,6 +29,7 @@ parser.add_argument('-a', '--axial', type=int, default=92, parser.add_argument('-d', '--depleted', action='store_true', help='Whether UO2 compositions should represent depleted fuel') parser.add_argument('-o', '--output-dir', type=Path, default=None) +parser.set_defaults(multipole=True) args = parser.parse_args() # Make directory for inputs diff --git a/smr/build-assembly-short.py b/smr/build-assembly-short.py index 64999fc..c06c67b 100644 --- a/smr/build-assembly-short.py +++ b/smr/build-assembly-short.py @@ -17,8 +17,10 @@ from smr.pins import pin_universes # Define command-line options parser = argparse.ArgumentParser() -parser.add_argument('-m', '--multipole', action='store_true', - help='Whether to use multipole cross sections') +parser.add_argument('--multipole', action='store_true', + help='Use multipole cross sections') +parser.add_argument('--no-multipole', action='store_false', + help='Do not use multipole cross sections') parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='mat', help='Whether to use distribmats or distribcells for tallies') parser.add_argument('-a', '--axial', type=int, default=10, @@ -26,6 +28,7 @@ parser.add_argument('-a', '--axial', type=int, default=10, parser.add_argument('-d', '--depleted', action='store_true', help='Whether UO2 compositions should represent depleted fuel') parser.add_argument('-o', '--output-dir', type=Path, default=None) +parser.set_defaults(multipole=True) args = parser.parse_args() # Make directory for inputs diff --git a/smr/build-assembly.py b/smr/build-assembly.py index e399b45..327f3bc 100644 --- a/smr/build-assembly.py +++ b/smr/build-assembly.py @@ -17,8 +17,10 @@ from smr import inlet_temperature # Define command-line options parser = argparse.ArgumentParser() -parser.add_argument('-m', '--multipole', action='store_true', - help='Whether to use multipole cross sections') +parser.add_argument('--multipole', action='store_true', + help='Use multipole cross sections') +parser.add_argument('--no-multipole', action='store_false', + help='Do not use multipole cross sections') parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='mat', help='Whether to use distribmats or distribcells for tallies') parser.add_argument('-r', '--rings', type=int, default=10, @@ -28,6 +30,7 @@ parser.add_argument('-a', '--axial', type=int, default=196, parser.add_argument('-d', '--depleted', action='store_true', help='Whether UO2 compositions should represent depleted fuel') parser.add_argument('-o', '--output-dir', type=Path, default=None) +parser.set_defaults(multipole=True) args = parser.parse_args() # Make directory for inputs diff --git a/smr/build-core-fresh.py b/smr/build-core-fresh.py index caba492..1272f3b 100644 --- a/smr/build-core-fresh.py +++ b/smr/build-core-fresh.py @@ -25,8 +25,10 @@ def clone(mat): # Define command-line options parser = argparse.ArgumentParser() -parser.add_argument('-m', '--multipole', action='store_true', - help='Whether to use multipole cross sections') +parser.add_argument('--multipole', action='store_true', + help='Use multipole cross sections') +parser.add_argument('--no-multipole', action='store_false', + help='Do not use multipole cross sections') parser.add_argument('-t', '--tallies', choices=('cell', 'mat'), default='cell', help='Whether to use distribmats or distribcells for tallies') parser.add_argument('-r', '--rings', type=int, default=10, @@ -36,6 +38,7 @@ parser.add_argument('-a', '--axial', type=int, default=196, parser.add_argument('-d', '--depleted', action='store_true', help='Whether UO2 compositions should represent depleted fuel') parser.add_argument('-o', '--output-dir', type=Path, default=None) +parser.set_defaults(multipole=True) args = parser.parse_args() # Make directory for inputs From b5b6f4611acef8183d7b9ce904cffa0f6d69117d Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 19 Dec 2019 14:59:44 -0600 Subject: [PATCH 8/8] Make sure model build doesn't fail if rings == 1 --- smr/build-assembly.py | 5 ++++- smr/build-core-fresh.py | 5 ++++- 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/smr/build-assembly.py b/smr/build-assembly.py index 327f3bc..392b114 100644 --- a/smr/build-assembly.py +++ b/smr/build-assembly.py @@ -44,7 +44,10 @@ else: directory.mkdir(exist_ok=True) # Define geometry with a single assembly -ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +if args.rings > 1: + ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +else: + ring_radii = None assembly = assembly_universes(ring_radii, args.axial, args.depleted) lattice_sides = openmc.model.get_rectangular_prism(lattice_pitch, lattice_pitch, boundary_type='reflective') diff --git a/smr/build-core-fresh.py b/smr/build-core-fresh.py index 1272f3b..a76ecda 100644 --- a/smr/build-core-fresh.py +++ b/smr/build-core-fresh.py @@ -51,7 +51,10 @@ else: directory = args.output_dir directory.mkdir(exist_ok=True) -ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +if args.rings > 1: + ring_radii = np.sqrt(np.arange(1, args.rings)*pellet_OR**2 / args.rings) +else: + ring_radii = None geometry = core_geometry(ring_radii, args.axial, args.depleted) #### "Differentiate" the geometry if using distribmats