From 968e9a95d4393490c0aa7a5f98c4638b7fee8819 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 14 Feb 2019 12:43:05 -0600 Subject: [PATCH] Convert particle restart to C++ --- CMakeLists.txt | 3 +- include/openmc/capi.h | 1 - include/openmc/hdf5_interface.h | 13 ++- include/openmc/particle_restart.h | 10 +++ include/openmc/simulation.h | 3 + src/constants.F90 | 12 --- src/main.cpp | 4 +- src/particle_restart.F90 | 142 ------------------------------ src/particle_restart.cpp | 109 +++++++++++++++++++++++ src/simulation.cpp | 2 - 10 files changed, 137 insertions(+), 162 deletions(-) create mode 100644 include/openmc/particle_restart.h delete mode 100644 src/particle_restart.F90 create mode 100644 src/particle_restart.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 2c625517a..43e250c34 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -319,10 +319,8 @@ add_library(libopenmc SHARED src/mgxs_interface.F90 src/nuclide_header.F90 src/particle_header.F90 - src/particle_restart.F90 src/photon_header.F90 src/pugixml/pugixml_f.F90 - src/random_lcg.F90 src/relaxng src/settings.F90 src/simulation_header.F90 @@ -370,6 +368,7 @@ add_library(libopenmc SHARED src/nuclide.cpp src/output.cpp src/particle.cpp + src/particle_restart.cpp src/photon.cpp src/physics.cpp src/physics_common.cpp diff --git a/include/openmc/capi.h b/include/openmc/capi.h index 80d431063..5ec852d2d 100644 --- a/include/openmc/capi.h +++ b/include/openmc/capi.h @@ -84,7 +84,6 @@ extern "C" { int openmc_meshsurface_filter_set_mesh(int32_t index, int32_t index_mesh); int openmc_next_batch(int* status); int openmc_nuclide_name(int index, const char** name); - int openmc_particle_restart(); int openmc_plot_geometry(); int openmc_reset(); int openmc_run(); diff --git a/include/openmc/hdf5_interface.h b/include/openmc/hdf5_interface.h index 1fc378bde..f55b14171 100644 --- a/include/openmc/hdf5_interface.h +++ b/include/openmc/hdf5_interface.h @@ -223,8 +223,9 @@ read_attribute(hid_t obj_id, const char* name, std::vector& vec) // Templates/overloads for read_dataset and related methods //============================================================================== -template -void read_dataset(hid_t obj_id, const char* name, T& buffer, bool indep=false) +template inline +std::enable_if_t>::value> +read_dataset(hid_t obj_id, const char* name, T& buffer, bool indep=false) { read_dataset(obj_id, name, H5TypeMap::type_id, &buffer, indep); } @@ -242,6 +243,14 @@ read_dataset(hid_t obj_id, const char* name, std::string& str, bool indep=false) str = std::string{buffer, n}; } +// array version +template inline void +read_dataset(hid_t dset, const char* name, std::array& buffer, bool indep=false) +{ + read_dataset(dset, name, H5TypeMap::type_id, buffer.data(), indep); +} + +// vector version template void read_dataset(hid_t dset, std::vector& vec, bool indep=false) { diff --git a/include/openmc/particle_restart.h b/include/openmc/particle_restart.h new file mode 100644 index 000000000..24ea237a4 --- /dev/null +++ b/include/openmc/particle_restart.h @@ -0,0 +1,10 @@ +#ifndef OPENMC_PARTICLE_RESTART_H +#define OPENMC_PARTICLE_RESTART_H + +namespace openmc { + +void run_particle_restart(); + +} // namespace openmc + +#endif // OPENMC_PARTICLE_RESTART_H diff --git a/include/openmc/simulation.h b/include/openmc/simulation.h index f9e1f7c96..af333f714 100644 --- a/include/openmc/simulation.h +++ b/include/openmc/simulation.h @@ -83,6 +83,9 @@ void finalize_generation(); //! Determine overall generation number extern "C" int overall_generation(); +extern "C" void simulation_init_f(); +extern "C" void simulation_finalize_f(); + #ifdef OPENMC_MPI void broadcast_results(); #endif diff --git a/src/constants.F90 b/src/constants.F90 index 9061bb800..c1566d776 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -147,18 +147,6 @@ module constants DIFF_NUCLIDE_DENSITY = 2, & DIFF_TEMPERATURE = 3 - ! ============================================================================ - ! RANDOM NUMBER STREAM CONSTANTS - - integer(C_INT), bind(C, name='N_STREAMS') :: N_STREAMS - integer(C_INT), bind(C, name='STREAM_TRACKING') :: STREAM_TRACKING - integer(C_INT), bind(C, name='STREAM_TALLIES') :: STREAM_TALLIES - integer(C_INT), bind(C, name='STREAM_SOURCE') :: STREAM_SOURCE - integer(C_INT), bind(C, name='STREAM_URR_PTABLE') :: STREAM_URR_PTABLE - integer(C_INT), bind(C, name='STREAM_VOLUME') :: STREAM_VOLUME - integer(C_INT), bind(C, name='STREAM_PHOTON') :: STREAM_PHOTON - integer(C_INT64_T), parameter :: DEFAULT_SEED = 1_8 - ! ============================================================================ ! MISCELLANEOUS CONSTANTS diff --git a/src/main.cpp b/src/main.cpp index b757c67ff..936c89b78 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -5,6 +5,7 @@ #include "openmc/constants.h" #include "openmc/error.h" #include "openmc/message_passing.h" +#include "openmc/particle_restart.h" #include "openmc/settings.h" @@ -36,7 +37,8 @@ int main(int argc, char* argv[]) { err = openmc_plot_geometry(); break; case RUN_MODE_PARTICLE: - if (mpi::master) err = openmc_particle_restart(); + if (mpi::master) run_particle_restart(); + err = 0; break; case RUN_MODE_VOLUME: err = openmc_calculate_volumes(); diff --git a/src/particle_restart.F90 b/src/particle_restart.F90 deleted file mode 100644 index c05f9f3df..000000000 --- a/src/particle_restart.F90 +++ /dev/null @@ -1,142 +0,0 @@ -module particle_restart - - use, intrinsic :: ISO_FORTRAN_ENV - - use bank_header, only: Bank - use constants - use error, only: write_message - use hdf5_interface, only: file_open, file_close, read_dataset, HID_T - use mgxs_interface, only: energy_bin_avg - use nuclide_header, only: micro_xs, nuclides_size - use particle_header - use photon_header, only: micro_photon_xs, n_elements - use random_lcg, only: set_particle_seed - use settings - use simulation_header - use tally_header, only: n_tallies - - implicit none - private - public :: openmc_particle_restart - -contains - -!=============================================================================== -! OPENMC_PARTICLE_RESTART is the main routine that runs the particle restart -!=============================================================================== - - function openmc_particle_restart() result(err) bind(C) - integer(C_INT) :: err - - integer(8) :: particle_seed - integer :: previous_run_mode - type(Particle) :: p - - interface - subroutine set_micro_xs() bind(C) - end subroutine - - subroutine print_particle(p) bind(C) - import Particle - type(Particle), intent(in) :: p - end subroutine - end interface - - err = 0 - - ! Set verbosity high - verbosity = 10 - - !$omp parallel - allocate(micro_xs(nuclides_size())) - allocate(micro_photon_xs(n_elements)) - !$omp end parallel - call set_micro_xs() - - ! Initialize the particle to be tracked - call particle_initialize(p) - - ! Read in the restart information - call read_particle_restart(p, previous_run_mode) - - ! Set all tallies to 0 for now (just tracking errors) - n_tallies = 0 - - ! Compute random number seed - select case (previous_run_mode) - case (MODE_EIGENVALUE) - particle_seed = (total_gen + overall_generation() - 1)*n_particles + p % id - case (MODE_FIXEDSOURCE) - particle_seed = p % id - end select - - call set_particle_seed(particle_seed) - - ! Transport neutron - ! FIXME - !call transport(p) - - ! Write output if particle made it - call print_particle(p) - - deallocate(micro_xs) - - end function openmc_particle_restart - -!=============================================================================== -! READ_PARTICLE_RESTART reads the particle restart file -!=============================================================================== - - subroutine read_particle_restart(p, previous_run_mode) - type(Particle), intent(inout) :: p - integer, intent(inout) :: previous_run_mode - - integer(HID_T) :: file_id - character(MAX_WORD_LEN) :: tempstr - - ! Write meessage - call write_message("Loading particle restart file " & - // trim(path_particle_restart) // "...", 5) - - ! Open file - file_id = file_open(path_particle_restart, 'r') - - ! Read data from file - call read_dataset(current_batch, file_id, 'current_batch') - call read_dataset(gen_per_batch, file_id, 'generations_per_batch') - call read_dataset(current_gen, file_id, 'current_generation') - call read_dataset(n_particles, file_id, 'n_particles') - call read_dataset(tempstr, file_id, 'run_mode') - select case (tempstr) - case ('eigenvalue') - previous_run_mode = MODE_EIGENVALUE - case ('fixed source') - previous_run_mode = MODE_FIXEDSOURCE - end select - call read_dataset(p % id, file_id, 'id') - call read_dataset(p % type, file_id, 'type') - call read_dataset(p % wgt, file_id, 'weight') - call read_dataset(p % E, file_id, 'energy') - call read_dataset(p % coord(1) % xyz, file_id, 'xyz') - call read_dataset(p % coord(1) % uvw, file_id, 'uvw') - - ! Set energy group and average energy in multi-group mode - if (.not. run_CE) then - p % g = int(p % E) - p % E = energy_bin_avg(p % g) - end if - - ! Set particle last attributes - p % last_wgt = p % wgt - p % last_xyz_current = p % coord(1)%xyz - p % last_xyz = p % coord(1)%xyz - p % last_uvw = p % coord(1)%uvw - p % last_E = p % E - p % last_g = p % g - - ! Close hdf5 file - call file_close(file_id) - - end subroutine read_particle_restart - -end module particle_restart diff --git a/src/particle_restart.cpp b/src/particle_restart.cpp new file mode 100644 index 000000000..89781bb6c --- /dev/null +++ b/src/particle_restart.cpp @@ -0,0 +1,109 @@ +#include "openmc/particle_restart.h" + +#include "openmc/constants.h" +#include "openmc/hdf5_interface.h" +#include "openmc/mgxs_interface.h" +#include "openmc/nuclide.h" +#include "openmc/output.h" +#include "openmc/particle.h" +#include "openmc/random_lcg.h" +#include "openmc/settings.h" +#include "openmc/simulation.h" +#include "openmc/tallies/tally.h" + +#include // for copy +#include +#include + +namespace openmc { + +void read_particle_restart(Particle& p, int& previous_run_mode) +{ + // Write meessage + write_message("Loading particle restart file " + + settings::path_particle_restart + "...", 5); + + // Open file + hid_t file_id = file_open(settings::path_particle_restart, 'r'); + + // Read data from file + read_dataset(file_id, "current_batch", simulation::current_batch); + read_dataset(file_id, "generations_per_batch", settings::gen_per_batch); + read_dataset(file_id, "current_generation", simulation::current_gen); + read_dataset(file_id, "n_particles", settings::n_particles); + std::string mode; + read_dataset(file_id, "run_mode", mode); + if (mode == "eigenvalue") { + previous_run_mode = RUN_MODE_EIGENVALUE; + } else if (mode == "fixed source") { + previous_run_mode = RUN_MODE_FIXEDSOURCE; + } + read_dataset(file_id, "id", p.id); + read_dataset(file_id, "type", p.type); + read_dataset(file_id, "weight", p.wgt); + read_dataset(file_id, "energy", p.E); + std::array x; + read_dataset(file_id, "xyz", x); + std::copy(x.data(), x.data() + 3, p.coord[0].xyz); + read_dataset(file_id, "uvw", x); + std::copy(x.data(), x.data() + 3, p.coord[0].uvw); + + // Set energy group and average energy in multi-group mode + if (!settings::run_CE) { + p.g = p.E; + p.E = data::energy_bin_avg[p.g - 1]; + } + + // Set particle last attributes + p.last_wgt = p.wgt; + std::copy(p.coord[0].xyz, p.coord[0].xyz + 3, p.last_xyz_current); + std::copy(p.coord[0].xyz, p.coord[0].xyz + 3, p.last_xyz); + std::copy(p.coord[0].uvw, p.coord[0].uvw + 3, p.last_uvw); + p.last_E = p.E; + p.last_g = p.g; + + // Close hdf5 file + file_close(file_id); +} + +void run_particle_restart() +{ + // Set verbosity high + settings::verbosity = 10; + + simulation_init_f(); + set_micro_xs(); + + // Initialize the particle to be tracked + Particle p; + p.initialize(); + + // Read in the restart information + int previous_run_mode; + read_particle_restart(p, previous_run_mode); + + // Set all tallies to 0 for now (just tracking errors) + model::tallies.clear(); + + // Compute random number seed + int64_t particle_seed; + switch (previous_run_mode) { + case RUN_MODE_EIGENVALUE: + particle_seed = (simulation::total_gen + overall_generation() - 1)*settings::n_particles + p.id; + break; + case RUN_MODE_FIXEDSOURCE: + particle_seed = p.id; + break; + } + set_particle_seed(particle_seed); + + // Transport neutron + p.transport(); + + // Write output if particle made it + print_particle(&p); + + simulation_finalize_f(); +} + +} // namespace openmc diff --git a/src/simulation.cpp b/src/simulation.cpp index acc488353..1fdda3061 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -30,8 +30,6 @@ namespace openmc { extern "C" void accumulate_tallies(); extern "C" void allocate_tally_results(); extern "C" void setup_active_tallies(); -extern "C" void simulation_init_f(); -extern "C" void simulation_finalize_f(); extern "C" void write_tallies(); } // namespace openmc