OpenMC Guide
Building a Reactor Model in OpenMC
Overview
This example builds a simple reactor model in OpenMC, from materials and pin cells through assembly lattices to a complete core configuration.
The model follows OpenMC's multi-level universe nesting strategy. A fuel pin cell is defined as its own universe, then tiled into a 17×17 rectangular lattice to form an assembly. That assembly lattice is wrapped in an assembly universe, which is then placed into a core-level lattice. Finally, the core lattice is embedded in the root universe that defines the outermost problem boundary. This hierarchical approach means each component — pin, assembly, core — is defined exactly once and reused wherever it appears, so modifying a single pin design automatically propagates to every assembly and the full core.
Step 1: Materials
Define UO2 fuel and light water moderator with thermal scattering treatment:
# Basic reactor model setup
import openmc
import numpy as np
# Define materials
fuel = openmc.Material(name='UO2 fuel')
fuel.add_nuclide('U235', 0.04)
fuel.add_nuclide('U238', 0.96)
fuel.add_element('O', 2.0)
fuel.set_density('g/cm3', 10.4)
water = openmc.Material(name='Water')
water.add_nuclide('H1', 2)
water.add_element('O', 1)
water.set_density('g/cm3', 1.0)
water.add_s_alpha_beta('c_H_in_H2O')
materials = openmc.Materials([fuel, water])
materials.export_to_xml()Tutorial snippet — no separate file in examples repo
For more complex material definitions, including burnable absorbers or control materials, see the Materials page.
Step 2: Pin Cell
Create a fuel pin cell — the fundamental building block — consisting of a fuel cylinder surrounded by moderator:
# Define pin cell geometry
fuel_radius = 0.4
pin_pitch = 1.26
# Create cylinder and reflective square boundary
fuel_surf = openmc.ZCylinder(r=fuel_radius)
square_region = openmc.model.RectangularPrism(
pin_pitch,
pin_pitch,
boundary_type='reflective'
)
# Create cells
fuel_cell = openmc.Cell(fill=fuel, region=-fuel_surf)
water_cell = openmc.Cell(fill=water, region=+fuel_surf & square_region)
# Create universe
pin_universe = openmc.Universe(cells=[fuel_cell, water_cell])Tutorial snippet — no separate file in examples repo
For detailed pin cell modeling, including gap and cladding, see thePin Cell Example.
Step 3: Fuel Assembly
Create a 17×17 fuel assembly using a lattice of pin cells with uniform loading:
# Create 17x17 assembly
assembly_pitch = 1.26
lattice = openmc.RectLattice()
lattice.lower_left = (-10.71, -10.71)
lattice.pitch = (assembly_pitch, assembly_pitch)
lattice.universes = [[pin_universe] * 17] * 17
# Reflective boundary around assembly
half_span = 17 * assembly_pitch / 2
xmin = openmc.XPlane(-half_span, boundary_type='reflective')
xmax = openmc.XPlane(half_span, boundary_type='reflective')
ymin = openmc.YPlane(-half_span, boundary_type='reflective')
ymax = openmc.YPlane(half_span, boundary_type='reflective')
zmin = openmc.ZPlane(-200.0, boundary_type='vacuum')
zmax = openmc.ZPlane(200.0, boundary_type='vacuum')
assembly_region = +xmin & -xmax & +ymin & -ymax & +zmin & -zmax
# Create assembly cell and universe
assembly_cell = openmc.Cell(fill=lattice, region=assembly_region)
assembly_universe = openmc.Universe(cells=[assembly_cell])Tutorial snippet — no separate file in examples repo
For advanced assembly configurations, including guide tubes and varying enrichments, see theAssembly Example.
Step 4: Core Configuration
Arrange assemblies into a 3×3 core layout with vacuum boundaries:
# Create 3x3 core
core_lattice = openmc.RectLattice()
core_lattice.lower_left = (-32.13, -32.13)
core_lattice.pitch = (21.42, 21.42)
core_lattice.universes = [[assembly_universe] * 3] * 3
# Vacuum boundary around core
core_half_span = 3 * 21.42 / 2
core_xmin = openmc.XPlane(-core_half_span, boundary_type='vacuum')
core_xmax = openmc.XPlane(core_half_span, boundary_type='vacuum')
core_ymin = openmc.YPlane(-core_half_span, boundary_type='vacuum')
core_ymax = openmc.YPlane(core_half_span, boundary_type='vacuum')
zmin = openmc.ZPlane(-200.0, boundary_type='vacuum')
zmax = openmc.ZPlane(200.0, boundary_type='vacuum')
core_region = +core_xmin & -core_xmax & +core_ymin & -core_ymax & +zmin & -zmax
# Create core cell and root universe
core_cell = openmc.Cell(fill=core_lattice, region=core_region)
root_universe = openmc.Universe(cells=[core_cell])
# Create and export geometry
geometry = openmc.Geometry(root_universe)
geometry.export_to_xml()Tutorial snippet — no separate file in examples repo
Vacuum boundaries terminate any particle that crosses them, effectively modeling a bare core with no surrounding reflector. This is appropriate for a simplified demonstration, but a production model would typically add a water or baffle reflector region outside the core to capture neutrons that would otherwise be lost — omitting the reflector underestimates k-effective and distorts the peripheral power distribution.
Step 5: Simulation Settings
Configure particle count, batches, and initial source distribution:
# Define simulation settings
settings = openmc.Settings()
settings.batches = 100
settings.inactive = 20
settings.particles = 10000
# Set initial fission source near the core center
lower_left = (-32.13, -32.13, -50.0)
upper_right = (32.13, 32.13, 50.0)
source = openmc.IndependentSource(space=openmc.stats.Box(lower_left, upper_right))
settings.source = source
settings.export_to_xml()Tutorial snippet — no separate file in examples repo
Step 6: Tallies
Score flux and fission rate distributions across the core:
# Define tallies
mesh = openmc.RegularMesh()
mesh.dimension = [34, 34]
mesh.lower_left = [-32.13, -32.13]
mesh.upper_right = [32.13, 32.13]
tally = openmc.Tally(name='core_mesh_flux')
tally.filters = [openmc.MeshFilter(mesh)]
tally.scores = ['flux', 'fission']
tallies = openmc.Tallies([tally])
tallies.export_to_xml()Tutorial snippet — no separate file in examples repo
For more advanced tally options, see theTallies page.
Next Steps
For more advanced topics, explore:
Step 7: Run the Model
Combine all components into a model and run:
# Build full model and run
model = openmc.Model(geometry, materials, settings, tallies)
sp = openmc.StatePoint(model.run())
print(f"k-effective: {sp.keff.nominal_value:.5f} ± {sp.keff.std_dev:.5f}")Tutorial snippet — no separate file in examples repo