Adding support for max particle events

docs/source/io_formats/settings.rst

@@ -262,6 +262,14 @@ may be more efficient computationally.
 
   *Default*: 100000
 
+----------------------------------------
+``<max_particle_events>`` Element
+----------------------------------------
+
+This element indicates the maximum number of events a particle can undergo.
+
+  *Default*: 100000
+
 ---------------------------
 ``<max_order>`` Element
 ---------------------------
@@ -280,7 +288,7 @@ then, OpenMC will only use up to the :math:`P_1` data.
 ``<max_splits>`` Element
 ---------------------------
 
-The ``<max_splits>`` element indicates the number of times a particle can split during a history. 
+The ``<max_splits>`` element indicates the number of times a particle can split during a history.
 
   *Default*: 1000
 

include/openmc/particle_data.h

@@ -27,9 +27,6 @@ constexpr int MAX_DELAYED_GROUPS {8};
 
 constexpr double CACHE_INVALID {-1.0};
 
-// Maximum number of collisions/crossings
-constexpr int MAX_EVENTS {1000000};
-
 //==========================================================================
 // Aliases and type definitions
 

include/openmc/settings.h

@@ -94,7 +94,8 @@ extern "C" int64_t n_particles;   //!< number of particles per generation
 
 extern int64_t
   max_particles_in_flight; //!< Max num. event-based particles in flight
-
+extern int
+  max_particle_events; //!< Maximum number of particle events
 extern ElectronTreatment
   electron_treatment; //!< how to treat secondary electrons
 extern array<double, 4>

openmc/settings.py

@@ -109,6 +109,10 @@ class Settings:
         parallelism.
 
         .. versionadded:: 0.12
+    max_particle_events : int
+        Maximum number of allowed particle events per particle history.
+
+        .. versionadded:: 0.14.1
     max_order : None or int
         Maximum scattering order to apply globally when in multi-group mode.
     max_splits : int
@@ -938,6 +942,16 @@ def max_particles_in_flight(self, value: int):
         cv.check_greater_than('max particles in flight', value, 0)
         self._max_particles_in_flight = value
 
+    @property
+    def max_particle_events(self) -> int:
+        return self._max_particle_events
+
+    @max_particle_events.setter
+    def max_particle_events(self, value: int):
+        cv.check_type('max particle events', value, Integral)
+        cv.check_greater_than('max particle events', value, 0)
+        self._max_particle_events = value
+
     @property
     def write_initial_source(self) -> bool:
         return self._write_initial_source
@@ -1341,6 +1355,11 @@ def _create_max_particles_in_flight_subelement(self, root):
             elem = ET.SubElement(root, "max_particles_in_flight")
             elem.text = str(self._max_particles_in_flight).lower()
 
+    def _create_max_events_subelement(self, root):
+        if self._max_particle_events is not None:
+            elem = ET.SubElement(root, "max_particle_events")
+            elem.text = str(self._max_particle_events).lower()
+
     def _create_material_cell_offsets_subelement(self, root):
         if self._material_cell_offsets is not None:
             elem = ET.SubElement(root, "material_cell_offsets")
@@ -1820,6 +1839,7 @@ def to_xml_element(self, mesh_memo=None):
         self._create_delayed_photon_scaling_subelement(element)
         self._create_event_based_subelement(element)
         self._create_max_particles_in_flight_subelement(element)
+        self._create_max_events_subelement(element)
         self._create_material_cell_offsets_subelement(element)
         self._create_log_grid_bins_subelement(element)
         self._create_write_initial_source_subelement(element)

src/particle.cpp

@@ -381,7 +381,7 @@ void Particle::event_revive_from_secondary()
 {
   // If particle has too many events, display warning and kill it
   ++n_event();
-  if (n_event() == MAX_EVENTS) {
+  if (n_event() == settings::max_events) {
     warning("Particle " + std::to_string(id()) +
             " underwent maximum number of events.");
     wgt() = 0.0;

src/settings.cpp

@@ -96,6 +96,7 @@ int32_t gen_per_batch {1};
 int64_t n_particles {-1};
 
 int64_t max_particles_in_flight {100000};
+int max_particle_events {1000000};
 
 ElectronTreatment electron_treatment {ElectronTreatment::TTB};
 array<double, 4> energy_cutoff {0.0, 1000.0, 0.0, 0.0};

tests/unit_tests/test_settings.py

@@ -27,8 +27,8 @@ 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, 'time_neutron': 1.0e-5, 
-                'time_photon': 1.0e-5, 'time_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}
     mesh = openmc.RegularMesh()
     mesh.lower_left = (-10., -10., -10.)
@@ -59,6 +59,8 @@ def test_export_to_xml(run_in_tmpdir):
     s.write_initial_source = True
     s.weight_window_checkpoints = {'surface': True, 'collision': False}
 
+    s.max_particle_events = 100
+
     # Make sure exporting XML works
     s.export_to_xml()
 
@@ -92,7 +94,7 @@ def test_export_to_xml(run_in_tmpdir):
     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,
-                        'time_neutron': 1.0e-5, 'time_photon': 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.]
@@ -128,3 +130,4 @@ def test_export_to_xml(run_in_tmpdir):
     assert vol.lower_left == (-10., -10., -10.)
     assert vol.upper_right == (10., 10., 10.)
     assert s.weight_window_checkpoints == {'surface': True, 'collision': False}
+    assert s.max_particle_events == 100