Hello. I'm attempting to model a system of 10 cylindrical segments stacked together via translational stacking. I've set up total reflection only on their top and bottom surfaces and placed them in the root universe. However, the code keeps giving errors

import openmc
UO2=openmc.Material(name=‘uo2’)
UO2.add_nuclide(‘U238’,6.91e-3,‘wo’)
UO2.add_nuclide(‘U235’, 0.22,‘wo’)
UO2.add_nuclide(‘O16’, 0.455,‘wo’)
UO2.set_density(‘g/cm3’,10.196)

clad_mat=openmc.Material(name=‘clad’)
clad_mat.add_nuclide(‘C0’,0.0009,‘wo’)
clad_mat.add_nuclide(‘P31’,0.0001,‘wo’)
clad_mat.add_nuclide(‘S32’,0.0002,‘wo’)
clad_mat.add_nuclide(‘Si28’,0.0048,‘wo’)
clad_mat.add_nuclide(‘Cr52’,0.0848,‘wo’)
clad_mat.add_nuclide(‘Mn55’,0.0053,‘wo’)
clad_mat.add_nuclide(‘Ni58’,0.0002,‘wo’)
clad_mat.add_nuclide(‘Fe56’,0.8911,‘wo’)
clad_mat.add_nuclide(‘Mo96’,0.0093,‘wo’)
clad_mat.add_nuclide(‘N14’,0.0004,‘wo’)
clad_mat.add_nuclide(‘Al27’,0.0001,‘wo’)
clad_mat.add_nuclide(‘Nb93’,0.0009,‘wo’)
clad_mat.set_density(‘g/cm3’,7.53)

H2O=openmc.Material(name=‘water’)
H2O.add_nuclide(‘H1’, 2.0)
H2O.add_nuclide(‘O16’, 1.0)
H2O.set_density(‘g/cm3’,1.0)

helium_mat=openmc.Material(name=‘Helium’)
helium_mat.add_nuclide(‘He4’,1.0,‘ao’)
helium_mat.set_density(‘g/cm3’, 0.0002)

mat_res=[UO2,clad_mat,helium_mat,H2O]
mat_xml_file=openmc.Materials(mat_res)
mat_xml_file.export_to_xml()
openmc.Materials.cross_sections=‘/home/zh/endfb71_hdf5/cross_sections.xml’

H=50
D=0.4
clad_thickness=0.2
gap=0.1
clad_D=D+(clad_thickness+gap)*2
out_D=3

axid=10

fuel_cyl=openmc.ZCylinder(r=D/2)
clad_in_cyl=openmc.ZCylinder(r=D/2+gap)
clad_out_cyl=openmc.ZCylinder(r=clad_D/2)
#out_cyl=openmc.ZCylinder(r=out_D, boundary_type=‘reflective’)
out_cyl=openmc.ZCylinder(r=out_D)

pz_list=
for i in range(axid+1):
if i==0 or i==axid:
res=openmc.ZPlane(z0=-H/2+iH/axid,boundary_type=‘reflective’)
else:
res=openmc.ZPlane(z0=-H/2+i
H/axid)
pz_list.append(res)

fuel_id=
fuel_universe=
rod_cells=
for i in range(axid):
clad_region=+clad_in_cyl &-clad_out_cyl & +pz_list[i] &-pz_list[i+1]
gap_region=-clad_in_cyl &+fuel_cyl & +pz_list[i] &-pz_list[i+1]
fuel_region=-fuel_cyl & +pz_list[i] &-pz_list[i+1]
water_region=-out_cyl &+clad_out_cyl & +pz_list[i] &-pz_list[i+1]

clad_cell=openmc.Cell(region=clad_region,fill=clad_mat)
gap_cell=openmc.Cell(region=gap_region,fill=helium_mat)
fuel_cell=openmc.Cell(region=fuel_region,fill=UO2)
fuel_id.append(fuel_cell.id)
water_cell=openmc.Cell(region=water_region,fill=H2O)

fuel_universe_res=openmc.Universe(cells=[fuel_cell,clad_cell,gap_cell,water_cell])
    


rod_cell=openmc.Cell(region=-out_cyl& +pz_list[i] &-pz_list[i+1])
rod_cell.fill=fuel_universe_res
rod_cell.translation=(0,0.001*i,0)

fuel_universe.append(fuel_universe_res)
rod_cells.append(rod_cell)

rod_unit=openmc.Universe(cells=rod_cells)

rod_unit.plot(basis=‘yz’,width=(10,60),pixels=(1000,1000),color_by=‘material’)

geomatry=openmc.Geometry()
geomatry.root_universe=rod_unit
geomatry.export_to_xml()

settings=openmc.Settings()
settings.batches=150
settings.inactive=50
settings.particles=1000

settings.export_to_xml()

fuel_cell_filter=openmc.CellFilter(fuel_id)
fuel_tally=openmc.Tally()
fuel_tally.filters.append(fuel_cell_filter)
fuel_tally.scores=[‘fission’,‘flux’,‘total’]

tallies=openmc.Tallies([fuel_tally])

tallies.export_to_xml()

openmc.run()

error:Maximum number of lost particles has been reached.

You just need to define an initial source distribution. May be something like this?

settings.source = openmc.IndependentSource(
space=openmc.stats.Point((0, 0, H / 2)),
angle=openmc.stats.Isotropic())

This define an isotropic point source. PS ( It doesn’t matter what initial source distribution you have as long as you are running sufficient inactive batches to achieve the true source distribution.)

Hey qvq,

I tried to run your code but I was getting too many errors in the construction of it (syntax and your -pz_list[i+1] was grabbing outside the index, but assuming it’s working properly on your end), just looking at it, it’s due to you not having a region outside of your geometry that has void boundary conditions. You didn’t apply boundary conditions to your out_cyl surface, therefore the particles are passing through it into no other cell and getting lost. Since you’re translating the cylinders in the y direction, there should be a bit of the z facing side of the cylinders not touching another z facing side (if i’m reading the geometry correct), which means you can’t just apply a boundary condition to the out_cyl, and you need to surround the entire geometry in a cell that has vaccuum conditions. I would do something thats just a slightly larger cylinder (taking in the fact that you are translating them) that ends at the same upper and lower z bounds, and with a boundary of vacuum. Also, I would always double check by plotting and seeing that it looks like what you are expecting, and also running in geometry debug mode with higher counting statistics (batchs/particles) and seeing what it tells you.