7. Execution Settings

Once you have created the materials and geometry for your simulation, the last step to have a complete model is to specify execution settings through the openmc.Settings class. At a minimum, you need to specify a source distribution and how many particles to run. Many other execution settings can be set using the openmc.Settings object, but they are generally optional.

7.1. Run Modes

The Settings.run_mode attribute controls what run mode is used when openmc is executed. There are five different run modes that can be specified:

‘eigenvalue’

Runs a \(k\) eigenvalue simulation. See Eigenvalue Calculations for a full description of eigenvalue calculations. In this mode, the Settings.source specifies a starting source that is only used for the first fission generation.

‘fixed source’

Runs a fixed-source calculation with a specified external source, specified in the Settings.source attribute.

‘volume’

Runs a stochastic volume calculation.

‘plot’

Generates slice or voxel plots (see Geometry Visualization).

‘particle restart’

Simulate a single source particle using a particle restart file.

So, for example, to specify that OpenMC should be run in fixed source mode, you would need to instantiate a openmc.Settings object and assign the Settings.run_mode attribute:

settings = openmc.Settings()
settings.run_mode = 'fixed source'

If you don’t specify a run mode, the default run mode is ‘eigenvalue’.

7.2. Run Strategy

For a fixed source simulation, the total number of source particle histories simulated is broken up into a number of batches, each corresponding to a realization of the tally random variables. Thus, you need to specify both the number of batches (Settings.batches) as well as the number of particles per batch (Settings.particles).

For a \(k\) eigenvalue simulation, particles are grouped into fission generations, as described in Eigenvalue Calculations. Successive fission generations can be combined into a batch for statistical purposes. By default, a batch will consist of only a single fission generation, but this can be changed with the Settings.generations_per_batch attribute. For problems with a high dominance ratio, using multiple generations per batch can help reduce underprediction of variance, thereby leading to more accurate confidence intervals. Tallies should not be scored to until the source distribution converges, as described in Method of Successive Generations, which may take many generations. To specify the number of batches that should be discarded before tallies begin to accumulate, use the Settings.inactive attribute.

The following example shows how one would simulate 10000 particles per generation, using 10 generations per batch, 150 total batches, and discarding 5 batches. Thus, a total of 145 active batches (or 1450 generations) will be used for accumulating tallies.

settings.particles = 10000
settings.generations_per_batch = 10
settings.batches = 150
settings.inactive = 5

7.2.1. Number of Batches

In general, the stochastic uncertainty in your simulation results is directly related to how many total active particles are simulated (the product of the number of active batches, number of generations per batch, and number of particles). At a minimum, you should use enough active batches so that the central limit theorem is satisfied (about 30). Otherwise, reducing the overall uncertainty in your simulation by a factor of 2 will require using 4 times as many batches (since the standard deviation decreases as \(1/\sqrt{N}\)).

7.2.2. Number of Inactive Batches

For \(k\) eigenvalue simulations, the source distribution is not known a priori. Thus, a “guess” of the source distribution is made and then iterated on, with the source evolving closer to the true distribution at each iteration. Once the source distribution has converged, it is then safe to start accumulating tallies. Consequently, a preset number of inactive batches are run before the active batches (where tallies are turned on) begin. The number of inactive batches necessary to reach a converged source depends on the spatial extent of the problem, its dominance ratio, what boundary conditions are used, and many other factors. For small problems, using 50–100 inactive batches is likely sufficient. For larger models, many hundreds of inactive batches may be necessary. Users are recommended to use the Shannon entropy diagnostic as a way of determining how many inactive batches are necessary.

Specifying the initial source used for the very first batch is described in below. Although the initial source is arbitrary in the sense that any source will eventually converge to the correct distribution, using a source guess that is closer to the actual converged source distribution will translate into needing fewer inactive batches (and hence less simulation time).

For fixed source simulations, the source distribution is known exactly, so no inactive batches are needed. In this case the Settings.inactive attribute can be omitted since it defaults to zero.

7.2.3. Number of Generations per Batch

The standard deviation of tally results is calculated assuming that all realizations (batches) are independent. However, in a \(k\) eigenvalue calculation, the source sites for each batch are produced from fissions in the preceding batch, resulting in a correlation between successive batches. This correlation can result in an underprediction of the variance. That is, the variance reported is actually less than the true variance. To mitigate this effect, OpenMC allows you to group together multiple fission generations into a single batch for statistical purposes, rather than having each fission generation be a separate batch, which is the default behavior.

7.2.4. Number of Particles per Generation

There are several considerations for choosing the number of particles per generation. As discussed in Number of Batches, the total number of active particles will determine the level of stochastic uncertainty in simulation results, so using a higher number of particles will result in less uncertainty. For parallel simulations that use OpenMP and/or MPI, the number of particles per generation should be large enough to ensure good load balancing between threads. For example, if you are running on a single processor with 32 cores, each core should have at least 100 particles or so (i.e., at least 3,200 particles per generation should be used). Using a larger number of particles per generation can also help reduce the cost of synchronization and communication between batches. For \(k\) eigenvalue calculations, experts recommend at least 10,000 particles per generation to avoid any bias in the estimate of \(k\) eigenvalue or tallies.

7.3. External Source Distributions

External source distributions can be specified through the Settings.source attribute. If you have a single external source, you can create an instance of openmc.Source and use it to set the Settings.source attribute. If you have multiple external sources with varying source strengths, Settings.source should be set to a list of openmc.Source objects.

The openmc.Source class has four main attributes that one can set: Source.space, which defines the spatial distribution, Source.angle, which defines the angular distribution, Source.energy, which defines the energy distribution, and Source.time, which defines the time distribution.

The spatial distribution can be set equal to a sub-class of openmc.stats.Spatial; common choices are openmc.stats.Point or openmc.stats.Box. To independently specify distributions in the \(x\), \(y\), and \(z\) coordinates, you can use openmc.stats.CartesianIndependent. To independently specify distributions using spherical or cylindrical coordinates, you can use openmc.stats.SphericalIndependent or openmc.stats.CylindricalIndependent, respectively.

The angular distribution can be set equal to a sub-class of openmc.stats.UnitSphere such as openmc.stats.Isotropic, openmc.stats.Monodirectional, or openmc.stats.PolarAzimuthal. By default, if no angular distribution is specified, an isotropic angular distribution is used. As an example of a non-trivial angular distribution, the following code would create a conical distribution with an aperture of 30 degrees pointed in the positive x direction:

from math import pi, cos
aperture = 30.0
mu = openmc.stats.Uniform(cos(aperture/2), 1.0)
phi = openmc.stats.Uniform(0.0, 2*pi)
angle = openmc.stats.PolarAzimuthal(mu, phi, reference_uvw=(1., 0., 0.))

The energy distribution can be set equal to any univariate probability distribution. This could be a probability mass function (openmc.stats.Discrete), a Watt fission spectrum (openmc.stats.Watt), or a tabular distribution (openmc.stats.Tabular). By default, if no energy distribution is specified, a Watt fission spectrum with \(a\) = 0.988 MeV and \(b\) = 2.249 MeV -1 is used.

The time distribution can be set equal to any univariate probability distribution. This could be a probability mass function (openmc.stats.Discrete), a uniform distribution (openmc.stats.Uniform), or a tabular distribution (openmc.stats.Tabular). By default, if no time distribution is specified, particles are started at \(t=0\).

As an example, to create an isotropic, 10 MeV monoenergetic source uniformly distributed over a cube centered at the origin with an edge length of 10 cm, and emitting a pulse of particles from 0 to 10 µs, one would run:

source = openmc.Source()
source.space = openmc.stats.Box((-5, -5, -5), (5, 5, 5))
source.angle = openmc.stats.Isotropic()
source.energy = openmc.stats.Discrete([10.0e6], [1.0])
source.time = openmc.stats.Uniform(0, 1e-6)
settings.source = source

The openmc.Source class also has a Source.strength attribute that indicates the relative strength of a source distribution if multiple are used. For example, to create two sources, one that should be sampled 70% of the time and another that should be sampled 30% of the time:

src1 = openmc.Source()
src1.strength = 0.7
...

src2 = openmc.Source()
src2.strength = 0.3
...

settings.source = [src1, src2]

Finally, the Source.particle attribute can be used to indicate the source should be composed of particles other than neutrons. For example, the following would generate a photon source:

source = openmc.Source()
source.particle = 'photon'
...

settings.source = source

For a full list of all classes related to statistical distributions, see openmc.stats – Statistics.

7.3.1. File-based Sources

OpenMC can use a pregenerated HDF5 source file by specifying the filename argument to openmc.Source:

settings.source = openmc.Source(filename='source.h5')

Statepoint and source files are generated automatically when a simulation is run and can be used as the starting source in a new simulation. Alternatively, a source file can be manually generated with the openmc.write_source_file() function. This is particularly useful for coupling OpenMC with another program that generates a source to be used in OpenMC.

A source file based on particles that cross one or more surfaces can be generated during a simulation using the Settings.surf_source_write attribute:

settings.surf_source_write = {
    'surfaces_ids': [1, 2, 3],
    'max_particles': 10000
}

In this example, at most 10,000 source particles are stored when particles cross surfaces with IDs of 1, 2, or 3.

7.3.2. Custom Sources

It is often the case that one may wish to simulate a complex source distribution that is not possible to represent with the classes described above. For these situations, it is possible to define a complex source class containing an externally defined source function that is loaded at runtime. A simple example source is shown below.

#include <memory> // for unique_ptr

#include "openmc/random_lcg.h"
#include "openmc/source.h"
#include "openmc/particle.h"

class CustomSource : public openmc::Source
{
  openmc::SourceSite sample(uint64_t* seed) const
  {
    openmc::SourceSite particle;
    // weight
    particle.particle = openmc::ParticleType::neutron;
    particle.wgt = 1.0;
    // position
    double angle = 2.0 * M_PI * openmc::prn(seed);
    double radius = 3.0;
    particle.r.x = radius * std::cos(angle);
    particle.r.y = radius * std::sin(angle);
    particle.r.z = 0.0;
    // angle
    particle.u = {1.0, 0.0, 0.0};
    particle.E = 14.08e6;
    particle.delayed_group = 0;
    return particle;
  }
};

extern "C" std::unique_ptr<CustomSource> openmc_create_source(std::string parameters)
{
  return std::make_unique<CustomSource>();
}

The above source creates monodirectional 14.08 MeV neutrons that are distributed in a ring with a 3 cm radius. This routine is not particularly complex, but should serve as an example upon which to build more complicated sources.

Note

The source class must inherit from openmc::Source and implement a sample() function.

Note

The openmc_create_source() function signature must be declared extern "C".

Note

You should only use the openmc::prn() random number generator.

In order to build your external source, you will need to link it against the OpenMC shared library. This can be done by writing a CMakeLists.txt file:

cmake_minimum_required(VERSION 3.3 FATAL_ERROR)
project(openmc_sources CXX)
add_library(source SHARED source_ring.cpp)
find_package(OpenMC REQUIRED HINTS <path to openmc>)
target_link_libraries(source OpenMC::libopenmc)

After running cmake and make, you will have a libsource.so (or .dylib) file in your build directory. Setting the openmc.Source.library attribute to the path of this shared library will indicate that it should be used for sampling source particles at runtime.

7.3.3. Custom Parameterized Sources

Some custom sources may have values (parameters) that can be changed between runs. This is supported by using the openmc_create_source() function to pass parameters defined in the openmc.Source.parameters attribute to the source class when it is created:

#include <memory> // for unique_ptr

#include "openmc/source.h"
#include "openmc/particle.h"

class CustomSource : public openmc::Source {
public:
  CustomSource(double energy) : energy_{energy} { }

  // Samples from an instance of this class.
  openmc::SourceSite sample(uint64_t* seed) const
  {
    openmc::SourceSite particle;
    // weight
    particle.particle = openmc::ParticleType::neutron;
    particle.wgt = 1.0;
    // position
    particle.r.x = 0.0;
    particle.r.y = 0.0;
    particle.r.z = 0.0;
    // angle
    particle.u = {1.0, 0.0, 0.0};
    particle.E = this->energy_;
    particle.delayed_group = 0;

    return particle;
  }

private:
  double energy_;
};

extern "C" std::unique_ptr<CustomSource> openmc_create_source(std::string parameter) {
  double energy = std::stod(parameter);
  return std::make_unique<CustomSource>(energy);
}

As with the basic custom source functionality, the custom source library location must be provided in the openmc.Source.library attribute.

7.4. Shannon Entropy

To assess convergence of the source distribution, the scalar Shannon entropy metric is often used in Monte Carlo codes. OpenMC also allows you to calculate Shannon entropy at each generation over a specified mesh, created using the openmc.RegularMesh class. After instantiating a RegularMesh, you need to specify the lower-left coordinates of the mesh (RegularMesh.lower_left), the number of mesh cells in each direction (RegularMesh.dimension) and either the upper-right coordinates of the mesh (RegularMesh.upper_right) or the width of each mesh cell (RegularMesh.width). Once you have a mesh, simply assign it to the Settings.entropy_mesh attribute.

entropy_mesh = openmc.RegularMesh()
entropy_mesh.lower_left = (-50, -50, -25)
entropy_mesh.upper_right = (50, 50, 25)
entropy_mesh.dimension = (8, 8, 8)

settings.entropy_mesh = entropy_mesh

If you’re unsure of what bounds to use for the entropy mesh, you can try getting a bounding box for the entire geometry using the Geometry.bounding_box property:

geom = openmc.Geometry()
...
m = openmc.RegularMesh()
m.lower_left, m.upper_right = geom.bounding_box
m.dimension = (8, 8, 8)

settings.entropy_mesh = m

7.5. Photon Transport

In addition to neutrons, OpenMC is also capable of simulating the passage of photons through matter. This allows the modeling of photon production from neutrons as well as pure photon calculations. The Settings.photon_transport attribute can be used to enable photon transport:

settings.photon_transport = True

The way in which OpenMC handles secondary charged particles can be specified with the Settings.electron_treatment attribute. By default, the thick-target bremsstrahlung (TTB) approximation is used to generate bremsstrahlung radiation emitted by electrons and positrons created in photon interactions. To neglect secondary bremsstrahlung photons and instead deposit all energy from electrons locally, the local energy deposition option can be selected:

settings.electron_treatment = 'led'

Note

Some features related to photon transport are not currently implemented, including:

  • Tallying photon energy deposition.

  • Generating a photon source from a neutron calculation that can be used for a later fixed source photon calculation.

  • Photoneutron reactions.

7.6. Generation of Output Files

A number of attributes of the openmc.Settings class can be used to control what files are output and how often. First, there is the Settings.output attribute which takes a dictionary having keys ‘summary’, ‘tallies’, and ‘path’. The first two keys controls whether a summary.h5 and tallies.out file are written, respectively (see Viewing and Analyzing Results for a description of those files). By default, output files are written to the current working directory; this can be changed by setting the ‘path’ key. For example, if you want to disable the tallies.out file and write the summary.h5 to a directory called ‘results’, you’d specify the Settings.output dictionary as:

settings.output = {
    'tallies': False,
    'path': 'results'
}

Generation of statepoint and source files is handled separately through the Settings.statepoint and Settings.sourcepoint attributes. Both of those attributes expect dictionaries and have a ‘batches’ key which indicates at which batches statepoints and source files should be written. Note that by default, the source is written as part of the statepoint file; this behavior can be changed by the ‘separate’ and ‘write’ keys of the Settings.sourcepoint dictionary, the first of which indicates whether the source should be written to a separate file and the second of which indicates whether the source should be written at all.

As an example, to write a statepoint file every five batches:

settings.batches = n
settings.statepoint = {'batches': range(5, n + 5, 5)}

7.6.1. Particle Track Files

OpenMC can generate a particle track file that contains track information (position, direction, energy, time, weight, cell ID, and material ID) for each state along a particle’s history. There are two ways to indicate which particles and/or how many particles should have their tracks written. First, you can identify specific source particles by their batch, generation, and particle ID numbers:

settings.tracks = [
  (1, 1, 50),
  (2, 1, 30),
  (5, 1, 75)
]

In this example, track information would be written for the 50th particle in the 1st generation of batch 1, the 30th particle in the first generation of batch 2, and the 75th particle in the 1st generation of batch 5. Unless you are using more than one generation per batch (see Run Strategy), the generation number should be 1. Alternatively, you can run OpenMC in a mode where track information is written for all particles, up to a user-specified limit:

openmc.run(tracks=True)

In this case, you can control the maximum number of source particles for which tracks will be written as follows:

settings.max_tracks = 1000

Particle track information is written to the tracks.h5 file, which can be analyzed using the Tracks class:

>>> tracks = openmc.Tracks('tracks.h5')
>>> tracks
[<Track (1, 1, 50): 151 particles>,
 <Track (2, 1, 30): 191 particles>,
 <Track (5, 1, 75): 81 particles>]

Each Track object stores a list of track information for every primary/secondary particle. In the above example, the first source particle produced 150 secondary particles for a total of 151 particles. Information for each primary/secondary particle can be accessed using the particle_tracks attribute:

>>> first_track = tracks[0]
>>> first_track.particle_tracks
[<ParticleTrack: neutron, 120 states>,
 <ParticleTrack: photon, 6 states>,
 <ParticleTrack: electron, 2 states>,
 <ParticleTrack: electron, 2 states>,
 <ParticleTrack: electron, 2 states>,
 ...
 <ParticleTrack: electron, 2 states>,
 <ParticleTrack: electron, 2 states>]
>>> photon = first_track.particle_tracks[1]

The ParticleTrack class is a named tuple indicating the particle type and then a NumPy array of the “states”. The states array is a compound type with a field for each physical quantity (position, direction, energy, time, weight, cell ID, and material ID). For example, to get the position for the above particle track:

>>> photon.states['r']
array([(-11.92987939, -12.28467295, 0.67837495),
       (-11.95213726, -12.2682    , 0.68783964),
       (-12.2682    , -12.03428339, 0.82223855),
       (-12.5913778 , -11.79510096, 0.95966298),
       (-12.6622572 , -11.74264344, 0.98980293),
       (-12.6907775 , -11.7215357 , 1.00193058)],
      dtype=[('x', '<f8'), ('y', '<f8'), ('z', '<f8')])

The full list of fields is as follows:

r

Position (each direction in [cm])

u

Direction

E

Energy in [eV]

time

Time in [s]

wgt

Weight

cell_id

Cell ID

cell_instance

Cell instance

material_id

Material ID

Both the Tracks and Track classes have a filter method that allows you to get a subset of tracks that meet a given criteria. For example, to get all tracks that involved a photon:

>>> tracks.filter(particle='photon')
[<Track (1, 1, 50): 151 particles>,
 <Track (2, 1, 30): 191 particles>,
 <Track (5, 1, 75): 81 particles>]

The openmc.Tracks.filter() method returns a new Tracks instance, whereas the openmc.Track.filter() method returns a new Track instance.

Note

If you are using an MPI-enabled install of OpenMC and run a simulation with more than one process, a separate track file will be written for each MPI process with the filename tracks_p#.h5 where # is the rank of the corresponding process. Multiple track files can be combined with the openmc-track-combine script:

openmc-track-combine tracks_p*.h5 --out tracks.h5