hello，everyone.
I’m modelling a gas-cooled microreactor with openmc and I’ve got problems running fuel assemble.
RuntimeError:No fission sites banked on MPI rank 0 Abort(-1) on node 0 (rank 0 in comm 0): application called MPI_Abort(MPI_COMM_WORLD, -1) - process 0
%matplotlib inline
from math import pi
import numpy as np
import matplotlib.pyplot as plt
import openmc
import openmc.model
import glob
import os
uo2 = openmc.Material(name = 'UO2')
uo2.add_nuclide('U235', 0.08599,'ao')
uo2.add_nuclide('U238', 0.91401,'ao')
uo2.add_nuclide('O16', 2.0,'ao')
uo2.set_density('g/cm3', 10.0)
print(uo2)
graphite = openmc.Material(name = 'Graphite')
graphite.add_element('C', 1.0)
graphite.set_density('g/cm3', 1.79)
print(graphite)
he = openmc.Material(name = 'helium')
he.add_element('He', 1.0)
he.set_density('g/cm3', 0.0001786)
print(he)
buffer = openmc.Material(name = 'Buffer')
buffer.add_element('C', 1.0)
buffer.set_density('g/cm3', 1.0)
buffer.add_s_alpha_beta('c_Graphite')
ipyc = openmc.Material(name = 'IPyC')
ipyc.add_element('C', 1.0)
ipyc.set_density('g/cm3', 1.9)
ipyc.add_s_alpha_beta('c_Graphite')
sic = openmc.Material(name = 'SiC')
sic.add_element('C', 0.3)
sic.add_element('Si', 0.7)
sic.set_density('g/cm3', 3.2)
opyc = openmc.Material(name = 'OPyC')
opyc.add_element('C', 1.0)
opyc.set_density('g/cm3', 1.9)
opyc.add_s_alpha_beta('c_Graphite')
spheres = [openmc.Sphere(r=1e-4*r)
for r in [250., 350., 390., 425.,465.]]
cells = [openmc.Cell(fill=uo2, region=-spheres[0]),
openmc.Cell(fill=buffer, region=+spheres[0] & -spheres[1]),
openmc.Cell(fill=ipyc, region=+spheres[1] & -spheres[2]),
openmc.Cell(fill=sic, region=+spheres[2] & -spheres[3]),
openmc.Cell(fill=opyc, region=+spheres[3]& -spheres[4])]
triso_univ = openmc.Universe(cells=cells)
materials = openmc.Materials([uo2,graphite,he,buffer,ipyc,sic,opyc])
materials.export_to_xml()
triso_univ.plot(width=(0.1,0.1))
plt.show()
fuel_rod_radius = 1.27
fuel_rod_height = 40
cylinder = openmc.ZCylinder(r=fuel_rod_radius,boundary_type='reflective')
z_min = openmc.ZPlane(z0=-fuel_rod_height/2,boundary_type='reflective')
z_max = openmc.ZPlane(z0=fuel_rod_height/2,boundary_type='reflective')
region = -cylinder & +z_min & -z_max
outer_radius = 465.*1e-4
centers = openmc.model.pack_spheres(radius=outer_radius, region=region, pf=0.3)
trisos = [openmc.model.TRISO(outer_radius, triso_univ, c) for c in centers]
print(trisos[0])
centers = np.vstack([t.center for t in trisos])
print(centers.min(axis=0))
print(centers.max(axis=0))
len(trisos)*4/3*pi*outer_radius**3
box = openmc.Cell(region=region)
lower_left, upper_right = box.region.bounding_box
shape = (3, 3, 3)
pitch = (upper_right - lower_left)/shape
lattice = openmc.model.create_triso_lattice(
trisos, lower_left, pitch, shape, sic)
box.fill = lattice
fuel_rod_univ = openmc.Universe(cells=[box])
fuel_rod_univ.plot(width = (3,3), colors={sic: 'yellow'}, color_by='material')
plt.show()
n_fuel_rods = 18
n_coolant_holes = 54
def generate_hex_positions(n_rings, fuel_pitch):
fuel_positions = []
for ring in range(n_rings):
radius = (ring + 1) * fuel_pitch
for i in range(6):
angle = 2*np.pi/6 * i
x = radius * np.cos(angle)
y = radius * np.sin(angle)
fuel_positions.append((x, y))
return fuel_positions
fuel_pitch = 4.6188
fuel_positions = generate_hex_positions(n_rings=2, fuel_pitch=fuel_pitch)
def generate_hex_positions2(fuel_positions):
coolant_positions = []
for (x, y) in fuel_positions:
for i in range(6):
angle =np.pi/6+ np.pi/3* (i)
dx = 2.6666 * np.cos(angle)
dy = 2.6666 * np.sin(angle)
coolant_positions.append((x + dx, y + dy))
return coolant_positions
coolant_positions = generate_hex_positions2(fuel_positions)
import numpy as np
from scipy.spatial import KDTree
tolerance = 1e-3
points = np.array(coolant_positions)
tree = KDTree(points)
neighbor_pairs = tree.query_pairs(tolerance)
keep_mask = np.ones(len(points), dtype=bool)
for i, j in neighbor_pairs:
if keep_mask[j]:
keep_mask[j] = False
unique_positions = points[keep_mask].tolist()
final_coolant_positions = [(round(x,4), round(y,4)) for x, y in unique_positions]
fuel_positions = fuel_positions +[(6.9281,3.6427),(0,7.999996269999128),(-6.9281,3.6427),(-6.9281,-3.6427),(0,-7.999996269999128),(6.9281,-3.6427)]
fuel_positions = [(round(x,4), round(y,4)) for x, y in fuel_positions]
final_coolant_positions = final_coolant_positions + [(9.2376,5.3333),(0,10.6666),(-9.2376,5.3333),(-9.2376,-5.3333),(0,-10.6666),(9.2376,-5.3333)]
final_coolant_positions = [(round(x,4), round(y,4)) for x, y in final_coolant_positions]
all_cells = []
cyl_air = openmc.ZCylinder(x0=0, y0=0, r=1.27,boundary_type='reflective')
region_air = -cyl_air & +z_min & -z_max
cell_air = openmc.Cell(
name=f'air hole',
fill=None,
region=region_air
)
all_cells.append(cell_air)
for i, (x, y) in enumerate(fuel_positions):
cyl_fuel = openmc.ZCylinder(x0=x, y0=y, r=1.27,boundary_type='reflective')
region_fuel = -cyl_fuel & +z_min & -z_max
cell_fuel = openmc.Cell(
name=f'Fuel Rod {i+1}',
fill=lattice,
region=region_fuel
)
all_cells.append(cell_fuel)
for i, (x, y) in enumerate(final_coolant_positions):
cyl_coolant = openmc.ZCylinder(x0=x, y0=y, r=0.4,boundary_type='reflective')
region_coolant = -cyl_coolant & +z_min & -z_max
cell_coolant = openmc.Cell(
name=f'Coolant Hole {i+1}',
fill=he,
region=region_coolant
)
all_cells.append(cell_coolant)
hex_prism = openmc.model.HexagonalPrism(edge_length=13.8564,orientation='x',origin=(0,0),boundary_type='reflective')
assembly_region =-hex_prism & +z_min & -z_max
from openmc import Union
fuel_regions = [c.region for c in all_cells if 'Fuel' in c.name]
coolant_regions = [c.region for c in all_cells if 'Coolant' in c.name]
air_region = [c.region for c in all_cells if 'air' in c.name]
all_holes = Union(fuel_regions + coolant_regions)
graphite_region1 = assembly_region & ~all_holes
graphite_region = graphite_region1 & ~region_air
graphite_cell = openmc.Cell(
name='Graphite',
fill=graphite,
region=graphite_region
)
all_cells.append(graphite_cell)
root_universe = openmc.Universe(cells=all_cells)
geometry = openmc.Geometry(root_universe)
materials.export_to_xml()
geometry.export_to_xml()
root_universe.plot(width = (30, 30),
color_by='material',
colors={he:'blue',
graphite:'gray',})
plt.show()
settings = openmc.Settings()
settings.run_mode = 'eigenvalue'
bounds = root_universe.bounding_box
lower_left = [3.0, -1.0, -1.0]
upper_right = [6.0, 1.0, 1.0]
src_space = openmc.stats.Box(lower_left, upper_right)
source = openmc.IndependentSource(space=src_space)
settings.source = source
settings.batches = 100
settings.inactive = 10
settings.particles = 1000
settings.export_to_xml()
fuel_cell = cells[0]
cell_filter = openmc.CellFilter(fuel_cell)
tally = openmc.Tally()
tally.filters = [cell_filter]
tally.nuclides = ['U235']
tally.scores = ['total', 'fission', 'absorption', '(n,gamma)']
tallies = openmc.Tallies([tally])
tallies.export_to_xml()
openmc.run()