crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-26 21:25:36 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Convert write_source_bank to C++, move some MPI initialization over too

Browse source
This commit is contained in:
Paul Romano 2018-04-20 11:03:40 -05:00
parent 0ee5ec4ecb
commit 44ff22dc14
11 changed files with 200 additions and 154 deletions
Download patch file Download diff file Expand all files Collapse all files

1
CMakeLists.txt
View file

@ -438,6 +438,7 @@ set(LIBOPENMC_CXX_SRC
src/message_passing.cpp
src/random_lcg.cpp
src/simulation.cpp
src/state_point.cpp
src/surface.cpp
src/xml_interface.cpp
src/pugixml/pugixml.cpp)

7
include/openmc.h
View file

@ -127,7 +127,6 @@ extern "C" {
extern char openmc_err_msg[256];
extern double openmc_keff;
extern double openmc_keff_std;
extern bool openmc_master;
extern int32_t n_batches;
extern int32_t n_cells;
extern int32_t n_filters;
@ -147,6 +146,12 @@ extern "C" {
extern bool openmc_simulation_initialized;
extern int openmc_verbosity;
// Variables that are shared by necessity (can be removed later)
extern bool openmc_master;
extern int openmc_n_procs;
extern int openmc_rank;
extern int64_t openmc_work;
// Run modes
constexpr int RUN_MODE_FIXEDSOURCE {1};
constexpr int RUN_MODE_EIGENVALUE {2};

14
openmc/capi/core.py
View file

@ -118,11 +118,15 @@ def init(intracomm=None):
if intracomm is not None:
# If an mpi4py communicator was passed, convert it to void* to be passed
# to openmc_init
from mpi4py import MPI
address = MPI._addressof(intracomm)
_dll.openmc_init(c_void_p(address))
else:
_dll.openmc_init(None)
try:
from mpi4py import MPI
except ImportError:
intracomm = None
else:
address = MPI._addressof(intracomm)
intracomm = c_void_p(address)
_dll.openmc_init(intracomm)
def iter_batches():

16
src/initialize.F90
View file

@ -131,27 +131,11 @@ contains
integer :: bank_types(5) ! Datatypes
#endif
integer(MPI_ADDRESS_KIND) :: bank_disp(5) ! Displacements
logical :: init_called
type(Bank) :: b
! Indicate that MPI is turned on
mpi_enabled = .true.
! Initialize MPI
call MPI_INITIALIZED(init_called, mpi_err)
if (.not. init_called) call MPI_INIT(mpi_err)
! Determine number of processors and rank of each processor
mpi_intracomm = intracomm
call MPI_COMM_SIZE(mpi_intracomm, n_procs, mpi_err)
call MPI_COMM_RANK(mpi_intracomm, rank, mpi_err)
! Determine master
if (rank == 0) then
master = .true.
else
master = .false.
end if
! ==========================================================================
! CREATE MPI_BANK TYPE

42
src/initialize.cpp
View file

@ -3,12 +3,14 @@
#include "message_passing.h"
int
openmc_init(const void* intracomm)
int openmc_init(const void* intracomm)
{
#ifdef OPENMC_MPI
openmc::mpi::intracomm = *static_cast<const MPI_Comm *>(intracomm);
// Initialize MPI for C++
MPI_Comm intracomm = *static_cast<const MPI_Comm *>(intracomm);
initialize_mpi(intracomm);
// Continue with rest of initialization
MPI_Fint fcomm = MPI_Comm_c2f(openmc::mpi::intracomm);
openmc_init_f(&fcomm);
#else
@ -16,3 +18,37 @@ openmc_init(const void* intracomm)
#endif
return 0;
}
#ifdef OPENMC_MPI
void initialize_mpi(MPI_Comm intracomm)
{
openmc::mpi::intracomm = intracomm;
// Initialize MPI
int flag;
if (!MPI_Initialized(&flag)) MPI_Init(nullptr, nullptr);
// Determine number of processes and rank for each
MPI_Comm_size(intracomm, openmc::mpi::n_procs);
MPI_Comm_rank(intracomm, openmc::mpi::rank);
// Set variable for Fortran side
openmc_n_procs = openmc::mpi::n_procs;
openmc_rank = openmc::mpi::rank;
openmc_master = (openmc::mpi::rank == 0);
// Create bank datatype
Bank b;
MPI_Aint disp[5];
disp[0] = &b.wgt - &b;
disp[1] = &b.xyz - &b;
disp[2] = &b.uvw - &b;
disp[3] = &b.E - &b;
disp[4] = &b.delayed_group - &b;
int blocks[] {1, 3, 3, 1, 1};
MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_INT};
MPI_Type_create_struct(5, blocks, disp, types, &openmc::mpi::bank);
MPI_Type_commit(&openmc::mpi::bank);
}
#endif

6
src/message_passing.F90
View file

@ -14,9 +14,9 @@ module message_passing
! mpi_enabled flag are for when MPI is not being used at all, i.e. a serial
! run. In this case, these variables are still used at times.
integer :: n_procs = 1 ! number of processes
integer :: rank = 0 ! rank of process
logical(C_BOOL), bind(C, name='openmc_master') :: master = .true. ! master process?
integer(C_INT), bind(C, name='openmc_n_procs') :: n_procs = 1 ! number of processes
integer(C_INT), bind(C, name='openmc_rank') :: rank = 0 ! rank of process
logical(C_BOOL), bind(C, name='openmc_master') :: master = .true. ! master process?
logical :: mpi_enabled = .false. ! is MPI in use and initialized?
#ifdef OPENMC_MPIF08
type(MPI_Datatype) :: MPI_BANK ! MPI datatype for fission bank

4
src/simulation_header.F90
View file

@ -32,8 +32,8 @@ module simulation_header
logical :: satisfy_triggers = .false. ! whether triggers are satisfied
integer(8) :: work ! number of particles per processor
integer(8), allocatable :: work_index(:) ! starting index in source bank for each process
integer(C_INT64_T), bind(C, name='openmc_work') :: work ! number of particles per processor
integer(C_INT64_T), allocatable :: work_index(:) ! starting index in source bank for each process
integer(8) :: current_work ! index in source bank of current history simulated
! ============================================================================

2
src/source.F90
View file

@ -88,7 +88,7 @@ contains
call write_message('Writing out initial source...', 5)
filename = trim(path_output) // 'initial_source.h5'
file_id = file_open(filename, 'w', parallel=.true.)
call write_source_bank(file_id)
call write_source_bank(file_id, work_index, source_bank)
call file_close(file_id)
end if

135
src/state_point.F90
View file

@ -15,6 +15,7 @@ module state_point
use hdf5
use bank_header, only: Bank
use cmfd_header
use constants
use eigenvalue, only: openmc_get_keff
@ -37,6 +38,15 @@ module state_point
implicit none
interface
subroutine write_source_bank(group_id, work_index, bank_) bind(C)
import HID_T, C_INT64_T, Bank
integer(HID_T), value :: group_id
integer(C_INT64_T), intent(in) :: work_index(*)
type(Bank), intent(in) :: bank_(*)
end subroutine write_source_bank
end interface
contains
!===============================================================================
@ -489,7 +499,7 @@ contains
end if
end if
call write_source_bank(file_id)
call write_source_bank(file_id, work_index, source_bank)
if (master .or. parallel) call file_close(file_id)
end if
@ -502,7 +512,7 @@ contains
call write_dataset(file_id, "filetype", 'source')
end if
call write_source_bank(file_id)
call write_source_bank(file_id, work_index, source_bank)
if (master .or. parallel) call file_close(file_id)
end if
@ -831,132 +841,11 @@ contains
end subroutine load_state_point
!===============================================================================
! WRITE_SOURCE_BANK writes OpenMC source_bank data
!===============================================================================
subroutine write_source_bank(group_id)
use bank_header, only: Bank
integer(HID_T), intent(in) :: group_id
integer :: hdf5_err
integer(HID_T) :: dset ! data set handle
integer(HID_T) :: dspace ! data or file space handle
integer(HID_T) :: memspace ! memory space handle
integer(HSIZE_T) :: offset(1) ! source data offset
integer(HSIZE_T) :: dims(1)
type(c_ptr) :: f_ptr
#ifdef PHDF5
integer(HID_T) :: plist ! property list
#else
integer :: i
#ifdef OPENMC_MPI
integer :: mpi_err ! MPI error code
type(Bank), allocatable, target :: temp_source(:)
#endif
#endif
#ifdef PHDF5
! Set size of total dataspace for all procs and rank
dims(1) = n_particles
call h5screate_simple_f(1, dims, dspace, hdf5_err)
call h5dcreate_f(group_id, "source_bank", hdf5_bank_t, dspace, dset, hdf5_err)
! Create another data space but for each proc individually
dims(1) = work
call h5screate_simple_f(1, dims, memspace, hdf5_err)
! Select hyperslab for this dataspace
offset(1) = work_index(rank)
call h5sselect_hyperslab_f(dspace, H5S_SELECT_SET_F, offset, dims, hdf5_err)
! Set up the property list for parallel writing
call h5pcreate_f(H5P_DATASET_XFER_F, plist, hdf5_err)
call h5pset_dxpl_mpio_f(plist, H5FD_MPIO_COLLECTIVE_F, hdf5_err)
! Set up pointer to data
f_ptr = c_loc(source_bank)
! Write data to file in parallel
call h5dwrite_f(dset, hdf5_bank_t, f_ptr, hdf5_err, &
file_space_id=dspace, mem_space_id=memspace, &
xfer_prp=plist)
! Close all ids
call h5sclose_f(dspace, hdf5_err)
call h5sclose_f(memspace, hdf5_err)
call h5dclose_f(dset, hdf5_err)
call h5pclose_f(plist, hdf5_err)
#else
if (master) then
! Create dataset big enough to hold all source sites
dims(1) = n_particles
call h5screate_simple_f(1, dims, dspace, hdf5_err)
call h5dcreate_f(group_id, "source_bank", hdf5_bank_t, &
dspace, dset, hdf5_err)
! Save source bank sites since the souce_bank array is overwritten below
#ifdef OPENMC_MPI
allocate(temp_source(work))
temp_source(:) = source_bank(:)
#endif
do i = 0, n_procs - 1
! Create memory space
dims(1) = work_index(i+1) - work_index(i)
call h5screate_simple_f(1, dims, memspace, hdf5_err)
#ifdef OPENMC_MPI
! Receive source sites from other processes
if (i > 0) then
call MPI_RECV(source_bank, int(dims(1)), MPI_BANK, i, i, &
mpi_intracomm, MPI_STATUS_IGNORE, mpi_err)
end if
#endif
! Select hyperslab for this dataspace
call h5dget_space_f(dset, dspace, hdf5_err)
offset(1) = work_index(i)
call h5sselect_hyperslab_f(dspace, H5S_SELECT_SET_F, offset, dims, hdf5_err)
! Set up pointer to data and write data to hyperslab
f_ptr = c_loc(source_bank)
call h5dwrite_f(dset, hdf5_bank_t, f_ptr, hdf5_err, &
file_space_id=dspace, mem_space_id=memspace)
call h5sclose_f(memspace, hdf5_err)
call h5sclose_f(dspace, hdf5_err)
end do
! Close all ids
call h5dclose_f(dset, hdf5_err)
! Restore state of source bank
#ifdef OPENMC_MPI
source_bank(:) = temp_source(:)
deallocate(temp_source)
#endif
else
#ifdef OPENMC_MPI
call MPI_SEND(source_bank, int(work), MPI_BANK, 0, rank, &
mpi_intracomm, mpi_err)
#endif
end if
#endif
end subroutine write_source_bank
!===============================================================================
! READ_SOURCE_BANK reads OpenMC source_bank data
!===============================================================================
subroutine read_source_bank(group_id)
use bank_header, only: Bank
integer(HID_T), intent(in) :: group_id
integer :: hdf5_err

112
src/state_point.cpp Normal file
View file

@ -0,0 +1,112 @@
#include "state_point.h"
#include <algorithm>
#include <vector>
#include "mpi.h"
#include "message_passing.h"
#include "openmc.h"
namespace openmc {
void
write_source_bank(hid_t group_id, int64_t* work_index, const Bank* source_bank)
{
// Create type for array of 3 reals
hsize_t dims[] {3};
hid_t triplet = H5Tarray_create(H5T_NATIVE_DOUBLE, 1, dims);
// Create bank datatype
hid_t banktype = H5Tcreate(H5T_COMPOUND, sizeof(struct Bank));
H5Tinsert(banktype, "wgt", HOFFSET(Bank, wgt), H5T_NATIVE_DOUBLE);
H5Tinsert(banktype, "xyz", HOFFSET(Bank, xyz), triplet);
H5Tinsert(banktype, "uvw", HOFFSET(Bank, uvw), triplet);
H5Tinsert(banktype, "E", HOFFSET(Bank, E), H5T_NATIVE_DOUBLE);
H5Tinsert(banktype, "delayed_group", HOFFSET(Bank, delayed_group), H5T_NATIVE_INT);
#ifdef PHDF5
// Set size of total dataspace for all procs and rank
dims[0] = n_particles;
hid_t dspace = H5Screate_simple(1, dims, H5P_DEFAULT);
hid_t dset = H5Dcreate(group_id, "source_bank", banktype, dspace,
H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
// Create another data space but for each proc individually
hsize_t count[] {openmc_work};
hid_t memspace = H5Screate_simple(1, count, H5P_DEFAULT);
// Select hyperslab for this dataspace
hsize_t start[] {work_index[openmc::mpi::rank]};
H5Sselect_hyperslab(dspace, H5S_SELECT_SET, start, nullptr, count, nullptr);
// Set up the property list for parallel writing
hid_t plist = H5Pcreate(H5P_DATASET_XFER);
H5Pset_dxpl_mpio_f(plist, H5FD_MPIO_COLLECTIVE);
// Write data to file in parallel
H5Dwrite(dset, banktype, memspace, dspace, memspace, plist, source_bank);
// Free resources
H5Sclose(dspace);
h5sclose(memspace);
H5Dclose(dset);
H5Pclose(plist);
#else
if (openmc_master) {
// Create dataset big enough to hold all source sites
hsize_t dims[] {static_cast<hsize_t>(n_particles)};
hid_t dspace = H5Screate_simple(1, dims, H5P_DEFAULT);
hid_t dset = H5Dcreate(group_id, "source_bank", banktype, dspace,
H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
// Save source bank sites since the souce_bank array is overwritten below
#ifdef OPENMC_MPI
std::vector<double> temp_source {source_bank, source_bank + openmc_work};
#endif
for (int i = 0; i < openmc::mpi::n_procs; ++i) {
// Create memory space
hsize_t count[] {static_cast<hsize_t>(work_index[i+1] - work_index[i])};
hid_t memspace = H5Screate_simple(1, count, H5P_DEFAULT);
#ifdef OPENMC_MPI
// Receive source sites from other processes
if (i > 0)
MPI_Recv(source_bank, count[0], openmc::mpi::bank, i, i,
openmc::mpi::intracomm, MPI_STATUS_IGNORE);
#endif
// Select hyperslab for this dataspace
dspace = H5Dget_space(dset);
hsize_t start[] {static_cast<hsize_t>(work_index[i])};
H5Sselect_hyperslab(dspace, H5S_SELECT_SET, start, nullptr, count, nullptr);
// Write data to hyperslab
H5Dwrite(dset, banktype, memspace, dspace, H5P_DEFAULT, source_bank);
H5Sclose(memspace);
H5Sclose(dspace);
}
// Close all ids
H5Dclose(dset);
#ifdef OPENMC_MPI
// Restore state of source bank
std::copy(temp_source.begin(), temp_source.end(), source_bank);
#endif
} else {
#ifdef OPENMC_MPI
MPI_Send(source_bank, openmc_work, openmc::mpi::bank, 0, openmc::mpi::rank,
openmc::mpi::mpi_intracomm);
#endif
}
#endif
H5Tclose(banktype);
H5Tclose(triplet);
}
} // namespace openmc

15
src/state_point.h Normal file
View file

@ -0,0 +1,15 @@
#ifndef STATE_POINT_H
#define STATE_POINT_H
#include <cstdint>
#include "hdf5.h"
#include "openmc.h"
namespace openmc {
extern "C" void write_source_bank(hid_t group_id, int64_t* work_index,
const struct Bank* bank);
} // namespace openmc
#endif // STATE_POINT_H