From 07fc0985bc4e05a7f8527fe8657f1e35334c125d Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Tue, 19 Feb 2019 22:06:20 -0600 Subject: [PATCH] Convert particle track writing to C++ --- CMakeLists.txt | 2 +- include/openmc/track_output.h | 18 +++++ src/constants.F90 | 6 -- src/particle.cpp | 14 ++-- src/track_output.F90 | 132 ---------------------------------- src/track_output.cpp | 80 +++++++++++++++++++++ 6 files changed, 104 insertions(+), 148 deletions(-) create mode 100644 include/openmc/track_output.h delete mode 100644 src/track_output.F90 create mode 100644 src/track_output.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 65292367d..3ab7cdb04 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -327,7 +327,6 @@ add_library(libopenmc SHARED src/state_point.F90 src/stl_vector.F90 src/string.F90 - src/track_output.F90 src/vector_header.F90 src/xml_interface.F90 src/tallies/tally.F90 @@ -419,6 +418,7 @@ add_library(libopenmc SHARED src/tallies/trigger.cpp src/timer.cpp src/thermal.cpp + src/track_output.cpp src/volume_calc.cpp src/wmp.cpp src/xml_interface.cpp diff --git a/include/openmc/track_output.h b/include/openmc/track_output.h new file mode 100644 index 000000000..f32529bc2 --- /dev/null +++ b/include/openmc/track_output.h @@ -0,0 +1,18 @@ +#ifndef OPENMC_TRACK_OUTPUT_H +#define OPENMC_TRACK_OUTPUT_H + +#include "openmc/particle.h" + +namespace openmc { + +//============================================================================== +// Non-member functions +//============================================================================== + +void add_particle_track(); +void write_particle_track(const Particle& p); +void finalize_particle_track(const Particle& p); + +} // namespace openmc + +#endif // OPENMC_TRACK_OUTPUT_H diff --git a/src/constants.F90 b/src/constants.F90 index c1566d776..1fd5d00e9 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -4,12 +4,6 @@ module constants implicit none - ! ============================================================================ - ! VERSIONING NUMBERS - - ! Version numbers for binary files - integer, parameter :: VERSION_TRACK(2) = [2, 0] - ! ============================================================================ ! ADJUSTABLE PARAMETERS diff --git a/src/particle.cpp b/src/particle.cpp index db2421d78..90d3d2835 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -24,6 +24,7 @@ #include "openmc/tallies/derivative.h" #include "openmc/tallies/tally.h" #include "openmc/tallies/tally_scoring.h" +#include "openmc/track_output.h" namespace openmc { @@ -132,11 +133,6 @@ Particle::from_source(const Bank* src) last_E = E; } -extern "C" void initialize_particle_track(); -extern "C" void write_particle_track(const Particle*); -extern "C" void add_particle_track(); -extern "C" void finalize_particle_track(const Particle*); - void Particle::transport() { @@ -160,7 +156,7 @@ Particle::transport() } // Prepare to write out particle track. - if (write_track) initialize_particle_track(); + if (write_track) add_particle_track(); // Every particle starts with no accumulated flux derivative. if (!model::active_tallies.empty()) zero_flux_derivs(); @@ -194,7 +190,7 @@ Particle::transport() } // Write particle track. - if (write_track) write_particle_track(this); + if (write_track) write_particle_track(*this); if (settings::check_overlaps) check_cell_overlap(this); @@ -396,8 +392,8 @@ Particle::transport() // Finish particle track output. if (write_track) { - write_particle_track(this); - finalize_particle_track(this); + write_particle_track(*this); + finalize_particle_track(*this); } } diff --git a/src/track_output.F90 b/src/track_output.F90 deleted file mode 100644 index cf1625ca3..000000000 --- a/src/track_output.F90 +++ /dev/null @@ -1,132 +0,0 @@ -!=============================================================================== -! TRACK_OUTPUT handles output of particle tracks (the paths taken by particles -! as they are transported through the geometry). -!=============================================================================== - -module track_output - - use, intrinsic :: ISO_C_BINDING - - use constants - use hdf5_interface - use particle_header, only: Particle - use settings, only: path_output - use simulation_header - use string, only: to_str - - implicit none - private - - type TrackCoordinates - real(8), allocatable :: coords(:,:) - end type TrackCoordinates - - type(TrackCoordinates), allocatable :: tracks(:) -!$omp threadprivate(tracks) - - public :: initialize_particle_track - public :: write_particle_track - public :: add_particle_track - public :: finalize_particle_track - -contains - -!=============================================================================== -! INITIALIZE_PARTICLE_TRACK allocates the array to store particle track -! information -!=============================================================================== - - subroutine initialize_particle_track() bind(C) - allocate(tracks(1)) - end subroutine initialize_particle_track - -!=============================================================================== -! WRITE_PARTICLE_TRACK copies particle position to an array. -!=============================================================================== - - subroutine write_particle_track(p) bind(C) - type(Particle), intent(in) :: p - real(8), allocatable :: new_coords(:, :) - - integer :: i - integer :: n_tracks - - ! Add another column to coords - i = size(tracks) - if (allocated(tracks(i)%coords)) then - n_tracks = size(tracks(i)%coords, 2) - allocate(new_coords(3, n_tracks + 1)) - new_coords(:, 1:n_tracks) = tracks(i)%coords - call move_alloc(FROM=new_coords, TO=tracks(i)%coords) - else - n_tracks = 0 - allocate(tracks(i)%coords(3, 1)) - end if - - ! Write current coordinates into the newest column. - n_tracks = n_tracks + 1 - tracks(i)%coords(:, n_tracks) = p%coord(1)%xyz - end subroutine write_particle_track - -!=============================================================================== -! ADD_PARTICLE_TRACK creates a new entry in the track coordinates for a -! secondary particle -!=============================================================================== - - subroutine add_particle_track() bind(C) - type(TrackCoordinates), allocatable :: new_tracks(:) - - integer :: i - - ! Determine current number of particle tracks - i = size(tracks) - - ! Create array one larger than current - allocate(new_tracks(i + 1)) - - ! Copy memory and move allocation - new_tracks(1:i) = tracks(i) - call move_alloc(FROM=new_tracks, TO=tracks) - - end subroutine add_particle_track - -!=============================================================================== -! FINALIZE_PARTICLE_TRACK writes the particle track array to disk. -!=============================================================================== - - subroutine finalize_particle_track(p) bind(C) - type(Particle), intent(in) :: p - - integer :: i - integer :: n_particle_tracks - integer(HID_T) :: file_id - character(MAX_FILE_LEN) :: fname - integer, allocatable :: n_coords(:) - - fname = trim(path_output) // 'track_' // trim(to_str(current_batch)) & - // '_' // trim(to_str(current_gen)) // '_' // trim(to_str(p%id)) & - // '.h5' - - ! Determine total number of particles and number of coordinates for each - n_particle_tracks = size(tracks) - allocate(n_coords(n_particle_tracks)) - do i = 1, n_particle_tracks - n_coords(i) = size(tracks(i)%coords, 2) - end do - -!$omp critical (FinalizeParticleTrack) - file_id = file_open(fname, 'w') - call write_attribute(file_id, 'filetype', 'track') - call write_attribute(file_id, 'version', VERSION_TRACK) - call write_attribute(file_id, 'n_particles', n_particle_tracks) - call write_attribute(file_id, 'n_coords', n_coords) - do i = 1, n_particle_tracks - call write_dataset(file_id, 'coordinates_' // trim(to_str(i)), & - tracks(i)%coords) - end do - call file_close(file_id) -!$omp end critical (FinalizeParticleTrack) - deallocate(tracks) - end subroutine finalize_particle_track - -end module track_output diff --git a/src/track_output.cpp b/src/track_output.cpp new file mode 100644 index 000000000..78548c76c --- /dev/null +++ b/src/track_output.cpp @@ -0,0 +1,80 @@ +#include "openmc/track_output.h" + +#include "openmc/constants.h" +#include "openmc/hdf5_interface.h" +#include "openmc/position.h" +#include "openmc/settings.h" +#include "openmc/simulation.h" + +#include "xtensor/xtensor.hpp" + +#include // for size_t +#include +#include +#include + +namespace openmc { + +//============================================================================== +// Global variables +//============================================================================== + +// Forward declaration needed in order to declare tracks as threadprivate +extern std::vector> tracks; +#pragma omp threadprivate(tracks) + +std::vector> tracks; + +//============================================================================== +// Non-member functions +//============================================================================== + +void add_particle_track() +{ + tracks.emplace_back(); +} + +void write_particle_track(const Particle& p) +{ + tracks.back().push_back({p.coord[0].xyz}); +} + +void finalize_particle_track(const Particle& p) +{ + std::stringstream filename; + filename << settings::path_output << "track_" << simulation::current_batch + << '_' << simulation::current_gen << '_' << p.id << ".h5"; + + // Determine number of coordinates for each particle + std::vector n_coords; + for (auto& coords : tracks) { + n_coords.push_back(coords.size()); + } + + #pragma omp critical (FinalizeParticleTrack) + { + hid_t file_id = file_open(filename.str().c_str(), 'w'); + write_attribute(file_id, "filetype", "track"); + write_attribute(file_id, "version", VERSION_TRACK); + write_attribute(file_id, "n_particles", tracks.size()); + write_attribute(file_id, "n_coords", n_coords); + for (int i = 1; i <= tracks.size(); ++i) { + const auto& t {tracks[i-1]}; + size_t n = t.size(); + xt::xtensor data({n,3}); + for (int j = 0; j < n; ++j) { + data(j, 0) = t[j].x; + data(j, 1) = t[j].y; + data(j, 2) = t[j].z; + } + std::string name = "coordinates_" + std::to_string(i); + write_dataset(file_id, name.c_str(), data); + } + file_close(file_id); + } + + // Clear particle tracks + tracks.clear(); +} + +} // namespace openmc