Move finalize_batch and write_source_point to C++

This commit is contained in:
Paul Romano 2018-10-17 10:06:16 -05:00
parent 8c64139c06
commit a4ddaf9ed6
21 changed files with 211 additions and 223 deletions

View file

@ -95,7 +95,7 @@ extern "C" {
int openmc_sphharm_filter_get_cosine(int32_t index, char cosine[]);
int openmc_sphharm_filter_set_order(int32_t index, int order);
int openmc_sphharm_filter_set_cosine(int32_t index, const char cosine[]);
int openmc_statepoint_write(const char filename[]);
int openmc_statepoint_write(const char filename[], bool* write_source);
int openmc_tally_allocate(int32_t index, const char* type);
int openmc_tally_get_active(int32_t index, bool* active);
int openmc_tally_get_estimator(int32_t index, int32_t* estimator);

View file

@ -0,0 +1,17 @@
#ifndef OPENMC_CONTAINER_UTIL_H
#define OPENMC_CONTAINER_UTIL_H
#include <algorithm> // for find
#include <iterator> // for begin, end
namespace openmc {
template<class C, class T>
inline bool contains(const C& v, const T& x)
{
return std::end(v) != std::find(std::begin(v), std::end(v), x);
}
}
#endif // OPENMC_CONTAINER_UTIL_H

View file

@ -7,6 +7,7 @@
#include <array>
#include <cstdint>
#include <string>
#include <unordered_set>
#include <vector>
#include "pugixml.hpp"
@ -76,6 +77,8 @@ extern "C" int res_scat_method; //!< resonance upscattering method
extern "C" double res_scat_energy_min; //!< Min energy in [eV] for res. upscattering
extern "C" double res_scat_energy_max; //!< Max energy in [eV] for res. upscattering
extern "C" int run_mode; //!< Run mode (eigenvalue, fixed src, etc.)
extern std::unordered_set<int> sourcepoint_batch; //!< Batches when source should be written
extern std::unordered_set<int> statepoint_batch; //!< Batches when state should be written
extern "C" int temperature_method; //!< method for choosing temperatures
extern "C" double temperature_tolerance; //!< Tolerance in [K] on choosing temperatures
extern "C" double temperature_default; //!< Default T in [K]

View file

@ -66,6 +66,13 @@ extern "C" void initialize_generation();
extern "C" void initialize_history(Particle* p, int64_t index_source);
//! Finalize a batch
//!
//! Handles synchronization and accumulation of tallies, calculation of Shannon
//! entropy, getting single-batch estimate of keff, and turning on tallies when
//! appropriate
extern "C" void finalize_batch();
//! Finalize a fission generation
extern "C" void finalize_generation();
@ -74,7 +81,6 @@ extern "C" int overall_generation();
#ifdef OPENMC_MPI
extern "C" void broadcast_results();
extern "C" void broadcast_triggers();
#endif
} // namespace openmc

View file

@ -9,6 +9,7 @@
namespace openmc {
void write_source_point(const char* filename);
extern "C" void write_source_bank(hid_t group_id, Bank* source_bank);
extern "C" void read_source_bank(hid_t group_id, Bank* source_bank);
extern "C" void write_tally_results_nr(hid_t file_id);

View file

@ -45,6 +45,7 @@ extern Timer time_bank_sendrecv;
extern Timer time_finalize;
extern Timer time_inactive;
extern Timer time_initialize;
extern Timer time_tallies;
extern Timer time_total;
} // namespace openmc

View file

@ -132,8 +132,6 @@ contains
n_batch_interval = 1
n_lost_particles = 0
n_particles = -1
n_source_points = 0
n_state_points = 0
n_tallies = 0
output_summary = .true.
output_tallies = .true.

View file

@ -29,7 +29,7 @@ contains
! EXECUTE_CMFD runs the CMFD calculation
!==============================================================================
subroutine execute_cmfd()
subroutine execute_cmfd() bind(C)
use cmfd_data, only: set_up_cmfd
use cmfd_solver, only: cmfd_solver_execute

View file

@ -143,7 +143,7 @@ module cmfd_header
logical, public :: cmfd_run_adjoint = .false.
! CMFD run logicals
logical, public :: cmfd_on = .false.
logical(C_BOOL), public, bind(C) :: cmfd_on = .false.
! CMFD display info
character(len=25), public :: cmfd_display = 'balance'

View file

@ -210,9 +210,7 @@ contains
integer :: i
integer :: n
integer, allocatable :: temp_int_array(:)
type(XMLNode) :: root
type(XMLNode) :: node_sp
type(XMLNode) :: node_res_scat
type(XMLNode) :: node_vol
type(XMLNode), allocatable :: node_vol_list(:)
@ -220,89 +218,6 @@ contains
! Get proper XMLNode type given pointer
root % ptr = root_ptr
! Check if the user has specified to write state points
if (check_for_node(root, "state_point")) then
! Get pointer to state_point node
node_sp = root % child("state_point")
! Determine number of batches at which to store state points
if (check_for_node(node_sp, "batches")) then
n_state_points = node_word_count(node_sp, "batches")
else
n_state_points = 0
end if
if (n_state_points > 0) then
! User gave specific batches to write state points
allocate(temp_int_array(n_state_points))
call get_node_array(node_sp, "batches", temp_int_array)
do i = 1, n_state_points
call statepoint_batch % add(temp_int_array(i))
end do
deallocate(temp_int_array)
else
! If neither were specified, write state point at last batch
n_state_points = 1
call statepoint_batch % add(n_batches)
end if
else
! If no <state_point> tag was present, by default write state point at
! last batch only
n_state_points = 1
call statepoint_batch % add(n_batches)
end if
! Check if the user has specified to write source points
if (check_for_node(root, "source_point")) then
! Get pointer to source_point node
node_sp = root % child("source_point")
! Determine number of batches at which to store source points
if (check_for_node(node_sp, "batches")) then
n_source_points = node_word_count(node_sp, "batches")
else
n_source_points = 0
end if
if (n_source_points > 0) then
! User gave specific batches to write source points
allocate(temp_int_array(n_source_points))
call get_node_array(node_sp, "batches", temp_int_array)
do i = 1, n_source_points
call sourcepoint_batch % add(temp_int_array(i))
end do
deallocate(temp_int_array)
else
! If neither were specified, write source points with state points
n_source_points = n_state_points
do i = 1, n_state_points
call sourcepoint_batch % add(statepoint_batch % get_item(i))
end do
end if
else
! If no <source_point> tag was present, by default we keep source bank in
! statepoint file and write it out at statepoints intervals
n_source_points = n_state_points
do i = 1, n_state_points
call sourcepoint_batch % add(statepoint_batch % get_item(i))
end do
end if
! If source is not seperate and is to be written out in the statepoint file,
! make sure that the sourcepoint batch numbers are contained in the
! statepoint list
if (.not. source_separate) then
do i = 1, n_source_points
if (.not. statepoint_batch % contains(sourcepoint_batch % &
get_item(i))) then
call fatal_error('Sourcepoint batches are not a subset&
& of statepoint batches.')
end if
end do
end if
! Resonance scattering parameters
if (check_for_node(root, "resonance_scattering")) then
node_res_scat = root % child("resonance_scattering")

View file

@ -360,7 +360,7 @@ contains
! multiplication factor as well as the average value if we're in active batches
!===============================================================================
subroutine print_batch_keff()
subroutine print_batch_keff() bind(C)
integer :: i ! overall generation
integer :: n ! number of active generations

View file

@ -109,14 +109,6 @@ module settings
! Whether create fission neutrons or not. Only applied for MODE_FIXEDSOURCE
logical(C_BOOL), bind(C) :: create_fission_neutrons
! Information about state points to be written
integer :: n_state_points = 0
type(SetInt) :: statepoint_batch
! Information about source points to be written
integer :: n_source_points = 0
type(SetInt) :: sourcepoint_batch
character(MAX_FILE_LEN) :: path_input ! Path to input file
character(MAX_FILE_LEN) :: path_cross_sections = '' ! Path to cross_sections.xml
character(MAX_FILE_LEN) :: path_multipole ! Path to wmp library
@ -149,10 +141,14 @@ contains
!===============================================================================
subroutine free_memory_settings()
interface
subroutine free_memory_settings_c() bind(C)
end subroutine
end interface
if (allocated(res_scat_nuclides)) deallocate(res_scat_nuclides)
call statepoint_batch % clear()
call sourcepoint_batch % clear()
call free_memory_settings_c()
end subroutine free_memory_settings
end module settings

View file

@ -9,6 +9,7 @@
#include "openmc/capi.h"
#include "openmc/constants.h"
#include "openmc/container_util.h"
#include "openmc/distribution.h"
#include "openmc/distribution_multi.h"
#include "openmc/distribution_spatial.h"
@ -87,6 +88,8 @@ int res_scat_method {RES_SCAT_ARES};
double res_scat_energy_min {0.01};
double res_scat_energy_max {1000.0};
int run_mode {-1};
std::unordered_set<int> sourcepoint_batch;
std::unordered_set<int> statepoint_batch;
int temperature_method {TEMPERATURE_NEAREST};
double temperature_tolerance {10.0};
double temperature_default {293.6};
@ -507,7 +510,7 @@ void read_settings_xml()
// Reshape into track_identifiers
int n_tracks = temp.size() / 3;
for (int i = 0; i < n_tracks; ++i) {
settings::track_identifiers.push_back({temp[3*i], temp[3*i + 1],
track_identifiers.push_back({temp[3*i], temp[3*i + 1],
temp[3*i + 2]});
}
}
@ -548,7 +551,7 @@ void read_settings_xml()
// If the user did not specify how many mesh cells are to be used in
// each direction, we automatically determine an appropriate number of
// cells
int n = std::ceil(std::pow(settings::n_particles / 20.0, 1.0/3.0));
int n = std::ceil(std::pow(n_particles / 20.0, 1.0/3.0));
m.shape_ = {n, n, n};
m.n_dimension_ = 3;
@ -557,7 +560,7 @@ void read_settings_xml()
}
// Turn on Shannon entropy calculation
settings::entropy_on = true;
entropy_on = true;
}
// Uniform fission source weighting mesh
@ -590,18 +593,48 @@ void read_settings_xml()
if (index_ufs_mesh >= 0) {
// Turn on uniform fission source weighting
settings::ufs_on = true;
ufs_on = true;
}
// TODO: Read <state_point>
// Check if the user has specified to write state points
if (check_for_node(root, "state_point")) {
// Get pointer to state_point node
auto node_sp = root.child("state_point");
// Determine number of batches at which to store state points
if (check_for_node(node_sp, "batches")) {
// User gave specific batches to write state points
auto temp = get_node_array<int>(node_sp, "batches");
for (const auto& b : temp) {
statepoint_batch.insert(b);
}
} else {
// If neither were specified, write state point at last batch
statepoint_batch.insert(n_batches);
}
} else {
// If no <state_point> tag was present, by default write state point at
// last batch only
statepoint_batch.insert(n_batches);
}
// Check if the user has specified to write source points
if (check_for_node(root, "source_point")) {
// Get source_point node
xml_node node_sp = root.child("source_point");
// TODO: Read source point batches
// Determine batches at which to store source points
if (check_for_node(node_sp, "batches")) {
// User gave specific batches to write source points
auto temp = get_node_array<int>(node_sp, "batches");
for (const auto& b : temp) {
sourcepoint_batch.insert(b);
}
} else {
// If neither were specified, write source points with state points
sourcepoint_batch = statepoint_batch;
}
// Check if the user has specified to write binary source file
if (check_for_node(node_sp, "separate")) {
@ -618,10 +651,19 @@ void read_settings_xml()
// If no <source_point> tag was present, by default we keep source bank in
// statepoint file and write it out at statepoints intervals
source_separate = false;
// TODO: add defaults
sourcepoint_batch = statepoint_batch;
}
// TODO: Check source points are subset
// If source is not seperate and is to be written out in the statepoint file,
// make sure that the sourcepoint batch numbers are contained in the
// statepoint list
if (!source_separate) {
for (const auto& b : sourcepoint_batch) {
if (!contains(statepoint_batch, b)) {
fatal_error("Sourcepoint batches are not a subset of statepoint batches.");
}
}
}
// Check if the user has specified to not reduce tallies at the end of every
// batch
@ -787,6 +829,11 @@ extern "C" {
const char* openmc_path_particle_restart() {
return settings::path_particle_restart.c_str();
}
void free_memory_settings_c() {
settings::statepoint_batch.clear();
settings::sourcepoint_batch.clear();
}
}
} // namespace openmc

View file

@ -26,7 +26,7 @@ module simulation
use random_lcg, only: set_particle_seed
use settings
use simulation_header
use state_point, only: openmc_statepoint_write, write_source_point, load_state_point
use state_point, only: openmc_statepoint_write, load_state_point
use string, only: to_str
use tally, only: accumulate_tallies, setup_active_tallies, &
init_tally_routines
@ -55,6 +55,9 @@ module simulation
subroutine finalize_generation() bind(C)
end subroutine finalize_generation
subroutine finalize_batch() bind(C)
end subroutine
subroutine initialize_history(p, index_source) bind(C)
import Particle, C_INT64_T
type(Particle), intent(inout) :: p
@ -141,73 +144,6 @@ contains
end function openmc_next_batch
!===============================================================================
! FINALIZE_BATCH handles synchronization and accumulation of tallies,
! calculation of Shannon entropy, getting single-batch estimate of keff, and
! turning on tallies when appropriate
!===============================================================================
subroutine finalize_batch()
integer(C_INT) :: err
character(MAX_FILE_LEN) :: filename
interface
subroutine broadcast_triggers() bind(C)
end subroutine broadcast_triggers
end interface
! Reduce tallies onto master process and accumulate
call time_tallies % start()
call accumulate_tallies()
call time_tallies % stop()
! Reset global tally results
if (current_batch <= n_inactive) then
global_tallies(:,:) = ZERO
n_realizations = 0
end if
if (run_mode == MODE_EIGENVALUE) then
! Perform CMFD calculation if on
if (cmfd_on) call execute_cmfd()
! Write batch output
if (master .and. verbosity >= 7) call print_batch_keff()
end if
! Check_triggers
if (master) call check_triggers()
#ifdef OPENMC_MPI
call broadcast_triggers()
#endif
if (satisfy_triggers .or. &
(trigger_on .and. current_batch == n_max_batches)) then
call statepoint_batch % add(current_batch)
end if
! Write out state point if it's been specified for this batch
if (statepoint_batch % contains(current_batch)) then
if (sourcepoint_batch % contains(current_batch) .and. source_write &
.and. .not. source_separate) then
err = openmc_statepoint_write(write_source=.true._C_BOOL)
else
err = openmc_statepoint_write(write_source=.false._C_BOOL)
end if
end if
! Write out a separate source point if it's been specified for this batch
if (sourcepoint_batch % contains(current_batch) .and. source_write &
.and. source_separate) call write_source_point()
! Write a continously-overwritten source point if requested.
if (source_latest) then
filename = trim(path_output) // 'source' // '.h5'
call write_source_point(filename)
end if
end subroutine finalize_batch
!===============================================================================
! INITIALIZE_SIMULATION
!===============================================================================

View file

@ -1,6 +1,7 @@
#include "openmc/simulation.h"
#include "openmc/capi.h"
#include "openmc/container_util.h"
#include "openmc/eigenvalue.h"
#include "openmc/error.h"
#include "openmc/message_passing.h"
@ -8,6 +9,7 @@
#include "openmc/random_lcg.h"
#include "openmc/settings.h"
#include "openmc/source.h"
#include "openmc/state_point.h"
#include "openmc/timer.h"
#include "openmc/tallies/tally.h"
@ -15,13 +17,19 @@
#include <string>
// data/functions from Fortran side
extern "C" bool cmfd_on;
extern "C" void accumulate_tallies();
extern "C" void allocate_banks();
extern "C" void check_triggers();
extern "C" void cmfd_init_batch();
extern "C" void cmfd_tally_init();
extern "C" void configure_tallies();
extern "C" void execute_cmfd();
extern "C" void init_tally_routines();
extern "C" void join_bank_from_threads();
extern "C" void load_state_point();
extern "C" void print_batch_keff();
extern "C" void print_columns();
extern "C" void print_generation();
extern "C" void print_results();
@ -234,6 +242,62 @@ void initialize_batch()
setup_active_tallies();
}
void finalize_batch()
{
// Reduce tallies onto master process and accumulate
time_tallies.start();
accumulate_tallies();
time_tallies.stop();
// Reset global tally results
if (simulation::current_batch <= settings::n_inactive) {
auto gt = global_tallies();
std::fill(gt.begin(), gt.end(), 0.0);
n_realizations = 0;
}
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
// Perform CMFD calculation if on
if (cmfd_on) execute_cmfd();
// Write batch output
if (mpi::master && settings::verbosity >= 7) print_batch_keff();
}
// Check_triggers
if (mpi::master) check_triggers();
#ifdef OPENMC_MPI
MPI_Bcast(&simulation::satisfy_triggers, 1, MPI_C_BOOL, 0, mpi::intracomm);
#endif
if (simulation::satisfy_triggers || (settings::trigger_on &&
simulation::current_batch == settings::n_max_batches)) {
settings::statepoint_batch.insert(simulation::current_batch);
}
// Write out state point if it's been specified for this batch
if (contains(settings::statepoint_batch, simulation::current_batch)) {
if (contains(settings::sourcepoint_batch, simulation::current_batch)
&& settings::source_write && !settings::source_separate) {
bool b = true;
openmc_statepoint_write(nullptr, &b);
} else {
bool b = false;
openmc_statepoint_write(nullptr, &b);
}
}
// Write out a separate source point if it's been specified for this batch
if (contains(settings::sourcepoint_batch, simulation::current_batch)
&& settings::source_write && settings::source_separate) {
write_source_point(nullptr);
}
// Write a continously-overwritten source point if requested.
if (settings::source_latest) {
auto filename = settings::path_output + "source.h5";
write_source_point(filename.c_str());
}
}
void initialize_generation()
{
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
@ -406,10 +470,6 @@ void broadcast_results() {
simulation::k_abs_tra = temp[2];
}
void broadcast_triggers()
{
MPI_Bcast(&simulation::satisfy_triggers, 1, MPI_C_BOOL, 0, mpi::intracomm);
}
#endif
//==============================================================================

View file

@ -443,44 +443,6 @@ contains
end if
end function openmc_statepoint_write
!===============================================================================
! WRITE_SOURCE_POINT
!===============================================================================
subroutine write_source_point(filename)
character(MAX_FILE_LEN), intent(in), optional :: filename
logical :: parallel
integer(HID_T) :: file_id
character(MAX_FILE_LEN) :: filename_
! When using parallel HDF5, the file is written to collectively by all
! processes. With MPI-only, the file is opened and written by the master
! (note that the call to write_source_bank is by all processes since slave
! processes need to send source bank data to the master.
#ifdef PHDF5
parallel = .true.
#else
parallel = .false.
#endif
if (present(filename)) then
filename_ = filename
else
filename_ = trim(path_output) // 'source.' // &
& zero_padded(current_batch, count_digits(n_max_batches))
filename_ = trim(filename_) // '.h5'
end if
if (master .or. parallel) then
file_id = file_open(filename_, 'w', parallel=.true.)
call write_attribute(file_id, "filetype", 'source')
end if
call write_source_bank(file_id, source_bank)
if (master .or. parallel) call file_close(file_id)
end subroutine write_source_point
!===============================================================================
! LOAD_STATE_POINT
!===============================================================================

View file

@ -1,6 +1,7 @@
#include "openmc/state_point.h"
#include <algorithm>
#include <iomanip> // for setfill, setw
#include <string>
#include <vector>
@ -19,7 +20,6 @@
namespace openmc {
hid_t h5banktype() {
// Create type for array of 3 reals
hsize_t dims[] {3};
@ -37,6 +37,46 @@ hid_t h5banktype() {
return banktype;
}
void
write_source_point(const char* filename)
{
// When using parallel HDF5, the file is written to collectively by all
// processes. With MPI-only, the file is opened and written by the master
// (note that the call to write_source_bank is by all processes since slave
// processes need to send source bank data to the master.
#ifdef PHDF5
bool parallel = true;
#else
bool parallel = false;
#endif
std::string filename_;
if (filename) {
filename_ = filename;
} else {
// Determine width for zero padding
int w = std::to_string(settings::n_max_batches).size();
std::stringstream s;
s << settings::path_output << "source." << std::setfill('0')
<< std::setw(w) << simulation::current_batch << ".h5";
filename_ = s.str();
}
hid_t file_id;
if (mpi::master || parallel) {
file_id = file_open(filename_, 'w', true);
write_attribute(file_id, "filetype", "source");
}
// Get pointer to source bank and write to file
Bank* source_bank;
int64_t n;
openmc_source_bank(&source_bank, &n);
write_source_bank(file_id, source_bank);
if (mpi::master || parallel) file_close(file_id);
}
void
write_source_bank(hid_t group_id, Bank* source_bank)

View file

@ -3769,7 +3769,7 @@ contains
! within the batch to a new random variable
!===============================================================================
subroutine accumulate_tallies()
subroutine accumulate_tallies() bind(C)
integer :: i
real(C_DOUBLE) :: k_col ! Copy of batch collision estimate of keff

View file

@ -29,7 +29,7 @@ contains
! and predicts the number of remainining batches to convergence.
!===============================================================================
subroutine check_triggers()
subroutine check_triggers() bind(C)
implicit none

View file

@ -45,6 +45,7 @@ Timer time_bank_sendrecv;
Timer time_finalize;
Timer time_inactive;
Timer time_initialize;
Timer time_tallies;
Timer time_total;
//==============================================================================
@ -58,6 +59,7 @@ extern "C" double time_bank_sendrecv_elapsed() { return time_bank_sendrecv.elaps
extern "C" double time_finalize_elapsed() { return time_finalize.elapsed(); }
extern "C" double time_inactive_elapsed() { return time_inactive.elapsed(); }
extern "C" double time_initialize_elapsed() { return time_initialize.elapsed(); }
extern "C" double time_tallies_elapsed() { return time_tallies.elapsed(); }
extern "C" double time_total_elapsed() { return time_total.elapsed(); }
extern "C" void reset_timers()

View file

@ -35,6 +35,10 @@ module timer_header
import C_DOUBLE
real(C_DOUBLE) :: t
end function
function time_tallies_elapsed() result(t) bind(C)
import C_DOUBLE
real(C_DOUBLE) :: t
end function
function time_total_elapsed() result(t) bind(C)
import C_DOUBLE
real(C_DOUBLE) :: t