I am currently conducting a geometric validation study to replicate the fuel log configuration of the TMSR-500 reactor, as proposed by ThorCon Indonesia. This structure consists of vertical fuel slabs embedded within a hexagonal prism filled with moderator, resembling the attached reference diagram. My goal is to reproduce this layout in OpenMC as accurately as possible to support a subsequent neutronics validation effort.
The attached script is my initial attempt to generate this geometry. It creates multiple vertical fuel slabs positioned within a hexagonal boundary and fills the remainder with moderator. However, the resulting geometry differs significantly from the reference—particularly in the number and arrangement of the fuel slabs, their spacing, and the interaction with the hexagonal moderator volume.
import openmc
# === PARAMETER GEOMETRI ===
apothem = 19.055 # cm
height_log = 17.0996 # cm
slot_thickness = 0.6 # cm
slab_thickness = 4.0 # cm
# === MATERIALS ===
#Material define in TMSR-500
#molten salts fuel NaF-BeF2-UF4-ThF4 with mole fraction (76% ; 12%; 9,5%; 2,5% )
T = 564 #inlet temperature as 564 celcius
fuel = openmc.Material(name='fuel_salt', temperature = T+273.2 )
fuel.add_nuclide('Na23', 0.306, 'ao') # Natrium-23
fuel.add_nuclide('Be9', 0.048, 'ao') # Berilium-9
fuel.add_nuclide('Th232', 0.0383, 'ao')# Thorium-232
fuel.add_nuclide('U233', 0.00101, 'ao')# Uranium-233 (10% from total uranium)
fuel.add_nuclide('U235', 0.00199, 'ao')# Uranium-235 (19.7% from total uranium)
fuel.add_nuclide('U238', 0.0071, 'ao') # Uranium-238 (70.3% from total uranium)
fuel.add_nuclide('F19', 0.597, 'ao') # Fluorine-19
#moderator core
moderator = openmc.Material(name='moderator', temperature = T+273.2)
moderator.add_nuclide('C12', 1.0)
moderator.set_density('g/cm3', 2.0)
moderator.add_s_alpha_beta('c_graphite')
#vessel core
vessel = openmc.Material(name='ss-16')
vessel.add_element('C',0.08,'wo')
vessel.add_element('N', 0.1, 'wo')
vessel.add_element('Si', 0.75,'wo')
vessel.add_element('P', 0.045,'wo')
vessel.add_element('S', 0.03,'wo')
vessel.add_element('Cr', 17,'wo')
vessel.add_element('Mn', 2, 'wo')
vessel.add_element('Fe', 65.497,'wo')
vessel.add_element('Ni', 12, 'wo')
vessel.add_element('Mo', 2.5,'wo')
vessel.set_density('g/cm3',8.00)
#reflector core
reflector = openmc.Material(name = 'reflector')
reflector.add_nuclide('C12', 1.0)
reflector.set_density('g/cm3', 2.0)
reflector.add_s_alpha_beta('c_graphite')
#shielding core
shielding = openmc.Material(name='shield')
shielding.add_nuclide('B10', 1.5920E-01)
shielding.add_nuclide('B11', 6480E-01)
shielding.add_nuclide('C12', 2.000E-01)
shielding.set_density('g/cm3', 2.57)
materials = openmc.Materials([fuel, moderator, vessel , reflector, shielding])
materials.export_to_xml()
# === HEXAGONAL FUEL LOG REGION ===
hex_prism = openmc.model.hexagonal_prism(edge_length=apothem, orientation='y', boundary_type='reflective')
bottom = openmc.ZPlane(z0=-height_log / 2, boundary_type='reflective')
top = openmc.ZPlane(z0=height_log / 2, boundary_type='reflective')
log_region = hex_prism & +bottom & -top
# === FUEL + SLOT SLAB REGION ===
slab_cells = []
slot_cells = []
#
num_slabs = int((2 * apothem) // (slab_thickness + slot_thickness))
x_start = -apothem + (slab_thickness / 2)
for i in range(num_slabs):
x_center = x_start + i * (slab_thickness + slot_thickness)
# boundary region of slab moderator
x_left = openmc.XPlane(x0=x_center - slab_thickness / 2)
x_right = openmc.XPlane(x0=x_center + slab_thickness / 2)
slab_region = +x_left & -x_right & log_region
slab_cells.append(openmc.Cell(name=f"Fuel Slab {i+1}", region=slab_region, fill=fuel))
if i < num_slabs - 1:
x_slot_left = openmc.XPlane(x0=x_center + slab_thickness / 2)
x_slot_right = openmc.XPlane(x0=x_center + slab_thickness / 2 + slot_thickness)
slot_region = +x_slot_left & -x_slot_right & log_region
slot_cells.append(openmc.Cell(name=f"Moderator Slot {i+1}", region=slot_region, fill=moderator))
occupied_region = openmc.Union([cell.region for cell in slab_cells + slot_cells])
remaining_moderator_region = log_region & ~occupied_region
remaining_moderator_cell = openmc.Cell(name="Remaining Moderator", region=remaining_moderator_region, fill=moderator)
# === formation FUEL LOG ===
fuel_log_universe = openmc.Universe(cells=slab_cells + slot_cells + [remaining_moderator_cell])
# === GEOMETRY ===
geometry = openmc.Geometry(fuel_log_universe)
# === EXPORT ===
materials.export_to_xml()
geometry.export_to_xml()
# === PLOTTING ===
geometry.plot(origin=(0, 0, 0), width=(2*apothem, 2*apothem), pixels=(600, 600), basis='xy')
geometry.plot(origin=(0, 0, 0), width=(2*apothem, height_log), pixels=(600, 600), basis='xz')
and the result, what i got is like this
