• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

openmc-dev / openmc / 27444237628

12 Jun 2026 09:31PM UTC coverage: 81.304% (+0.02%) from 81.281%
27444237628

Pull #3971

github

web-flow
Merge 9487507b3 into 02eb999af
Pull Request #3971: Delta tracking

18480 of 26797 branches covered (68.96%)

Branch coverage included in aggregate %.

592 of 644 new or added lines in 20 files covered. (91.93%)

31 existing lines in 3 files now uncovered.

59819 of 69507 relevant lines covered (86.06%)

49514105.74 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

76.46
/src/settings.cpp
1
#include "openmc/settings.h"
2
#include "openmc/random_ray/flat_source_domain.h"
3

4
#include <cmath>  // for ceil, pow
5
#include <limits> // for numeric_limits
6
#include <string>
7

8
#include <fmt/core.h>
9
#ifdef _OPENMP
10
#include <omp.h>
11
#endif
12

13
#include "openmc/capi.h"
14
#include "openmc/collision_track.h"
15
#include "openmc/constants.h"
16
#include "openmc/container_util.h"
17
#include "openmc/distribution.h"
18
#include "openmc/distribution_multi.h"
19
#include "openmc/distribution_spatial.h"
20
#include "openmc/eigenvalue.h"
21
#include "openmc/error.h"
22
#include "openmc/file_utils.h"
23
#include "openmc/majorant.h"
24
#include "openmc/mcpl_interface.h"
25
#include "openmc/mesh.h"
26
#include "openmc/message_passing.h"
27
#include "openmc/output.h"
28
#include "openmc/plot.h"
29
#include "openmc/random_lcg.h"
30
#include "openmc/random_ray/random_ray.h"
31
#include "openmc/reaction.h"
32
#include "openmc/simulation.h"
33
#include "openmc/source.h"
34
#include "openmc/string_utils.h"
35
#include "openmc/tallies/trigger.h"
36
#include "openmc/volume_calc.h"
37
#include "openmc/weight_windows.h"
38
#include "openmc/xml_interface.h"
39

40
namespace openmc {
41

42
//==============================================================================
43
// Global variables
44
//==============================================================================
45

46
namespace settings {
47

48
// Default values for boolean flags
49
bool assume_separate {false};
50
bool check_overlaps {false};
51
bool collision_track {false};
52
bool cmfd_run {false};
53
bool confidence_intervals {false};
54
bool create_delayed_neutrons {true};
55
bool create_fission_neutrons {true};
56
bool delayed_photon_scaling {true};
57
bool delta_tracking {false};
58
bool entropy_on {false};
59
bool event_based {false};
60
bool ifp_on {false};
61
bool legendre_to_tabular {true};
62
bool material_cell_offsets {true};
63
bool output_summary {true};
64
bool output_tallies {true};
65
bool particle_restart_run {false};
66
bool photon_transport {false};
67
bool atomic_relaxation {true};
68
bool reduce_tallies {true};
69
bool res_scat_on {false};
70
bool restart_run {false};
71
bool run_CE {true};
72
bool source_latest {false};
73
bool source_separate {false};
74
bool source_write {true};
75
bool source_mcpl_write {false};
76
bool surf_source_write {false};
77
bool surf_mcpl_write {false};
78
bool surf_source_read {false};
79
bool survival_biasing {false};
80
bool survival_normalization {false};
81
bool temperature_multipole {false};
82
bool trigger_on {false};
83
bool trigger_predict {false};
84
bool uniform_source_sampling {false};
85
bool ufs_on {false};
86
bool urr_ptables_on {true};
87
bool use_decay_photons {false};
88
bool use_shared_secondary_bank {false};
89
bool weight_windows_on {false};
90
bool weight_window_checkpoint_surface {false};
91
bool weight_window_checkpoint_collision {true};
92
bool write_all_tracks {false};
93
bool write_initial_source {false};
94

95
std::string path_cross_sections;
96
std::string path_input;
97
std::string path_output;
98
std::string path_particle_restart;
99
std::string path_sourcepoint;
100
std::string path_statepoint;
101
const char* path_statepoint_c {path_statepoint.c_str()};
102
std::string weight_windows_file;
103
std::string properties_file;
104

105
int32_t n_inactive {0};
106
int32_t max_lost_particles {10};
107
double rel_max_lost_particles {1.0e-6};
108
int32_t max_write_lost_particles {-1};
109
int32_t gen_per_batch {1};
110
int64_t n_particles {-1};
111

112
int64_t max_particles_in_flight {100000};
113
int max_particle_events {1000000};
114

115
ElectronTreatment electron_treatment {ElectronTreatment::TTB};
116
array<double, 4> energy_cutoff {0.0, 1000.0, 0.0, 0.0};
117
array<double, 4> time_cutoff {INFTY, INFTY, INFTY, INFTY};
118
int ifp_n_generation {-1};
119
IFPParameter ifp_parameter {IFPParameter::None};
120
int legendre_to_tabular_points {C_NONE};
121
int max_order {0};
122
int n_log_bins {8000};
123
int n_batches;
124
int n_max_batches;
125
int max_secondaries {10000};
126
int max_history_splits {10'000'000};
127
int max_tracks {1000};
128
ResScatMethod res_scat_method {ResScatMethod::rvs};
129
double res_scat_energy_min {0.01};
130
double res_scat_energy_max {1000.0};
131
vector<std::string> res_scat_nuclides;
132
RunMode run_mode {RunMode::UNSET};
133
SolverType solver_type {SolverType::MONTE_CARLO};
134
std::unordered_set<int> sourcepoint_batch;
135
std::unordered_set<int> statepoint_batch;
136
double source_rejection_fraction {0.05};
137
double free_gas_threshold {400.0};
138
std::unordered_set<int> source_write_surf_id;
139
CollisionTrackConfig collision_track_config {};
140
int64_t ssw_max_particles;
141
int64_t ssw_max_files;
142
int64_t ssw_cell_id {C_NONE};
143
SSWCellType ssw_cell_type {SSWCellType::None};
144
double surface_grazing_cutoff {0.001};
145
double surface_grazing_ratio {0.5};
146
TemperatureMethod temperature_method {TemperatureMethod::NEAREST};
147
double temperature_tolerance {10.0};
148
double temperature_default {293.6};
149
array<double, 2> temperature_range {0.0, 0.0};
150
int trace_batch;
151
int trace_gen;
152
int64_t trace_particle;
153
vector<array<int, 3>> track_identifiers;
154
int trigger_batch_interval {1};
155
int verbosity {-1};
156
double weight_cutoff {0.25};
157
double weight_survive {1.0};
158

159
} // namespace settings
160

161
//==============================================================================
162
// Functions
163
//==============================================================================
164

165
void get_run_parameters(pugi::xml_node node_base)
7,969✔
166
{
167
  using namespace settings;
7,969✔
168
  using namespace pugi;
7,969✔
169

170
  // Check number of particles
171
  if (!check_for_node(node_base, "particles")) {
7,969!
172
    fatal_error("Need to specify number of particles.");
×
173
  }
174

175
  // Get number of particles if it wasn't specified as a command-line argument
176
  if (n_particles == -1) {
7,969✔
177
    n_particles = std::stoll(get_node_value(node_base, "particles"));
7,958✔
178
  }
179

180
  // Get maximum number of in flight particles for event-based mode
181
  if (check_for_node(node_base, "max_particles_in_flight")) {
7,969!
182
    max_particles_in_flight =
×
183
      std::stoll(get_node_value(node_base, "max_particles_in_flight"));
×
184
  }
185

186
  // Get maximum number of events allowed per particle
187
  if (check_for_node(node_base, "max_particle_events")) {
7,969!
188
    max_particle_events =
×
189
      std::stoll(get_node_value(node_base, "max_particle_events"));
×
190
  }
191

192
  // Get number of basic batches
193
  if (check_for_node(node_base, "batches")) {
7,969!
194
    n_batches = std::stoi(get_node_value(node_base, "batches"));
7,969✔
195
  }
196
  if (!trigger_on)
7,969✔
197
    n_max_batches = n_batches;
7,828✔
198

199
  // Get max number of lost particles
200
  if (check_for_node(node_base, "max_lost_particles")) {
7,969✔
201
    max_lost_particles =
92✔
202
      std::stoi(get_node_value(node_base, "max_lost_particles"));
46✔
203
  }
204

205
  // Get relative number of lost particles
206
  if (check_for_node(node_base, "rel_max_lost_particles")) {
7,969!
207
    rel_max_lost_particles =
×
208
      std::stod(get_node_value(node_base, "rel_max_lost_particles"));
×
209
  }
210

211
  // Get relative number of lost particles
212
  if (check_for_node(node_base, "max_write_lost_particles")) {
7,969✔
213
    max_write_lost_particles =
30✔
214
      std::stoi(get_node_value(node_base, "max_write_lost_particles"));
15✔
215
  }
216

217
  // Get number of inactive batches
218
  if (run_mode == RunMode::EIGENVALUE ||
7,969✔
219
      solver_type == SolverType::RANDOM_RAY) {
3,188✔
220
    if (check_for_node(node_base, "inactive")) {
5,217✔
221
      n_inactive = std::stoi(get_node_value(node_base, "inactive"));
4,995✔
222
    }
223
    if (check_for_node(node_base, "generations_per_batch")) {
5,217✔
224
      gen_per_batch =
30✔
225
        std::stoi(get_node_value(node_base, "generations_per_batch"));
15✔
226
    }
227

228
    // Preallocate space for keff and entropy by generation
229
    int m = settings::n_max_batches * settings::gen_per_batch;
5,217✔
230
    simulation::k_generation.reserve(m);
5,217✔
231
    simulation::entropy.reserve(m);
5,217✔
232

233
    // Get the trigger information for keff
234
    if (check_for_node(node_base, "keff_trigger")) {
5,217✔
235
      xml_node node_keff_trigger = node_base.child("keff_trigger");
101✔
236

237
      if (check_for_node(node_keff_trigger, "type")) {
101!
238
        auto temp = get_node_value(node_keff_trigger, "type", true, true);
101✔
239
        if (temp == "std_dev") {
101!
240
          keff_trigger.metric = TriggerMetric::standard_deviation;
101✔
241
        } else if (temp == "variance") {
×
242
          keff_trigger.metric = TriggerMetric::variance;
×
243
        } else if (temp == "rel_err") {
×
244
          keff_trigger.metric = TriggerMetric::relative_error;
×
245
        } else {
246
          fatal_error("Unrecognized keff trigger type " + temp);
×
247
        }
248
      } else {
×
249
        fatal_error("Specify keff trigger type in settings XML");
×
250
      }
251

252
      if (check_for_node(node_keff_trigger, "threshold")) {
101!
253
        keff_trigger.threshold =
202✔
254
          std::stod(get_node_value(node_keff_trigger, "threshold"));
202✔
255
        if (keff_trigger.threshold <= 0) {
101!
256
          fatal_error("keff trigger threshold must be positive");
×
257
        }
258
      } else {
259
        fatal_error("Specify keff trigger threshold in settings XML");
×
260
      }
261
    }
262
  }
263

264
  // Random ray variables
265
  if (solver_type == SolverType::RANDOM_RAY) {
7,969✔
266
    xml_node random_ray_node = node_base.child("random_ray");
792✔
267
    if (check_for_node(random_ray_node, "distance_active")) {
792!
268
      RandomRay::distance_active_ =
1,584✔
269
        std::stod(get_node_value(random_ray_node, "distance_active"));
1,584✔
270
      if (RandomRay::distance_active_ <= 0.0) {
792!
271
        fatal_error("Random ray active distance must be greater than 0");
×
272
      }
273
    } else {
274
      fatal_error("Specify random ray active distance in settings XML");
×
275
    }
276
    if (check_for_node(random_ray_node, "distance_inactive")) {
792!
277
      RandomRay::distance_inactive_ =
1,584✔
278
        std::stod(get_node_value(random_ray_node, "distance_inactive"));
1,584✔
279
      if (RandomRay::distance_inactive_ < 0) {
792!
280
        fatal_error(
×
281
          "Random ray inactive distance must be greater than or equal to 0");
282
      }
283
    } else {
284
      fatal_error("Specify random ray inactive distance in settings XML");
×
285
    }
286
    if (check_for_node(random_ray_node, "ray_source")) {
792!
287
      xml_node ray_source_node = random_ray_node.child("ray_source");
792✔
288
      xml_node source_node = ray_source_node.child("source");
792✔
289
      // Get point to list of <source> elements and make sure there is at least
290
      // one
291
      RandomRay::ray_source_ = Source::create(source_node);
1,584✔
292
    } else {
293
      fatal_error("Specify random ray source in settings XML");
×
294
    }
295
    if (check_for_node(random_ray_node, "volume_estimator")) {
792✔
296
      std::string temp_str =
151✔
297
        get_node_value(random_ray_node, "volume_estimator", true, true);
151✔
298
      if (temp_str == "simulation_averaged") {
151✔
299
        FlatSourceDomain::volume_estimator_ =
30✔
300
          RandomRayVolumeEstimator::SIMULATION_AVERAGED;
301
      } else if (temp_str == "naive") {
121✔
302
        FlatSourceDomain::volume_estimator_ = RandomRayVolumeEstimator::NAIVE;
91✔
303
      } else if (temp_str == "hybrid") {
30!
304
        FlatSourceDomain::volume_estimator_ = RandomRayVolumeEstimator::HYBRID;
30✔
305
      } else {
306
        fatal_error("Unrecognized volume estimator: " + temp_str);
×
307
      }
308
    }
151✔
309
    if (check_for_node(random_ray_node, "source_shape")) {
792✔
310
      std::string temp_str =
450✔
311
        get_node_value(random_ray_node, "source_shape", true, true);
450✔
312
      if (temp_str == "flat") {
450✔
313
        RandomRay::source_shape_ = RandomRaySourceShape::FLAT;
75✔
314
      } else if (temp_str == "linear") {
375✔
315
        RandomRay::source_shape_ = RandomRaySourceShape::LINEAR;
330✔
316
      } else if (temp_str == "linear_xy") {
45!
317
        RandomRay::source_shape_ = RandomRaySourceShape::LINEAR_XY;
45✔
318
      } else {
319
        fatal_error("Unrecognized source shape: " + temp_str);
×
320
      }
321
    }
450✔
322
    if (check_for_node(random_ray_node, "volume_normalized_flux_tallies")) {
792✔
323
      FlatSourceDomain::volume_normalized_flux_tallies_ =
551✔
324
        get_node_value_bool(random_ray_node, "volume_normalized_flux_tallies");
551✔
325
    }
326
    if (check_for_node(random_ray_node, "adjoint")) {
792✔
327
      FlatSourceDomain::adjoint_ =
45✔
328
        get_node_value_bool(random_ray_node, "adjoint");
45✔
329
    }
330
    if (check_for_node(random_ray_node, "sample_method")) {
792✔
331
      std::string temp_str =
30✔
332
        get_node_value(random_ray_node, "sample_method", true, true);
30✔
333
      if (temp_str == "prng") {
30!
334
        RandomRay::sample_method_ = RandomRaySampleMethod::PRNG;
×
335
      } else if (temp_str == "halton") {
30✔
336
        RandomRay::sample_method_ = RandomRaySampleMethod::HALTON;
15✔
337
      } else if (temp_str == "s2") {
15!
338
        RandomRay::sample_method_ = RandomRaySampleMethod::S2;
15✔
339
      } else {
340
        fatal_error("Unrecognized sample method: " + temp_str);
×
341
      }
342
    }
30✔
343
    if (check_for_node(random_ray_node, "source_region_meshes")) {
792✔
344
      pugi::xml_node node_source_region_meshes =
346✔
345
        random_ray_node.child("source_region_meshes");
346✔
346
      for (pugi::xml_node node_mesh :
752✔
347
        node_source_region_meshes.children("mesh")) {
752✔
348
        int mesh_id = std::stoi(node_mesh.attribute("id").value());
812✔
349
        for (pugi::xml_node node_domain : node_mesh.children("domain")) {
812✔
350
          int domain_id = std::stoi(node_domain.attribute("id").value());
812✔
351
          std::string domain_type = node_domain.attribute("type").value();
406✔
352
          Source::DomainType type;
406✔
353
          if (domain_type == "material") {
406✔
354
            type = Source::DomainType::MATERIAL;
30✔
355
          } else if (domain_type == "cell") {
376✔
356
            type = Source::DomainType::CELL;
30✔
357
          } else if (domain_type == "universe") {
346!
358
            type = Source::DomainType::UNIVERSE;
346✔
359
          } else {
360
            throw std::runtime_error("Unknown domain type: " + domain_type);
×
361
          }
362
          FlatSourceDomain::mesh_domain_map_[mesh_id].emplace_back(
406✔
363
            type, domain_id);
364
        }
406✔
365
      }
366
    }
367
    if (check_for_node(random_ray_node, "diagonal_stabilization_rho")) {
792✔
368
      FlatSourceDomain::diagonal_stabilization_rho_ = std::stod(
15✔
369
        get_node_value(random_ray_node, "diagonal_stabilization_rho"));
15✔
370
      if (FlatSourceDomain::diagonal_stabilization_rho_ < 0.0 ||
15!
371
          FlatSourceDomain::diagonal_stabilization_rho_ > 1.0) {
372
        fatal_error("Random ray diagonal stabilization rho factor must be "
×
373
                    "between 0 and 1");
374
      }
375
    }
376
    if (check_for_node(random_ray_node, "adjoint_source")) {
792✔
377
      pugi::xml_node adj_source_node = random_ray_node.child("adjoint_source");
15✔
378
      for (pugi::xml_node source_node : adj_source_node.children("source")) {
30✔
379
        // Find any local adjoint sources
380
        model::adjoint_sources.push_back(Source::create(source_node));
30✔
381
      }
382
    }
383
  }
384
}
7,969✔
385

386
void read_settings_xml()
1,409✔
387
{
388
  using namespace settings;
1,409✔
389
  using namespace pugi;
1,409✔
390
  // Check if settings.xml exists
391
  std::string filename = settings::path_input + "settings.xml";
1,409✔
392
  if (!file_exists(filename)) {
1,409✔
393
    if (run_mode != RunMode::PLOTTING) {
22!
394
      fatal_error("Could not find any XML input files! In order to run OpenMC, "
×
395
                  "you first need a set of input files; at a minimum, this "
396
                  "includes settings.xml, geometry.xml, and materials.xml or a "
397
                  "single model XML file. Please consult the user's guide at "
398
                  "https://docs.openmc.org for further information.");
399
    } else {
400
      // The settings.xml file is optional if we just want to make a plot.
401
      return;
22✔
402
    }
403
  }
404

405
  // Parse settings.xml file
406
  xml_document doc;
1,387✔
407
  auto result = doc.load_file(filename.c_str());
1,387✔
408
  if (!result) {
1,387!
409
    fatal_error("Error processing settings.xml file.");
×
410
  }
411

412
  // Get root element
413
  xml_node root = doc.document_element();
1,387✔
414

415
  // Verbosity
416
  if (check_for_node(root, "verbosity") && verbosity == -1) {
1,387!
417
    verbosity = std::stoi(get_node_value(root, "verbosity"));
428✔
418
  } else if (verbosity == -1) {
1,173!
419
    verbosity = 7;
1,173✔
420
  }
421

422
  // To this point, we haven't displayed any output since we didn't know what
423
  // the verbosity is. Now that we checked for it, show the title if necessary
424
  if (mpi::master) {
1,387✔
425
    if (verbosity >= 2)
1,199✔
426
      title();
993✔
427
  }
428

429
  write_message("Reading settings XML file...", 5);
1,387✔
430

431
  read_settings_xml(root);
1,387✔
432
}
1,399✔
433

434
void read_settings_xml(pugi::xml_node root)
8,920✔
435
{
436
  using namespace settings;
8,920✔
437
  using namespace pugi;
8,920✔
438

439
  // Find if a multi-group or continuous-energy simulation is desired
440
  if (check_for_node(root, "energy_mode")) {
8,920✔
441
    std::string temp_str = get_node_value(root, "energy_mode", true, true);
1,357✔
442
    if (temp_str == "mg" || temp_str == "multi-group") {
2,714!
443
      run_CE = false;
1,357✔
444
    } else if (temp_str == "ce" || temp_str == "continuous-energy") {
×
445
      run_CE = true;
×
446
    }
447
  }
1,357✔
448

449
  // Check for user meshes and allocate
450
  read_meshes(root);
8,920✔
451

452
  // Look for deprecated cross_sections.xml file in settings.xml
453
  if (check_for_node(root, "cross_sections")) {
8,920!
454
    warning(
×
455
      "Setting cross_sections in settings.xml has been deprecated."
456
      " The cross_sections are now set in materials.xml and the "
457
      "cross_sections input to materials.xml and the OPENMC_CROSS_SECTIONS"
458
      " environment variable will take precendent over setting "
459
      "cross_sections in settings.xml.");
460
    path_cross_sections = get_node_value(root, "cross_sections");
×
461
  }
462

463
  if (!run_CE) {
8,920✔
464
    // Scattering Treatments
465
    if (check_for_node(root, "max_order")) {
1,357✔
466
      max_order = std::stoi(get_node_value(root, "max_order"));
30✔
467
    } else {
468
      // Set to default of largest int - 1, which means to use whatever is
469
      // contained in library. This is largest int - 1 because for legendre
470
      // scattering, a value of 1 is added to the order; adding 1 to the largest
471
      // int gets you the largest negative integer, which is not what we want.
472
      max_order = std::numeric_limits<int>::max() - 1;
1,342✔
473
    }
474
  }
475

476
  // Check for a trigger node and get trigger information
477
  if (check_for_node(root, "trigger")) {
8,920✔
478
    xml_node node_trigger = root.child("trigger");
156✔
479

480
    // Check if trigger(s) are to be turned on
481
    trigger_on = get_node_value_bool(node_trigger, "active");
156✔
482

483
    if (trigger_on) {
156✔
484
      if (check_for_node(node_trigger, "max_batches")) {
141!
485
        n_max_batches = std::stoi(get_node_value(node_trigger, "max_batches"));
282✔
486
      } else {
487
        fatal_error("<max_batches> must be specified with triggers");
×
488
      }
489

490
      // Get the batch interval to check triggers
491
      if (!check_for_node(node_trigger, "batch_interval")) {
141✔
492
        trigger_predict = true;
15✔
493
      } else {
494
        trigger_batch_interval =
252✔
495
          std::stoi(get_node_value(node_trigger, "batch_interval"));
252✔
496
        if (trigger_batch_interval <= 0) {
126!
497
          fatal_error("Trigger batch interval must be greater than zero");
×
498
        }
499
      }
500
    }
501
  }
502

503
  // Check run mode if it hasn't been set from the command line
504
  xml_node node_mode;
8,920✔
505
  if (run_mode == RunMode::UNSET) {
8,920✔
506
    if (check_for_node(root, "run_mode")) {
8,001✔
507
      std::string temp_str = get_node_value(root, "run_mode", true, true);
7,971✔
508
      if (temp_str == "eigenvalue") {
7,971✔
509
        run_mode = RunMode::EIGENVALUE;
4,751✔
510
      } else if (temp_str == "fixed source") {
3,220✔
511
        run_mode = RunMode::FIXED_SOURCE;
3,188✔
512
      } else if (temp_str == "plot") {
32!
513
        run_mode = RunMode::PLOTTING;
×
514
      } else if (temp_str == "particle restart") {
32!
515
        run_mode = RunMode::PARTICLE;
×
516
      } else if (temp_str == "volume") {
32!
517
        run_mode = RunMode::VOLUME;
32✔
518
      } else {
519
        fatal_error("Unrecognized run mode: " + temp_str);
×
520
      }
521

522
      // Assume XML specifies <particles>, <batches>, etc. directly
523
      node_mode = root;
7,971✔
524
    } else {
7,971✔
525
      warning("<run_mode> should be specified.");
30✔
526

527
      // Make sure that either eigenvalue or fixed source was specified
528
      node_mode = root.child("eigenvalue");
30✔
529
      if (node_mode) {
30!
530
        run_mode = RunMode::EIGENVALUE;
30✔
531
      } else {
532
        node_mode = root.child("fixed_source");
×
533
        if (node_mode) {
×
534
          run_mode = RunMode::FIXED_SOURCE;
×
535
        } else {
536
          fatal_error("<eigenvalue> or <fixed_source> not specified.");
×
537
        }
538
      }
539
    }
540
  }
541

542
  // Check solver type
543
  if (check_for_node(root, "random_ray")) {
8,920✔
544
    solver_type = SolverType::RANDOM_RAY;
792✔
545
    if (run_CE)
792!
546
      fatal_error("multi-group energy mode must be specified in settings XML "
×
547
                  "when using the random ray solver.");
548
  }
549

550
  if (run_mode == RunMode::EIGENVALUE || run_mode == RunMode::FIXED_SOURCE) {
8,920✔
551
    // Read run parameters
552
    get_run_parameters(node_mode);
7,969✔
553

554
    // Check number of active batches, inactive batches, max lost particles and
555
    // particles
556
    if (n_batches <= n_inactive) {
7,969!
557
      fatal_error("Number of active batches must be greater than zero.");
×
558
    } else if (n_inactive < 0) {
7,969!
559
      fatal_error("Number of inactive batches must be non-negative.");
×
560
    } else if (n_particles <= 0) {
7,969!
561
      fatal_error("Number of particles must be greater than zero.");
×
562
    } else if (max_lost_particles <= 0) {
7,969!
563
      fatal_error("Number of max lost particles must be greater than zero.");
×
564
    } else if (rel_max_lost_particles <= 0.0 || rel_max_lost_particles >= 1.0) {
7,969!
565
      fatal_error("Relative max lost particles must be between zero and one.");
×
566
    }
567

568
    // Check for user value for the number of generation of the Iterated Fission
569
    // Probability (IFP) method
570
    if (check_for_node(root, "ifp_n_generation")) {
7,969✔
571
      ifp_n_generation = std::stoi(get_node_value(root, "ifp_n_generation"));
166✔
572
      if (ifp_n_generation <= 0) {
83!
573
        fatal_error("'ifp_n_generation' must be greater than 0.");
×
574
      }
575
      // Avoid tallying 0 if IFP logs are not complete when active cycles start
576
      if (ifp_n_generation > n_inactive) {
83✔
577
        fatal_error("'ifp_n_generation' must be lower than or equal to the "
9✔
578
                    "number of inactive cycles.");
579
      }
580
    }
581
  }
582

583
  // Copy plotting random number seed if specified
584
  if (check_for_node(root, "plot_seed")) {
8,911!
585
    auto seed = std::stoll(get_node_value(root, "plot_seed"));
×
586
    model::plotter_seed = seed;
×
587
  }
588

589
  // Copy random number seed if specified
590
  if (check_for_node(root, "seed")) {
8,911✔
591
    auto seed = std::stoll(get_node_value(root, "seed"));
1,178✔
592
    openmc_set_seed(seed);
589✔
593
  }
594

595
  // Copy random number stride if specified
596
  if (check_for_node(root, "stride")) {
8,911✔
597
    auto stride = std::stoull(get_node_value(root, "stride"));
30✔
598
    openmc_set_stride(stride);
15✔
599
  }
600

601
  // Check for electron treatment
602
  if (check_for_node(root, "electron_treatment")) {
8,911✔
603
    auto temp_str = get_node_value(root, "electron_treatment", true, true);
82✔
604
    if (temp_str == "led") {
82✔
605
      electron_treatment = ElectronTreatment::LED;
26✔
606
    } else if (temp_str == "ttb") {
56!
607
      electron_treatment = ElectronTreatment::TTB;
56✔
608
    } else {
609
      fatal_error("Unrecognized electron treatment: " + temp_str + ".");
×
610
    }
611
  }
82✔
612

613
  // Check for photon transport
614
  if (check_for_node(root, "photon_transport")) {
8,911✔
615
    photon_transport = get_node_value_bool(root, "photon_transport");
522✔
616

617
    if (!run_CE && photon_transport) {
522!
618
      fatal_error("Photon transport is not currently supported in "
×
619
                  "multigroup mode");
620
    }
621
  }
622

623
  // Check for atomic relaxation
624
  if (check_for_node(root, "atomic_relaxation")) {
8,911✔
625
    atomic_relaxation = get_node_value_bool(root, "atomic_relaxation");
15✔
626
  }
627

628
  // Number of bins for logarithmic grid
629
  if (check_for_node(root, "log_grid_bins")) {
8,911✔
630
    n_log_bins = std::stoi(get_node_value(root, "log_grid_bins"));
30✔
631
    if (n_log_bins < 1) {
15!
632
      fatal_error("Number of bins for logarithmic grid must be greater "
×
633
                  "than zero.");
634
    }
635
  }
636

637
  // Number of OpenMP threads
638
  if (check_for_node(root, "threads")) {
8,911!
639
    if (mpi::master)
×
640
      warning("The <threads> element has been deprecated. Use "
×
641
              "the OMP_NUM_THREADS environment variable to set the number of "
642
              "threads.");
643
  }
644

645
  // ==========================================================================
646
  // EXTERNAL SOURCE
647

648
  // Get point to list of <source> elements and make sure there is at least one
649
  for (pugi::xml_node node : root.children("source")) {
17,096✔
650
    model::external_sources.push_back(Source::create(node));
16,380✔
651
  }
652

653
  // Check if the user has specified to read surface source
654
  if (check_for_node(root, "surf_source_read")) {
8,901✔
655
    surf_source_read = true;
30✔
656
    // Get surface source read node
657
    xml_node node_ssr = root.child("surf_source_read");
30✔
658

659
    std::string path = "surface_source.h5";
30✔
660
    // Check if the user has specified different file for surface source reading
661
    if (check_for_node(node_ssr, "path")) {
30!
662
      path = get_node_value(node_ssr, "path", false, true);
30✔
663
    }
664
    model::external_sources.push_back(make_unique<FileSource>(path));
30✔
665
  }
30✔
666

667
  // If no source specified, default to isotropic point source at origin with
668
  // Watt spectrum. No default source is needed in random ray mode.
669
  if (model::external_sources.empty() &&
8,901✔
670
      settings::solver_type != SolverType::RANDOM_RAY) {
2,409✔
671
    double T[] {0.0};
2,263✔
672
    double p[] {1.0};
2,263✔
673
    model::external_sources.push_back(make_unique<IndependentSource>(
2,263✔
674
      UPtrSpace {new SpatialPoint({0.0, 0.0, 0.0})},
4,526✔
675
      UPtrAngle {new Isotropic()}, UPtrDist {new Watt(0.988e6, 2.249e-6)},
4,526✔
676
      UPtrDist {new Discrete(T, p, 1)}));
4,526✔
677
  }
678

679
  // Build probability mass function for sampling external sources
680
  vector<double> source_strengths;
8,901✔
681
  for (auto& s : model::external_sources) {
19,379✔
682
    source_strengths.push_back(s->strength());
10,478✔
683
  }
684
  model::external_sources_probability.assign(source_strengths);
8,901✔
685

686
  // Check if we want to write out source
687
  if (check_for_node(root, "write_initial_source")) {
8,901!
688
    write_initial_source = get_node_value_bool(root, "write_initial_source");
×
689
  }
690

691
  // Get relative number of lost particles
692
  if (check_for_node(root, "source_rejection_fraction")) {
8,901✔
693
    source_rejection_fraction =
14✔
694
      std::stod(get_node_value(root, "source_rejection_fraction"));
14!
695
  }
696

697
  if (check_for_node(root, "free_gas_threshold")) {
8,901!
698
    free_gas_threshold = std::stod(get_node_value(root, "free_gas_threshold"));
×
699
  }
700

701
  // Surface grazing
702
  if (check_for_node(root, "surface_grazing_cutoff"))
8,901!
703
    surface_grazing_cutoff =
×
704
      std::stod(get_node_value(root, "surface_grazing_cutoff"));
×
705
  if (check_for_node(root, "surface_grazing_ratio"))
8,901!
706
    surface_grazing_ratio =
×
707
      std::stod(get_node_value(root, "surface_grazing_ratio"));
×
708

709
  // Survival biasing
710
  if (check_for_node(root, "survival_biasing")) {
8,901✔
711
    survival_biasing = get_node_value_bool(root, "survival_biasing");
205✔
712
  }
713

714
  // Probability tables
715
  if (check_for_node(root, "ptables")) {
8,901✔
716
    urr_ptables_on = get_node_value_bool(root, "ptables");
15✔
717
  }
718

719
  // Cutoffs
720
  if (check_for_node(root, "cutoff")) {
8,901✔
721
    xml_node node_cutoff = root.child("cutoff");
138✔
722
    if (check_for_node(node_cutoff, "weight")) {
138✔
723
      weight_cutoff = std::stod(get_node_value(node_cutoff, "weight"));
30✔
724
    }
725
    if (check_for_node(node_cutoff, "weight_avg")) {
138✔
726
      weight_survive = std::stod(get_node_value(node_cutoff, "weight_avg"));
30✔
727
    }
728
    if (check_for_node(node_cutoff, "survival_normalization")) {
138!
729
      survival_normalization =
×
730
        get_node_value_bool(node_cutoff, "survival_normalization");
×
731
    }
732
    if (check_for_node(node_cutoff, "energy_neutron")) {
138✔
733
      energy_cutoff[0] =
15✔
734
        std::stod(get_node_value(node_cutoff, "energy_neutron"));
30✔
735
    } else if (check_for_node(node_cutoff, "energy")) {
123!
736
      warning("The use of an <energy> cutoff is deprecated and should "
×
737
              "be replaced by <energy_neutron>.");
738
      energy_cutoff[0] = std::stod(get_node_value(node_cutoff, "energy"));
×
739
    }
740
    if (check_for_node(node_cutoff, "energy_photon")) {
138✔
741
      energy_cutoff[1] =
82✔
742
        std::stod(get_node_value(node_cutoff, "energy_photon"));
164✔
743
    }
744
    if (check_for_node(node_cutoff, "energy_electron")) {
138!
745
      energy_cutoff[2] =
×
746
        std::stof(get_node_value(node_cutoff, "energy_electron"));
×
747
    }
748
    if (check_for_node(node_cutoff, "energy_positron")) {
138!
749
      energy_cutoff[3] =
×
750
        std::stod(get_node_value(node_cutoff, "energy_positron"));
×
751
    }
752
    if (check_for_node(node_cutoff, "time_neutron")) {
138✔
753
      time_cutoff[0] = std::stod(get_node_value(node_cutoff, "time_neutron"));
26✔
754
    }
755
    if (check_for_node(node_cutoff, "time_photon")) {
138!
756
      time_cutoff[1] = std::stod(get_node_value(node_cutoff, "time_photon"));
×
757
    }
758
    if (check_for_node(node_cutoff, "time_electron")) {
138!
759
      time_cutoff[2] = std::stod(get_node_value(node_cutoff, "time_electron"));
×
760
    }
761
    if (check_for_node(node_cutoff, "time_positron")) {
138!
762
      time_cutoff[3] = std::stod(get_node_value(node_cutoff, "time_positron"));
×
763
    }
764
  }
765

766
  // read properties from file
767
  if (check_for_node(root, "properties_file")) {
8,901✔
768
    properties_file = get_node_value(root, "properties_file");
11✔
769
    if (!file_exists(properties_file)) {
11!
770
      fatal_error(fmt::format("File '{}' does not exist.", properties_file));
×
771
    }
772
  }
773

774
  // Particle trace
775
  if (check_for_node(root, "trace")) {
8,901✔
776
    auto temp = get_node_array<int64_t>(root, "trace");
15✔
777
    if (temp.size() != 3) {
15!
778
      fatal_error("Must provide 3 integers for <trace> that specify the "
×
779
                  "batch, generation, and particle number.");
780
    }
781
    trace_batch = temp.at(0);
15✔
782
    trace_gen = temp.at(1);
15✔
783
    trace_particle = temp.at(2);
15✔
784
  }
15✔
785

786
  // Particle tracks
787
  if (check_for_node(root, "track")) {
8,901✔
788
    // Get values and make sure there are three per particle
789
    auto temp = get_node_array<int>(root, "track");
45✔
790
    if (temp.size() % 3 != 0) {
45!
791
      fatal_error(
×
792
        "Number of integers specified in 'track' is not "
793
        "divisible by 3.  Please provide 3 integers per particle to be "
794
        "tracked.");
795
    }
796

797
    // Reshape into track_identifiers
798
    int n_tracks = temp.size() / 3;
45✔
799
    for (int i = 0; i < n_tracks; ++i) {
180✔
800
      track_identifiers.push_back(
135✔
801
        {temp[3 * i], temp[3 * i + 1], temp[3 * i + 2]});
135✔
802
    }
803
  }
45✔
804

805
  // Shannon entropy
806
  if (solver_type == SolverType::RANDOM_RAY) {
8,901✔
807
    if (check_for_node(root, "entropy_mesh")) {
792!
808
      fatal_error("Random ray uses FSRs to compute the Shannon entropy. "
×
809
                  "No user-defined entropy mesh is supported.");
810
    }
811
    entropy_on = true;
792✔
812
  } else if (solver_type == SolverType::MONTE_CARLO) {
8,109!
813
    if (check_for_node(root, "entropy_mesh")) {
8,109✔
814
      int temp = std::stoi(get_node_value(root, "entropy_mesh"));
668✔
815
      if (model::mesh_map.find(temp) == model::mesh_map.end()) {
334!
816
        fatal_error(fmt::format(
×
817
          "Mesh {} specified for Shannon entropy does not exist.", temp));
818
      }
819

820
      auto* m = dynamic_cast<RegularMesh*>(
334!
821
        model::meshes[model::mesh_map.at(temp)].get());
334!
822
      if (!m)
334!
823
        fatal_error("Only regular meshes can be used as an entropy mesh");
×
824
      simulation::entropy_mesh = m;
334✔
825

826
      // Turn on Shannon entropy calculation
827
      entropy_on = true;
334✔
828

829
    } else if (check_for_node(root, "entropy")) {
7,775!
830
      fatal_error(
×
831
        "Specifying a Shannon entropy mesh via the <entropy> element "
832
        "is deprecated. Please create a mesh using <mesh> and then reference "
833
        "it by specifying its ID in an <entropy_mesh> element.");
834
    }
835
  }
836
  // Uniform fission source weighting mesh
837
  if (check_for_node(root, "ufs_mesh")) {
8,901✔
838
    auto temp = std::stoi(get_node_value(root, "ufs_mesh"));
30✔
839
    if (model::mesh_map.find(temp) == model::mesh_map.end()) {
15!
840
      fatal_error(fmt::format("Mesh {} specified for uniform fission site "
×
841
                              "method does not exist.",
842
        temp));
843
    }
844

845
    auto* m =
15✔
846
      dynamic_cast<RegularMesh*>(model::meshes[model::mesh_map.at(temp)].get());
15!
847
    if (!m)
15!
848
      fatal_error("Only regular meshes can be used as a UFS mesh");
×
849
    simulation::ufs_mesh = m;
15✔
850

851
    // Turn on uniform fission source weighting
852
    ufs_on = true;
15✔
853

854
  } else if (check_for_node(root, "uniform_fs")) {
8,886!
855
    fatal_error(
×
856
      "Specifying a UFS mesh via the <uniform_fs> element "
857
      "is deprecated. Please create a mesh using <mesh> and then reference "
858
      "it by specifying its ID in a <ufs_mesh> element.");
859
  }
860

861
  // Check if the user has specified to write state points
862
  if (check_for_node(root, "state_point")) {
8,901✔
863

864
    // Get pointer to state_point node
865
    auto node_sp = root.child("state_point");
160✔
866

867
    // Determine number of batches at which to store state points
868
    if (check_for_node(node_sp, "batches")) {
160!
869
      // User gave specific batches to write state points
870
      auto temp = get_node_array<int>(node_sp, "batches");
160✔
871
      for (const auto& b : temp) {
491✔
872
        statepoint_batch.insert(b);
331✔
873
      }
874
    } else {
160✔
875
      // If neither were specified, write state point at last batch
876
      statepoint_batch.insert(n_batches);
×
877
    }
878
  } else {
879
    // If no <state_point> tag was present, by default write state point at
880
    // last batch only
881
    statepoint_batch.insert(n_batches);
8,741✔
882
  }
883

884
  // Check if the user has specified to write source points
885
  if (check_for_node(root, "source_point")) {
8,901✔
886
    // Get source_point node
887
    xml_node node_sp = root.child("source_point");
101✔
888

889
    // Determine batches at which to store source points
890
    if (check_for_node(node_sp, "batches")) {
101✔
891
      // User gave specific batches to write source points
892
      auto temp = get_node_array<int>(node_sp, "batches");
45✔
893
      for (const auto& b : temp) {
120✔
894
        sourcepoint_batch.insert(b);
75✔
895
      }
896
    } else {
45✔
897
      // If neither were specified, write source points with state points
898
      sourcepoint_batch = statepoint_batch;
56!
899
    }
900

901
    // Check if the user has specified to write binary source file
902
    if (check_for_node(node_sp, "separate")) {
101✔
903
      source_separate = get_node_value_bool(node_sp, "separate");
71✔
904
    }
905
    if (check_for_node(node_sp, "write")) {
101!
906
      source_write = get_node_value_bool(node_sp, "write");
×
907
    }
908
    if (check_for_node(node_sp, "mcpl")) {
101✔
909
      source_mcpl_write = get_node_value_bool(node_sp, "mcpl");
26✔
910
    }
911
    if (check_for_node(node_sp, "overwrite_latest")) {
101✔
912
      source_latest = get_node_value_bool(node_sp, "overwrite_latest");
15✔
913
      source_separate = source_latest;
15✔
914
    }
915
  } else {
916
    // If no <source_point> tag was present, by default we keep source bank in
917
    // statepoint file and write it out at statepoints intervals
918
    source_separate = false;
8,800✔
919
    sourcepoint_batch = statepoint_batch;
8,800!
920
  }
921

922
  // Check is the user specified to convert strength to statistical weight
923
  if (check_for_node(root, "uniform_source_sampling")) {
8,901✔
924
    uniform_source_sampling =
55✔
925
      get_node_value_bool(root, "uniform_source_sampling");
55✔
926
  }
927

928
  // Check if the user has specified to write surface source
929
  if (check_for_node(root, "surf_source_write")) {
8,901✔
930
    surf_source_write = true;
412✔
931
    // Get surface source write node
932
    xml_node node_ssw = root.child("surf_source_write");
412✔
933

934
    // Determine surface ids at which crossing particles are to be banked.
935
    // If no surfaces are specified, all surfaces in the model will be used
936
    // to bank source points.
937
    if (check_for_node(node_ssw, "surface_ids")) {
412✔
938
      auto temp = get_node_array<int>(node_ssw, "surface_ids");
202✔
939
      for (const auto& b : temp) {
994✔
940
        source_write_surf_id.insert(b);
792✔
941
      }
942
    }
202✔
943

944
    // Get maximum number of particles to be banked per surface
945
    if (check_for_node(node_ssw, "max_particles")) {
412✔
946
      ssw_max_particles = std::stoll(get_node_value(node_ssw, "max_particles"));
806✔
947
    } else {
948
      fatal_error("A maximum number of particles needs to be specified "
9✔
949
                  "using the 'max_particles' parameter to store surface "
950
                  "source points.");
951
    }
952

953
    // Get maximum number of surface source files to be created
954
    if (check_for_node(node_ssw, "max_source_files")) {
403✔
955
      ssw_max_files = std::stoll(get_node_value(node_ssw, "max_source_files"));
66✔
956
    } else {
957
      ssw_max_files = 1;
370✔
958
    }
959

960
    if (check_for_node(node_ssw, "mcpl")) {
403✔
961
      surf_mcpl_write = get_node_value_bool(node_ssw, "mcpl");
11✔
962
    }
963
    // Get cell information
964
    if (check_for_node(node_ssw, "cell")) {
403✔
965
      ssw_cell_id = std::stoll(get_node_value(node_ssw, "cell"));
208✔
966
      ssw_cell_type = SSWCellType::Both;
104✔
967
    }
968
    if (check_for_node(node_ssw, "cellfrom")) {
403✔
969
      if (ssw_cell_id != C_NONE) {
90✔
970
        fatal_error(
18✔
971
          "'cell', 'cellfrom' and 'cellto' cannot be used at the same time.");
972
      }
973
      ssw_cell_id = std::stoll(get_node_value(node_ssw, "cellfrom"));
144✔
974
      ssw_cell_type = SSWCellType::From;
72✔
975
    }
976
    if (check_for_node(node_ssw, "cellto")) {
385✔
977
      if (ssw_cell_id != C_NONE) {
71✔
978
        fatal_error(
18✔
979
          "'cell', 'cellfrom' and 'cellto' cannot be used at the same time.");
980
      }
981
      ssw_cell_id = std::stoll(get_node_value(node_ssw, "cellto"));
106✔
982
      ssw_cell_type = SSWCellType::To;
53✔
983
    }
984
  }
985

986
  // Check if the user has specified to write specific collisions
987
  if (check_for_node(root, "collision_track")) {
8,856✔
988
    settings::collision_track = true;
160✔
989
    // Get collision track node
990
    xml_node node_ct = root.child("collision_track");
160✔
991
    collision_track_config = CollisionTrackConfig {};
160✔
992

993
    // Determine cell ids at which crossing particles are to be banked
994
    if (check_for_node(node_ct, "cell_ids")) {
160✔
995
      auto temp = get_node_array<int>(node_ct, "cell_ids");
89✔
996
      for (const auto& b : temp) {
237✔
997
        collision_track_config.cell_ids.insert(b);
148✔
998
      }
999
    }
89✔
1000
    if (check_for_node(node_ct, "reactions")) {
160✔
1001
      auto temp = get_node_array<std::string>(node_ct, "reactions");
63✔
1002
      for (const auto& b : temp) {
171✔
1003
        int reaction_int = reaction_mt(b);
108✔
1004
        if (reaction_int > 0) {
108!
1005
          collision_track_config.mt_numbers.insert(reaction_int);
108✔
1006
        }
1007
      }
1008
    }
63✔
1009
    if (check_for_node(node_ct, "universe_ids")) {
160✔
1010
      auto temp = get_node_array<int>(node_ct, "universe_ids");
30✔
1011
      for (const auto& b : temp) {
60✔
1012
        collision_track_config.universe_ids.insert(b);
30✔
1013
      }
1014
    }
30✔
1015
    if (check_for_node(node_ct, "material_ids")) {
160✔
1016
      auto temp = get_node_array<int>(node_ct, "material_ids");
30✔
1017
      for (const auto& b : temp) {
75✔
1018
        collision_track_config.material_ids.insert(b);
45✔
1019
      }
1020
    }
30✔
1021
    if (check_for_node(node_ct, "nuclides")) {
160✔
1022
      auto temp = get_node_array<std::string>(node_ct, "nuclides");
30✔
1023
      for (const auto& b : temp) {
120✔
1024
        collision_track_config.nuclides.insert(b);
90✔
1025
      }
1026
    }
30✔
1027
    if (check_for_node(node_ct, "deposited_E_threshold")) {
160✔
1028
      collision_track_config.deposited_energy_threshold =
60✔
1029
        std::stod(get_node_value(node_ct, "deposited_E_threshold"));
60✔
1030
    }
1031
    // Get maximum number of particles to be banked per collision
1032
    if (check_for_node(node_ct, "max_collisions")) {
160!
1033
      collision_track_config.max_collisions =
320✔
1034
        std::stoll(get_node_value(node_ct, "max_collisions"));
320✔
1035
    } else {
1036
      warning("A maximum number of collisions needs to be specified. "
×
1037
              "By default the code sets 'max_collisions' parameter equals to "
1038
              "1000.");
1039
    }
1040
    // Get maximum number of collision_track files to be created
1041
    if (check_for_node(node_ct, "max_collision_track_files")) {
160!
1042
      collision_track_config.max_files =
×
1043
        std::stoll(get_node_value(node_ct, "max_collision_track_files"));
×
1044
    }
1045
    if (check_for_node(node_ct, "mcpl")) {
160✔
1046
      collision_track_config.mcpl_write = get_node_value_bool(node_ct, "mcpl");
22✔
1047
    }
1048
  }
1049

1050
  // If source is not separate and is to be written out in the statepoint
1051
  // file, make sure that the sourcepoint batch numbers are contained in the
1052
  // statepoint list
1053
  if (!source_separate) {
8,856✔
1054
    for (const auto& b : sourcepoint_batch) {
17,666✔
1055
      if (!contains(statepoint_batch, b)) {
17,792!
1056
        fatal_error(
×
1057
          "Sourcepoint batches are not a subset of statepoint batches.");
1058
      }
1059
    }
1060
  }
1061

1062
  // Check if the user has specified to not reduce tallies at the end of every
1063
  // batch
1064
  if (check_for_node(root, "no_reduce")) {
8,856✔
1065
    reduce_tallies = !get_node_value_bool(root, "no_reduce");
30✔
1066
  }
1067

1068
  // Check if the user has specified to use confidence intervals for
1069
  // uncertainties rather than standard deviations
1070
  if (check_for_node(root, "confidence_intervals")) {
8,856✔
1071
    confidence_intervals = get_node_value_bool(root, "confidence_intervals");
15✔
1072
  }
1073

1074
  // Check for output options
1075
  if (check_for_node(root, "output")) {
8,856✔
1076
    // Get pointer to output node
1077
    pugi::xml_node node_output = root.child("output");
790✔
1078

1079
    // Check for summary option
1080
    if (check_for_node(node_output, "summary")) {
790✔
1081
      output_summary = get_node_value_bool(node_output, "summary");
764✔
1082
    }
1083

1084
    // Check for ASCII tallies output option
1085
    if (check_for_node(node_output, "tallies")) {
790✔
1086
      output_tallies = get_node_value_bool(node_output, "tallies");
349✔
1087
    }
1088

1089
    // Set output directory if a path has been specified
1090
    if (check_for_node(node_output, "path")) {
790!
1091
      path_output = get_node_value(node_output, "path");
×
1092
      if (!ends_with(path_output, "/")) {
×
1093
        path_output += "/";
790!
1094
      }
1095
    }
1096
  }
1097

1098
  // Resonance scattering parameters
1099
  if (check_for_node(root, "resonance_scattering")) {
8,856✔
1100
    xml_node node_res_scat = root.child("resonance_scattering");
15✔
1101

1102
    // See if resonance scattering is enabled
1103
    if (check_for_node(node_res_scat, "enable")) {
15!
1104
      res_scat_on = get_node_value_bool(node_res_scat, "enable");
15✔
1105
    } else {
1106
      res_scat_on = true;
×
1107
    }
1108

1109
    // Determine what method is used
1110
    if (check_for_node(node_res_scat, "method")) {
15!
1111
      auto temp = get_node_value(node_res_scat, "method", true, true);
15✔
1112
      if (temp == "rvs") {
15!
1113
        res_scat_method = ResScatMethod::rvs;
15✔
1114
      } else if (temp == "dbrc") {
×
1115
        res_scat_method = ResScatMethod::dbrc;
×
1116
      } else {
1117
        fatal_error(
×
1118
          "Unrecognized resonance elastic scattering method: " + temp + ".");
×
1119
      }
1120
    }
15✔
1121

1122
    // Minimum energy for resonance scattering
1123
    if (check_for_node(node_res_scat, "energy_min")) {
15!
1124
      res_scat_energy_min =
30✔
1125
        std::stod(get_node_value(node_res_scat, "energy_min"));
30✔
1126
    }
1127
    if (res_scat_energy_min < 0.0) {
15!
1128
      fatal_error("Lower resonance scattering energy bound is negative");
×
1129
    }
1130

1131
    // Maximum energy for resonance scattering
1132
    if (check_for_node(node_res_scat, "energy_max")) {
15!
1133
      res_scat_energy_max =
30✔
1134
        std::stod(get_node_value(node_res_scat, "energy_max"));
30✔
1135
    }
1136
    if (res_scat_energy_max < res_scat_energy_min) {
15!
1137
      fatal_error("Upper resonance scattering energy bound is below the "
×
1138
                  "lower resonance scattering energy bound.");
1139
    }
1140

1141
    // Get resonance scattering nuclides
1142
    if (check_for_node(node_res_scat, "nuclides")) {
15!
1143
      res_scat_nuclides =
15✔
1144
        get_node_array<std::string>(node_res_scat, "nuclides");
30✔
1145
    }
1146
  }
1147

1148
  // Get volume calculations
1149
  for (pugi::xml_node node_vol : root.children("volume_calc")) {
9,164✔
1150
    model::volume_calcs.emplace_back(node_vol);
308✔
1151
  }
1152

1153
  // Get temperature settings
1154
  if (check_for_node(root, "temperature_default")) {
8,856✔
1155
    temperature_default =
342✔
1156
      std::stod(get_node_value(root, "temperature_default"));
342✔
1157
  }
1158
  if (check_for_node(root, "temperature_method")) {
8,856✔
1159
    auto temp = get_node_value(root, "temperature_method", true, true);
485✔
1160
    if (temp == "nearest") {
485✔
1161
      temperature_method = TemperatureMethod::NEAREST;
304✔
1162
    } else if (temp == "interpolation") {
181!
1163
      temperature_method = TemperatureMethod::INTERPOLATION;
181✔
1164
    } else {
1165
      fatal_error("Unknown temperature method: " + temp);
×
1166
    }
1167
  }
485✔
1168
  if (check_for_node(root, "temperature_tolerance")) {
8,856✔
1169
    temperature_tolerance =
680✔
1170
      std::stod(get_node_value(root, "temperature_tolerance"));
680✔
1171
  }
1172
  if (check_for_node(root, "temperature_multipole")) {
8,856✔
1173
    temperature_multipole = get_node_value_bool(root, "temperature_multipole");
185✔
1174

1175
    // Multipole currently doesn't work with photon transport
1176
    if (temperature_multipole && photon_transport) {
185!
1177
      fatal_error("Multipole data cannot currently be used in conjunction with "
×
1178
                  "photon transport.");
1179
    }
1180
  }
1181
  if (check_for_node(root, "temperature_range")) {
8,856✔
1182
    auto range = get_node_array<double>(root, "temperature_range");
170✔
1183
    temperature_range[0] = range.at(0);
170✔
1184
    temperature_range[1] = range.at(1);
170✔
1185
  }
170✔
1186

1187
  // Check for tabular_legendre options
1188
  if (check_for_node(root, "tabular_legendre")) {
8,856✔
1189
    // Get pointer to tabular_legendre node
1190
    xml_node node_tab_leg = root.child("tabular_legendre");
105✔
1191

1192
    // Check for enable option
1193
    if (check_for_node(node_tab_leg, "enable")) {
105!
1194
      legendre_to_tabular = get_node_value_bool(node_tab_leg, "enable");
105✔
1195
    }
1196

1197
    // Check for the number of points
1198
    if (check_for_node(node_tab_leg, "num_points")) {
105!
1199
      legendre_to_tabular_points =
×
1200
        std::stoi(get_node_value(node_tab_leg, "num_points"));
×
1201
      if (legendre_to_tabular_points <= 1 && !run_CE) {
×
1202
        fatal_error(
×
1203
          "The 'num_points' subelement/attribute of the "
1204
          "<tabular_legendre> element must contain a value greater than 1");
1205
      }
1206
    }
1207
  }
1208

1209
  // Check whether create delayed neutrons in fission
1210
  if (check_for_node(root, "create_delayed_neutrons")) {
8,856!
1211
    create_delayed_neutrons =
×
1212
      get_node_value_bool(root, "create_delayed_neutrons");
×
1213
  }
1214

1215
  // Check whether create fission sites
1216
  if (run_mode == RunMode::FIXED_SOURCE) {
8,856✔
1217
    if (check_for_node(root, "create_fission_neutrons")) {
3,142✔
1218
      create_fission_neutrons =
325✔
1219
        get_node_value_bool(root, "create_fission_neutrons");
325✔
1220
    }
1221
  }
1222

1223
  // Check whether to scale fission photon yields
1224
  if (check_for_node(root, "delayed_photon_scaling")) {
8,856!
1225
    delayed_photon_scaling =
×
1226
      get_node_value_bool(root, "delayed_photon_scaling");
×
1227
  }
1228

1229
  // Check whether to use event-based parallelism
1230
  if (check_for_node(root, "event_based")) {
8,856✔
1231
    event_based = get_node_value_bool(root, "event_based");
30✔
1232
  }
1233

1234
  // Check whether or not to use delta tracking
1235
  if (check_for_node(root, "delta_tracking")) {
8,856✔
1236
    delta_tracking = get_node_value_bool(root, "delta_tracking");
150✔
1237

1238
    if (temperature_multipole && delta_tracking) {
150!
NEW
1239
      fatal_error("Delta tracking cannot be used with a windowed multipole "
×
1240
                  "temperature treatment.");
1241
    }
1242

1243
    if (!run_CE && delta_tracking) {
150!
NEW
1244
      fatal_error("At present, delta tracking can only be used in continuous "
×
1245
                  "energy simulations.");
1246
    }
1247
  }
1248

1249
  // Check whether material cell offsets should be generated
1250
  if (check_for_node(root, "material_cell_offsets")) {
8,856!
1251
    material_cell_offsets = get_node_value_bool(root, "material_cell_offsets");
×
1252
  }
1253

1254
  // Weight window information
1255
  for (pugi::xml_node node_ww : root.children("weight_windows")) {
9,073✔
1256
    variance_reduction::weight_windows.emplace_back(
217✔
1257
      std::make_unique<WeightWindows>(node_ww));
434✔
1258
  }
1259

1260
  // Enable weight windows by default if one or more are present
1261
  if (variance_reduction::weight_windows.size() > 0)
8,856✔
1262
    settings::weight_windows_on = true;
165✔
1263

1264
  // read weight windows from file
1265
  if (check_for_node(root, "weight_windows_file")) {
8,856!
1266
    weight_windows_file = get_node_value(root, "weight_windows_file");
×
1267
  }
1268

1269
  // read settings for weight windows value, this will override
1270
  // the automatic setting even if weight windows are present
1271
  if (check_for_node(root, "weight_windows_on")) {
8,856✔
1272
    weight_windows_on = get_node_value_bool(root, "weight_windows_on");
60✔
1273
  }
1274

1275
  if (check_for_node(root, "max_secondaries")) {
8,856!
1276
    settings::max_secondaries =
×
1277
      std::stoi(get_node_value(root, "max_secondaries"));
×
1278
  }
1279

1280
  if (check_for_node(root, "max_history_splits")) {
8,856✔
1281
    settings::max_history_splits =
592✔
1282
      std::stoi(get_node_value(root, "max_history_splits"));
592✔
1283
  }
1284

1285
  if (check_for_node(root, "max_tracks")) {
8,856✔
1286
    settings::max_tracks = std::stoi(get_node_value(root, "max_tracks"));
90✔
1287
  }
1288

1289
  // Create weight window generator objects
1290
  if (check_for_node(root, "weight_window_generators")) {
8,856✔
1291
    auto wwgs_node = root.child("weight_window_generators");
105✔
1292
    for (pugi::xml_node node_wwg :
210✔
1293
      wwgs_node.children("weight_windows_generator")) {
210✔
1294
      variance_reduction::weight_windows_generators.emplace_back(
105✔
1295
        std::make_unique<WeightWindowsGenerator>(node_wwg));
210✔
1296
    }
1297
    // if any of the weight windows are intended to be generated otf, make
1298
    // sure they're applied
1299
    for (const auto& wwg : variance_reduction::weight_windows_generators) {
105!
1300
      if (wwg->on_the_fly_) {
105!
1301
        settings::weight_windows_on = true;
105✔
1302
        break;
105✔
1303
      }
1304
    }
1305
    // If any weight window generators have local FW-CADIS target tallies,
1306
    // user-defined adjoint sources cannot be used at the same time.
1307
    if (!model::adjoint_sources.empty()) {
105!
1308
      for (const auto& wwg : variance_reduction::weight_windows_generators) {
×
1309
        if (!wwg->targets_.empty()) {
×
1310
          fatal_error("Cannot use both user-defined adjoint sources and "
×
1311
                      "FW-CADIS target tallies at the same time.");
1312
        }
1313
      }
1314
    }
1315
  }
1316

1317
  // Set up weight window checkpoints
1318
  if (check_for_node(root, "weight_window_checkpoints")) {
8,856✔
1319
    xml_node ww_checkpoints = root.child("weight_window_checkpoints");
92✔
1320
    if (check_for_node(ww_checkpoints, "collision")) {
92!
1321
      weight_window_checkpoint_collision =
92✔
1322
        get_node_value_bool(ww_checkpoints, "collision");
92✔
1323
    }
1324
    if (check_for_node(ww_checkpoints, "surface")) {
92!
1325
      weight_window_checkpoint_surface =
92✔
1326
        get_node_value_bool(ww_checkpoints, "surface");
92✔
1327
    }
1328
  }
1329

1330
  if (weight_windows_on) {
8,856✔
1331
    if (!weight_window_checkpoint_surface &&
270✔
1332
        !weight_window_checkpoint_collision)
178!
1333
      fatal_error(
×
1334
        "Weight Windows are enabled but there are no valid checkpoints.");
1335
  }
1336

1337
  if (check_for_node(root, "use_decay_photons")) {
8,856✔
1338
    settings::use_decay_photons =
11✔
1339
      get_node_value_bool(root, "use_decay_photons");
11✔
1340
  }
1341

1342
  // If weight windows are on, also enable shared secondary bank (unless
1343
  // explicitly disabled by user).
1344
  if (check_for_node(root, "shared_secondary_bank")) {
8,856✔
1345
    bool val = get_node_value_bool(root, "shared_secondary_bank");
291✔
1346
    if (val && run_mode == RunMode::EIGENVALUE) {
291!
1347
      warning(
×
1348
        "Shared secondary bank is not supported in eigenvalue calculations. "
1349
        "Setting will be ignored.");
1350
    } else {
1351
      settings::use_shared_secondary_bank = val;
291✔
1352
    }
1353
  } else if (settings::weight_windows_on) {
8,565✔
1354
    if (run_mode == RunMode::EIGENVALUE) {
121✔
1355
      warning(
22✔
1356
        "Shared secondary bank is not supported in eigenvalue calculations. "
1357
        "Particle local secondary banks will be used instead.");
1358
    } else if (run_mode == RunMode::FIXED_SOURCE) {
110!
1359
      settings::use_shared_secondary_bank = true;
110✔
1360
    }
1361
  }
1362
}
8,856✔
1363

1364
void free_memory_settings()
9,004✔
1365
{
1366
  settings::statepoint_batch.clear();
9,004✔
1367
  settings::sourcepoint_batch.clear();
9,004✔
1368
  settings::source_write_surf_id.clear();
9,004✔
1369
  settings::res_scat_nuclides.clear();
9,004✔
1370
}
9,004✔
1371

1372
//==============================================================================
1373
// C API functions
1374
//==============================================================================
1375

1376
extern "C" int openmc_set_n_batches(
220✔
1377
  int32_t n_batches, bool set_max_batches, bool add_statepoint_batch)
1378
{
1379
  if (settings::n_inactive >= n_batches) {
220✔
1380
    set_errmsg("Number of active batches must be greater than zero.");
11✔
1381
    return OPENMC_E_INVALID_ARGUMENT;
11✔
1382
  }
1383

1384
  if (!settings::trigger_on) {
209✔
1385
    // Set n_batches and n_max_batches to same value
1386
    settings::n_batches = n_batches;
187✔
1387
    settings::n_max_batches = n_batches;
187✔
1388
  } else {
1389
    // Set n_batches and n_max_batches based on value of set_max_batches
1390
    if (set_max_batches) {
22✔
1391
      settings::n_max_batches = n_batches;
11✔
1392
    } else {
1393
      settings::n_batches = n_batches;
11✔
1394
    }
1395
  }
1396

1397
  // Update size of k_generation and entropy
1398
  int m = settings::n_max_batches * settings::gen_per_batch;
209✔
1399
  simulation::k_generation.reserve(m);
209✔
1400
  simulation::entropy.reserve(m);
209✔
1401

1402
  // Add value of n_batches to statepoint_batch
1403
  if (add_statepoint_batch &&
209✔
1404
      !(contains(settings::statepoint_batch, n_batches)))
198✔
1405
    settings::statepoint_batch.insert(n_batches);
33✔
1406

1407
  return 0;
1408
}
1409

1410
extern "C" int openmc_get_n_batches(int* n_batches, bool get_max_batches)
2,530✔
1411
{
1412
  *n_batches = get_max_batches ? settings::n_max_batches : settings::n_batches;
2,530✔
1413

1414
  return 0;
2,530✔
1415
}
1416

1417
} // namespace openmc
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2026 Coveralls, Inc