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

openmc-dev / openmc / 21043561037

15 Jan 2026 07:23PM UTC coverage: 81.388% (-0.7%) from 82.044%
21043561037

Pull #3734

github

web-flow
Merge 5e79b7b6f into 179048b80
Pull Request #3734: Specify temperature from a field (structured mesh only)

16703 of 22995 branches covered (72.64%)

Branch coverage included in aggregate %.

156 of 179 new or added lines in 12 files covered. (87.15%)

843 existing lines in 36 files now uncovered.

54487 of 64475 relevant lines covered (84.51%)

27592585.27 hits per line

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

75.65
/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)
3,578✔
46
{
47
  using namespace openmc;
48

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

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

62
  // Parse command-line arguments
63
  int err = parse_command_line(argc, argv);
3,578✔
64
  if (err)
3,563✔
65
    return err;
5✔
66

67
#ifdef OPENMC_LIBMESH_ENABLED
68
  const int n_threads = num_threads();
1,463✔
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,463!
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 =
78
      make_unique<libMesh::LibMeshInit>(argc, argv, comm, n_threads);
843✔
79
#else
80
    // pass command line args, empty MPI communicator, and number of threads
81
    settings::libmesh_init =
82
      make_unique<libMesh::LibMeshInit>(argc, argv, 0, n_threads);
620✔
83
#endif
84

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

88
#endif
89

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

94
#ifdef _OPENMP
95
  // If OMP_SCHEDULE is not set, default to a static schedule
96
  char* envvar = std::getenv("OMP_SCHEDULE");
1,463✔
97
  if (!envvar) {
1,463!
98
    omp_set_schedule(omp_sched_static, 0);
1,463✔
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);
3,558✔
105
  openmc::openmc_set_stride(DEFAULT_STRIDE);
3,558✔
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);
3,558✔
113
  if (std::setlocale(LC_ALL, "C") == NULL) {
3,558!
114
    fatal_error("Cannot set locale to C.");
×
115
  }
116

117
  // Read XML input files
118
  if (!read_model_xml())
3,558✔
119
    read_separate_xml_files();
610✔
120

121
  // Reset locale to previous state
122
  if (std::setlocale(LC_ALL, prev_locale.c_str()) == NULL) {
3,528!
123
    fatal_error("Cannot reset locale.");
×
124
  }
125

126
  // Write some initial output under the header if needed
127
  initial_output();
3,528✔
128

129
  // Check for particle restart run
130
  if (settings::particle_restart_run)
3,528✔
131
    settings::run_mode = RunMode::PARTICLE;
19✔
132

133
  // Stop initialization timer
134
  simulation::time_initialize.stop();
3,528✔
135
  simulation::time_total.stop();
3,528✔
136

137
  return 0;
3,528✔
138
}
3,528✔
139

140
namespace openmc {
141

142
#ifdef OPENMC_MPI
143
void initialize_mpi(MPI_Comm intracomm)
1,694✔
144
{
145
  mpi::intracomm = intracomm;
1,694✔
146

147
  // Initialize MPI
148
  int flag;
149
  MPI_Initialized(&flag);
1,694✔
150
  if (!flag)
1,694✔
151
    MPI_Init(nullptr, nullptr);
1,391✔
152

153
  // Determine number of processes and rank for each
154
  MPI_Comm_size(intracomm, &mpi::n_procs);
1,694✔
155
  MPI_Comm_rank(intracomm, &mpi::rank);
1,694✔
156
  mpi::master = (mpi::rank == 0);
1,694✔
157

158
  // Create bank datatype
159
  SourceSite b;
1,694✔
160
  MPI_Aint disp[11];
161
  MPI_Get_address(&b.r, &disp[0]);
1,694✔
162
  MPI_Get_address(&b.u, &disp[1]);
1,694✔
163
  MPI_Get_address(&b.E, &disp[2]);
1,694✔
164
  MPI_Get_address(&b.time, &disp[3]);
1,694✔
165
  MPI_Get_address(&b.wgt, &disp[4]);
1,694✔
166
  MPI_Get_address(&b.delayed_group, &disp[5]);
1,694✔
167
  MPI_Get_address(&b.surf_id, &disp[6]);
1,694✔
168
  MPI_Get_address(&b.particle, &disp[7]);
1,694✔
169
  MPI_Get_address(&b.parent_nuclide, &disp[8]);
1,694✔
170
  MPI_Get_address(&b.parent_id, &disp[9]);
1,694✔
171
  MPI_Get_address(&b.progeny_id, &disp[10]);
1,694✔
172
  for (int i = 10; i >= 0; --i) {
20,328✔
173
    disp[i] -= disp[0];
18,634✔
174
  }
175

176
  int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1};
1,694✔
177
  MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
1,694✔
178
    MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, MPI_LONG};
179
  MPI_Type_create_struct(11, blocks, disp, types, &mpi::source_site);
1,694✔
180
  MPI_Type_commit(&mpi::source_site);
1,694✔
181

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

204
  int blocksc[] = {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1};
1,694✔
205
  MPI_Datatype typesc[] = {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
1,694✔
206
    MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_INT,
207
    MPI_INT, MPI_INT, MPI_INT, MPI_INT64_T, MPI_INT64_T};
208

209
  MPI_Type_create_struct(
1,694✔
210
    16, blocksc, dispc, typesc, &mpi::collision_track_site);
211
  MPI_Type_commit(&mpi::collision_track_site);
1,694✔
212
}
1,694✔
213
#endif // OPENMC_MPI
214

215
int parse_command_line(int argc, char* argv[])
3,578✔
216
{
217
  int last_flag = 0;
3,578✔
218
  for (int i = 1; i < argc; ++i) {
4,049✔
219
    std::string arg {argv[i]};
476✔
220
    if (arg[0] == '-') {
476✔
221
      if (arg == "-p" || arg == "--plot") {
431!
222
        settings::run_mode = RunMode::PLOTTING;
53✔
223
        settings::check_overlaps = true;
53✔
224

225
      } else if (arg == "-n" || arg == "--particles") {
378!
226
        i += 1;
×
227
        settings::n_particles = std::stoll(argv[i]);
×
228

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

238
      } else if (arg == "-e" || arg == "--event") {
378!
UNCOV
239
        settings::event_based = true;
×
240
      } else if (arg == "-r" || arg == "--restart") {
378!
241
        i += 1;
48✔
242
        // Check what type of file this is
243
        hid_t file_id = file_open(argv[i], 'r', true);
48✔
244
        std::string filetype;
48✔
245
        read_attribute(file_id, "filetype", filetype);
48✔
246
        file_close(file_id);
48✔
247

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

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

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

279
            // It is a source file
280
            settings::path_sourcepoint = argv[i + 1];
×
281
            i += 1;
×
282

283
          } else {
284
            // Source is in statepoint file
285
            settings::path_sourcepoint = settings::path_statepoint;
×
286
          }
287

288
        } else {
×
289
          // Source is assumed to be in statepoint file
290
          settings::path_sourcepoint = settings::path_statepoint;
48✔
291
        }
292

293
      } else if (arg == "-g" || arg == "--geometry-debug") {
378!
294
        settings::check_overlaps = true;
×
295
      } else if (arg == "-c" || arg == "--volume") {
330✔
296
        settings::run_mode = RunMode::VOLUME;
289✔
297
      } else if (arg == "-s" || arg == "--threads") {
41!
298
        // Read number of threads
299
        i += 1;
10✔
300

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

316
      } else if (arg == "-?" || arg == "-h" || arg == "--help") {
31!
317
        print_usage();
×
318
        return OPENMC_E_UNASSIGNED;
×
319

320
      } else if (arg == "-v" || arg == "--version") {
31!
321
        print_version();
5✔
322
        print_build_info();
5✔
323
        return OPENMC_E_UNASSIGNED;
5✔
324

325
      } else if (arg == "-t" || arg == "--track") {
26!
326
        settings::write_all_tracks = true;
26✔
327

328
      } else {
329
        fmt::print(stderr, "Unknown option: {}\n", argv[i]);
×
330
        print_usage();
×
331
        return OPENMC_E_UNASSIGNED;
×
332
      }
333

334
      last_flag = i;
426✔
335
    }
336
  }
476✔
337

338
  // Determine directory where XML input files are
339
  if (argc > 1 && last_flag < argc - 1) {
3,573✔
340
    settings::path_input = std::string(argv[last_flag + 1]);
45✔
341

342
    // check that the path is either a valid directory or file
343
    if (!dir_exists(settings::path_input) &&
85✔
344
        !file_exists(settings::path_input)) {
40✔
345
      fatal_error(fmt::format(
15✔
346
        "The path specified to the OpenMC executable '{}' does not exist.",
347
        settings::path_input));
348
    }
349

350
    // Add slash at end of directory if it isn't there
351
    if (!ends_with(settings::path_input, "/") &&
60!
352
        dir_exists(settings::path_input)) {
30✔
353
      settings::path_input += "/";
5✔
354
    }
355
  }
356

357
  return 0;
3,558✔
358
}
359

360
bool read_model_xml()
3,558✔
361
{
362
  std::string model_filename = settings::path_input;
3,558✔
363

364
  // if the current filename is a directory, append the default model filename
365
  if (model_filename.empty() || dir_exists(model_filename))
3,558✔
366
    model_filename += "model.xml";
3,533✔
367

368
  // if this file doesn't exist, stop here
369
  if (!file_exists(model_filename))
3,558✔
370
    return false;
610✔
371

372
  // try to process the path input as an XML file
373
  pugi::xml_document doc;
2,948✔
374
  if (!doc.load_file(model_filename.c_str())) {
2,948!
375
    fatal_error(fmt::format(
×
376
      "Error reading from single XML input file '{}'", model_filename));
377
  }
378

379
  pugi::xml_node root = doc.document_element();
2,948✔
380

381
  // Read settings
382
  if (!check_for_node(root, "settings")) {
2,948!
383
    fatal_error("No <settings> node present in the model.xml file.");
×
384
  }
385
  auto settings_root = root.child("settings");
2,948✔
386

387
  // Verbosity
388
  if (check_for_node(settings_root, "verbosity") && settings::verbosity == -1) {
2,948!
389
    settings::verbosity = std::stoi(get_node_value(settings_root, "verbosity"));
14✔
390
  } else if (settings::verbosity == -1) {
2,934!
391
    settings::verbosity = 7;
2,934✔
392
  }
393

394
  // To this point, we haven't displayed any output since we didn't know what
395
  // the verbosity is. Now that we checked for it, show the title if necessary
396
  if (mpi::master) {
2,948✔
397
    if (settings::verbosity >= 2)
2,593✔
398
      title();
2,583✔
399
  }
400

401
  write_message(
2,948✔
402
    fmt::format("Reading model XML file '{}' ...", model_filename), 5);
4,687✔
403

404
  read_settings_xml(settings_root);
2,948✔
405

406
  // If other XML files are present, display warning
407
  // that they will be ignored
408
  auto other_inputs = {"materials.xml", "geometry.xml", "settings.xml",
2,930✔
409
    "tallies.xml", "plots.xml"};
2,930✔
410
  for (const auto& input : other_inputs) {
17,222✔
411
    if (file_exists(settings::path_input + input)) {
14,378✔
412
      warning((fmt::format("Other XML file input(s) are present. These files "
86✔
413
                           "may be ignored in favor of the {} file.",
414
        model_filename)));
415
      break;
86✔
416
    }
417
  }
418

419
  // Read data from chain file
420
  read_chain_file_xml();
2,930✔
421

422
  // Read materials and cross sections
423
  if (!check_for_node(root, "materials")) {
2,930!
424
    fatal_error(fmt::format(
×
425
      "No <materials> node present in the {} file.", model_filename));
426
  }
427

428
  if (settings::run_mode != RunMode::PLOTTING) {
2,930✔
429
    read_cross_sections_xml(root.child("materials"));
2,892✔
430
  }
431
  read_materials_xml(root.child("materials"));
2,930✔
432

433
  // Read geometry
434
  if (!check_for_node(root, "geometry")) {
2,930!
435
    fatal_error(fmt::format(
×
436
      "No <geometry> node present in the {} file.", model_filename));
437
  }
438
  read_geometry_xml(root.child("geometry"));
2,930✔
439

440
  // Final geometry setup and assign temperatures
441
  finalize_geometry();
2,930✔
442

443
  // Finalize cross sections having assigned temperatures
444
  finalize_cross_sections();
2,930✔
445

446
  // Compute cell density multipliers now that material densities
447
  // have been finalized (from geometry_aux.h)
448
  finalize_cell_densities();
2,930✔
449

450
  if (check_for_node(root, "tallies"))
2,930✔
451
    read_tallies_xml(root.child("tallies"));
1,653✔
452

453
  // Initialize distribcell_filters
454
  prepare_distribcell();
2,924✔
455

456
  if (check_for_node(root, "plots")) {
2,924✔
457
    read_plots_xml(root.child("plots"));
140✔
458
  } else {
459
    // When no <plots> element is present in the model.xml file, check for a
460
    // regular plots.xml file
461
    std::string filename = settings::path_input + "plots.xml";
2,784✔
462
    if (file_exists(filename)) {
2,784!
463
      read_plots_xml();
×
464
    }
465
  }
2,784✔
466

467
  finalize_variance_reduction();
2,921✔
468

469
  return true;
2,921✔
470
}
3,531✔
471

472
void read_separate_xml_files()
610✔
473
{
474
  read_settings_xml();
610✔
475
  if (settings::run_mode != RunMode::PLOTTING) {
607✔
476
    read_cross_sections_xml();
592✔
477
  }
478

479
  // Read data from chain file
480
  read_chain_file_xml();
607✔
481

482
  read_materials_xml();
607✔
483
  read_geometry_xml();
607✔
484

485
  // Final geometry setup and assign temperatures
486
  finalize_geometry();
607✔
487

488
  // Finalize cross sections having assigned temperatures
489
  finalize_cross_sections();
607✔
490

491
  // Compute cell density multipliers now that material densities
492
  // have been finalized (from geometry_aux.h)
493
  finalize_cell_densities();
607✔
494

495
  read_tallies_xml();
607✔
496

497
  // Initialize distribcell_filters
498
  prepare_distribcell();
607✔
499

500
  // Read the plots.xml regardless of plot mode in case plots are requested
501
  // via the API
502
  read_plots_xml();
607✔
503

504
  finalize_variance_reduction();
607✔
505
}
607✔
506

507
void initial_output()
3,528✔
508
{
509
  // write initial output
510
  if (settings::run_mode == RunMode::PLOTTING) {
3,528✔
511
    // Read plots.xml if it exists
512
    if (mpi::master && settings::verbosity >= 5)
50!
513
      print_plot();
40✔
514

515
  } else {
516
    // Write summary information
517
    if (mpi::master && settings::output_summary)
3,478✔
518
      write_summary();
2,888✔
519

520
    // Warn if overlap checking is on
521
    if (mpi::master && settings::check_overlaps) {
3,478!
522
      warning("Cell overlap checking is ON.");
×
523
    }
524
  }
525
}
3,528✔
526

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