crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 14:15:42 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Move remainder of MCPL-related code to mcpl_interface.h/cpp

Browse source
This commit is contained in:
Paul Romano 2022-12-22 12:27:56 -06:00
parent a3065b479c
commit 361b3486ef
7 changed files with 201 additions and 179 deletions
Download patch file Download diff file Expand all files Collapse all files

1
CMakeLists.txt
View file

@ -164,6 +164,7 @@ endif()
#===============================================================================
if (OPENMC_USE_MCPL)
find_package(MCPL REQUIRED)
message(STATUS "Found MCPL: ${MCPL_DIR} (found version \"${MCPL_VERSION}\")")
endif()
#===============================================================================

19
include/openmc/mcpl_interface.h
View file

@ -6,20 +6,31 @@
#include <string>
#include <vector>
#ifdef OPENMC_MCPL
#include <mcpl.h>
#endif
namespace openmc {
//==============================================================================
// Constants
//==============================================================================
extern "C" const bool MCPL_ENABLED;
//==============================================================================
// Functions
//==============================================================================
//! Get a vector of source sites from an MCPL file
//
//! \param[in] path Path to MCPL file
//! \return Vector of source sites
vector<SourceSite> mcpl_source_sites(std::string path);
//! Write an MCPL source file
//
//! \param[in] filename Path to MCPL file
//! \param[in] surf_source_bank Whether to use the surface source bank
void write_mcpl_source_point(
const char* filename, bool surf_source_bank = false);
} // namespace openmc
#endif // OPENMC_MCPL_INTERFACE_H

10
include/openmc/state_point.h
View file

@ -9,10 +9,6 @@
#include "openmc/particle.h"
#include "openmc/vector.h"
#ifdef OPENMC_MCPL
#include <mcpl.h>
#endif
namespace openmc {
void load_state_point();
@ -25,11 +21,5 @@ void write_tally_results_nr(hid_t file_id);
void restart_set_keff();
void write_unstructured_mesh_results();
#ifdef OPENMC_MCPL
void write_mcpl_source_point(
const char* filename, bool surf_source_bank = false);
void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank);
#endif
} // namespace openmc
#endif // OPENMC_STATE_POINT_H

185
src/mcpl_interface.cpp
View file

@ -1,6 +1,13 @@
#include "openmc/mcpl_interface.h"
#include "openmc/bank.h"
#include "openmc/error.h"
#include "openmc/message_passing.h"
#include "openmc/settings.h"
#include "openmc/simulation.h"
#include "openmc/state_point.h"
#include <fmt/core.h>
#ifdef OPENMC_MCPL
#include <mcpl.h>
@ -8,14 +15,22 @@
namespace openmc {
//==============================================================================
// Constants
//==============================================================================
#ifdef OPENMC_MCPL
const bool MCPL_ENABLED = true;
#else
const bool MCPL_ENABLED = false;
#endif
//==============================================================================
// Functions
//==============================================================================
#ifdef OPENMC_MCPL
SourceSite mcpl_particle_to_site(const mcpl_particle* particle)
SourceSite mcpl_particle_to_site(const mcpl_particle_t* particle)
{
SourceSite site;
@ -35,40 +50,42 @@ SourceSite mcpl_particle_to_site(const mcpl_particle* particle)
}
// Copy position and direction
site.r.x = mcpl_particle->position[0];
site.r.y = mcpl_particle->position[1];
site.r.z = mcpl_particle->position[2];
site.u.x = mcpl_particle->direction[0];
site.u.y = mcpl_particle->direction[1];
site.u.z = mcpl_particle->direction[2];
site.r.x = particle->position[0];
site.r.y = particle->position[1];
site.r.z = particle->position[2];
site.u.x = particle->direction[0];
site.u.y = particle->direction[1];
site.u.z = particle->direction[2];
// mcpl stores kinetic energy in MeV
site.E = mcpl_particle->ekin * 1e6;
site.E = particle->ekin * 1e6;
// mcpl stores time in ms
site.time = mcpl_particle->time * 1e-3;
site.wgt = mcpl_particle->weight;
site.time = particle->time * 1e-3;
site.wgt = particle->weight;
return site;
}
#endif
//==============================================================================
vector<SourceSite> mcpl_source_sites(std::string path)
{
vector<SourceSite> sites;
#ifdef OPENMC_MCPL
// Open MCPL file and determine number of particles
auto mcpl_file = mcpl_open_file(path.c_str());
size_t n_sites = mcpl_hdr_nparticles(mcpl_file);
for (int i = 0; i < n_sites; i++) {
SourceSite site;
// Extract particle from mcpl-file, checking if it is a neutron, photon,
// electron, or positron. Otherwise skip.
const mcpl_particle_t* particle;
int pdg = 0;
while (pdg != 2112 && pdg != 22 && pdg != 11 && pdg != -11) {
particle = mcpl_read(mcpl_file);
pdg = mcpl_particle->pdgcode;
pdg = particle->pdgcode;
}
// Convert to source site and add to vector
@ -76,7 +93,7 @@ vector<SourceSite> mcpl_source_sites(std::string path)
}
// Check that some sites were read
if (sites_.empty()) {
if (sites.empty()) {
fatal_error("MCPL file contained no neutron, photon, electron, or positron "
"source particles.");
}
@ -90,4 +107,144 @@ vector<SourceSite> mcpl_source_sites(std::string path)
return sites;
}
//==============================================================================
#ifdef OPENMC_MCPL
void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank)
{
int64_t dims_size = settings::n_particles;
int64_t count_size = simulation::work_per_rank;
// Set vectors for source bank and starting bank index of each process
vector<int64_t>* bank_index = &simulation::work_index;
vector<SourceSite>* source_bank = &simulation::source_bank;
vector<int64_t> surf_source_index_vector;
vector<SourceSite> surf_source_bank_vector;
if (surf_source_bank) {
surf_source_index_vector = calculate_surf_source_size();
dims_size = surf_source_index_vector[mpi::n_procs];
count_size = simulation::surf_source_bank.size();
bank_index = &surf_source_index_vector;
// Copy data in a SharedArray into a vector.
surf_source_bank_vector.resize(count_size);
surf_source_bank_vector.assign(simulation::surf_source_bank.data(),
simulation::surf_source_bank.data() + count_size);
source_bank = &surf_source_bank_vector;
}
if (mpi::master) {
// Particles are writeen to disk from the master node only
// Save source bank sites since the array is overwritten below
#ifdef OPENMC_MPI
vector<SourceSite> temp_source {source_bank->begin(), source_bank->end()};
#endif
// loop over the other nodes and receive data - then write those.
for (int i = 0; i < mpi::n_procs; ++i) {
// number of particles for node node i
size_t count[] {
static_cast<size_t>((*bank_index)[i + 1] - (*bank_index)[i])};
#ifdef OPENMC_MPI
if (i > 0)
MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i,
mpi::intracomm, MPI_STATUS_IGNORE);
#endif
// now write the source_bank data again.
for (vector<SourceSite>::iterator _site = source_bank->begin();
_site != source_bank->end(); _site++) {
// particle is now at the iterator
// write it to the mcpl-file
mcpl_particle_t p;
p.position[0] = _site->r.x;
p.position[1] = _site->r.y;
p.position[2] = _site->r.z;
// mcpl requires that the direction vector is unit length
// which is also the case in openmc
p.direction[0] = _site->u.x;
p.direction[1] = _site->u.y;
p.direction[2] = _site->u.z;
// mcpl stores kinetic energy in MeV
p.ekin = _site->E * 1e-6;
p.time = _site->time * 1e3;
p.weight = _site->wgt;
switch (_site->particle) {
case ParticleType::neutron:
p.pdgcode = 2112;
break;
case ParticleType::photon:
p.pdgcode = 22;
break;
case ParticleType::electron:
p.pdgcode = 11;
break;
case ParticleType::positron:
p.pdgcode = -11;
break;
}
mcpl_add_particle(file_id, &p);
}
}
#ifdef OPENMC_MPI
// Restore state of source bank
std::copy(temp_source.begin(), temp_source.end(), source_bank->begin());
#endif
} else {
#ifdef OPENMC_MPI
MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank,
mpi::intracomm);
#endif
}
}
#endif
//==============================================================================
void write_mcpl_source_point(const char* filename, bool surf_source_bank)
{
std::string filename_;
if (filename) {
filename_ = filename;
} else {
// Determine width for zero padding
int w = std::to_string(settings::n_max_batches).size();
filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output,
simulation::current_batch, w);
}
#ifdef OPENMC_MCPL
mcpl_outfile_t file_id;
std::string line;
if (mpi::master) {
file_id = mcpl_create_outfile(filename_.c_str());
if (VERSION_DEV) {
line = fmt::format("OpenMC {0}.{1}.{2}-development", VERSION_MAJOR,
VERSION_MINOR, VERSION_RELEASE);
} else {
line = fmt::format(
"OpenMC {0}.{1}.{2}", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE);
}
mcpl_hdr_set_srcname(file_id, line.c_str());
}
write_mcpl_source_bank(file_id, surf_source_bank);
if (mpi::master) {
mcpl_close_outfile(file_id);
}
#endif
}
} // namespace openmc

18
src/settings.cpp
View file

@ -433,10 +433,8 @@ void read_settings_xml()
if (check_for_node(node, "file")) {
auto path = get_node_value(node, "file", false, true);
if (ends_with(path, ".mcpl") || ends_with(path, ".mcpl.gz")) {
#ifdef OPENMC_MCPL
auto sites = mcpl_source_sites(path);
model::external_sources.push_back(make_unique<FileSource>(sites));
#endif
} else {
model::external_sources.push_back(make_unique<FileSource>(path));
}
@ -659,6 +657,12 @@ void read_settings_xml()
}
if (check_for_node(node_sp, "mcpl")) {
source_mcpl_write = get_node_value_bool(node_sp, "mcpl");
// Make sure MCPL support is enabled
if (source_mcpl_write && !MCPL_ENABLED) {
fatal_error(
"Your build of OpenMC does not support writing MCPL source files.");
}
}
if (check_for_node(node_sp, "overwrite_latest")) {
source_latest = get_node_value_bool(node_sp, "overwrite_latest");
@ -690,11 +694,15 @@ void read_settings_xml()
max_surface_particles =
std::stoll(get_node_value(node_ssw, "max_particles"));
}
#ifdef OPENMC_MCPL
if (check_for_node(node_ssw, "mcpl")) {
surf_mcpl_write = true;
surf_mcpl_write = get_node_value_bool(node_ssw, "mcpl");
// Make sure MCPL support is enabled
if (surf_mcpl_write && !MCPL_ENABLED) {
fatal_error("Your build of OpenMC does not support writing MCPL "
"surface source files.");
}
}
#endif
}
// If source is not separate and is to be written out in the statepoint file,

7
src/simulation.cpp
View file

@ -8,6 +8,7 @@
#include "openmc/event.h"
#include "openmc/geometry_aux.h"
#include "openmc/material.h"
#include "openmc/mcpl_interface.h"
#include "openmc/message_passing.h"
#include "openmc/nuclide.h"
#include "openmc/output.h"
@ -389,9 +390,7 @@ void finalize_batch()
if (settings::run_mode == RunMode::EIGENVALUE) {
// Write out a separate source point if it's been specified for this batch
#ifdef OPENMC_MCPL
if (!settings::source_mcpl_write) {
#endif
if (contains(settings::sourcepoint_batch, simulation::current_batch) &&
settings::source_write && settings::source_separate) {
write_source_point(nullptr);
@ -402,7 +401,6 @@ void finalize_batch()
auto filename = settings::path_output + "source.h5";
write_source_point(filename.c_str());
}
#ifdef OPENMC_MCPL
} else {
if (contains(settings::sourcepoint_batch, simulation::current_batch) &&
settings::source_mcpl_write && settings::source_separate) {
@ -415,7 +413,6 @@ void finalize_batch()
write_mcpl_source_point(filename.c_str());
}
}
#endif
}
// Write out surface source if requested.
@ -424,13 +421,11 @@ void finalize_batch()
auto filename = settings::path_output + "surface_source.h5";
write_source_point(filename.c_str(), true);
}
#ifdef OPENMC_MCPL
if (settings::surf_mcpl_write &&
simulation::current_batch == settings::n_batches) {
auto filename = settings::path_output + "surface_source.mcpl";
write_mcpl_source_point(filename.c_str(), true);
}
#endif
}
void initialize_generation()

140
src/state_point.cpp
View file

@ -28,10 +28,6 @@
#include "openmc/timer.h"
#include "openmc/vector.h"
#ifdef OPENMC_MCPL
#include <mcpl.h>
#endif
namespace openmc {
extern "C" int openmc_statepoint_write(const char* filename, bool* write_source)
@ -603,43 +599,6 @@ void write_source_point(const char* filename, bool surf_source_bank)
file_close(file_id);
}
#ifdef OPENMC_MCPL
void write_mcpl_source_point(const char* filename, bool surf_source_bank)
{
std::string filename_;
if (filename) {
filename_ = filename;
} else {
// Determine width for zero padding
int w = std::to_string(settings::n_max_batches).size();
filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output,
simulation::current_batch, w);
}
mcpl_outfile_t file_id;
std::string line;
if (mpi::master) {
file_id = mcpl_create_outfile(filename_.c_str());
if (VERSION_DEV) {
line = fmt::format("OpenMC {0}.{1}.{2}-development", VERSION_MAJOR,
VERSION_MINOR, VERSION_RELEASE);
} else {
line = fmt::format(
"OpenMC {0}.{1}.{2}", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE);
}
mcpl_hdr_set_srcname(file_id, line.c_str());
}
write_mcpl_source_bank(file_id, surf_source_bank);
if (mpi::master) {
mcpl_close_outfile(file_id);
}
}
#endif
void write_source_bank(hid_t group_id, bool surf_source_bank)
{
hid_t banktype = h5banktype();
@ -756,105 +715,6 @@ void write_source_bank(hid_t group_id, bool surf_source_bank)
H5Tclose(banktype);
}
#ifdef OPENMC_MCPL
void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank)
{
int64_t dims_size = settings::n_particles;
int64_t count_size = simulation::work_per_rank;
// Set vectors for source bank and starting bank index of each process
vector<int64_t>* bank_index = &simulation::work_index;
vector<SourceSite>* source_bank = &simulation::source_bank;
vector<int64_t> surf_source_index_vector;
vector<SourceSite> surf_source_bank_vector;
if (surf_source_bank) {
surf_source_index_vector = calculate_surf_source_size();
dims_size = surf_source_index_vector[mpi::n_procs];
count_size = simulation::surf_source_bank.size();
bank_index = &surf_source_index_vector;
// Copy data in a SharedArray into a vector.
surf_source_bank_vector.resize(count_size);
surf_source_bank_vector.assign(simulation::surf_source_bank.data(),
simulation::surf_source_bank.data() + count_size);
source_bank = &surf_source_bank_vector;
}
if (mpi::master) {
// Particles are writeen to disk from the master node only
// Save source bank sites since the array is overwritten below
#ifdef OPENMC_MPI
vector<SourceSite> temp_source {source_bank->begin(), source_bank->end()};
#endif
// loop over the other nodes and receive data - then write those.
for (int i = 0; i < mpi::n_procs; ++i) {
// number of particles for node node i
size_t count[] {
static_cast<size_t>((*bank_index)[i + 1] - (*bank_index)[i])};
#ifdef OPENMC_MPI
if (i > 0)
MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i,
mpi::intracomm, MPI_STATUS_IGNORE);
#endif
// now write the source_bank data again.
for (vector<SourceSite>::iterator _site = source_bank->begin();
_site != source_bank->end(); _site++) {
// particle is now at the iterator
// write it to the mcpl-file
mcpl_particle_t p;
p.position[0] = _site->r.x;
p.position[1] = _site->r.y;
p.position[2] = _site->r.z;
// mcpl requires that the direction vector is unit length
// which is also the case in openmc
p.direction[0] = _site->u.x;
p.direction[1] = _site->u.y;
p.direction[2] = _site->u.z;
// mcpl stores kinetic energy in MeV
p.ekin = _site->E * 1e-6;
p.time = _site->time * 1e3;
p.weight = _site->wgt;
switch (_site->particle) {
case ParticleType::neutron:
p.pdgcode = 2112;
break;
case ParticleType::photon:
p.pdgcode = 22;
break;
case ParticleType::electron:
p.pdgcode = 11;
break;
case ParticleType::positron:
p.pdgcode = -11;
break;
}
mcpl_add_particle(file_id, &p);
}
}
#ifdef OPENMC_MPI
// Restore state of source bank
std::copy(temp_source.begin(), temp_source.end(), source_bank->begin());
#endif
} else {
#ifdef OPENMC_MPI
MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank,
mpi::intracomm);
#endif
}
}
#endif
// Determine member names of a compound HDF5 datatype
std::string dtype_member_names(hid_t dtype_id)
{