Photon Flux Reading Zero

User Support
1

I am brand new to openMC and am trying to solve a relatively simple photon shielding problem. The problem consists of a 2mCi Co-60 point source centered in a 100cm tall, 10cm radius iron tank that’s filled with water. The tank has 1 cm thick iron walls. I am looking to find the photon flux on the outside of the tank at the point closest to the source (along centerline height of the tank).

In order to do this I set up my geometry with a spherical source cell centered in the tank, the tank and water, and a “detector” cell located at the centerline height on the outside of the tank. I believe my geometry is correct based on the geometry plots below.

geometry_plot_xy
geometry_plot_xy800×800 1.11 KB

geometry_plot_xz
geometry_plot_xz800×800 902 Bytes

However I have been unsuccessful running this and getting a flux that isn’t zero. This makes me think something is incorrect in my tally or source parts of my code. I am attempting to tally with energy bins so I can manually convert those fluxes to exposure using tabulated response functions. My code is shown below:

import openmc

# Constants
mCi_to_Bq = 3.7e7  # 1 mCi = 3.7 × 10^7 Bq
activity_mCi = 2.0  # Source activity in mCi
source_strength = activity_mCi * mCi_to_Bq  # Particles per second

# Define materials
# Water inside the tank
water = openmc.Material(name="Water")
water.add_element("H", 2)
water.add_element("O", 1)
water.set_density("g/cm3", 1.0)

# Iron for the tank walls
iron = openmc.Material(name="Iron")
iron.add_element("Fe", 1)
iron.set_density("g/cm3", 7.87)

# Create a materials collection
materials = openmc.Materials([water, iron])
materials.export_to_xml()

# Define geometry
# Cylindrical water-filled tank
tank_radius = 10.0  # cm
tank_height = 100.0  # cm
wall_thickness = 1.0  # cm
detector_radius = 5 # cm
# Define surfaces
cylinder_inner = openmc.ZCylinder(r=tank_radius)
cylinder_outer = openmc.ZCylinder(r=tank_radius + wall_thickness)
bottom_surface = openmc.ZPlane(z0=-tank_height / 2, boundary_type= 'vacuum')
top_surface = openmc.ZPlane(z0=tank_height / 2, boundary_type= 'vacuum')
void_surface = openmc.Sphere(r=500, boundary_type= 'vacuum')
detector_surface = openmc.Sphere(x0=tank_radius + wall_thickness + detector_radius+0.1, y0=0, z0=0, r=detector_radius)
source_surface =openmc.Sphere(x0=0,y0=0, z0=0, r=.5)
# Define cells
source_region = -source_surface
water_region = +source_surface & -cylinder_inner & +bottom_surface & -top_surface 
iron_region = +cylinder_inner & -cylinder_outer & +bottom_surface & -top_surface
detector_region = -detector_surface

water_cell = openmc.Cell(fill=water, region=water_region, name="Water inside the tank")
iron_cell = openmc.Cell(fill=iron, region=iron_region, name="Iron tank wall")
detector_cell = openmc.Cell(fill=water, region=detector_region, name="detector")
source_cell = openmc.Cell(fill=None, region=source_region, name='source')

# Outer void region
outer_void_region = +cylinder_outer & +bottom_surface & -top_surface & -void_surface
outer_void_cell = openmc.Cell(region=outer_void_region, name="Outer void")

# Universe
universe = openmc.Universe(cells=[water_cell, iron_cell, outer_void_cell, detector_cell, source_cell])

# Geometry
geometry = openmc.Geometry(universe)
geometry.export_to_xml()

# Define Co-60 source
energies = [1.17e6, 1.33e6]  # Gamma-ray energies in eV
weights = [0.502, 0.498]  # Relative probabilities
source = openmc.Source(space=openmc.stats.Point((0,0,0)), angle=openmc.stats.Isotropic(),energy=openmc.stats.Discrete(energies, weights),particle='photon',strength=source_strength)

# Settings
settings = openmc.Settings()
settings.source = source
settings.batches = 15
settings.particles = 30000000
settings.inactive = 5
settings.run_mode = "fixed source"
settings.photon_transport = True
settings.export_to_xml()

# Define energy bins from the provided table (convert MeV to eV)
energy_bins = [
    0.0,        # Start from 0 eV
    0.01e6,     # 10 keV
    0.08e6,     # 80 keV
    0.1e6,      # 100 keV
    0.15e6,     # 150 keV
    0.2e6,      # 200 keV
    0.3e6,      # 300 keV
    0.4e6,      # 400 keV
    0.5e6,      # 500 keV
    0.6e6,      # 600 keV
    0.8e6,      # 800 keV
    1.0e6,      # 1.0 MeV
    1.25e6,     # 1.25 MeV
    1.5e6       # 1.5 MeV
]
energy_filter = openmc.EnergyFilter(energy_bins)

# Add a cell filter for the spherical region
point_cell_filter = openmc.CellFilter([detector_cell])

# Add a particle filter for photons
photon_filter = openmc.ParticleFilter(["photon"])

# Define a tally for the flux at the spherical region
point_tally = openmc.Tally(name="point_flux")
point_tally.filters = [point_cell_filter, energy_filter, photon_filter]
point_tally.scores = ["flux"]

# Append tally to the tallies object
tallies = openmc.Tallies()
tallies.append(point_tally)
tallies.export_to_xml()

# Run the simulation
openmc.run()

Any assistance would be appreciated. I haven’t succesfully been alble to apply anything I’ve seen on other posts.

Ryan

2

I ended up being able to fix this by defining a (1,1,1) mesh cell that was a 1 cm^2 box centered around the point closest to the centerline of the tank. It seems that in order to do a flux tally you must use a mesh? I was able to do a current tally across a surface or geometry cell for a different shielding problem. But I only received zero as a flux when doing a flux tally across a geometry cell. Am I right in thinking that a flux tally requires a meshed volume as opposed to say a geometry cell?