Time Cutoff

docs/source/io_formats/settings.rst

@@ -60,13 +60,16 @@ fission.
 ``<cutoff>`` Element
 --------------------
 
-The ``<cutoff>`` element indicates two kinds of cutoffs. The first is the weight
-cutoff used below which particles undergo Russian roulette. Surviving particles
-are assigned a user-determined weight. Note that weight cutoffs and Russian
-rouletting are not turned on by default. The second is the energy cutoff which
-is used to kill particles under certain energy. The energy cutoff should not be
-used unless you know particles under the energy are of no importance to results
-you care. This element has the following attributes/sub-elements:
+The ``<cutoff>`` element indicates three kinds of cutoffs. The first is the
+weight cutoff used below which particles undergo Russian roulette. Surviving
+particles are assigned a user-determined weight. Note that weight cutoffs and
+Russian rouletting are not turned on by default. The second is the energy cutoff
+which is used to kill particles under certain energy. The energy cutoff should
+not be used unless you know particles under the energy are of no importance to
+results you care. The third is the time cutoff used to kill particles whose time
+exceeds a specific cutoff.
+
+This element has the following attributes/sub-elements:
 
   :weight:
     The weight below which particles undergo Russian roulette.
@@ -99,6 +102,26 @@ you care. This element has the following attributes/sub-elements:
 
     *Default*: 0.0
 
+  :time_neutrons 
+    The time above which neutrons will be killed.
+
+    *Default*: 1.7976931348623157E+308
+
+  :time_photons
+    The time above which photons will be killed.
+
+    *Default*: 1.7976931348623157E+308
+
+  :time_electrons 
+    The time above which electrons will be killed.
+
+    *Default*: 1.7976931348623157E+308
+
+  :time_positrons
+    The time above which positorns will be killed.
+
+    *Default*: 1.7976931348623157E+308
+
 ----------------------------
 ``<delayed_photon_scaling>``
 ----------------------------

include/openmc/settings.h

@@ -27,14 +27,14 @@ extern bool assume_separate;      //!< assume tallies are spatially separate?
 extern bool check_overlaps;       //!< check overlaps in geometry?
 extern bool confidence_intervals; //!< use confidence intervals for results?
 extern bool
-  create_fission_neutrons; //!< create fission neutrons (fixed source)?
+  create_fission_neutrons;        //!< create fission neutrons (fixed source)?
 extern bool create_delayed_neutrons; //!< create delayed fission neutrons?
-extern "C" bool cmfd_run;  //!< is a CMFD run?
+extern "C" bool cmfd_run;            //!< is a CMFD run?
 extern bool
   delayed_photon_scaling;   //!< Scale fission photon yield to include delayed
 extern "C" bool entropy_on; //!< calculate Shannon entropy?
 extern "C" bool
-  event_based; //!< use event-based mode (instead of history-based)
+  event_based;              //!< use event-based mode (instead of history-based)
 extern bool legendre_to_tabular; //!< convert Legendre distributions to tabular?
 extern bool material_cell_offsets; //!< create material cells offsets?
 extern "C" bool output_summary;    //!< write summary.h5?
@@ -70,8 +70,8 @@ extern std::string path_particle_restart; //!< path to a particle restart file
 extern std::string path_sourcepoint;      //!< path to a source file
 extern "C" std::string path_statepoint;   //!< path to a statepoint file
 
-extern "C" int32_t n_inactive;         //!< number of inactive batches
-extern "C" int32_t max_lost_particles; //!< maximum number of lost particles
+extern "C" int32_t n_inactive;            //!< number of inactive batches
+extern "C" int32_t max_lost_particles;    //!< maximum number of lost particles
 extern double
   rel_max_lost_particles; //!< maximum number of lost particles, relative to the
                           //!< total number of particles
@@ -82,9 +82,11 @@ extern int64_t
   max_particles_in_flight; //!< Max num. event-based particles in flight
 
 extern ElectronTreatment
-  electron_treatment; //!< how to treat secondary electrons
+  electron_treatment;         //!< how to treat secondary electrons
 extern array<double, 4>
-  energy_cutoff; //!< Energy cutoff in [eV] for each particle type
+  energy_cutoff;              //!< Energy cutoff in [eV] for each particle type
+extern array<double, 4>
+  time_cutoff;                //!< Time cutoff in [s] for each particle type
 extern int
   legendre_to_tabular_points; //!< number of points to convert Legendres
 extern int max_order;         //!< Maximum Legendre order for multigroup data
@@ -99,24 +101,24 @@ extern vector<std::string>
   res_scat_nuclides;     //!< Nuclides using res. upscattering treatment
 extern RunMode run_mode; //!< Run mode (eigenvalue, fixed src, etc.)
 extern std::unordered_set<int>
-  sourcepoint_batch; //!< Batches when source should be written
+  sourcepoint_batch;     //!< Batches when source should be written
 extern std::unordered_set<int>
-  statepoint_batch; //!< Batches when state should be written
+  statepoint_batch;      //!< Batches when state should be written
 extern std::unordered_set<int>
-  source_write_surf_id; //!< Surface ids where sources will be written
+  source_write_surf_id;  //!< Surface ids where sources will be written
 extern int max_splits; //!< maximum number of particle splits for weight windows
 extern int64_t max_surface_particles; //!< maximum number of particles to be
                                       //!< banked on surfaces per process
 extern TemperatureMethod
-  temperature_method; //!< method for choosing temperatures
+  temperature_method;                 //!< method for choosing temperatures
 extern double
   temperature_tolerance; //!< Tolerance in [K] on choosing temperatures
 extern double temperature_default; //!< Default T in [K]
 extern array<double, 2>
-  temperature_range;           //!< Min/max T in [K] over which to load xs
-extern int trace_batch;        //!< Batch to trace particle on
-extern int trace_gen;          //!< Generation to trace particle on
-extern int64_t trace_particle; //!< Particle ID to enable trace on
+  temperature_range;               //!< Min/max T in [K] over which to load xs
+extern int trace_batch;            //!< Batch to trace particle on
+extern int trace_gen;              //!< Generation to trace particle on
+extern int64_t trace_particle;     //!< Particle ID to enable trace on
 extern vector<array<int, 3>>
   track_identifiers;               //!< Particle numbers for writing tracks
 extern int trigger_batch_interval; //!< Batch interval for triggers

openmc/settings.py

@@ -49,15 +49,17 @@ class Settings:
         deviation.
     create_fission_neutrons : bool
         Indicate whether fission neutrons should be created or not.
-    cutoff : dict
-        Dictionary defining weight cutoff and energy cutoff. The dictionary may
-        have six keys, 'weight', 'weight_avg', 'energy_neutron', 'energy_photon',
-        'energy_electron', and 'energy_positron'. Value for 'weight'
+    cutoff : dict 
+        Dictionary defining weight cutoff, energy cutoff and time cutoff. The
+        dictionary may have ten keys, 'weight', 'weight_avg', 'energy_neutron',
+        'energy_photon', 'energy_electron', 'energy_positron', 'time_neutron',
+        'time_photon', 'time_electron', and 'time_positron'. Value for 'weight'
         should be a float indicating weight cutoff below which particle undergo
         Russian roulette. Value for 'weight_avg' should be a float indicating
-        weight assigned to particles that are not killed after Russian
-        roulette. Value of energy should be a float indicating energy in eV
-        below which particle type will be killed.
+        weight assigned to particles that are not killed after Russian roulette.
+        Value of energy should be a float indicating energy in eV below which
+        particle type will be killed. Value of time should be a float indicating
+        time in s above which particle type will be killed.
     delayed_photon_scaling : bool
         Indicate whether to scale the fission photon yield by (EGP + EGD)/EGP
         where EGP is the energy release of prompt photons and EGD is the energy
@@ -736,6 +738,10 @@ def cutoff(self, cutoff: dict):
                          'energy_positron']:
                 cv.check_type('energy cutoff', cutoff[key], Real)
                 cv.check_greater_than('energy cutoff', cutoff[key], 0.0)
+            elif key in ['time_neutron', 'time_photon', 'time_electron',
+                         'time_positron']:
+                cv.check_type('time cutoff', cutoff[key], Real)
+                cv.check_greater_than('time cutoff', cutoff[key], 0.0)
             else:
                 msg = f'Unable to set cutoff to "{key}" which is unsupported ' \
                       'by OpenMC'
@@ -1457,7 +1463,8 @@ def _cutoff_from_xml_element(self, root):
         if elem is not None:
             self.cutoff = {}
             for key in ('energy_neutron', 'energy_photon', 'energy_electron',
-                        'energy_positron', 'weight', 'weight_avg'):
+                        'energy_positron', 'weight', 'weight_avg', 'time_neutron',
+                        'time_photon', 'time_electron', 'time_positron'):
                 value = get_text(elem, key)
                 if value is not None:
                     self.cutoff[key] = float(value)

src/finalize.cpp

@@ -28,6 +28,7 @@
 #include "openmc/weight_windows.h"
 
 #include "xtensor/xview.hpp"
+#include <limits> // for numeric_limits
 
 namespace openmc {
 
@@ -78,6 +79,9 @@ int openmc_finalize()
   settings::electron_treatment = ElectronTreatment::LED;
   settings::delayed_photon_scaling = true;
   settings::energy_cutoff = {0.0, 1000.0, 0.0, 0.0};
+  settings::time_cutoff = {std::numeric_limits<double>::max(),
+    std::numeric_limits<double>::max(), std::numeric_limits<double>::max(),
+    std::numeric_limits<double>::max()};
   settings::entropy_on = false;
   settings::event_based = false;
   settings::gen_per_batch = 1;

src/particle.cpp

@@ -197,6 +197,12 @@ void Particle::event_advance()
   // Select smaller of the two distances
   double distance = std::min(boundary().distance, collision_distance());
 
+  // Kill particle if its time exceeds the cutoff
+  if (time() + (distance / this->speed()) >
+      settings::time_cutoff[static_cast<int>(type())]) {
+    wgt() = 0.0;
+  }
+
   // Advance particle in space and time
   for (int j = 0; j < n_coord(); ++j) {
     coord(j).r += distance * coord(j).u;

src/settings.cpp

@@ -93,6 +93,9 @@ int64_t max_particles_in_flight {100000};
 
 ElectronTreatment electron_treatment {ElectronTreatment::TTB};
 array<double, 4> energy_cutoff {0.0, 1000.0, 0.0, 0.0};
+array<double, 4> time_cutoff {std::numeric_limits<double>::max(),
+  std::numeric_limits<double>::max(), std::numeric_limits<double>::max(),
+  std::numeric_limits<double>::max()};
 int legendre_to_tabular_points {C_NONE};
 int max_order {0};
 int n_log_bins {8000};
@@ -536,6 +539,18 @@ void read_settings_xml(pugi::xml_node root)
       energy_cutoff[3] =
         std::stod(get_node_value(node_cutoff, "energy_positron"));
     }
+    if (check_for_node(node_cutoff, "time_neutron")) {
+      time_cutoff[0] = std::stod(get_node_value(node_cutoff, "time_neutron"));
+    }
+    if (check_for_node(node_cutoff, "time_photon")) {
+      time_cutoff[1] = std::stod(get_node_value(node_cutoff, "time_photon"));
+    }
+    if (check_for_node(node_cutoff, "time_electron")) {
+      time_cutoff[2] = std::stod(get_node_value(node_cutoff, "time_electron"));
+    }
+    if (check_for_node(node_cutoff, "time_positron")) {
+      time_cutoff[3] = std::stod(get_node_value(node_cutoff, "time_positron"));
+    }
   }
 
   // Particle trace

tests/unit_tests/test_settings.py

@@ -26,7 +26,9 @@ def test_export_to_xml(run_in_tmpdir):
     s.survival_biasing = True
     s.cutoff = {'weight': 0.25, 'weight_avg': 0.5, 'energy_neutron': 1.0e-5,
                 'energy_photon': 1000.0, 'energy_electron': 1.0e-5,
-                'energy_positron': 1.0e-5}
+                'energy_positron': 1.0e-5, 'time_neutron': 1.0e-5, 
+                'time_photon': 1.0e-5, 'time_electron': 1.0e-5, 
+                'time_positron': 1.0e-5}
     mesh = openmc.RegularMesh()
     mesh.lower_left = (-10., -10., -10.)
     mesh.upper_right = (10., 10., 10.)
@@ -86,7 +88,9 @@ def test_export_to_xml(run_in_tmpdir):
     assert s.survival_biasing
     assert s.cutoff == {'weight': 0.25, 'weight_avg': 0.5,
                         'energy_neutron': 1.0e-5, 'energy_photon': 1000.0,
-                        'energy_electron': 1.0e-5, 'energy_positron': 1.0e-5}
+                        'energy_electron': 1.0e-5, 'energy_positron': 1.0e-5,
+                        'time_neutron': 1.0e-5, 'time_photon': 1.0e-5, 
+                        'time_electron': 1.0e-5, 'time_positron': 1.0e-5}
     assert isinstance(s.entropy_mesh, openmc.RegularMesh)
     assert s.entropy_mesh.lower_left == [-10., -10., -10.]
     assert s.entropy_mesh.upper_right == [10., 10., 10.]