mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 05:35:49 -04:00
Translate most of load_state_point to C++
This commit is contained in:
Paul Romano
parent
102aa934ae
commit
ab0b9eb0bc
8 changed files with 190 additions and 220 deletions
|
|
@ -83,6 +83,10 @@ extern "C" double ufs_get_weight(const Particle* p);
|
|||
//! \param[in] group HDF5 group
|
||||
void write_eigenvalue_hdf5(hid_t group);
|
||||
|
||||
//! Read data related to k-eigenvalue from statepoint
|
||||
//! \param[in] group HDF5 group
|
||||
void read_eigenvalue_hdf5(hid_t group);
|
||||
|
||||
} // namespace openmc
|
||||
|
||||
#endif // OPENMC_EIGENVALUE_H
|
||||
|
|
|
|||
|
|
@ -229,6 +229,19 @@ void read_dataset(hid_t obj_id, const char* name, T& buffer, bool indep=false)
|
|||
read_dataset(obj_id, name, H5TypeMap<T>::type_id, &buffer, indep);
|
||||
}
|
||||
|
||||
// overload for std::string
|
||||
inline void
|
||||
read_dataset(hid_t obj_id, const char* name, std::string& str, bool indep=false)
|
||||
{
|
||||
// Create buffer to read data into
|
||||
auto n = attribute_typesize(obj_id, name);
|
||||
char buffer[n];
|
||||
|
||||
// Read attribute and set string
|
||||
read_string(obj_id, name, n, buffer, indep);
|
||||
str = std::string{buffer, n};
|
||||
}
|
||||
|
||||
template <typename T>
|
||||
void read_dataset(hid_t dset, std::vector<T>& vec, bool indep=false)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -9,11 +9,12 @@
|
|||
|
||||
namespace openmc {
|
||||
|
||||
void load_state_point();
|
||||
void write_source_point(const char* filename);
|
||||
extern "C" void write_source_bank(hid_t group_id);
|
||||
extern "C" void read_source_bank(hid_t group_id);
|
||||
extern "C" void write_tally_results_nr(hid_t file_id);
|
||||
extern "C" void restart_set_keff();
|
||||
void write_source_bank(hid_t group_id);
|
||||
void read_source_bank(hid_t group_id);
|
||||
void write_tally_results_nr(hid_t file_id);
|
||||
void restart_set_keff();
|
||||
|
||||
} // namespace openmc
|
||||
#endif // OPENMC_STATE_POINT_H
|
||||
|
|
|
|||
|
|
@ -638,7 +638,7 @@ void write_eigenvalue_hdf5(hid_t group)
|
|||
write_dataset(group, "k_combined", k_combined);
|
||||
}
|
||||
|
||||
extern "C" void read_eigenvalue_hdf5(hid_t group)
|
||||
void read_eigenvalue_hdf5(hid_t group)
|
||||
{
|
||||
read_dataset(group, "generations_per_batch", settings::gen_per_batch);
|
||||
int n = simulation::restart_batch*settings::gen_per_batch;
|
||||
|
|
|
|||
|
|
@ -30,7 +30,6 @@ extern "C" void accumulate_tallies();
|
|||
extern "C" void allocate_tally_results();
|
||||
extern "C" void check_triggers();
|
||||
extern "C" void init_tally_routines();
|
||||
extern "C" void load_state_point();
|
||||
extern "C" void setup_active_tallies();
|
||||
extern "C" void simulation_init_f();
|
||||
extern "C" void simulation_finalize_f();
|
||||
|
|
|
|||
|
|
@ -13,7 +13,6 @@ module state_point
|
|||
|
||||
use, intrinsic :: ISO_C_BINDING
|
||||
|
||||
use bank_header, only: Bank
|
||||
use constants
|
||||
use endf, only: reaction_name
|
||||
use error, only: fatal_error, warning, write_message
|
||||
|
|
@ -21,11 +20,9 @@ module state_point
|
|||
use message_passing
|
||||
use mgxs_interface
|
||||
use nuclide_header, only: nuclides
|
||||
use output, only: time_stamp
|
||||
use random_lcg, only: openmc_get_seed, openmc_set_seed
|
||||
use settings
|
||||
use simulation_header
|
||||
use string, only: to_str, count_digits, zero_padded, to_f_string
|
||||
use string, only: to_str, zero_padded
|
||||
use tally_header
|
||||
use tally_filter_header
|
||||
use tally_derivative_header, only: tally_derivs
|
||||
|
|
@ -33,18 +30,6 @@ module state_point
|
|||
|
||||
implicit none
|
||||
|
||||
interface
|
||||
subroutine write_source_bank(group_id) bind(C)
|
||||
import HID_T
|
||||
integer(HID_T), value :: group_id
|
||||
end subroutine write_source_bank
|
||||
|
||||
subroutine read_source_bank(group_id) bind(C)
|
||||
import HID_T
|
||||
integer(HID_T), value :: group_id
|
||||
end subroutine read_source_bank
|
||||
end interface
|
||||
|
||||
contains
|
||||
|
||||
!===============================================================================
|
||||
|
|
@ -248,181 +233,39 @@ contains
|
|||
! LOAD_STATE_POINT
|
||||
!===============================================================================
|
||||
|
||||
subroutine load_state_point() bind(C)
|
||||
subroutine load_state_point_f(file_id) bind(C)
|
||||
integer(HID_T), value :: file_id
|
||||
|
||||
integer :: i
|
||||
integer :: int_array(3)
|
||||
integer, allocatable :: array(:)
|
||||
integer(C_INT64_T) :: seed
|
||||
integer(HID_T) :: file_id
|
||||
integer :: temp
|
||||
integer(HID_T) :: tallies_group
|
||||
integer(HID_T) :: tally_group
|
||||
logical :: source_present
|
||||
character(MAX_WORD_LEN) :: word
|
||||
|
||||
interface
|
||||
subroutine read_eigenvalue_hdf5(group) bind(C)
|
||||
import HID_T
|
||||
integer(HID_T), value :: group
|
||||
end subroutine
|
||||
! Read global tally data
|
||||
call read_dataset(global_tallies, file_id, "global_tallies")
|
||||
|
||||
subroutine restart_set_keff() bind(C)
|
||||
end subroutine
|
||||
end interface
|
||||
! Check if tally results are present
|
||||
call read_attribute(temp, file_id, "tallies_present")
|
||||
|
||||
! Write message
|
||||
call write_message("Loading state point " // trim(path_state_point) &
|
||||
// "...", 5)
|
||||
! Read in sum and sum squared
|
||||
if (temp == 1) then
|
||||
tallies_group = open_group(file_id, "tallies")
|
||||
|
||||
! Open file for reading
|
||||
file_id = file_open(path_state_point, 'r', parallel=.true.)
|
||||
TALLY_RESULTS: do i = 1, n_tallies
|
||||
associate (t => tallies(i) % obj)
|
||||
! Read sum, sum_sq, and N for each bin
|
||||
tally_group = open_group(tallies_group, "tally " // &
|
||||
trim(to_str(t % id)))
|
||||
call t % read_results_hdf5(tally_group)
|
||||
call read_dataset(t % n_realizations, tally_group, &
|
||||
"n_realizations")
|
||||
call close_group(tally_group)
|
||||
end associate
|
||||
end do TALLY_RESULTS
|
||||
|
||||
! Read filetype
|
||||
call read_attribute(word, file_id, "filetype")
|
||||
if (word /= 'statepoint') then
|
||||
call fatal_error("OpenMC tried to restart from a non-statepoint file.")
|
||||
call close_group(tallies_group)
|
||||
end if
|
||||
|
||||
! Read revision number for state point file and make sure it matches with
|
||||
! current version
|
||||
call read_attribute(array, file_id, "version")
|
||||
if (any(array /= VERSION_STATEPOINT)) then
|
||||
call fatal_error("State point version does not match current version &
|
||||
&in OpenMC.")
|
||||
end if
|
||||
|
||||
! Read and overwrite random number seed
|
||||
call read_dataset(seed, file_id, "seed")
|
||||
call openmc_set_seed(seed)
|
||||
|
||||
! It is not impossible for a state point to be generated from a CE run but
|
||||
! to be loaded in to an MG run (or vice versa), check to prevent that.
|
||||
call read_dataset(word, file_id, "energy_mode")
|
||||
if (word == "multi-group" .and. run_CE) then
|
||||
call fatal_error("State point file is from multi-group run but &
|
||||
& current run is continous-energy!")
|
||||
else if (word == "continuous-energy" .and. .not. run_CE) then
|
||||
call fatal_error("State point file is from continuous-energy run but &
|
||||
& current run is multi-group!")
|
||||
end if
|
||||
|
||||
! Read and overwrite run information except number of batches
|
||||
call read_dataset(word, file_id, "run_mode")
|
||||
select case(word)
|
||||
case ('fixed source')
|
||||
run_mode = MODE_FIXEDSOURCE
|
||||
case ('eigenvalue')
|
||||
run_mode = MODE_EIGENVALUE
|
||||
end select
|
||||
call read_attribute(int_array(1), file_id, "photon_transport")
|
||||
if (int_array(1) == 1) then
|
||||
photon_transport = .true.
|
||||
else
|
||||
photon_transport = .false.
|
||||
end if
|
||||
call read_dataset(n_particles, file_id, "n_particles")
|
||||
call read_dataset(int_array(1), file_id, "n_batches")
|
||||
|
||||
! Take maximum of statepoint n_batches and input n_batches
|
||||
n_batches = max(n_batches, int_array(1))
|
||||
|
||||
! Read batch number to restart at
|
||||
call read_dataset(restart_batch, file_id, "current_batch")
|
||||
|
||||
! Check for source in statepoint if needed
|
||||
call read_attribute(int_array(1), file_id, "source_present")
|
||||
if (int_array(1) == 1) then
|
||||
source_present = .true.
|
||||
else
|
||||
source_present = .false.
|
||||
end if
|
||||
|
||||
if (restart_batch > n_batches) then
|
||||
call fatal_error("The number batches specified in settings.xml is fewer &
|
||||
& than the number of batches in the given statepoint file.")
|
||||
end if
|
||||
|
||||
! Read information specific to eigenvalue run
|
||||
if (run_mode == MODE_EIGENVALUE) then
|
||||
call read_dataset(int_array(1), file_id, "n_inactive")
|
||||
call read_eigenvalue_hdf5(file_id)
|
||||
|
||||
! Take maximum of statepoint n_inactive and input n_inactive
|
||||
n_inactive = max(n_inactive, int_array(1))
|
||||
end if
|
||||
|
||||
! Read number of realizations for global tallies
|
||||
call read_dataset(n_realizations, file_id, "n_realizations", indep=.true.)
|
||||
|
||||
! Set k_sum, keff, and current_batch based on whether restart file is part
|
||||
! of active cycle or inactive cycle
|
||||
call restart_set_keff()
|
||||
current_batch = restart_batch
|
||||
|
||||
! Check to make sure source bank is present
|
||||
if (path_source_point == path_state_point .and. .not. source_present) then
|
||||
call fatal_error("Source bank must be contained in statepoint restart &
|
||||
&file")
|
||||
end if
|
||||
|
||||
! Read tallies to master. If we are using Parallel HDF5, all processes
|
||||
! need to be included in the HDF5 calls.
|
||||
#ifdef PHDF5
|
||||
if (.true.) then
|
||||
#else
|
||||
if (master) then
|
||||
#endif
|
||||
! Read global tally data
|
||||
call read_dataset(global_tallies, file_id, "global_tallies")
|
||||
|
||||
! Check if tally results are present
|
||||
call read_attribute(int_array(1), file_id, "tallies_present")
|
||||
|
||||
! Read in sum and sum squared
|
||||
if (int_array(1) == 1) then
|
||||
tallies_group = open_group(file_id, "tallies")
|
||||
|
||||
TALLY_RESULTS: do i = 1, n_tallies
|
||||
associate (t => tallies(i) % obj)
|
||||
! Read sum, sum_sq, and N for each bin
|
||||
tally_group = open_group(tallies_group, "tally " // &
|
||||
trim(to_str(t % id)))
|
||||
call t % read_results_hdf5(tally_group)
|
||||
call read_dataset(t % n_realizations, tally_group, &
|
||||
"n_realizations")
|
||||
call close_group(tally_group)
|
||||
end associate
|
||||
end do TALLY_RESULTS
|
||||
|
||||
call close_group(tallies_group)
|
||||
end if
|
||||
end if
|
||||
|
||||
! Read source if in eigenvalue mode
|
||||
if (run_mode == MODE_EIGENVALUE) then
|
||||
|
||||
! Check if source was written out separately
|
||||
if (.not. source_present) then
|
||||
|
||||
! Close statepoint file
|
||||
call file_close(file_id)
|
||||
|
||||
! Write message
|
||||
call write_message("Loading source file " // trim(path_source_point) &
|
||||
// "...", 5)
|
||||
|
||||
! Open source file
|
||||
file_id = file_open(path_source_point, 'r', parallel=.true.)
|
||||
end if
|
||||
|
||||
! Read source
|
||||
call read_source_bank(file_id)
|
||||
|
||||
end if
|
||||
|
||||
! Close file
|
||||
call file_close(file_id)
|
||||
|
||||
end subroutine load_state_point
|
||||
end subroutine load_state_point_f
|
||||
|
||||
end module state_point
|
||||
|
|
|
|||
|
|
@ -1,6 +1,7 @@
|
|||
#include "openmc/state_point.h"
|
||||
|
||||
#include <algorithm>
|
||||
#include <cstdint> // for int64_t
|
||||
#include <iomanip> // for setfill, setw
|
||||
#include <string>
|
||||
#include <vector>
|
||||
|
|
@ -26,6 +27,7 @@
|
|||
namespace openmc {
|
||||
|
||||
extern "C" void statepoint_write_f(hid_t file_id);
|
||||
extern "C" void load_state_point_f(hid_t file_id);
|
||||
|
||||
extern "C" int
|
||||
openmc_statepoint_write(const char* filename, bool* write_source)
|
||||
|
|
@ -160,6 +162,147 @@ openmc_statepoint_write(const char* filename, bool* write_source)
|
|||
return 0;
|
||||
}
|
||||
|
||||
void restart_set_keff()
|
||||
{
|
||||
if (simulation::restart_batch > settings::n_inactive) {
|
||||
for (int i = settings::n_inactive; i < simulation::restart_batch; ++i) {
|
||||
simulation::k_sum[0] += simulation::k_generation[i];
|
||||
simulation::k_sum[1] += std::pow(simulation::k_generation[i], 2);
|
||||
}
|
||||
int n = settings::gen_per_batch*n_realizations;
|
||||
simulation::keff = simulation::k_sum[0] / n;
|
||||
} else {
|
||||
simulation::keff = simulation::k_generation.back();
|
||||
}
|
||||
}
|
||||
|
||||
void load_state_point()
|
||||
{
|
||||
// Write message
|
||||
write_message("Loading state point " + settings::path_statepoint + "...", 5);
|
||||
|
||||
// Open file for reading
|
||||
hid_t file_id = file_open(settings::path_statepoint.c_str(), 'r', true);
|
||||
|
||||
// Read filetype
|
||||
std::string word;
|
||||
read_attribute(file_id, "filetype", word);
|
||||
if (word != "statepoint") {
|
||||
fatal_error("OpenMC tried to restart from a non-statepoint file.");
|
||||
}
|
||||
|
||||
// Read revision number for state point file and make sure it matches with
|
||||
// current version
|
||||
std::array<int, 2> array;
|
||||
read_attribute(file_id, "version", array);
|
||||
if (array != VERSION_STATEPOINT) {
|
||||
fatal_error("State point version does not match current version in OpenMC.");
|
||||
}
|
||||
|
||||
// Read and overwrite random number seed
|
||||
int64_t seed;
|
||||
read_dataset(file_id, "seed", seed);
|
||||
openmc_set_seed(seed);
|
||||
|
||||
// It is not impossible for a state point to be generated from a CE run but
|
||||
// to be loaded in to an MG run (or vice versa), check to prevent that.
|
||||
read_dataset(file_id, "energy_mode", word);
|
||||
if (word == "multi-group" && settings::run_CE) {
|
||||
fatal_error("State point file is from multigroup run but current run is "
|
||||
"continous energy.");
|
||||
} else if (word == "continuous-energy" && !settings::run_CE) {
|
||||
fatal_error("State point file is from continuous-energy run but current "
|
||||
"run is multigroup!");
|
||||
}
|
||||
|
||||
// Read and overwrite run information except number of batches
|
||||
read_dataset(file_id, "run_mode", word);
|
||||
if (word == "fixed source") {
|
||||
settings::run_mode = RUN_MODE_FIXEDSOURCE;
|
||||
} else if (word == "eigenvalue") {
|
||||
settings::run_mode = RUN_MODE_EIGENVALUE;
|
||||
}
|
||||
read_attribute(file_id, "photon_transport", settings::photon_transport);
|
||||
read_dataset(file_id, "n_particles", settings::n_particles);
|
||||
int temp;
|
||||
read_dataset(file_id, "n_batches", temp);
|
||||
|
||||
// Take maximum of statepoint n_batches and input n_batches
|
||||
settings::n_batches = std::max(settings::n_batches, temp);
|
||||
|
||||
// Read batch number to restart at
|
||||
read_dataset(file_id, "current_batch", simulation::restart_batch);
|
||||
|
||||
// Check for source in statepoint if needed
|
||||
bool source_present;
|
||||
read_attribute(file_id, "source_present", source_present);
|
||||
|
||||
if (simulation::restart_batch > settings::n_batches) {
|
||||
fatal_error("The number batches specified in settings.xml is fewer "
|
||||
" than the number of batches in the given statepoint file.");
|
||||
}
|
||||
|
||||
// Read information specific to eigenvalue run
|
||||
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
|
||||
read_dataset(file_id, "n_inactive", temp);
|
||||
read_eigenvalue_hdf5(file_id);
|
||||
|
||||
// Take maximum of statepoint n_inactive and input n_inactive
|
||||
settings::n_inactive = std::max(settings::n_inactive, temp);
|
||||
}
|
||||
|
||||
// Read number of realizations for global tallies
|
||||
read_dataset(file_id, "n_realizations", n_realizations);
|
||||
|
||||
// Set k_sum, keff, and current_batch based on whether restart file is part
|
||||
// of active cycle or inactive cycle
|
||||
restart_set_keff();
|
||||
simulation::current_batch = simulation::restart_batch;
|
||||
|
||||
// Check to make sure source bank is present
|
||||
if (settings::path_sourcepoint == settings::path_statepoint &&
|
||||
!source_present) {
|
||||
fatal_error("Source bank must be contained in statepoint restart file");
|
||||
}
|
||||
|
||||
// Read tallies to master. If we are using Parallel HDF5, all processes
|
||||
// need to be included in the HDF5 calls.
|
||||
#ifdef PHDF5
|
||||
if (true) {
|
||||
#else
|
||||
if (mpi::master) {
|
||||
#endif
|
||||
// Read tally results
|
||||
load_state_point_f(file_id);
|
||||
}
|
||||
|
||||
// Read source if in eigenvalue mode
|
||||
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
|
||||
|
||||
// Check if source was written out separately
|
||||
if (!source_present) {
|
||||
|
||||
// Close statepoint file
|
||||
file_close(file_id);
|
||||
|
||||
// Write message
|
||||
write_message("Loading source file " + settings::path_sourcepoint
|
||||
+ "...", 5);
|
||||
|
||||
// Open source file
|
||||
file_id = file_open(settings::path_source.c_str(), 'r', true);
|
||||
}
|
||||
|
||||
// Read source
|
||||
read_source_bank(file_id);
|
||||
|
||||
}
|
||||
|
||||
// Close file
|
||||
file_close(file_id);
|
||||
}
|
||||
|
||||
|
||||
hid_t h5banktype() {
|
||||
// Create type for array of 3 reals
|
||||
hsize_t dims[] {3};
|
||||
|
|
@ -455,18 +598,4 @@ void write_tally_results_nr(hid_t file_id)
|
|||
}
|
||||
}
|
||||
|
||||
void restart_set_keff()
|
||||
{
|
||||
if (simulation::restart_batch > settings::n_inactive) {
|
||||
for (int i = settings::n_inactive; i < simulation::restart_batch; ++i) {
|
||||
simulation::k_sum[0] += simulation::k_generation[i];
|
||||
simulation::k_sum[1] += std::pow(simulation::k_generation[i], 2);
|
||||
}
|
||||
int n = settings::gen_per_batch*n_realizations;
|
||||
simulation::keff = simulation::k_sum[0] / n;
|
||||
} else {
|
||||
simulation::keff = simulation::k_generation.back();
|
||||
}
|
||||
}
|
||||
|
||||
} // namespace openmc
|
||||
|
|
|
|||
|
|
@ -247,25 +247,6 @@ contains
|
|||
|
||||
end function ends_with
|
||||
|
||||
!===============================================================================
|
||||
! COUNT_DIGITS returns the number of digits needed to represent the input
|
||||
! integer.
|
||||
!===============================================================================
|
||||
|
||||
pure function count_digits(num) result(n_digits)
|
||||
integer, intent(in) :: num
|
||||
integer :: n_digits
|
||||
|
||||
n_digits = 1
|
||||
do while (num / 10**(n_digits) /= 0 .and. abs(num / 10 **(n_digits-1)) /= 1&
|
||||
&.and. n_digits /= 10)
|
||||
! Note that 10 digits is the maximum needed to represent an integer(4) so
|
||||
! the loop automatically exits when n_digits = 10.
|
||||
n_digits = n_digits + 1
|
||||
end do
|
||||
|
||||
end function count_digits
|
||||
|
||||
!===============================================================================
|
||||
! INT4_TO_STR converts an integer(4) to a string.
|
||||
!===============================================================================
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue