OpenMC Guide

❮ Example: Pin CellExample: Full Reactor ❯

Example: Fuel Assembly

Building an Assembly

This example creates a 17×17 PWR fuel assembly using a hierarchical approach: define pin universes, arrange them in a lattice, then place the lattice in an assembly cell. This pattern scales to full core models.

Materials Definition

python
import openmc
import numpy as np

# =============================================================================
# MATERIALS
# =============================================================================

# UO2 fuel - 4.5% enriched
fuel = openmc.Material(name='UO2 Fuel')
fuel.set_density('g/cm3', 10.4)
fuel.add_nuclide('U235', 0.045)
fuel.add_nuclide('U238', 0.955)
fuel.add_element('O', 2.0)

# Zircaloy-4 cladding
clad = openmc.Material(name='Zircaloy-4')
clad.set_density('g/cm3', 6.55)
clad.add_element('Zr', 0.982)
clad.add_element('Sn', 0.015)
clad.add_element('Fe', 0.002)
clad.add_element('Cr', 0.001)

# Light water moderator
water = openmc.Material(name='Light Water')
water.set_density('g/cm3', 1.0)
water.add_nuclide('H1', 2.0)
water.add_element('O', 1.0)
water.add_s_alpha_beta('c_H_in_H2O')

# Stainless steel for guide tubes
steel = openmc.Material(name='Stainless Steel')
steel.set_density('g/cm3', 8.0)
steel.add_element('Fe', 0.68)
steel.add_element('Cr', 0.20)
steel.add_element('Ni', 0.12)

materials = openmc.Materials([fuel, clad, water, steel])
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Pin Universe Definitions

python
# =============================================================================
# PIN UNIVERSES
# =============================================================================

def create_fuel_pin():
    """Standard fuel pin with cladding"""
    fuel_radius = 0.4096
    clad_inner = 0.4178
    clad_outer = 0.4750
    
    # Surfaces
    fuel_surf = openmc.ZCylinder(r=fuel_radius)
    clad_inner_surf = openmc.ZCylinder(r=clad_inner)
    clad_outer_surf = openmc.ZCylinder(r=clad_outer)
    
    # Regions
    fuel_region = -fuel_surf
    gap_region = +fuel_surf & -clad_inner_surf
    clad_region = +clad_inner_surf & -clad_outer_surf
    water_region = +clad_outer_surf
    
    # Cells
    fuel_cell = openmc.Cell(fill=fuel, region=fuel_region)
    gap_cell = openmc.Cell(fill=None, region=gap_region)
    clad_cell = openmc.Cell(fill=clad, region=clad_region)
    water_cell = openmc.Cell(fill=water, region=water_region)
    
    return openmc.Universe(cells=[fuel_cell, gap_cell, clad_cell, water_cell])

def create_guide_tube():
    """Guide tube for control rods"""
    guide_inner = 0.5610
    guide_outer = 0.6020
    
    # Surfaces
    guide_inner_surf = openmc.ZCylinder(r=guide_inner)
    guide_outer_surf = openmc.ZCylinder(r=guide_outer)
    
    # Regions
    water_inner = -guide_inner_surf
    tube_region = +guide_inner_surf & -guide_outer_surf
    water_outer = +guide_outer_surf
    
    # Cells
    water_inner_cell = openmc.Cell(fill=water, region=water_inner)
    tube_cell = openmc.Cell(fill=steel, region=tube_region)
    water_outer_cell = openmc.Cell(fill=water, region=water_outer)
    
    return openmc.Universe(cells=[water_inner_cell, tube_cell, water_outer_cell])

def create_instrumentation_tube():
    """Instrumentation tube (water-filled)"""
    tube_inner = 0.5590
    tube_outer = 0.6050
    
    # Surfaces
    tube_inner_surf = openmc.ZCylinder(r=tube_inner)
    tube_outer_surf = openmc.ZCylinder(r=tube_outer)
    
    # Regions
    water_inner = -tube_inner_surf
    tube_region = +tube_inner_surf & -tube_outer_surf
    water_outer = +tube_outer_surf
    
    # Cells
    water_inner_cell = openmc.Cell(fill=water, region=water_inner)
    tube_cell = openmc.Cell(fill=steel, region=tube_region)
    water_outer_cell = openmc.Cell(fill=water, region=water_outer)
    
    return openmc.Universe(cells=[water_inner_cell, tube_cell, water_outer_cell])

# Create the pin universes
fuel_pin = create_fuel_pin()
guide_tube = create_guide_tube()
instr_tube = create_instrumentation_tube()
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Assembly Lattice

python
# =============================================================================
# ASSEMBLY LATTICE
# =============================================================================

# Create 17x17 lattice
lattice = openmc.RectLattice(name='Assembly Lattice')
lattice.lower_left = [-10.71, -10.71]  # Half of 17 * 1.26 cm pitch
lattice.pitch = [1.26, 1.26]

# Define assembly pattern
# Standard Westinghouse 17x17 assembly layout
# F = fuel pin, G = guide tube, I = instrumentation tube

assembly_pattern = [
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','F','F','G','F','F','G','F','F','G','F','F','F','F','F'],
    ['F','F','F','G','F','F','F','F','F','F','F','F','F','G','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','G','F','F','G','F','F','G','F','F','G','F','F','G','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','G','F','F','G','F','F','I','F','F','G','F','F','G','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','G','F','F','G','F','F','G','F','F','G','F','F','G','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','G','F','F','F','F','F','F','F','F','F','G','F','F','F'],
    ['F','F','F','F','F','G','F','F','G','F','F','G','F','F','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F'],
    ['F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F','F']
]

# Convert pattern to universe array
universe_map = {'F': fuel_pin, 'G': guide_tube, 'I': instr_tube}
lattice_universes = []
for row in assembly_pattern:
    universe_row = [universe_map[symbol] for symbol in row]
    lattice_universes.append(universe_row)

lattice.universes = lattice_universes
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Assembly Cell and Geometry

python
# =============================================================================
# ASSEMBLY GEOMETRY
# =============================================================================

# Assembly boundary (21.42 cm square)
pitch = 1.26
assembly_size = 17 * pitch

# Create assembly boundary
assembly_left = openmc.XPlane(-assembly_size/2, boundary_type='reflective')
assembly_right = openmc.XPlane(assembly_size/2, boundary_type='reflective')
assembly_bottom = openmc.YPlane(-assembly_size/2, boundary_type='reflective')
assembly_top = openmc.YPlane(assembly_size/2, boundary_type='reflective')
assembly_lower = openmc.ZPlane(-200, boundary_type='vacuum')
assembly_upper = openmc.ZPlane(200, boundary_type='vacuum')

# Assembly region
assembly_region = (+assembly_left & -assembly_right & 
                  +assembly_bottom & -assembly_top &
                  +assembly_lower & -assembly_upper)

# Assembly cell
assembly_cell = openmc.Cell(fill=lattice, region=assembly_region)

# Root universe and geometry
root_universe = openmc.Universe(cells=[assembly_cell])
geometry = openmc.Geometry(root_universe)
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Simulation Settings

python
# =============================================================================
# SIMULATION SETTINGS
# =============================================================================

settings = openmc.Settings()
settings.particles = 10000
settings.batches = 120
settings.inactive = 20

# Source distribution across the assembly
source_box = openmc.stats.Box(
    [-assembly_size/2, -assembly_size/2, -100],
    [assembly_size/2, assembly_size/2, 100]
)
settings.source = openmc.IndependentSource(space=source_box)

# Set temperature on cells (temperature is a Cell property, not Material)
# In the fuel pin function above, you would set:
#   fuel_cell.temperature = 900   # K
#   gap_cell.temperature = 600    # K
#   clad_cell.temperature = 600   # K
#   water_cell.temperature = 574  # K
# Here we apply it after geometry is built:
for cell in geometry.get_all_cells().values():
    if cell.fill == fuel:
        cell.temperature = 900   # K
    elif cell.fill == clad or cell.fill == steel:
        cell.temperature = 600   # K
    elif cell.fill == water:
        cell.temperature = 574   # K
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Tallies

python
# =============================================================================
# TALLIES
# =============================================================================

tallies = openmc.Tallies()

# Assembly-averaged flux and reaction rates
assembly_tally = openmc.Tally(name='assembly_avg')
assembly_tally.filters = [openmc.CellFilter(assembly_cell)]
assembly_tally.scores = ['flux', 'fission', 'absorption', 'nu-fission']
tallies.append(assembly_tally)

# Pin-by-pin power distribution
mesh = openmc.RegularMesh()
mesh.lower_left = [-assembly_size/2, -assembly_size/2]
mesh.upper_right = [assembly_size/2, assembly_size/2]
mesh.dimension = [17, 17]  # One bin per pin

mesh_tally = openmc.Tally(name='pin_powers')
mesh_tally.filters = [openmc.MeshFilter(mesh)]
mesh_tally.scores = ['kappa-fission']  # Energy released per fission
tallies.append(mesh_tally)

# Energy spectrum in fuel regions
energy_bins = np.logspace(-3, 7, 50)  # 1 meV to 10 MeV
energy_filter = openmc.EnergyFilter(energy_bins)

spectrum_tally = openmc.Tally(name='fuel_spectrum')
spectrum_tally.filters = [openmc.MaterialFilter(fuel), energy_filter]
spectrum_tally.scores = ['flux']
tallies.append(spectrum_tally)
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

Running the Model

python
# =============================================================================
# RUN SIMULATION
# =============================================================================

# Create complete model
model = openmc.Model(geometry, materials, settings, tallies)

# Run simulation (model.run() returns path to statepoint)
statepoint_path = model.run()
sp = openmc.StatePoint(statepoint_path)
keff = sp.keff
print(f"k-effective: {keff.nominal_value:.5f}")
print(f"Standard deviation: {keff.std_dev:.5f}")

# =============================================================================
# BASIC ANALYSIS
# =============================================================================

# Load results
sp = openmc.StatePoint('statepoint.120.h5')

# Get pin power distribution
pin_powers = sp.get_tally(name='pin_powers')
power_array = pin_powers.mean.reshape((17, 17))

# Simple power peaking analysis
max_power = np.max(power_array)
avg_power = np.mean(power_array)
peaking_factor = max_power / avg_power

print(f"Peak-to-average power ratio: {peaking_factor:.3f}")

# Find location of peak power
peak_location = np.unravel_index(np.argmax(power_array), power_array.shape)
print(f"Peak power at pin ({peak_location[0]+1}, {peak_location[1]+1})")
caalh/ReactorMCBrowse filesopenmc-examples/example_assembly/build_model.py

Section of the same example script in examples repo

The same hierarchical pattern extends to full core models, where assemblies become the building blocks.

Related Topics

❮ Example: Pin CellExample: Full Reactor ❯