Neutron flux and units#
This example creates a simple sphere of water and tallies neutron flux averaged across a cell.
First import OpenMC and configure the nuclear data path
import openmc
from pathlib import Path
# Setting the cross section path to the correct location in the docker image.
# If you are running this outside the docker image you will have to change this path to your local cross section path.
openmc.config['cross_sections'] = Path.home() / 'nuclear_data' / 'cross_sections.xml'
This section creates a simple material, geometry and settings. This model is used in both the neutron current tally and the neutron flux tally.
# MATERIALS
# Due to the hydrogen content water is a very good neutron moderator
my_material = openmc.Material()
my_material.add_element('H', 1, percent_type='ao')
my_material.add_element('O', 2, percent_type='ao')
my_material.set_density('g/cm3', 1)
my_materials = openmc.Materials([my_material])
# GEOMETRY
# surfaces
outer_surface = openmc.Sphere(r=500, boundary_type='vacuum')
# cells
cell_1 = openmc.Cell(region=-outer_surface)
cell_1.fill = my_material
# we will need this volume later to convert the flux units
cell_1.volume = 5.24e8 # using (4/3)pi * r^3
my_geometry = openmc.Geometry([cell_1])
# SIMULATION SETTINGS
# Instantiate a Settings object
my_settings = openmc.Settings()
my_settings.batches = 10
my_settings.particles = 1000
my_settings.run_mode = 'fixed source'
# Create a DT point source
my_source = openmc.IndependentSource()
my_source.space = openmc.stats.Point((0, 0, 0))
my_source.angle = openmc.stats.Isotropic()
my_source.energy = openmc.stats.Discrete([14e6], [1])
my_settings.source = my_source
This section section adds a tally for the average neutron flux across a cell.
# sets up filters for the tallies
neutron_particle_filter = openmc.ParticleFilter(['neutron'])
# setup the filters for the cell tally
cell_filter = openmc.CellFilter(cell_1)
# create the tally
cell_spectra_tally = openmc.Tally(name='cell_flux_tally')
cell_spectra_tally.scores = ['flux']
cell_spectra_tally.filters = [cell_filter, neutron_particle_filter]
my_tallies = openmc.Tallies([cell_spectra_tally])
This section adds two surface current tallies - one on the inner sphere surface and one on the outer sphere surface.
This section runs the simulation.
# combine all the required parts to make a model
model = openmc.model.Model(my_geometry, my_materials, my_settings, my_tallies)
# remove old files and runs OpenMC
!rm *.h5 || true
results_filename = model.run()
%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%
############### %%%%%%%%%%%%%%%%%%%%%%%%
################## %%%%%%%%%%%%%%%%%%%%%%%
################### %%%%%%%%%%%%%%%%%%%%%%%
#################### %%%%%%%%%%%%%%%%%%%%%%
##################### %%%%%%%%%%%%%%%%%%%%%
###################### %%%%%%%%%%%%%%%%%%%%
####################### %%%%%%%%%%%%%%%%%%
####################### %%%%%%%%%%%%%%%%%
###################### %%%%%%%%%%%%%%%%%
#################### %%%%%%%%%%%%%%%%%
################# %%%%%%%%%%%%%%%%%
############### %%%%%%%%%%%%%%%%
############ %%%%%%%%%%%%%%%
######## %%%%%%%%%%%%%%
%%%%%%%%%%%
| The OpenMC Monte Carlo Code
Copyright | 2011-2025 MIT, UChicago Argonne LLC, and contributors
License | https://docs.openmc.org/en/latest/license.html
Version | 0.15.3-dev144
Commit Hash | 14ab4d3b6f3d77b73e8bb480782bbe1bdd56483d
Date/Time | 2025-06-04 19:46:55
OpenMP Threads | 16
Reading model XML file 'model.xml' ...
Reading cross sections XML file...
Reading H1 from /home/jon/nuclear_data/neutron/H1.h5
Reading H2 from /home/jon/nuclear_data/neutron/H2.h5
Reading O16 from /home/jon/nuclear_data/neutron/O16.h5
Reading O17 from /home/jon/nuclear_data/neutron/O17.h5
Reading O18 from /home/jon/nuclear_data/neutron/O18.h5
Minimum neutron data temperature: 294.0 K
Maximum neutron data temperature: 294.0 K
Preparing distributed cell instances...
Writing summary.h5 file...
Maximum neutron transport energy: 20000000.0 eV for H1
===============> FIXED SOURCE TRANSPORT SIMULATION <===============
Simulating batch 1
Simulating batch 2
Simulating batch 3
Simulating batch 4
Simulating batch 5
Simulating batch 6
Simulating batch 7
Simulating batch 8
Simulating batch 9
Simulating batch 10
Creating state point statepoint.10.h5...
=======================> TIMING STATISTICS <=======================
Total time for initialization = 1.8620e-01 seconds
Reading cross sections = 1.7641e-01 seconds
Total time in simulation = 5.8975e-01 seconds
Time in transport only = 5.1639e-01 seconds
Time in active batches = 5.8975e-01 seconds
Time accumulating tallies = 3.9983e-02 seconds
Time writing statepoints = 3.1614e-02 seconds
Total time for finalization = 6.2546e-04 seconds
Total time elapsed = 7.7765e-01 seconds
Calculation Rate (active) = 16956.5 particles/second
============================> RESULTS <============================
Leakage Fraction = 0.00000 +/- 0.00000
This section extracts the cell tally data from the results file and plots neutron flux across the cell. Selecting log-log scale will allow you to see a distribution of thermal neutrons.
# open the results file
results = openmc.StatePoint(results_filename)
#extracts the tally values from the simulation results
cell_tally = results.get_tally(name='cell_flux_tally')
# flattens the ndarray into a 1d array
openmc_flux = cell_tally.mean.flatten()
Discussion on results and units of flux
Openmc like most of other neutronics codes accumulates track lengths within cell volumes, i.e the length that a neutron travels in a material.
A track length has units of centimeters (cm).
Neutronics codes typically make use of cm instead of the SI base unit for length of meters, this is partly historical and partly due to the format of nuclear data files.
A flux score on a cell in OpenMC therefore returns the average length that neutrons travel through the cell.
As we have simulated many neutrons we can get an average result and as we have batches we can get a standard deviation on that result.
OpenMC returns a flux tally in units of “neutron cm per source neutron”
To convert this into more common units of flux “neutrons per cm2 per second” we must first divide by the volume.
As this is a cell tally we divide by the cell volume, if this was a mesh tally we would divide by the voxel volume.
This gives us units of “neutrons per cm2 per source neutron”.
In a fixed source simulation (not an fission eigenvalue simulation) we can then scale the result by the number of neutrons per second that the source emits.
To find the number of neutrons you would typically know the power of the fusion reactor in watts.
Assuming we are using DT fuel then we know that each neutron resulted in 17.6MeV or 2.8e-18Joules of energy.
Therefore the neutrons per second is power in watts / 2.8e-18. For a 500MW fusion power reactor we would get 500e6/2.8e-18=1.785e+26 neutrons per second
So our source strength is 1.785e+26 and we multiply out flux by this to get units of “neutrons per cm2 per second”
For ICF (Inertial Confinement Fusion) you might use units of per shot instead of per second
volume_of_cell = 5.24e8 # in units of cm3
reactor_power = 500e6 # in units of watts
energy_of_each_fusion_reaction = 17.5e6 * 1.60218e-19 # eV converted to Joules
neutrons_per_second = reactor_power / energy_of_each_fusion_reaction
flux = (openmc_flux / volume_of_cell) * neutrons_per_second # divide by cell volume and then multiply by source strength
print(f'neutron flux = {flux} neutrons per cm2 per second')
neutron flux = [7.74281414e+13] neutrons per cm2 per second
Learning Outcomes for Part 1:
Neutron flux in found in units of “neutron cm per source neutron” but can be converted to “neutrons per cm2 per second”