#include "macros.h" ! fpp types module shared_array use iso_c_binding, only: c_sizeof use constants, only: int32, int64, dp, MPIArg use shared_memory_mpi, only: shared_allocate_mpi, shared_deallocate_mpi use mpi_f08, only: MPI_Barrier, MPI_Win_Sync, MPI_Win use MPI_wrapper, only: mpi_comm_intra use MemoryManager, only: LogMemALloc, LogMemDealloc, TagIntType use util_mod, only: ByteSize_t better_implicit_none private ! Shared memory array types are defined here public :: shared_array_real_t type :: shared_array_real_t ! They contain a ptr (access to the array) ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF real(dp), pointer :: ptr(:) => null() ! and an MPI window type(MPI_Win) :: win ! Tag for the NECI memory manager integer(TagIntType) :: tag = 0 contains ! allocation and deallocation routines procedure :: shared_alloc => safe_shared_memory_alloc_real procedure :: shared_dealloc => safe_shared_memory_dealloc_real procedure :: sync => sync_real procedure :: get_memory_demand => get_memory_demand_real end type shared_array_real_t public :: shared_array_int64_t type :: shared_array_int64_t ! They contain a ptr (access to the array) ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF integer(int64), pointer :: ptr(:) => null() ! and an MPI window type(MPI_Win) :: win ! Tag for the NECI memory manager integer(TagIntType) :: tag = 0 contains ! allocation and deallocation routines procedure :: shared_alloc => safe_shared_memory_alloc_int64 procedure :: shared_dealloc => safe_shared_memory_dealloc_int64 procedure :: sync => sync_int64 procedure :: get_memory_demand => get_memory_demand_int64 end type shared_array_int64_t public :: shared_array_int32_t type :: shared_array_int32_t ! They contain a ptr (access to the array) ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF integer(int32), pointer :: ptr(:) => null() ! and an MPI window type(MPI_Win) :: win ! Tag for the NECI memory manager integer(TagIntType) :: tag = 0 contains ! allocation and deallocation routines procedure :: shared_alloc => safe_shared_memory_alloc_int32 procedure :: shared_dealloc => safe_shared_memory_dealloc_int32 procedure :: sync => sync_int32 procedure :: get_memory_demand => get_memory_demand_int32 end type shared_array_int32_t public :: shared_array_cmplx_t type :: shared_array_cmplx_t ! They contain a ptr (access to the array) ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF complex(dp), pointer :: ptr(:) => null() ! and an MPI window type(MPI_Win) :: win ! Tag for the NECI memory manager integer(TagIntType) :: tag = 0 contains ! allocation and deallocation routines procedure :: shared_alloc => safe_shared_memory_alloc_cmplx procedure :: shared_dealloc => safe_shared_memory_dealloc_cmplx procedure :: sync => sync_cmplx procedure :: get_memory_demand => get_memory_demand_cmplx end type shared_array_cmplx_t public :: shared_array_bool_t type :: shared_array_bool_t ! They contain a ptr (access to the array) ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF logical, pointer :: ptr(:) => null() ! and an MPI window type(MPI_Win) :: win ! Tag for the NECI memory manager integer(TagIntType) :: tag = 0 contains ! allocation and deallocation routines procedure :: shared_alloc => safe_shared_memory_alloc_bool procedure :: shared_dealloc => safe_shared_memory_dealloc_bool procedure :: sync => sync_bool procedure :: get_memory_demand => get_memory_demand_bool end type shared_array_bool_t contains !------------------------------------------------------------------------------------------! ! Auxiliary functions to prevent code duplication !------------------------------------------------------------------------------------------! !> Wrapper for shared_allocate_mpi that tests if the pointer is associated !> @param[out] win MPI shared memory window for internal MPI usage !> @param[out] ptr pointer to be allocated, on return points to a shared memory segment of given size !> @param[in] size size of the memory segment to be allocated subroutine safe_shared_memory_alloc_real (this, size, name) class(shared_array_real_t), intent(inout) :: this integer(int64), intent(in) :: size character(*), intent(in), optional :: name character(*), parameter :: t_r = "shared_alloc" ! if pointer was allocated prior, re-allocate the probabilities ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF call safe_shared_memory_dealloc_real (this) call shared_allocate_mpi(this%win, this%ptr, (/size/)) ! If a name is given, log the allocation if (associated(this%ptr) .and. present(name)) & call LogMemAlloc(name, size, sizeof(this%ptr(1)), t_r, this%tag) end subroutine safe_shared_memory_alloc_real !------------------------------------------------------------------------------------------! !> wrapper for shared_deallocate_mpi that tests if the pointer is associated !> @param[in,out] win MPI shared memory window for internal MPI usage !> @param[in,out] ptr pointer to be deallocated (if associated) ! WARNING: THIS ASSUMES THAT IF ptr IS ASSOCIATED, IT POINTS TO AN MPI SHARED MEMORY ! WINDOW win subroutine safe_shared_memory_dealloc_real (this) class(shared_array_real_t), intent(inout) :: this character(*), parameter :: t_r = "shared_dealloc" ! assume that if ptr is associated, it points to mpi shared memory if (associated(this%ptr)) call shared_deallocate_mpi(this%win, this%ptr) ! First, check if we have to log the deallocation if (this%tag /= 0) call LogMemDealloc(t_r, this%tag) end subroutine safe_shared_memory_dealloc_real !> callls MPI_Win_Sync on the array's shared memory window to sync rma !! This has to be called between read/write epochs to ensure all tasks of a node are !! looking at the same shared data subroutine sync_real (this) class(shared_array_real_t), intent(inout) :: this integer(MPIArg) :: ierr call MPI_Win_Sync(this%win, ierr) call MPI_Barrier(mpi_comm_intra, ierr) end subroutine sync_real elemental function get_memory_demand_real(this) result(res) class(shared_array_real_t), intent(in) :: this type(ByteSize_t) :: res if (associated(this%ptr)) then res = ByteSize_t(size(this%ptr) * int(c_sizeof(this%ptr(1)), int64)) else res = ByteSize_t(0_int64) end if end function !> Wrapper for shared_allocate_mpi that tests if the pointer is associated !> @param[out] win MPI shared memory window for internal MPI usage !> @param[out] ptr pointer to be allocated, on return points to a shared memory segment of given size !> @param[in] size size of the memory segment to be allocated subroutine safe_shared_memory_alloc_int64 (this, size, name) class(shared_array_int64_t), intent(inout) :: this integer(int64), intent(in) :: size character(*), intent(in), optional :: name character(*), parameter :: t_r = "shared_alloc" ! if pointer was allocated prior, re-allocate the probabilities ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF call safe_shared_memory_dealloc_int64 (this) call shared_allocate_mpi(this%win, this%ptr, (/size/)) ! If a name is given, log the allocation if (associated(this%ptr) .and. present(name)) & call LogMemAlloc(name, size, sizeof(this%ptr(1)), t_r, this%tag) end subroutine safe_shared_memory_alloc_int64 !------------------------------------------------------------------------------------------! !> wrapper for shared_deallocate_mpi that tests if the pointer is associated !> @param[in,out] win MPI shared memory window for internal MPI usage !> @param[in,out] ptr pointer to be deallocated (if associated) ! WARNING: THIS ASSUMES THAT IF ptr IS ASSOCIATED, IT POINTS TO AN MPI SHARED MEMORY ! WINDOW win subroutine safe_shared_memory_dealloc_int64 (this) class(shared_array_int64_t), intent(inout) :: this character(*), parameter :: t_r = "shared_dealloc" ! assume that if ptr is associated, it points to mpi shared memory if (associated(this%ptr)) call shared_deallocate_mpi(this%win, this%ptr) ! First, check if we have to log the deallocation if (this%tag /= 0) call LogMemDealloc(t_r, this%tag) end subroutine safe_shared_memory_dealloc_int64 !> callls MPI_Win_Sync on the array's shared memory window to sync rma !! This has to be called between read/write epochs to ensure all tasks of a node are !! looking at the same shared data subroutine sync_int64 (this) class(shared_array_int64_t), intent(inout) :: this integer(MPIArg) :: ierr call MPI_Win_Sync(this%win, ierr) call MPI_Barrier(mpi_comm_intra, ierr) end subroutine sync_int64 elemental function get_memory_demand_int64(this) result(res) class(shared_array_int64_t), intent(in) :: this type(ByteSize_t) :: res if (associated(this%ptr)) then res = ByteSize_t(size(this%ptr) * int(c_sizeof(this%ptr(1)), int64)) else res = ByteSize_t(0_int64) end if end function !> Wrapper for shared_allocate_mpi that tests if the pointer is associated !> @param[out] win MPI shared memory window for internal MPI usage !> @param[out] ptr pointer to be allocated, on return points to a shared memory segment of given size !> @param[in] size size of the memory segment to be allocated subroutine safe_shared_memory_alloc_int32 (this, size, name) class(shared_array_int32_t), intent(inout) :: this integer(int64), intent(in) :: size character(*), intent(in), optional :: name character(*), parameter :: t_r = "shared_alloc" ! if pointer was allocated prior, re-allocate the probabilities ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF call safe_shared_memory_dealloc_int32 (this) call shared_allocate_mpi(this%win, this%ptr, (/size/)) ! If a name is given, log the allocation if (associated(this%ptr) .and. present(name)) & call LogMemAlloc(name, size, sizeof(this%ptr(1)), t_r, this%tag) end subroutine safe_shared_memory_alloc_int32 !------------------------------------------------------------------------------------------! !> wrapper for shared_deallocate_mpi that tests if the pointer is associated !> @param[in,out] win MPI shared memory window for internal MPI usage !> @param[in,out] ptr pointer to be deallocated (if associated) ! WARNING: THIS ASSUMES THAT IF ptr IS ASSOCIATED, IT POINTS TO AN MPI SHARED MEMORY ! WINDOW win subroutine safe_shared_memory_dealloc_int32 (this) class(shared_array_int32_t), intent(inout) :: this character(*), parameter :: t_r = "shared_dealloc" ! assume that if ptr is associated, it points to mpi shared memory if (associated(this%ptr)) call shared_deallocate_mpi(this%win, this%ptr) ! First, check if we have to log the deallocation if (this%tag /= 0) call LogMemDealloc(t_r, this%tag) end subroutine safe_shared_memory_dealloc_int32 !> callls MPI_Win_Sync on the array's shared memory window to sync rma !! This has to be called between read/write epochs to ensure all tasks of a node are !! looking at the same shared data subroutine sync_int32 (this) class(shared_array_int32_t), intent(inout) :: this integer(MPIArg) :: ierr call MPI_Win_Sync(this%win, ierr) call MPI_Barrier(mpi_comm_intra, ierr) end subroutine sync_int32 elemental function get_memory_demand_int32(this) result(res) class(shared_array_int32_t), intent(in) :: this type(ByteSize_t) :: res if (associated(this%ptr)) then res = ByteSize_t(size(this%ptr) * int(c_sizeof(this%ptr(1)), int64)) else res = ByteSize_t(0_int64) end if end function !> Wrapper for shared_allocate_mpi that tests if the pointer is associated !> @param[out] win MPI shared memory window for internal MPI usage !> @param[out] ptr pointer to be allocated, on return points to a shared memory segment of given size !> @param[in] size size of the memory segment to be allocated subroutine safe_shared_memory_alloc_cmplx (this, size, name) class(shared_array_cmplx_t), intent(inout) :: this integer(int64), intent(in) :: size character(*), intent(in), optional :: name character(*), parameter :: t_r = "shared_alloc" ! if pointer was allocated prior, re-allocate the probabilities ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF call safe_shared_memory_dealloc_cmplx (this) call shared_allocate_mpi(this%win, this%ptr, (/size/)) ! If a name is given, log the allocation if (associated(this%ptr) .and. present(name)) & call LogMemAlloc(name, size, sizeof(this%ptr(1)), t_r, this%tag) end subroutine safe_shared_memory_alloc_cmplx !------------------------------------------------------------------------------------------! !> wrapper for shared_deallocate_mpi that tests if the pointer is associated !> @param[in,out] win MPI shared memory window for internal MPI usage !> @param[in,out] ptr pointer to be deallocated (if associated) ! WARNING: THIS ASSUMES THAT IF ptr IS ASSOCIATED, IT POINTS TO AN MPI SHARED MEMORY ! WINDOW win subroutine safe_shared_memory_dealloc_cmplx (this) class(shared_array_cmplx_t), intent(inout) :: this character(*), parameter :: t_r = "shared_dealloc" ! assume that if ptr is associated, it points to mpi shared memory if (associated(this%ptr)) call shared_deallocate_mpi(this%win, this%ptr) ! First, check if we have to log the deallocation if (this%tag /= 0) call LogMemDealloc(t_r, this%tag) end subroutine safe_shared_memory_dealloc_cmplx !> callls MPI_Win_Sync on the array's shared memory window to sync rma !! This has to be called between read/write epochs to ensure all tasks of a node are !! looking at the same shared data subroutine sync_cmplx (this) class(shared_array_cmplx_t), intent(inout) :: this integer(MPIArg) :: ierr call MPI_Win_Sync(this%win, ierr) call MPI_Barrier(mpi_comm_intra, ierr) end subroutine sync_cmplx elemental function get_memory_demand_cmplx(this) result(res) class(shared_array_cmplx_t), intent(in) :: this type(ByteSize_t) :: res if (associated(this%ptr)) then res = ByteSize_t(size(this%ptr) * int(c_sizeof(this%ptr(1)), int64)) else res = ByteSize_t(0_int64) end if end function !> Wrapper for shared_allocate_mpi that tests if the pointer is associated !> @param[out] win MPI shared memory window for internal MPI usage !> @param[out] ptr pointer to be allocated, on return points to a shared memory segment of given size !> @param[in] size size of the memory segment to be allocated subroutine safe_shared_memory_alloc_bool (this, size, name) class(shared_array_bool_t), intent(inout) :: this integer(int64), intent(in) :: size character(*), intent(in), optional :: name character(*), parameter :: t_r = "shared_alloc" ! if pointer was allocated prior, re-allocate the probabilities ! WARNING: DO NOT MANUALLY RE-ASSIGN ptr, THIS WILL MOST LIKELY BREAK STUFF call safe_shared_memory_dealloc_bool (this) call shared_allocate_mpi(this%win, this%ptr, (/size/)) ! If a name is given, log the allocation if (associated(this%ptr) .and. present(name)) & call LogMemAlloc(name, size, sizeof(this%ptr(1)), t_r, this%tag) end subroutine safe_shared_memory_alloc_bool !------------------------------------------------------------------------------------------! !> wrapper for shared_deallocate_mpi that tests if the pointer is associated !> @param[in,out] win MPI shared memory window for internal MPI usage !> @param[in,out] ptr pointer to be deallocated (if associated) ! WARNING: THIS ASSUMES THAT IF ptr IS ASSOCIATED, IT POINTS TO AN MPI SHARED MEMORY ! WINDOW win subroutine safe_shared_memory_dealloc_bool (this) class(shared_array_bool_t), intent(inout) :: this character(*), parameter :: t_r = "shared_dealloc" ! assume that if ptr is associated, it points to mpi shared memory if (associated(this%ptr)) call shared_deallocate_mpi(this%win, this%ptr) ! First, check if we have to log the deallocation if (this%tag /= 0) call LogMemDealloc(t_r, this%tag) end subroutine safe_shared_memory_dealloc_bool !> callls MPI_Win_Sync on the array's shared memory window to sync rma !! This has to be called between read/write epochs to ensure all tasks of a node are !! looking at the same shared data subroutine sync_bool (this) class(shared_array_bool_t), intent(inout) :: this integer(MPIArg) :: ierr call MPI_Win_Sync(this%win, ierr) call MPI_Barrier(mpi_comm_intra, ierr) end subroutine sync_bool elemental function get_memory_demand_bool(this) result(res) class(shared_array_bool_t), intent(in) :: this type(ByteSize_t) :: res if (associated(this%ptr)) then res = ByteSize_t(size(this%ptr) * int(c_sizeof(this%ptr(1)), int64)) else res = ByteSize_t(0_int64) end if end function end module shared_array