16898768927

Committed 12 Aug 2025 04:13AM UTC coverage: 85.139% (-0.07%) from 85.205%

Build # 16898768927

Build Type

Pull #3528

github

Committed by

web-flow

Commit Message

Merge ed24df12a into a34365396

Pull Request Pull Request #3528: feature/split plot classes

Run Details

225 of 293 new or added lines in 2 files covered. (76.79%)

51 existing lines in 2 files now uncovered.

53131 of 62405 relevant lines covered (85.14%)

37146686.55 hits per line

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

64.52

/src/mcpl_interface.cpp

#include "openmc/mcpl_interface.h"

#include "openmc/bank.h"
#include "openmc/error.h"
#include "openmc/file_utils.h"
#include "openmc/message_passing.h"
#include "openmc/settings.h"
#include "openmc/simulation.h"
#include "openmc/state_point.h"
#include "openmc/vector.h"

#include <fmt/core.h>

#include <cctype>
#include <cstdint>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <memory>
#include <mutex>
#include <stdexcept>
#include <string>

#ifdef _WIN32
#define WIN32_LEAN_AND_MEAN
#include <windows.h>
#else
#include <dlfcn.h>
#endif

// WARNING: These declarations MUST EXACTLY MATCH the structure and function
// signatures of the libmcpl being loaded at runtime. Any discrepancy will
// likely lead to crashes or incorrect behavior. This is a maintenance risk.
// MCPL 2.2.0

#pragma pack(push, 1)
struct mcpl_particle_repr_t {
  double ekin;
  double polarisation[3];
  double position[3];
  double direction[3];
  double time;
  double weight;
  int32_t pdgcode;
  uint32_t userflags;
};
#pragma pack(pop)

// Opaque struct definitions replicating the MCPL C-API to ensure ABI
// compatibility without including mcpl.h. These must be kept in sync.
struct mcpl_file_t {
  void* internal;
};
struct mcpl_outfile_t {
  void* internal;
};

// Function pointer types for the dynamically loaded MCPL library
using mcpl_open_file_fpt = mcpl_file_t* (*)(const char* filename);
using mcpl_hdr_nparticles_fpt = uint64_t (*)(mcpl_file_t* file_handle);
using mcpl_read_fpt = const mcpl_particle_repr_t* (*)(mcpl_file_t *
                                                      file_handle);
using mcpl_close_file_fpt = void (*)(mcpl_file_t* file_handle);

using mcpl_create_outfile_fpt = mcpl_outfile_t* (*)(const char* filename);
using mcpl_hdr_set_srcname_fpt = void (*)(
  mcpl_outfile_t* outfile_handle, const char* srcname);
using mcpl_add_particle_fpt = void (*)(
  mcpl_outfile_t* outfile_handle, const mcpl_particle_repr_t* particle);
using mcpl_close_outfile_fpt = void (*)(mcpl_outfile_t* outfile_handle);
using mcpl_hdr_add_data_fpt = void (*)(mcpl_outfile_t* outfile_handle,
  const char* key, uint32_t datalength, const char* data);
using mcpl_hdr_add_stat_sum_fpt = void (*)(
  mcpl_outfile_t* outfile_handle, const char* key, double value);

namespace openmc {

#ifdef _WIN32
using LibraryHandleType = HMODULE;
#else
using LibraryHandleType = void*;
#endif

std::string get_last_library_error()
{
#ifdef _WIN32
  DWORD error_code = GetLastError();
  if (error_code == 0)
    return "No error reported by system."; // More accurate than "No error."
  LPSTR message_buffer = nullptr;
  size_t size =
    FormatMessageA(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM |
                     FORMAT_MESSAGE_IGNORE_INSERTS,
      NULL, error_code, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT),
      (LPSTR)&message_buffer, 0, NULL);
  std::string message(message_buffer, size);
  LocalFree(message_buffer);
  while (
    !message.empty() && (message.back() == '\n' || message.back() == '\r')) {
    message.pop_back();
  }
  return message;
#else
  const char* err = dlerror();
  return err ? std::string(err) : "No error reported by dlerror.";
#endif
}

struct McplApi {
  mcpl_open_file_fpt open_file;
  mcpl_hdr_nparticles_fpt hdr_nparticles;
  mcpl_read_fpt read;
  mcpl_close_file_fpt close_file;
  mcpl_create_outfile_fpt create_outfile;
  mcpl_hdr_set_srcname_fpt hdr_set_srcname;
  mcpl_add_particle_fpt add_particle;
  mcpl_close_outfile_fpt close_outfile;
  mcpl_hdr_add_data_fpt hdr_add_data;
  mcpl_hdr_add_stat_sum_fpt hdr_add_stat_sum;

  explicit McplApi(LibraryHandleType lib_handle)
  {
    if (!lib_handle)
      throw std::runtime_error(
        "MCPL library handle is null during API binding.");

    auto load_symbol_platform = [lib_handle](const char* name) {
      void* sym = nullptr;
#ifdef _WIN32
      sym = (void*)GetProcAddress(lib_handle, name);
#else
      sym = dlsym(lib_handle, name);
#endif
      if (!sym) {
        throw std::runtime_error(
          fmt::format("Failed to load MCPL symbol '{}': {}", name,
            get_last_library_error()));
      }
      return sym;
    };

    open_file = reinterpret_cast<mcpl_open_file_fpt>(
      load_symbol_platform("mcpl_open_file"));
    hdr_nparticles = reinterpret_cast<mcpl_hdr_nparticles_fpt>(
      load_symbol_platform("mcpl_hdr_nparticles"));
    read = reinterpret_cast<mcpl_read_fpt>(load_symbol_platform("mcpl_read"));
    close_file = reinterpret_cast<mcpl_close_file_fpt>(
      load_symbol_platform("mcpl_close_file"));
    create_outfile = reinterpret_cast<mcpl_create_outfile_fpt>(
      load_symbol_platform("mcpl_create_outfile"));
    hdr_set_srcname = reinterpret_cast<mcpl_hdr_set_srcname_fpt>(
      load_symbol_platform("mcpl_hdr_set_srcname"));
    add_particle = reinterpret_cast<mcpl_add_particle_fpt>(
      load_symbol_platform("mcpl_add_particle"));
    close_outfile = reinterpret_cast<mcpl_close_outfile_fpt>(
      load_symbol_platform("mcpl_close_outfile"));

    // Try to load mcpl_hdr_add_data (available in MCPL >= 2.1.0)
    // Set to nullptr if not available for graceful fallback
    try {
      hdr_add_data = reinterpret_cast<mcpl_hdr_add_data_fpt>(
        load_symbol_platform("mcpl_hdr_add_data"));
    } catch (const std::runtime_error&) {
      hdr_add_data = nullptr;
    }

    // Try to load mcpl_hdr_add_stat_sum (available in MCPL >= 2.1.0)
    // Set to nullptr if not available for graceful fallback
    try {
      hdr_add_stat_sum = reinterpret_cast<mcpl_hdr_add_stat_sum_fpt>(
        load_symbol_platform("mcpl_hdr_add_stat_sum"));
    } catch (const std::runtime_error&) {
      hdr_add_stat_sum = nullptr;
    }
  }
};

static LibraryHandleType g_mcpl_lib_handle = nullptr;
static std::unique_ptr<McplApi> g_mcpl_api;
static bool g_mcpl_init_attempted = false;
static bool g_mcpl_successfully_loaded = false;
static std::string g_mcpl_load_error_msg;
static std::once_flag g_mcpl_init_flag;

void append_error(std::string& existing_msg, const std::string& new_error)
{
  if (!existing_msg.empty()) {
    existing_msg += "; ";
  }
  existing_msg += new_error;
}

void initialize_mcpl_interface_impl()
{
  g_mcpl_init_attempted = true;
  g_mcpl_load_error_msg.clear();

  // Try mcpl-config
  if (!g_mcpl_lib_handle) {
    FILE* pipe = nullptr;
#ifdef _WIN32
    pipe = _popen("mcpl-config --show libpath", "r");
#else
    pipe = popen("mcpl-config --show libpath 2>/dev/null", "r");
#endif
    if (pipe) {
      char buffer[512];
      if (fgets(buffer, sizeof(buffer), pipe) != nullptr) {
        std::string shlibpath = buffer;
        // Remove trailing whitespace
        while (!shlibpath.empty() &&
               std::isspace(static_cast<unsigned char>(shlibpath.back()))) {
          shlibpath.pop_back();
        }

        if (!shlibpath.empty()) {
#ifdef _WIN32
          g_mcpl_lib_handle = LoadLibraryA(shlibpath.c_str());
#else
          g_mcpl_lib_handle = dlopen(shlibpath.c_str(), RTLD_LAZY);
#endif
          if (!g_mcpl_lib_handle) {
            append_error(
              g_mcpl_load_error_msg, fmt::format("From mcpl-config ({}): {}",
                                       shlibpath, get_last_library_error()));
          }
        }
      }
#ifdef _WIN32
      _pclose(pipe);
#else
      pclose(pipe);
#endif
    } else { // pipe failed to open
      append_error(g_mcpl_load_error_msg,
        "mcpl-config command not found or failed to execute");
    }
  }

  // Try standard library names
  if (!g_mcpl_lib_handle) {
#ifdef _WIN32
    const char* standard_names[] = {"mcpl.dll", "libmcpl.dll"};
#else
    const char* standard_names[] = {"libmcpl.so", "libmcpl.dylib"};
#endif
    for (const char* name : standard_names) {
#ifdef _WIN32
      g_mcpl_lib_handle = LoadLibraryA(name);
#else
      g_mcpl_lib_handle = dlopen(name, RTLD_LAZY);
#endif
      if (g_mcpl_lib_handle)
        break;
    }
    if (!g_mcpl_lib_handle) {
      append_error(
        g_mcpl_load_error_msg, fmt::format("Using standard names (e.g. {}): {}",
                                 standard_names[0], get_last_library_error()));
    }
  }

  if (!g_mcpl_lib_handle) {
    if (mpi::master) {
      warning(fmt::format("MCPL library could not be loaded. MCPL-dependent "
                          "features will be unavailable. Load attempts: {}",
        g_mcpl_load_error_msg.empty()
          ? "No specific error during load attempts."
          : g_mcpl_load_error_msg));
    }
    g_mcpl_successfully_loaded = false;
    return;
  }

  try {
    g_mcpl_api = std::make_unique<McplApi>(g_mcpl_lib_handle);
    g_mcpl_successfully_loaded = true;
    // Do not call dlclose/FreeLibrary at exit. Leaking the handle is safer
    // and standard practice for libraries used for the application's lifetime.
  } catch (const std::runtime_error& e) {
    append_error(g_mcpl_load_error_msg,
      fmt::format(
        "MCPL library loaded, but failed to bind symbols: {}", e.what()));
    if (mpi::master) {
      warning(g_mcpl_load_error_msg);
    }
#ifdef _WIN32
    FreeLibrary(g_mcpl_lib_handle);
#else
    dlclose(g_mcpl_lib_handle);
#endif
    g_mcpl_lib_handle = nullptr;
    g_mcpl_successfully_loaded = false;
  }
}

void initialize_mcpl_interface_if_needed()
{
  std::call_once(g_mcpl_init_flag, initialize_mcpl_interface_impl);
}

bool is_mcpl_interface_available()
{
  initialize_mcpl_interface_if_needed();
  return g_mcpl_successfully_loaded;
}

inline void ensure_mcpl_ready_or_fatal()
{
  initialize_mcpl_interface_if_needed();
  if (!g_mcpl_successfully_loaded) {
    fatal_error("MCPL functionality is required, but the MCPL library is not "
                "available or failed to initialize. Please ensure MCPL is "
                "installed and its library can be found (e.g., via PATH on "
                "Windows, LD_LIBRARY_PATH on Linux, or DYLD_LIBRARY_PATH on "
                "macOS). You can often install MCPL with 'pip install mcpl' or "
                "'conda install mcpl'.");
  }
}

SourceSite mcpl_particle_to_site(const mcpl_particle_repr_t* particle_repr)
{
  SourceSite site;
  switch (particle_repr->pdgcode) {
  case 2112:
    site.particle = ParticleType::neutron;
    break;
  case 22:
    site.particle = ParticleType::photon;
    break;
  case 11:
    site.particle = ParticleType::electron;
    break;
  case -11:
    site.particle = ParticleType::positron;
    break;
  default:
    fatal_error(fmt::format(
      "MCPL: Encountered unexpected PDG code {} when converting to SourceSite.",
      particle_repr->pdgcode));
    break;
  }

  // Copy position and direction
  site.r.x = particle_repr->position[0];
  site.r.y = particle_repr->position[1];
  site.r.z = particle_repr->position[2];
  site.u.x = particle_repr->direction[0];
  site.u.y = particle_repr->direction[1];
  site.u.z = particle_repr->direction[2];
  // MCPL stores kinetic energy in [MeV], time in [ms]
  site.E = particle_repr->ekin * 1e6;
  site.time = particle_repr->time * 1e-3;
  site.wgt = particle_repr->weight;
  return site;
}

vector<SourceSite> mcpl_source_sites(std::string path)
{
  ensure_mcpl_ready_or_fatal();
  vector<SourceSite> sites;

  mcpl_file_t* mcpl_file = g_mcpl_api->open_file(path.c_str());
  if (!mcpl_file) {
    fatal_error(fmt::format("MCPL: Could not open file '{}'. It might be "
                            "missing, inaccessible, or not a valid MCPL file.",
      path));
  }

  size_t n_particles_in_file = g_mcpl_api->hdr_nparticles(mcpl_file);
  size_t n_skipped = 0;
  if (n_particles_in_file > 0) {
    sites.reserve(n_particles_in_file);
  }

  for (size_t i = 0; i < n_particles_in_file; ++i) {
    const mcpl_particle_repr_t* p_repr = g_mcpl_api->read(mcpl_file);
    if (!p_repr) {
      warning(fmt::format("MCPL: Read error or unexpected end of file '{}' "
                          "after reading {} of {} expected particles.",
        path, sites.size(), n_particles_in_file));
      break;
    }
    if (p_repr->pdgcode == 2112 || p_repr->pdgcode == 22 ||
        p_repr->pdgcode == 11 || p_repr->pdgcode == -11) {
      sites.push_back(mcpl_particle_to_site(p_repr));
    } else {
      n_skipped++;
    }
  }

  g_mcpl_api->close_file(mcpl_file);

  if (n_skipped > 0 && n_particles_in_file > 0) {
    double percent_skipped =
      100.0 * static_cast<double>(n_skipped) / n_particles_in_file;
    warning(fmt::format(
      "MCPL: Skipped {} of {} total particles ({:.1f}%) in file '{}' because "
      "their type is not supported by OpenMC.",
      n_skipped, n_particles_in_file, percent_skipped, path));
  }

  if (sites.empty()) {
    if (n_particles_in_file > 0) {
      fatal_error(fmt::format(
        "MCPL file '{}' contained {} particles, but none were of the supported "
        "types (neutron, photon, electron, positron). OpenMC cannot proceed "
        "without source particles.",
        path, n_particles_in_file));
    } else {
      fatal_error(fmt::format(
        "MCPL file '{}' is empty or contains no particle data.", path));
    }
  }
  return sites;
}

void write_mcpl_source_bank_internal(mcpl_outfile_t* file_id,
  span<SourceSite> local_source_bank,
  const vector<int64_t>& bank_index_all_ranks)
{
  if (mpi::master) {
    if (!file_id) {
      fatal_error("MCPL: Internal error - master rank called "
                  "write_mcpl_source_bank_internal with null file_id.");
    }
    vector<SourceSite> receive_buffer;

    for (int rank_idx = 0; rank_idx < mpi::n_procs; ++rank_idx) {
      size_t num_sites_on_rank = static_cast<size_t>(
        bank_index_all_ranks[rank_idx + 1] - bank_index_all_ranks[rank_idx]);
      if (num_sites_on_rank == 0)
        continue;

      span<const SourceSite> sites_to_write;
#ifdef OPENMC_MPI
      if (rank_idx == mpi::rank) {
        sites_to_write = openmc::span<const SourceSite>(
          local_source_bank.data(), num_sites_on_rank);
      } else {
        if (receive_buffer.size() < num_sites_on_rank) {
          receive_buffer.resize(num_sites_on_rank);
        }
        MPI_Recv(receive_buffer.data(), num_sites_on_rank, mpi::source_site,
          rank_idx, rank_idx, mpi::intracomm, MPI_STATUS_IGNORE);
        sites_to_write = openmc::span<const SourceSite>(
          receive_buffer.data(), num_sites_on_rank);
      }
#else
      sites_to_write = openmc::span<const SourceSite>(
        local_source_bank.data(), num_sites_on_rank);
#endif
      for (const auto& site : sites_to_write) {
        mcpl_particle_repr_t p_repr {};
        p_repr.position[0] = site.r.x;
        p_repr.position[1] = site.r.y;
        p_repr.position[2] = site.r.z;
        p_repr.direction[0] = site.u.x;
        p_repr.direction[1] = site.u.y;
        p_repr.direction[2] = site.u.z;
        p_repr.ekin = site.E * 1e-6;
        p_repr.time = site.time * 1e3;
        p_repr.weight = site.wgt;
        switch (site.particle) {
        case ParticleType::neutron:
          p_repr.pdgcode = 2112;
          break;
        case ParticleType::photon:
          p_repr.pdgcode = 22;
          break;
        case ParticleType::electron:
          p_repr.pdgcode = 11;
          break;
        case ParticleType::positron:
          p_repr.pdgcode = -11;
          break;
        default:
          continue;
        }
        g_mcpl_api->add_particle(file_id, &p_repr);
      }
    }
  } else {
#ifdef OPENMC_MPI
    if (!local_source_bank.empty()) {
      MPI_Send(local_source_bank.data(), local_source_bank.size(),
        mpi::source_site, 0, mpi::rank, mpi::intracomm);
    }
#endif
  }
}

void write_mcpl_source_point(const char* filename, span<SourceSite> source_bank,
  const vector<int64_t>& bank_index)
{
  ensure_mcpl_ready_or_fatal();

  std::string filename_(filename);
  const auto extension = get_file_extension(filename_);
  if (extension.empty()) {
    filename_.append(".mcpl");
  } else if (extension != "mcpl") {
    warning(fmt::format("Specified filename '{}' has an extension '.{}', but "
                        "an MCPL file (.mcpl) will be written using this name.",
      filename, extension));
  }

  mcpl_outfile_t* file_id = nullptr;

  if (mpi::master) {
    file_id = g_mcpl_api->create_outfile(filename_.c_str());
    if (!file_id) {
      fatal_error(fmt::format(
        "MCPL: Failed to create output file '{}'. Check permissions and path.",
        filename_));
    }
    std::string src_line;
    if (VERSION_DEV) {
      src_line = fmt::format("OpenMC {}.{}.{}-dev{}", VERSION_MAJOR,
        VERSION_MINOR, VERSION_RELEASE, VERSION_COMMIT_COUNT);
    } else {
      src_line = fmt::format(
        "OpenMC {}.{}.{}", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE);
    }
    g_mcpl_api->hdr_set_srcname(file_id, src_line.c_str());

    // Initialize stat:sum with -1 to indicate incomplete file (issue #3514)
    // This follows MCPL >= 2.1.0 convention for tracking simulation statistics
    // The -1 value indicates "not available" if file creation is interrupted
    if (g_mcpl_api->hdr_add_stat_sum) {
      // Using key "openmc_np1" following tkittel's recommendation
      // Initial value of -1 prevents misleading values in case of crashes
      g_mcpl_api->hdr_add_stat_sum(file_id, "openmc_np1", -1.0);
    }
  }

  write_mcpl_source_bank_internal(file_id, source_bank, bank_index);

  if (mpi::master) {
    if (file_id) {
      // Update stat:sum with actual particle count before closing (issue #3514)
      // This represents the original number of source particles in the
      // simulation (not the number of particles in the file)
      if (g_mcpl_api->hdr_add_stat_sum) {
        // Calculate total source particles from active batches
        // Per issue #3514: this should be the original number of source
        // particles, not the number written to the file
        int64_t total_source_particles =
          static_cast<int64_t>(settings::n_batches - settings::n_inactive) *
          settings::gen_per_batch * settings::n_particles;
        // Update with actual count - this overwrites the initial -1 value
        g_mcpl_api->hdr_add_stat_sum(
          file_id, "openmc_np1", static_cast<double>(total_source_particles));
      }

      g_mcpl_api->close_outfile(file_id);
    }
  }
}

} // namespace openmc

1	#include "openmc/mcpl_interface.h"
2
3	#include "openmc/bank.h"
4	#include "openmc/error.h"
5	#include "openmc/file_utils.h"
6	#include "openmc/message_passing.h"
7	#include "openmc/settings.h"
8	#include "openmc/simulation.h"
9	#include "openmc/state_point.h"
10	#include "openmc/vector.h"
11
12	#include <fmt/core.h>
13
14	#include <cctype>
15	#include <cstdint>
16	#include <cstdio>
17	#include <cstdlib>
18	#include <cstring>
19	#include <memory>
20	#include <mutex>
21	#include <stdexcept>
22	#include <string>
23
24	#ifdef _WIN32
25	#define WIN32_LEAN_AND_MEAN
26	#include <windows.h>
27	#else
28	#include <dlfcn.h>
29	#endif
30
31	// WARNING: These declarations MUST EXACTLY MATCH the structure and function
32	// signatures of the libmcpl being loaded at runtime. Any discrepancy will
33	// likely lead to crashes or incorrect behavior. This is a maintenance risk.
34	// MCPL 2.2.0
35
36	#pragma pack(push, 1)
37	struct mcpl_particle_repr_t {
38	double ekin;
39	double polarisation[3];
40	double position[3];
41	double direction[3];
42	double time;
43	double weight;
44	int32_t pdgcode;
45	uint32_t userflags;
46	};
47	#pragma pack(pop)
48
49	// Opaque struct definitions replicating the MCPL C-API to ensure ABI
50	// compatibility without including mcpl.h. These must be kept in sync.
51	struct mcpl_file_t {
52	void* internal;
53	};
54	struct mcpl_outfile_t {
55	void* internal;
56	};
57
58	// Function pointer types for the dynamically loaded MCPL library
59	using mcpl_open_file_fpt = mcpl_file_t* ()(const char filename);
60	using mcpl_hdr_nparticles_fpt = uint64_t ()(mcpl_file_t file_handle);
61	using mcpl_read_fpt = const mcpl_particle_repr_t* ()(mcpl_file_t
62	file_handle);
63	using mcpl_close_file_fpt = void ()(mcpl_file_t file_handle);
64
65	using mcpl_create_outfile_fpt = mcpl_outfile_t* ()(const char filename);
66	using mcpl_hdr_set_srcname_fpt = void (*)(
67	mcpl_outfile_t* outfile_handle, const char* srcname);
68	using mcpl_add_particle_fpt = void (*)(
69	mcpl_outfile_t* outfile_handle, const mcpl_particle_repr_t* particle);
70	using mcpl_close_outfile_fpt = void ()(mcpl_outfile_t outfile_handle);
71	using mcpl_hdr_add_data_fpt = void ()(mcpl_outfile_t outfile_handle,
72	const char* key, uint32_t datalength, const char* data);
73	using mcpl_hdr_add_stat_sum_fpt = void (*)(
74	mcpl_outfile_t* outfile_handle, const char* key, double value);
75
76	namespace openmc {
77
78	#ifdef _WIN32
79	using LibraryHandleType = HMODULE;
80	#else
81	using LibraryHandleType = void*;
82	#endif
83
UNCOV 84	std::string get_last_library_error()	×
85	{
86	#ifdef _WIN32
87	DWORD error_code = GetLastError();
88	if (error_code == 0)
89	return "No error reported by system."; // More accurate than "No error."
90	LPSTR message_buffer = nullptr;
91	size_t size =
92	FormatMessageA(FORMAT_MESSAGE_ALLOCATE_BUFFER \| FORMAT_MESSAGE_FROM_SYSTEM \|
93	FORMAT_MESSAGE_IGNORE_INSERTS,
94	NULL, error_code, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT),
95	(LPSTR)&message_buffer, 0, NULL);
96	std::string message(message_buffer, size);
97	LocalFree(message_buffer);
98	while (
99	!message.empty() && (message.back() == '\n' \|\| message.back() == '\r')) {
100	message.pop_back();
101	}
102	return message;
103	#else
104	const char* err = dlerror();	×
UNCOV 105	return err ? std::string(err) : "No error reported by dlerror.";	×
106	#endif
107	}
108
109	struct McplApi {
110	mcpl_open_file_fpt open_file;
111	mcpl_hdr_nparticles_fpt hdr_nparticles;
112	mcpl_read_fpt read;
113	mcpl_close_file_fpt close_file;
114	mcpl_create_outfile_fpt create_outfile;
115	mcpl_hdr_set_srcname_fpt hdr_set_srcname;
116	mcpl_add_particle_fpt add_particle;
117	mcpl_close_outfile_fpt close_outfile;
118	mcpl_hdr_add_data_fpt hdr_add_data;
119	mcpl_hdr_add_stat_sum_fpt hdr_add_stat_sum;
120
121	explicit McplApi(LibraryHandleType lib_handle)	81✔
122	{	81✔
123	if (!lib_handle)	81✔
124	throw std::runtime_error(	×
UNCOV 125	"MCPL library handle is null during API binding.");	×
126
127	auto load_symbol_platform = [lib_handle](const char* name) {	810✔
128	void* sym = nullptr;	810✔
129	#ifdef _WIN32
130	sym = (void*)GetProcAddress(lib_handle, name);
131	#else
132	sym = dlsym(lib_handle, name);	810✔
133	#endif
134	if (!sym) {	810✔
135	throw std::runtime_error(	×
136	fmt::format("Failed to load MCPL symbol '{}': {}", name,	×
UNCOV 137	get_last_library_error()));	×
138	}
139	return sym;	810✔
140	};	81✔
141
142	open_file = reinterpret_cast<mcpl_open_file_fpt>(	81✔
143	load_symbol_platform("mcpl_open_file"));	81✔
144	hdr_nparticles = reinterpret_cast<mcpl_hdr_nparticles_fpt>(	81✔
145	load_symbol_platform("mcpl_hdr_nparticles"));	81✔
146	read = reinterpret_cast<mcpl_read_fpt>(load_symbol_platform("mcpl_read"));	81✔
147	close_file = reinterpret_cast<mcpl_close_file_fpt>(	81✔
148	load_symbol_platform("mcpl_close_file"));	81✔
149	create_outfile = reinterpret_cast<mcpl_create_outfile_fpt>(	81✔
150	load_symbol_platform("mcpl_create_outfile"));	81✔
151	hdr_set_srcname = reinterpret_cast<mcpl_hdr_set_srcname_fpt>(	81✔
152	load_symbol_platform("mcpl_hdr_set_srcname"));	81✔
153	add_particle = reinterpret_cast<mcpl_add_particle_fpt>(	81✔
154	load_symbol_platform("mcpl_add_particle"));	81✔
155	close_outfile = reinterpret_cast<mcpl_close_outfile_fpt>(	81✔
156	load_symbol_platform("mcpl_close_outfile"));	81✔
157
158	// Try to load mcpl_hdr_add_data (available in MCPL >= 2.1.0)
159	// Set to nullptr if not available for graceful fallback
160	try {
161	hdr_add_data = reinterpret_cast<mcpl_hdr_add_data_fpt>(	81✔
162	load_symbol_platform("mcpl_hdr_add_data"));	81✔
NEW 163	} catch (const std::runtime_error&) {	×
NEW 164	hdr_add_data = nullptr;	×
NEW 165	}	×
166
167	// Try to load mcpl_hdr_add_stat_sum (available in MCPL >= 2.1.0)
168	// Set to nullptr if not available for graceful fallback
169	try {
170	hdr_add_stat_sum = reinterpret_cast<mcpl_hdr_add_stat_sum_fpt>(	81✔
171	load_symbol_platform("mcpl_hdr_add_stat_sum"));	81✔
172	} catch (const std::runtime_error&) {	×
173	hdr_add_stat_sum = nullptr;	×
UNCOV 174	}	×
175	}	81✔
176	};
177
178	static LibraryHandleType g_mcpl_lib_handle = nullptr;
179	static std::unique_ptr<McplApi> g_mcpl_api;
180	static bool g_mcpl_init_attempted = false;
181	static bool g_mcpl_successfully_loaded = false;
182	static std::string g_mcpl_load_error_msg;
183	static std::once_flag g_mcpl_init_flag;
184
UNCOV 185	void append_error(std::string& existing_msg, const std::string& new_error)	×
186	{
187	if (!existing_msg.empty()) {	×
UNCOV 188	existing_msg += "; ";	×
189	}
UNCOV 190	existing_msg += new_error;	×
191	}
192
193	void initialize_mcpl_interface_impl()	81✔
194	{
195	g_mcpl_init_attempted = true;	81✔
196	g_mcpl_load_error_msg.clear();	81✔
197
198	// Try mcpl-config
199	if (!g_mcpl_lib_handle) {	81✔
200	FILE* pipe = nullptr;	81✔
201	#ifdef _WIN32
202	pipe = _popen("mcpl-config --show libpath", "r");
203	#else
204	pipe = popen("mcpl-config --show libpath 2>/dev/null", "r");	81✔
205	#endif
206	if (pipe) {	81✔
207	char buffer[512];
208	if (fgets(buffer, sizeof(buffer), pipe) != nullptr) {	81✔
209	std::string shlibpath = buffer;	81✔
210	// Remove trailing whitespace
211	while (!shlibpath.empty() &&	324✔
212	std::isspace(static_cast<unsigned char>(shlibpath.back()))) {	162✔
213	shlibpath.pop_back();	81✔
214	}
215
216	if (!shlibpath.empty()) {	81✔
217	#ifdef _WIN32
218	g_mcpl_lib_handle = LoadLibraryA(shlibpath.c_str());
219	#else
220	g_mcpl_lib_handle = dlopen(shlibpath.c_str(), RTLD_LAZY);	81✔
221	#endif
222	if (!g_mcpl_lib_handle) {	81✔
223	append_error(	×
224	g_mcpl_load_error_msg, fmt::format("From mcpl-config ({}): {}",	×
UNCOV 225	shlibpath, get_last_library_error()));	×
226	}
227	}
228	}	81✔
229	#ifdef _WIN32
230	_pclose(pipe);
231	#else
232	pclose(pipe);	81✔
233	#endif
234	} else { // pipe failed to open
UNCOV 235	append_error(g_mcpl_load_error_msg,	×
236	"mcpl-config command not found or failed to execute");
237	}
238	}
239
240	// Try standard library names
241	if (!g_mcpl_lib_handle) {	81✔
242	#ifdef _WIN32
243	const char* standard_names[] = {"mcpl.dll", "libmcpl.dll"};
244	#else
UNCOV 245	const char* standard_names[] = {"libmcpl.so", "libmcpl.dylib"};	×
246	#endif
UNCOV 247	for (const char* name : standard_names) {	×
248	#ifdef _WIN32
249	g_mcpl_lib_handle = LoadLibraryA(name);
250	#else
UNCOV 251	g_mcpl_lib_handle = dlopen(name, RTLD_LAZY);	×
252	#endif
253	if (g_mcpl_lib_handle)	×
UNCOV 254	break;	×
255	}
256	if (!g_mcpl_lib_handle) {	×
257	append_error(	×
258	g_mcpl_load_error_msg, fmt::format("Using standard names (e.g. {}): {}",	×
UNCOV 259	standard_names[0], get_last_library_error()));	×
260	}
261	}
262
263	if (!g_mcpl_lib_handle) {	81✔
264	if (mpi::master) {	×
UNCOV 265	warning(fmt::format("MCPL library could not be loaded. MCPL-dependent "	×
266	"features will be unavailable. Load attempts: {}",
267	g_mcpl_load_error_msg.empty()	×
UNCOV 268	? "No specific error during load attempts."	×
269	: g_mcpl_load_error_msg));
270	}
271	g_mcpl_successfully_loaded = false;	×
UNCOV 272	return;	×
273	}
274
275	try {
276	g_mcpl_api = std::make_unique<McplApi>(g_mcpl_lib_handle);	81✔
277	g_mcpl_successfully_loaded = true;	81✔
278	// Do not call dlclose/FreeLibrary at exit. Leaking the handle is safer
279	// and standard practice for libraries used for the application's lifetime.
280	} catch (const std::runtime_error& e) {	×
281	append_error(g_mcpl_load_error_msg,	×
282	fmt::format(	×
283	"MCPL library loaded, but failed to bind symbols: {}", e.what()));	×
284	if (mpi::master) {	×
UNCOV 285	warning(g_mcpl_load_error_msg);	×
286	}
287	#ifdef _WIN32
288	FreeLibrary(g_mcpl_lib_handle);
289	#else
UNCOV 290	dlclose(g_mcpl_lib_handle);	×
291	#endif
292	g_mcpl_lib_handle = nullptr;	×
293	g_mcpl_successfully_loaded = false;	×
UNCOV 294	}	×
295	}
296
297	void initialize_mcpl_interface_if_needed()	169✔
298	{
299	std::call_once(g_mcpl_init_flag, initialize_mcpl_interface_impl);	169✔
300	}	169✔
301
302	bool is_mcpl_interface_available()	33✔
303	{
304	initialize_mcpl_interface_if_needed();	33✔
305	return g_mcpl_successfully_loaded;	33✔
306	}
307
308	inline void ensure_mcpl_ready_or_fatal()	136✔
309	{
310	initialize_mcpl_interface_if_needed();	136✔
311	if (!g_mcpl_successfully_loaded) {	136✔
UNCOV 312	fatal_error("MCPL functionality is required, but the MCPL library is not "	×
313	"available or failed to initialize. Please ensure MCPL is "
314	"installed and its library can be found (e.g., via PATH on "
315	"Windows, LD_LIBRARY_PATH on Linux, or DYLD_LIBRARY_PATH on "
316	"macOS). You can often install MCPL with 'pip install mcpl' or "
317	"'conda install mcpl'.");
318	}
319	}	136✔
320
321	SourceSite mcpl_particle_to_site(const mcpl_particle_repr_t* particle_repr)	32,000✔
322	{
323	SourceSite site;	32,000✔
324	switch (particle_repr->pdgcode) {	32,000✔
325	case 2112:	32,000✔
326	site.particle = ParticleType::neutron;	32,000✔
327	break;	32,000✔
328	case 22:	×
329	site.particle = ParticleType::photon;	×
330	break;	×
331	case 11:	×
332	site.particle = ParticleType::electron;	×
333	break;	×
334	case -11:	×
335	site.particle = ParticleType::positron;	×
336	break;	×
337	default:	×
UNCOV 338	fatal_error(fmt::format(	×
339	"MCPL: Encountered unexpected PDG code {} when converting to SourceSite.",
UNCOV 340	particle_repr->pdgcode));	×
341	break;
342	}
343
344	// Copy position and direction
345	site.r.x = particle_repr->position[0];	32,000✔
346	site.r.y = particle_repr->position[1];	32,000✔
347	site.r.z = particle_repr->position[2];	32,000✔
348	site.u.x = particle_repr->direction[0];	32,000✔
349	site.u.y = particle_repr->direction[1];	32,000✔
350	site.u.z = particle_repr->direction[2];	32,000✔
351	// MCPL stores kinetic energy in [MeV], time in [ms]
352	site.E = particle_repr->ekin * 1e6;	32,000✔
353	site.time = particle_repr->time * 1e-3;	32,000✔
354	site.wgt = particle_repr->weight;	32,000✔
355	return site;	32,000✔
356	}
357
358	vector<SourceSite> mcpl_source_sites(std::string path)	32✔
359	{
360	ensure_mcpl_ready_or_fatal();	32✔
361	vector<SourceSite> sites;	32✔
362
363	mcpl_file_t* mcpl_file = g_mcpl_api->open_file(path.c_str());	32✔
364	if (!mcpl_file) {	32✔
UNCOV 365	fatal_error(fmt::format("MCPL: Could not open file '{}'. It might be "	×
366	"missing, inaccessible, or not a valid MCPL file.",
367	path));
368	}
369
370	size_t n_particles_in_file = g_mcpl_api->hdr_nparticles(mcpl_file);	32✔
371	size_t n_skipped = 0;	32✔
372	if (n_particles_in_file > 0) {	32✔
373	sites.reserve(n_particles_in_file);	32✔
374	}
375
376	for (size_t i = 0; i < n_particles_in_file; ++i) {	32,032✔
377	const mcpl_particle_repr_t* p_repr = g_mcpl_api->read(mcpl_file);	32,000✔
378	if (!p_repr) {	32,000✔
UNCOV 379	warning(fmt::format("MCPL: Read error or unexpected end of file '{}' "	×
380	"after reading {} of {} expected particles.",
381	path, sites.size(), n_particles_in_file));	×
UNCOV 382	break;	×
383	}
384	if (p_repr->pdgcode == 2112 \|\| p_repr->pdgcode == 22 \|\|	32,000✔
UNCOV 385	p_repr->pdgcode == 11 \|\| p_repr->pdgcode == -11) {	×
386	sites.push_back(mcpl_particle_to_site(p_repr));	32,000✔
387	} else {
UNCOV 388	n_skipped++;	×
389	}
390	}
391
392	g_mcpl_api->close_file(mcpl_file);	32✔
393
394	if (n_skipped > 0 && n_particles_in_file > 0) {	32✔
395	double percent_skipped =	×
396	100.0 * static_cast<double>(n_skipped) / n_particles_in_file;	×
UNCOV 397	warning(fmt::format(	×
398	"MCPL: Skipped {} of {} total particles ({:.1f}%) in file '{}' because "
399	"their type is not supported by OpenMC.",
400	n_skipped, n_particles_in_file, percent_skipped, path));
401	}
402
403	if (sites.empty()) {	32✔
404	if (n_particles_in_file > 0) {	×
UNCOV 405	fatal_error(fmt::format(	×
406	"MCPL file '{}' contained {} particles, but none were of the supported "
407	"types (neutron, photon, electron, positron). OpenMC cannot proceed "
408	"without source particles.",
409	path, n_particles_in_file));
410	} else {
UNCOV 411	fatal_error(fmt::format(	×
412	"MCPL file '{}' is empty or contains no particle data.", path));
413	}
414	}
415	return sites;	64✔
UNCOV 416	}	×
417
418	void write_mcpl_source_bank_internal(mcpl_outfile_t* file_id,	104✔
419	span<SourceSite> local_source_bank,
420	const vector<int64_t>& bank_index_all_ranks)
421	{
422	if (mpi::master) {	104✔
423	if (!file_id) {	99✔
UNCOV 424	fatal_error("MCPL: Internal error - master rank called "	×
425	"write_mcpl_source_bank_internal with null file_id.");
426	}
427	vector<SourceSite> receive_buffer;	99✔
428
429	for (int rank_idx = 0; rank_idx < mpi::n_procs; ++rank_idx) {	203✔
430	size_t num_sites_on_rank = static_cast<size_t>(
431	bank_index_all_ranks[rank_idx + 1] - bank_index_all_ranks[rank_idx]);	104✔
432	if (num_sites_on_rank == 0)	104✔
433	continue;	11✔
434
435	span<const SourceSite> sites_to_write;	93✔
436	#ifdef OPENMC_MPI
437	if (rank_idx == mpi::rank) {	45✔
438	sites_to_write = openmc::span<const SourceSite>(	40✔
439	local_source_bank.data(), num_sites_on_rank);	40✔
440	} else {
441	if (receive_buffer.size() < num_sites_on_rank) {	5✔
442	receive_buffer.resize(num_sites_on_rank);	5✔
443	}
444	MPI_Recv(receive_buffer.data(), num_sites_on_rank, mpi::source_site,	5✔
445	rank_idx, rank_idx, mpi::intracomm, MPI_STATUS_IGNORE);
446	sites_to_write = openmc::span<const SourceSite>(	5✔
447	receive_buffer.data(), num_sites_on_rank);	5✔
448	}
449	#else
450	sites_to_write = openmc::span<const SourceSite>(	48✔
451	local_source_bank.data(), num_sites_on_rank);	48✔
452	#endif
453	for (const auto& site : sites_to_write) {	46,414✔
454	mcpl_particle_repr_t p_repr {};	46,321✔
455	p_repr.position[0] = site.r.x;	46,321✔
456	p_repr.position[1] = site.r.y;	46,321✔
457	p_repr.position[2] = site.r.z;	46,321✔
458	p_repr.direction[0] = site.u.x;	46,321✔
459	p_repr.direction[1] = site.u.y;	46,321✔
460	p_repr.direction[2] = site.u.z;	46,321✔
461	p_repr.ekin = site.E * 1e-6;	46,321✔
462	p_repr.time = site.time * 1e3;	46,321✔
463	p_repr.weight = site.wgt;	46,321✔
464	switch (site.particle) {	46,321✔
465	case ParticleType::neutron:	46,321✔
466	p_repr.pdgcode = 2112;	46,321✔
467	break;	46,321✔
468	case ParticleType::photon:	×
469	p_repr.pdgcode = 22;	×
470	break;	×
471	case ParticleType::electron:	×
472	p_repr.pdgcode = 11;	×
473	break;	×
474	case ParticleType::positron:	×
475	p_repr.pdgcode = -11;	×
476	break;	×
477	default:	×
UNCOV 478	continue;	×
479	}
480	g_mcpl_api->add_particle(file_id, &p_repr);	46,321✔
481	}
482	}
483	} else {	99✔
484	#ifdef OPENMC_MPI
485	if (!local_source_bank.empty()) {	5✔
486	MPI_Send(local_source_bank.data(), local_source_bank.size(),	5✔
487	mpi::source_site, 0, mpi::rank, mpi::intracomm);
488	}
489	#endif
490	}
491	}	104✔
492
493	void write_mcpl_source_point(const char* filename, span<SourceSite> source_bank,	104✔
494	const vector<int64_t>& bank_index)
495	{
496	ensure_mcpl_ready_or_fatal();	104✔
497
498	std::string filename_(filename);	104✔
499	const auto extension = get_file_extension(filename_);	104✔
500	if (extension.empty()) {	104✔
UNCOV 501	filename_.append(".mcpl");	×
502	} else if (extension != "mcpl") {	104✔
UNCOV 503	warning(fmt::format("Specified filename '{}' has an extension '.{}', but "	×
504	"an MCPL file (.mcpl) will be written using this name.",
505	filename, extension));
506	}
507
508	mcpl_outfile_t* file_id = nullptr;	104✔
509
510	if (mpi::master) {	104✔
511	file_id = g_mcpl_api->create_outfile(filename_.c_str());	99✔
512	if (!file_id) {	99✔
UNCOV 513	fatal_error(fmt::format(	×
514	"MCPL: Failed to create output file '{}'. Check permissions and path.",
515	filename_));
516	}
517	std::string src_line;	99✔
518	if (VERSION_DEV) {
519	src_line = fmt::format("OpenMC {}.{}.{}-dev{}", VERSION_MAJOR,	180✔
520	VERSION_MINOR, VERSION_RELEASE, VERSION_COMMIT_COUNT);	99✔
521	} else {
522	src_line = fmt::format(
523	"OpenMC {}.{}.{}", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE);
524	}
525	g_mcpl_api->hdr_set_srcname(file_id, src_line.c_str());	99✔
526
527	// Initialize stat:sum with -1 to indicate incomplete file (issue #3514)
528	// This follows MCPL >= 2.1.0 convention for tracking simulation statistics
529	// The -1 value indicates "not available" if file creation is interrupted
530	if (g_mcpl_api->hdr_add_stat_sum) {	99✔
531	// Using key "openmc_np1" following tkittel's recommendation
532	// Initial value of -1 prevents misleading values in case of crashes
533	g_mcpl_api->hdr_add_stat_sum(file_id, "openmc_np1", -1.0);	99✔
534	}
535	}	99✔
536
537	write_mcpl_source_bank_internal(file_id, source_bank, bank_index);	104✔
538
539	if (mpi::master) {	104✔
540	if (file_id) {	99✔
541	// Update stat:sum with actual particle count before closing (issue #3514)
542	// This represents the original number of source particles in the
543	// simulation (not the number of particles in the file)
544	if (g_mcpl_api->hdr_add_stat_sum) {	99✔
545	// Calculate total source particles from active batches
546	// Per issue #3514: this should be the original number of source
547	// particles, not the number written to the file
548	int64_t total_source_particles =	99✔
549	static_cast<int64_t>(settings::n_batches - settings::n_inactive) *	99✔
550	settings::gen_per_batch * settings::n_particles;	99✔
551	// Update with actual count - this overwrites the initial -1 value
552	g_mcpl_api->hdr_add_stat_sum(	99✔
553	file_id, "openmc_np1", static_cast<double>(total_source_particles));
554	}
555
556	g_mcpl_api->close_outfile(file_id);	99✔
557	}
558	}
559	}	104✔
560
561	} // namespace openmc

openmc-dev / openmc / 16898768927

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

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