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

openmc-dev / openmc / 25930342666

15 May 2026 04:56PM UTC coverage: 81.414% (+0.09%) from 81.324%
25930342666

Pull #3734

github

web-flow
Merge 0af48acaf into d56cda254
Pull Request #3734: Specify temperature from a field (structured mesh only)

17919 of 25855 branches covered (69.31%)

Branch coverage included in aggregate %.

275 of 305 new or added lines in 18 files covered. (90.16%)

1321 existing lines in 29 files now uncovered.

58986 of 68607 relevant lines covered (85.98%)

47884537.9 hits per line

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

75.67
/src/initialize.cpp
1
#include "openmc/initialize.h"
2

3
#include <clocale>
4
#include <cstddef>
5
#include <cstdlib> // for getenv
6
#include <cstring>
7
#include <string>
8

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

14
#include "openmc/capi.h"
15
#include "openmc/chain.h"
16
#include "openmc/constants.h"
17
#include "openmc/cross_sections.h"
18
#include "openmc/error.h"
19
#include "openmc/file_utils.h"
20
#include "openmc/geometry_aux.h"
21
#include "openmc/hdf5_interface.h"
22
#include "openmc/material.h"
23
#include "openmc/memory.h"
24
#include "openmc/message_passing.h"
25
#include "openmc/mgxs_interface.h"
26
#include "openmc/nuclide.h"
27
#include "openmc/openmp_interface.h"
28
#include "openmc/output.h"
29
#include "openmc/plot.h"
30
#include "openmc/random_lcg.h"
31
#include "openmc/settings.h"
32
#include "openmc/simulation.h"
33
#include "openmc/string_utils.h"
34
#include "openmc/summary.h"
35
#include "openmc/tallies/tally.h"
36
#include "openmc/thermal.h"
37
#include "openmc/timer.h"
38
#include "openmc/vector.h"
39
#include "openmc/weight_windows.h"
40

41
#ifdef OPENMC_LIBMESH_ENABLED
42
#include "libmesh/libmesh.h"
43
#endif
44

45
int openmc_init(int argc, char* argv[], const void* intracomm)
8,540✔
46
{
47
  using namespace openmc;
8,540✔
48

49
#ifdef OPENMC_MPI
50
  // Check if intracomm was passed
51
  MPI_Comm comm;
3,759✔
52
  if (intracomm) {
3,759✔
53
    comm = *static_cast<const MPI_Comm*>(intracomm);
3,561✔
54
  } else {
55
    comm = MPI_COMM_WORLD;
47✔
56
  }
57

58
  // Initialize MPI for C++
59
  initialize_mpi(comm);
3,759✔
60
#endif
61

62
  // Parse command-line arguments
63
  int err = parse_command_line(argc, argv);
8,540✔
64
  if (err)
8,507✔
65
    return err;
66

67
#ifdef OPENMC_LIBMESH_ENABLED
68
  const int n_threads = num_threads();
1,595✔
69
  // initialize libMesh if it hasn't been initialized already
70
  // (if initialized externally, the libmesh_init object needs to be provided
71
  // also)
72
  if (!settings::libmesh_init && !libMesh::initialized()) {
1,595!
73
#ifdef OPENMC_MPI
74
    // pass command line args, empty MPI communicator, and number of threads.
75
    // Because libMesh was not initialized, we assume that OpenMC is the primary
76
    // application and that its main MPI comm should be used.
77
    settings::libmesh_init =
918✔
78
      make_unique<libMesh::LibMeshInit>(argc, argv, comm, n_threads);
1,836✔
79
#else
80
    // pass command line args, empty MPI communicator, and number of threads
81
    settings::libmesh_init =
677✔
82
      make_unique<libMesh::LibMeshInit>(argc, argv, 0, n_threads);
1,354✔
83
#endif
84

85
    settings::libmesh_comm = &(settings::libmesh_init->comm());
1,595✔
86
  }
87

88
#endif
89

90
  // Start total and initialization timer
91
  simulation::time_total.start();
8,496✔
92
  simulation::time_initialize.start();
8,496✔
93

94
#ifdef _OPENMP
95
  // If OMP_SCHEDULE is not set, default to a static schedule
96
  char* envvar = std::getenv("OMP_SCHEDULE");
4,846✔
97
  if (!envvar) {
4,846!
98
    omp_set_schedule(omp_sched_static, 0);
4,846✔
99
  }
100
#endif
101

102
  // Initialize random number generator -- if the user specifies a seed and/or
103
  // stride, it will be re-initialized later
104
  openmc::openmc_set_seed(DEFAULT_SEED);
8,496✔
105
  openmc::openmc_set_stride(DEFAULT_STRIDE);
8,496✔
106

107
  // Copy previous locale and set locale to C. This is a workaround for an issue
108
  // whereby when openmc_init is called from the plotter, the Qt application
109
  // framework first calls std::setlocale, which affects how pugixml reads
110
  // floating point numbers due to a bug:
111
  // https://github.com/zeux/pugixml/issues/469
112
  std::string prev_locale = std::setlocale(LC_ALL, nullptr);
8,496✔
113
  if (std::setlocale(LC_ALL, "C") == NULL) {
8,496!
114
    fatal_error("Cannot set locale to C.");
×
115
  }
116

117
  // Read XML input files
118
  if (!read_model_xml())
8,496✔
119
    read_separate_xml_files();
1,382✔
120

121
  if (!settings::properties_file.empty()) {
8,403✔
122
    openmc_properties_import(settings::properties_file.c_str());
11✔
123
  }
124

125
  // Reset locale to previous state
126
  if (std::setlocale(LC_ALL, prev_locale.c_str()) == NULL) {
8,403!
127
    fatal_error("Cannot reset locale.");
×
128
  }
129

130
  // Write some initial output under the header if needed
131
  initial_output();
8,403✔
132

133
  // Check for particle restart run
134
  if (settings::particle_restart_run)
8,403✔
135
    settings::run_mode = RunMode::PARTICLE;
41✔
136

137
  // Stop initialization timer
138
  simulation::time_initialize.stop();
8,403✔
139
  simulation::time_total.stop();
8,403✔
140

141
  return 0;
8,403✔
142
}
8,414✔
143

144
namespace openmc {
145

146
#ifdef OPENMC_MPI
147
void initialize_mpi(MPI_Comm intracomm)
3,759✔
148
{
149
  mpi::intracomm = intracomm;
3,759✔
150

151
  // Initialize MPI
152
  int flag;
3,759✔
153
  MPI_Initialized(&flag);
3,759✔
154
  if (!flag)
3,759✔
155
    MPI_Init(nullptr, nullptr);
3,075✔
156

157
  // Determine number of processes and rank for each
158
  MPI_Comm_size(intracomm, &mpi::n_procs);
3,759✔
159
  MPI_Comm_rank(intracomm, &mpi::rank);
3,759✔
160
  mpi::master = (mpi::rank == 0);
3,759✔
161

162
  // Create bank datatype
163
  SourceSite b;
3,759✔
164
  MPI_Aint disp[11];
3,759✔
165
  MPI_Get_address(&b.r, &disp[0]);
3,759✔
166
  MPI_Get_address(&b.u, &disp[1]);
3,759✔
167
  MPI_Get_address(&b.E, &disp[2]);
3,759✔
168
  MPI_Get_address(&b.time, &disp[3]);
3,759✔
169
  MPI_Get_address(&b.wgt, &disp[4]);
3,759✔
170
  MPI_Get_address(&b.delayed_group, &disp[5]);
3,759✔
171
  MPI_Get_address(&b.surf_id, &disp[6]);
3,759✔
172
  MPI_Get_address(&b.particle, &disp[7]);
3,759✔
173
  MPI_Get_address(&b.parent_nuclide, &disp[8]);
3,759✔
174
  MPI_Get_address(&b.parent_id, &disp[9]);
3,759✔
175
  MPI_Get_address(&b.progeny_id, &disp[10]);
3,759✔
176
  for (int i = 10; i >= 0; --i) {
45,108✔
177
    disp[i] -= disp[0];
41,349✔
178
  }
179

180
  int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1};
3,759✔
181
  MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
3,759✔
182
    MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, MPI_LONG};
183
  MPI_Type_create_struct(11, blocks, disp, types, &mpi::source_site);
3,759✔
184
  MPI_Type_commit(&mpi::source_site);
3,759✔
185

186
  CollisionTrackSite bc;
3,759✔
187
  MPI_Aint dispc[16];
3,759✔
188
  MPI_Get_address(&bc.r, &dispc[0]);             // double
3,759✔
189
  MPI_Get_address(&bc.u, &dispc[1]);             // double
3,759✔
190
  MPI_Get_address(&bc.E, &dispc[2]);             // double
3,759✔
191
  MPI_Get_address(&bc.dE, &dispc[3]);            // double
3,759✔
192
  MPI_Get_address(&bc.time, &dispc[4]);          // double
3,759✔
193
  MPI_Get_address(&bc.wgt, &dispc[5]);           // double
3,759✔
194
  MPI_Get_address(&bc.event_mt, &dispc[6]);      // int
3,759✔
195
  MPI_Get_address(&bc.delayed_group, &dispc[7]); // int
3,759✔
196
  MPI_Get_address(&bc.cell_id, &dispc[8]);       // int
3,759✔
197
  MPI_Get_address(&bc.nuclide_id, &dispc[9]);    // int
3,759✔
198
  MPI_Get_address(&bc.material_id, &dispc[10]);  // int
3,759✔
199
  MPI_Get_address(&bc.universe_id, &dispc[11]);  // int
3,759✔
200
  MPI_Get_address(&bc.n_collision, &dispc[12]);  // int
3,759✔
201
  MPI_Get_address(&bc.particle, &dispc[13]);     // int
3,759✔
202
  MPI_Get_address(&bc.parent_id, &dispc[14]);    // int64_t
3,759✔
203
  MPI_Get_address(&bc.progeny_id, &dispc[15]);   // int64_t
3,759✔
204
  for (int i = 15; i >= 0; --i) {
63,903✔
205
    dispc[i] -= dispc[0];
60,144✔
206
  }
207

208
  int blocksc[] = {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
3,759✔
209
  MPI_Datatype typesc[] = {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
3,759✔
210
    MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_INT,
211
    MPI_INT, MPI_INT, MPI_INT, MPI_INT64_T, MPI_INT64_T};
212

213
  MPI_Type_create_struct(
3,759✔
214
    16, blocksc, dispc, typesc, &mpi::collision_track_site);
215
  MPI_Type_commit(&mpi::collision_track_site);
3,759✔
216
}
3,759✔
217
#endif // OPENMC_MPI
218

219
int parse_command_line(int argc, char* argv[])
8,540✔
220
{
221
  int last_flag = 0;
8,540✔
222
  for (int i = 1; i < argc; ++i) {
9,860✔
223
    std::string arg {argv[i]};
1,331✔
224
    if (arg[0] == '-') {
1,331✔
225
      if (arg == "-p" || arg == "--plot") {
2,337!
226
        settings::run_mode = RunMode::PLOTTING;
119✔
227
        settings::check_overlaps = true;
119✔
228

229
      } else if (arg == "-n" || arg == "--particles") {
2,218!
230
        i += 1;
×
231
        settings::n_particles = std::stoll(argv[i]);
×
232

233
      } else if (arg == "-q" || arg == "--verbosity") {
2,218!
234
        i += 1;
×
235
        settings::verbosity = std::stoi(argv[i]);
×
236
        if (settings::verbosity > 10 || settings::verbosity < 1) {
×
237
          auto msg = fmt::format("Invalid verbosity: {}.", settings::verbosity);
×
238
          strcpy(openmc_err_msg, msg.c_str());
×
239
          return OPENMC_E_INVALID_ARGUMENT;
×
240
        }
×
241

242
      } else if (arg == "-e" || arg == "--event") {
2,005!
243
        settings::event_based = true;
213✔
244
      } else if (arg == "-r" || arg == "--restart") {
1,688!
245
        i += 1;
104✔
246
        // Check what type of file this is
247
        hid_t file_id = file_open(argv[i], 'r', true);
104✔
248
        std::string filetype;
104✔
249
        read_attribute(file_id, "filetype", filetype);
104✔
250
        file_close(file_id);
104✔
251

252
        // Set path and flag for type of run
253
        if (filetype == "statepoint") {
104✔
254
          settings::path_statepoint = argv[i];
63✔
255
          settings::path_statepoint_c = settings::path_statepoint.c_str();
63✔
256
          settings::restart_run = true;
63✔
257
        } else if (filetype == "particle restart") {
41!
258
          settings::path_particle_restart = argv[i];
41✔
259
          settings::particle_restart_run = true;
41✔
260
        } else {
261
          auto msg =
×
262
            fmt::format("Unrecognized file after restart flag: {}.", filetype);
×
263
          strcpy(openmc_err_msg, msg.c_str());
×
264
          return OPENMC_E_INVALID_ARGUMENT;
×
265
        }
×
266

267
        // If its a restart run check for additional source file
268
        if (settings::restart_run && i + 1 < argc) {
104!
269
          // Check if it has extension we can read
270
          if (ends_with(argv[i + 1], ".h5")) {
×
271

272
            // Check file type is a source file
273
            file_id = file_open(argv[i + 1], 'r', true);
×
274
            read_attribute(file_id, "filetype", filetype);
×
275
            file_close(file_id);
×
276
            if (filetype != "source") {
×
277
              std::string msg {
×
278
                "Second file after restart flag must be a source file"};
×
279
              strcpy(openmc_err_msg, msg.c_str());
×
280
              return OPENMC_E_INVALID_ARGUMENT;
×
281
            }
×
282

283
            // It is a source file
284
            settings::path_sourcepoint = argv[i + 1];
104!
285
            i += 1;
286

287
          } else {
288
            // Source is in statepoint file
289
            settings::path_sourcepoint = settings::path_statepoint;
×
290
          }
291

292
        } else {
293
          // Source is assumed to be in statepoint file
294
          settings::path_sourcepoint = settings::path_statepoint;
208✔
295
        }
296

297
      } else if (arg == "-g" || arg == "--geometry-debug") {
1,688!
298
        settings::check_overlaps = true;
×
299
      } else if (arg == "-c" || arg == "--volume") {
965✔
300
        settings::run_mode = RunMode::VOLUME;
701✔
301
      } else if (arg == "-s" || arg == "--threads") {
158!
302
        // Read number of threads
303
        i += 1;
24✔
304

305
#ifdef _OPENMP
306
        // Read and set number of OpenMP threads
307
        int n_threads = std::stoi(argv[i]);
26✔
308
        if (n_threads < 1) {
13!
309
          std::string msg {"Number of threads must be positive."};
×
310
          strcpy(openmc_err_msg, msg.c_str());
311
          return OPENMC_E_INVALID_ARGUMENT;
312
        }
313
        omp_set_num_threads(n_threads);
13✔
314
#else
315
        if (mpi::master) {
11✔
316
          warning("Ignoring number of threads specified on command line.");
20✔
317
        }
318
#endif
319

320
      } else if (arg == "-?" || arg == "-h" || arg == "--help") {
201!
321
        print_usage();
×
322
        return OPENMC_E_UNASSIGNED;
×
323

324
      } else if (arg == "-v" || arg == "--version") {
123!
325
        print_version();
11✔
326
        print_build_info();
11✔
327
        return OPENMC_E_UNASSIGNED;
11✔
328

329
      } else if (arg == "-t" || arg == "--track") {
56!
330
        settings::write_all_tracks = true;
56✔
331

332
      } else {
333
        fmt::print(stderr, "Unknown option: {}\n", argv[i]);
×
334
        print_usage();
×
335
        return OPENMC_E_UNASSIGNED;
×
336
      }
337

338
      last_flag = i;
339
    }
340
  }
1,331✔
341

342
  // Determine directory where XML input files are
343
  if (argc > 1 && last_flag < argc - 1) {
8,529✔
344
    settings::path_input = std::string(argv[last_flag + 1]);
103✔
345

346
    // check that the path is either a valid directory or file
347
    if (!dir_exists(settings::path_input) &&
193✔
348
        !file_exists(settings::path_input)) {
90✔
349
      fatal_error(fmt::format(
39✔
350
        "The path specified to the OpenMC executable '{}' does not exist.",
351
        settings::path_input));
352
    }
353

354
    // Add slash at end of directory if it isn't there
355
    if (!ends_with(settings::path_input, "/") &&
210!
356
        dir_exists(settings::path_input)) {
70✔
357
      settings::path_input += "/";
8,507✔
358
    }
359
  }
360

361
  return 0;
362
}
363

364
bool read_model_xml()
8,496✔
365
{
366
  std::string model_filename = settings::path_input;
8,496✔
367

368
  // if the current filename is a directory, append the default model filename
369
  if (model_filename.empty() || dir_exists(model_filename))
8,496✔
370
    model_filename += "model.xml";
8,439✔
371

372
  // if this file doesn't exist, stop here
373
  if (!file_exists(model_filename))
8,496✔
374
    return false;
375

376
  // try to process the path input as an XML file
377
  pugi::xml_document doc;
7,114✔
378
  if (!doc.load_file(model_filename.c_str())) {
7,114!
379
    fatal_error(fmt::format(
×
380
      "Error reading from single XML input file '{}'", model_filename));
381
  }
382

383
  pugi::xml_node root = doc.document_element();
7,114✔
384

385
  // Read settings
386
  if (!check_for_node(root, "settings")) {
7,114!
387
    fatal_error("No <settings> node present in the model.xml file.");
×
388
  }
389
  auto settings_root = root.child("settings");
7,114✔
390

391
  // Verbosity
392
  if (check_for_node(settings_root, "verbosity") && settings::verbosity == -1) {
7,114!
393
    settings::verbosity = std::stoi(get_node_value(settings_root, "verbosity"));
60✔
394
  } else if (settings::verbosity == -1) {
7,084!
395
    settings::verbosity = 7;
7,084✔
396
  }
397

398
  // To this point, we haven't displayed any output since we didn't know what
399
  // the verbosity is. Now that we checked for it, show the title if necessary
400
  if (mpi::master) {
7,114✔
401
    if (settings::verbosity >= 2)
6,322✔
402
      title();
6,300✔
403
  }
404

405
  write_message(
7,114✔
406
    fmt::format("Reading model XML file '{}' ...", model_filename), 5);
8,447✔
407

408
  // Read chain data before settings so DecaySpectrum source distributions can
409
  // resolve nuclides while sources are constructed.
410
  read_chain_file_xml();
7,114✔
411

412
  read_settings_xml(settings_root);
7,114✔
413

414
  // If other XML files are present, display warning
415
  // that they will be ignored
416
  auto other_inputs = {"materials.xml", "geometry.xml", "settings.xml",
7,060✔
417
    "tallies.xml", "plots.xml"};
7,060✔
418
  for (const auto& input : other_inputs) {
41,566✔
419
    if (file_exists(settings::path_input + input)) {
34,696✔
420
      warning((fmt::format("Other XML file input(s) are present. These files "
225✔
421
                           "may be ignored in favor of the {} file.",
422
        model_filename)));
423
      break;
190✔
424
    }
425
  }
426

427
  // Read materials and cross sections
428
  if (!check_for_node(root, "materials")) {
7,060!
UNCOV
429
    fatal_error(fmt::format(
×
430
      "No <materials> node present in the {} file.", model_filename));
431
  }
432

433
  if (settings::run_mode != RunMode::PLOTTING) {
7,060✔
434
    read_cross_sections_xml(root.child("materials"));
6,974✔
435
  }
436
  read_materials_xml(root.child("materials"));
7,060✔
437

438
  // Read geometry
439
  if (!check_for_node(root, "geometry")) {
7,060!
UNCOV
440
    fatal_error(fmt::format(
×
441
      "No <geometry> node present in the {} file.", model_filename));
442
  }
443
  read_geometry_xml(root.child("geometry"));
7,060✔
444

445
  // Final geometry setup and assign temperatures
446
  finalize_geometry();
7,058✔
447

448
  // Finalize cross sections having assigned temperatures
449
  finalize_cross_sections();
7,058✔
450

451
  // Compute cell density multipliers now that material densities
452
  // have been finalized (from geometry_aux.h)
453
  finalize_cell_densities();
7,058✔
454

455
  if (check_for_node(root, "tallies"))
7,058✔
456
    read_tallies_xml(root.child("tallies"));
4,168✔
457

458
  // Initialize distribcell_filters
459
  prepare_distribcell();
7,040✔
460

461
  if (check_for_node(root, "plots")) {
7,040✔
462
    read_plots_xml(root.child("plots"));
365✔
463
  } else {
464
    // When no <plots> element is present in the model.xml file, check for a
465
    // regular plots.xml file
466
    std::string filename = settings::path_input + "plots.xml";
6,675✔
467
    if (file_exists(filename)) {
6,675!
UNCOV
468
      read_plots_xml();
×
469
    }
470
  }
6,675✔
471

472
  finalize_variance_reduction();
7,031✔
473

474
  return true;
7,031✔
475
}
15,444✔
476

477
void read_separate_xml_files()
1,382✔
478
{
479
  // Read chain data before settings so DecaySpectrum source distributions can
480
  // resolve nuclides while sources are constructed.
481
  read_chain_file_xml();
1,382✔
482

483
  read_settings_xml();
1,382✔
484
  if (settings::run_mode != RunMode::PLOTTING) {
1,372✔
485
    read_cross_sections_xml();
1,339✔
486
  }
487

488
  read_materials_xml();
1,372✔
489
  read_geometry_xml();
1,372✔
490

491
  // Final geometry setup and assign temperatures
492
  finalize_geometry();
1,372✔
493

494
  // Finalize cross sections having assigned temperatures
495
  finalize_cross_sections();
1,372✔
496

497
  // Compute cell density multipliers now that material densities
498
  // have been finalized (from geometry_aux.h)
499
  finalize_cell_densities();
1,372✔
500

501
  read_tallies_xml();
1,372✔
502

503
  // Initialize distribcell_filters
504
  prepare_distribcell();
1,372✔
505

506
  // Read the plots.xml regardless of plot mode in case plots are requested
507
  // via the API
508
  read_plots_xml();
1,372✔
509

510
  finalize_variance_reduction();
1,372✔
511
}
1,372✔
512

513
void initial_output()
8,403✔
514
{
515
  // write initial output
516
  if (settings::run_mode == RunMode::PLOTTING) {
8,403✔
517
    // Read plots.xml if it exists
518
    if (mpi::master && settings::verbosity >= 5)
110!
519
      print_plot();
88✔
520

521
  } else {
522
    // Write summary information
523
    if (mpi::master && settings::output_summary)
8,293✔
524
      write_summary();
6,977✔
525

526
    // Warn if overlap checking is on
527
    if (mpi::master && settings::check_overlaps) {
8,293!
UNCOV
528
      warning("Cell overlap checking is ON.");
×
529
    }
530
  }
531
}
8,403✔
532

533
} // 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